Самый быстрый способ расчета инверсии матрицы — Производительность
yingqiuz
1
Привет всем
Мне нужно вычислить обратную положительно определенную матрицу H вида H = (X’ * X + Диагональ(d)). X — плоская матрица, имеющая (преимущественно) больше столбцов, чем строк. d состоит из положительных больших значений.
Максимум, что я мог сделать, это превратить задачу в нахождение обратной величины (I + X * Diagonal(1 ./ d) * X), используя формулу Вудбери. Интересно, есть ли какой-то особый способ ускорить обращение матрицы, учитывая, что она имеет специфическую структуру? холецкого!(I + X * Диагональ(1 ./ d) * X) или inv(I + X * Диагональ(1 ./ d) * X), или что-нибудь более осмысленное?
Большое спасибо. Любые предложения будут высоко ценится!
Стивендж
2
Инцюз:
Мне нужно вычислить обратную положительно определенную матрицу
Вы уверены, что вам нужно обратное? Почему нельзя использовать факторизацию Холецкого как есть?
1 Нравится
yingqiuz
3
Привет, спасибо за ответ. Вы предполагаете, что я могу использовать C = холецкий! (I + X * Диагональ (1 ./ d) * X) так, чтобы его инверсия была C.U \ (C.L \ I)?
Мне нужна обратная сторона (предпочтительнее) или хотя бы диагонали обратной…
stevengj
4
Инцюз:
Вы предполагаете, что я могу использовать C = холецкий!(I + X * Диагональ(1 ./ d) * X) так, чтобы его инверсия была C.U \ (C.L \ I)?
Нет, я предлагал вам подумать, можно ли использовать факторы Холецкого без явного вычисления обратной величины. То есть, учитывая C = cholesky!(X'X + Diagonal(d)) , вы можете быстро решить линейную систему для любой заданной правой части, поэтому во многих случаях вам не нужна обратная матрица в явном виде. .
Если вам действительно нужна вся обратная матрица, я бы предложил
LinearAlgebra.inv!(холецкий!(X'X + Diagonal(d)))
(не вижу смысла в ваших Диагональ(1 ./d) переформулировка.)
Мне нужна инверсия (предпочтительнее) или хотя бы диагонали инверсии…
Что вы в конечном итоге пытаетесь вычислить и каковы типичные размеры вашей матрицы?
9 лайков
yingqiuz
5
Большое спасибо за ответ. Я предполагаю, что мне нужна инверсия — это промежуточный шаг в моих итерациях, и мне нужна диагональ инверсии и ее определитель, которые нужно передать в следующую итерацию.
Вычисление детерминанта не требует обращения, но я не смог найти, как извлечь диагональ обратной матрицы без кубической сложности…
Учитывая разложение Холецкого, есть ли более быстрый способ только вычислить диагональ обратной? К сожалению, интересующая матрица огромна, по крайней мере, 30000 x 30000… но она становится меньше по мере итераций.
Большое спасибо!
Стивендж
6
Инцюз:
Полагаю, мне нужна обратная функция — это промежуточный шаг в моих итерациях, и мне нужна диагональ обратной и ее определитель, которые нужно передать в следующую итерацию.
Что пытается сделать ваша итерация? Может быть какой-то другой алгоритм для вашей конечной цели…
1 Нравится
yingqiuz
7
Большое спасибо! Моя итерация пытается найти оптимальное значение параметра, которое минимизирует функцию, например, argmin_{w} f(X, w, a). Мне нужно отрицательное значение, обратное гессиану, чтобы аппроксимировать w (апостериорное) распределение, , по крайней мере, на последней итерации. В промежуточных итерациях мне (только) нужна диагональ отрицательного обратного гессиана, что необходимо для оптимизации другого параметра, a.
Гессиан функции w имеет особую структуру (X’X + диагональ(d)). К сожалению, размеры X и w огромны — X 30 000 на 2 000 000 и w 2 000 000 на 1, что делает невозможным вычисление диагонали отрицательного обратного Гессеана.
В последней итерации мне нужен весь отрицательный обратный гессиан, но он будет гораздо меньшего размера.
Большое спасибо, что нашли время, чтобы прочитать это — есть ли лучший способ сделать это? Я рассматривал BFGS, но хранить фактический (аппроксимированный) обратный гессиан на каждом шаге также очень дорого.
Стивендж
8
Инцюз:
В промежуточных итерациях мне (только) нужна диагональ отрицательного обратного гессиана, что необходимо для оптимизации другого параметра, a.
Почему? Вы так и не сказали, какую функцию вы минимизируете.
Инцюз:
Я думал о BFGS, но хранить фактический (аппроксимированный) обратный гессиан на каждом шаге также очень дорого.
Рассматривали ли вы L-BFGS?
Инцюз
9
Спасибо! Моя итерация состоит из двух шагов. Предположим, что X имеет размерность n на d, w и a имеют размерность d на 1.решить argmin_{w}f(X, t, w, a) = y’y + 0,5 * w’ * диагональ(a) * w, сохраняя и фиксированными.
Здесь y = t .* log.(сигмоид.(X * w)) .+ (1 .- t) .* log(1 .- сигмоид.(X * w)). t n на 1.после нахождения оптимального значения w, сохранения его фиксированным и решения argmin_{a}f2(X, t, w, a) = y’y + 0,5 * w’ * Diagonal(a) * w + 0,5sum(log .
(а)) + 0,5logdet(Сигма). Здесь Sigma — отрицательный обратный гессиан для w при его оптимальном значении, а гессиан для w имеет вид (X’ * Diagonal(y .* (1 .-y)) * X + Diagonal(a))
(а)) + 0,5logdet(Сигма). Здесь Sigma — отрицательный обратный гессиан для w при его оптимальном значении, а гессиан для w имеет вид (X’ * Diagonal(y .* (1 .-y)) * X + Diagonal(a))Это проблема вектора релевантности.
На самом деле я уже давно не понимаю этой проблемы. f(X, t, w, a) — апостериорное значение w, а f2(X, t, w, a) — функция с проинтегрированным w (т. е. вероятность второго типа). Поскольку точное интегрирование трудноразрешимо, я использовал приближение Лапласа, как это обычно делают люди, вычислив f2(X, t, w, a) как приближение вероятности типа 2. Мои друзья предлагают выполнить SGD для f(X, t, w, a), чтобы найти w и a вместе, чтобы избежать f2, но я не исследовал, даст ли это хорошие результаты. …
Насколько близок приблизительный обратный гессиан к реальному обратному гессиану, полученному с помощью такой программы, как BFGS/L-BFGS?
Еще раз большое спасибо за ваше время.
Стивендж
10
Инцюз:
- 0,5логдет(Сигма). Здесь Sigma — отрицательный обратный гессиан для w при его оптимальном значении, а гессиан для w имеет вид (X’ * Diagonal(y .* (1 .-y)) * X + Diagonal(a))
Хорошо, так вот почему вы думаете, что вам нужна полная обратная матрица — вы вычисляете логарифмический определитель? Существует множество методов оценки детерминанта журнала, например. с помощью стохастических квадратур Ланцоша + оценок следов Хатчинсона. Более того, поскольку оценки трасс Хатчинсона представляют собой именно форму вычисления ожидаемого значения, они совместимы со стохастическим градиентным спуском (поэтому вам не нужно вычислять точный логарифмический детерминант на каждом шаге, чтобы минимизировать его). Это итерационные методы, поэтому вам нужны только произведения матрицы на вектор, что полезно, если ваша матрица большая и разреженная (или структурированная).