Как вычислить scipy разреженный матричный определитель, не превращая его в плотный?

Я пытаюсь найти самый быстрый способ найти определитель разреженных симметричных и вещественных матриц в python. используя модуль scipy sparse, но на самом деле удивлен, что нет детерминантной функции. Я знаю, что могу использовать LU-факторизацию для вычисления детерминанта, но не вижу простого способа сделать это, потому что возвращение scipy.sparse.linalg.splu является объектом и создание экземпляра плотной матрицы L и U не стоит того, что я могу сделать sp.linalg.det(A.todense()) где A - моя скудная разреженная матрица.

Я также немного удивлен, почему другие не столкнулись с проблемой эффективного детерминантного вычисления в scipy. Как использовать splu для вычисления детерминанта?

Я просмотрел pySparse и scikits.sparse.chlmod. Последнее для меня нецелесообразно - требует установки пакетов, а также не уверен, насколько быстро код, прежде чем я войду во все проблемы. Любые решения? Заранее спасибо.

Ответы

Ответ 1

Вот несколько ссылок, которые я представил как часть ответа здесь. Я думаю, что они решают проблему, которую вы пытаетесь решить:

  • примечания для реализации в библиотеке Shogun
  • Эрленд Ауне, Даниэль П. Симпсон: оценка параметров в высокоразмерных гауссовских распределениях, в частности раздел 2.1 (arxiv: 1105.5256)
  • Ilse C.F. Ipsen, Dean J. Lee: Определяющие аппроксимации (arxiv: 1105.0437)
  • Арнольд Реускен: Приближение детерминанта больших разреженных симметричных положительно определенных матриц (arxiv: hep-lat/0008007)

Цитата из примечаний по сёгуну:

Обычная методика вычисления логарифминантного члена в выражении правдоподобия опирается на факторизацию Cholesky матрицы, т.е. Σ = LLT, (L - нижний треугольный фактор Холецкого), а затем используя диагональные элементы фактора для вычисления журнал (Det (Σ)) = 2Σni = 1log (ЛИИ). Однако для разреженных матриц, как правило, являются ковариационные матрицы, факторы Холецкого часто страдают от феноменов заполнения - они оказываются не такими скудными. Поэтому для больших размеров эта техника становится невозможной из-за огромной потребности в памяти для хранения всех этих нерелевантных недиагональных коэффи- циентов фактора. В то время как были разработаны методы упорядочения, чтобы предварительно переставить строки и столбцы, чтобы уменьшить заполнение, например. ориентировочная минимальная степень (AMD), эти методы во многом зависят от шаблона разреженности и, следовательно, не гарантируют лучшего результата.

Недавние исследования показывают, что, используя ряд методов комплексного анализа, числовую линейную алгебру и жадную градиентную окраску, мы можем, однако, приблизить логарифмический определитель до произвольной точности [Aune et. al., 2012]. Основной трюк лежит в наблюдении, что мы можем написать log (det (Σ)) как след (log (Σ)), где log (Σ) - матричный логарифм.

Ответ 2

"Стандартный" способ решить эту проблему - с холесковой декомпозицией, но если вы не используете какой-либо новый скомпилированный код, вам не повезло. Лучшей редкой холеской реализацией является Тим Дэвис CHOLMOD, который лицензируется под LGPL и, следовательно, не доступен в scipy (Scipy - BSD).