Эффективно вычислить диагональ обратной разреженной матрицы ⇐ Python

Программы на Python
Anonymous
Эффективно вычислить диагональ обратной разреженной матрицы

Сообщение Anonymous »

В Python имеется большая разреженная матрица NxN, представляющая матрицу точности Q (обратную ковариационную матрицу) многомерного распределения Гаусса. Меня интересует получение предельных отклонений модели. Это будут элементы на диагонали обратной матрицы точности (диагональные элементы ковариационной матрицы).
Поскольку Q очень велико, мне следует избегать инвертирования всей матрицы. Решение линейной системы для получения результата также происходит недостаточно быстро. Однако известно, что диагональные элементы можно эффективно вычислять с использованием разреженной факторизации Холецкого. Это описано в уравнении 12.23 в книге Гельфанда «Справочник по пространственной статистике» (2010). Книга Рю и Хелда (2005) «Гауссовы марковские случайные поля» также содержит некоторую информацию по этому поводу.
В R функция inla.qinv из R-INLA< Библиотека /code> делает именно это. Как я могу получить функциональность, аналогичную функции inla.qinv в Python? В качестве альтернативы, есть ли какие-либо другие способы эффективного вычисления или оценки диагональных элементов?
На данный момент я использую пакет sksparse.cholmod для разреженного разложения Холецкого, но открыт для другие методы. Просто решив обратную задачу,

Код: Выделить всё

import scipy.sparse as sparse
from sksparse.cholmod import cholesky
inverse_diagonal = cholesky(Q).solve_A(sparse.eye_array(Q.shape[0])).diagonal()
недостаточно быстро для моих нужд.
Одной альтернативой является использование пакета rpy2 в Python для запуска R-функций. Однако этого осложнения я хочу избежать. Кроме того, в ходе тестирования я обнаружил, что простой запуск функции qinv с использованием rpy2 не работает должным образом.
Если существуют какие-либо C-библиотеки, которые выполняют желаемое вычисление, безусловно, возможно! qinv использует C-библиотеку GMRFLib под капотом, так что, возможно, стоит попробовать?

Подробнее здесь: https://stackoverflow.com/questions/790 ... rse-matrix

Вернуться в «Python»