Выборка из многомерного гауссова с графом лапласовой (обратной) ковариации

12

Мы знаем, например, из Koutis-Miller-Peng (на основе работы Spielman & Teng), что мы можем очень быстро решить линейные системы Ax=b для матриц A которые представляют собой матрицу Лапласа графа для некоторого разреженного графа с неотрицательными весами ребер ,

Теперь (первый вопрос) рассмотрим использование одной из этих графов лапласовских матриц A в качестве ковариационной или (второй вопрос) обратной ковариационной матрицы многомерного нормального распределения с нулевым средним или . Для каждого из этих случаев у меня есть два вопроса:N ( 0 , A - 1 )N(0,A)N(0,A1)

A. Насколько эффективно мы можем извлечь образец из этого распределения? (Как правило, чтобы нарисовать образец, мы вычисляем разложение Холецкого , нарисуем стандартную нормаль y \ sim \ mathcal {N} (\ boldsymbol {0}, I) , затем вычисляем образец как x = L ^ { -1} у ). у ~ N ( 0 , я ) х = L - 1 гA=LLTyN(0,I)x=L1y

B. Насколько эффективно мы можем вычислить определитель A ?

Обратите внимание, что оба из них могут быть легко решены с помощью разложения Холецкого, но я не сразу вижу, как извлечь L более эффективно, чем просто с помощью стандартного разреженного алгоритма Холецкого, который не будет использовать методы, представленные в вышеупомянутой ссылке работает, и который будет иметь кубическую сложность для графов разреженной, но высокой ширины.

dan_x
источник
Я думаю, что это может помочь немного конкретизировать то, что вы считаете «эффективным» в обоих случаях. Является ли «эффективность» такой же, как «не зависит от разложения Холецкого»?
Суреш Венкат
Спасибо за предложение. Возможно, ответ на все вопросы таков: «вам нужно вычислить разложение Холецкого, и нет структуры, которую можно было бы использовать за пределами редкости матрицы». Мне было бы интересно узнать, правда ли это (но я надеюсь, что это не так). Что касается «эффективно» в последнем абзаце, да, я в основном имею в виду более эффективно, чем стандартные разреженные алгоритмы Холецкого. Хотя, если бы был способ использовать методы вышеупомянутой работы, чтобы вычислить число Холецки так же быстро, как это можно сделать с помощью других средств, это также было бы интересно.
dan_x
Если вы хотите произвести выборку из , вы можете использовать это , где - матрица заболеваемости на графике. Таким образом, вы можете попробовать от стандартной гауссовой на ( являются ребра) и применить линейное преобразование . Я не знаю, как это соотносится с предложениями ниже, но вам не нужно вычислять разложение Холецкого. A = B T B B R E E BN(0,A)A=BTBBREEB
Лоренцо Найт

Ответы:

3

Здесь есть два отдельных вопроса.

  1. Как использовать эффективные решатель для для применения A 1 / 2 б .Ax=bA1/2b
  2. Как вычислить определитель.

Краткие ответы: 1) использовать приближения рациональных матричных функций и 2) нет, но в любом случае вам это не нужно. Я рассматриваю оба эти вопроса ниже.

Матрица квадратного корня

Идея состоит в том, чтобы преобразовать приближение рациональной функции для скалярных функций в приближение рациональной функции для матричных функций.

Мы знаем, что существуют рациональные функции, которые могут очень хорошо аппроксимировать функцию квадратного корня, для положительногоbi. Действительно, чтобы получить высокую точность на интервале[m,M], вам нужноO(logM

xr(x):=a1x+b1+a2x+b2++aNx+bN,
bi[m,M]условия в серии. Чтобы получить подходящие веса (ai) и полюса (-bi), просто найдите рациональное приближение функции онлайн или в книге.O(logMm)aibi

Теперь рассмотрим применение этой рациональной функции к вашей матрице:

r(A)=a1(A+b1I)1+a2(A+b2I)1++aN(A+bNI)1.

AA=UΣUA

||A1/2r(A)||2=||U(Σ1/2r(Σ))U||2,=maxi|σir(σi)|
A=UΣUA

Обозначив число условия через , мы можем применить к любому желаемому допуску, выполнив положительно смещенный граф лапласовых решений вида Κ 1 / 2 б O ( журнал κ ) ( + б I ) х = Ь .AκA1/2bO(logκ)

(A+bI)x=b.

Эти решения могут быть сделаны с вашим любимым графическим лапласовским решателем - я предпочитаю методы многосеточного типа, но тот, который вы цитируете, тоже подойдет. Дополнительный только помогает сходимости решателя.bI

Отличная статья, в которой обсуждается это, а также более общие методы комплексного анализа, применимые к несимметричным матрицам, см. разделе « Вычисление , и связанных матричных функций с помощью контурных интегралов» лог ( )Aαlog(A) , Hale, Higham и Trefethen (2008). ).

Определитель "вычисления"

Детерминант сложнее вычислить. Насколько я знаю, лучший способ вычислить разложение Шура с помощью QR - алгоритма, а затем считывать собственные от диагонали верхней треугольной матрицы . Это занимает времени, где - количество узлов в графе.A=QUQUO(n3)n

Однако вычисление детерминантов является изначально плохо обусловленной проблемой, поэтому, если вы когда-либо читали статью, которая основывается на вычислении детерминантов большой матрицы, вы должны очень скептически относиться к этому методу.

К счастью, вам, вероятно, не нужен определитель. Например,

  • Чтобы нарисовать выборки из одного гауссовского распределения , константа нормализации одинакова во всех точках, поэтому вам никогда не придется ее вычислять.N(0,A1)
  • Если ваша матрица Лапласа представляет обратную ковариацию локального гауссова приближения в точке к негауссову распределению, то определитель действительно меняется от точки к точке. Тем не менее, в каждой эффективной схеме выборки, которую я знаю (включая цепочку Маркова, Монте-Карло, выборку по важности и т. Д.), Вам действительно нужно соотношение детерминантов , где - текущая точка, а - предлагаемый следующий образец.A=Axx
    det(Ax01Axp),
    x0xp

Мы можем рассматривать как низкоуровневое обновление идентификатора, где эффективный числовой ранг низшего ранга является локальной мерой того, насколько негауссово истинное распределение; обычно это намного ниже, чем полный ранг матрицы. В самом деле, если велико, то истинное распределение локально настолько негауссово, что следует подвергнуть сомнению всю стратегию попытки выборки этого распределения с использованием локальных гауссовых приближений.Ax01Axp

Ax01Axp=I+QDQ,
rr

Факторы низкого ранга и можно найти с помощью рандомизированных SVD или Lanczos, применяя матрицу к различным векторам, для каждого применения которых требуется один граф Лапласово решение. Таким образом, общая работа по получению этих факторов низкого ранга составляет .QD

Ax01AxpI
O(r)O(rmax(n,E))

Зная , отношение определителей будет тогда det ( A - 1 x 0 A x p ) = det ( I + Q D Q ) = exp ( r i = 1 log d i ) .D=diag(d1,d2,,dr)

det(Ax01Axp)=det(I+QDQ)=exp(i=1rlogdi).

Эти методы расчета коэффициента детерминанта низкого ранга можно найти в методе стохастического Ньютона MCMC для крупномасштабных статистических обратных задач с применением к сейсмической инверсии , Martin, et al. (2012). В этой статье он применяется к задачам континуума, поэтому «граф» - это сетка в трехмерном пространстве, а граф Лапласа - это действительная матрица Лапласа. Однако все техники применимы к лапласианам общего графа. Вероятно, есть и другие статьи, в которых этот метод применяется к общим графам (расширение тривиально и, в основном, то, что я только что написал)

Ник Алджер
источник