Какой самый быстрый способ вычислить наибольшее собственное значение общей матрицы?

27

РЕДАКТИРОВАТЬ: я проверяю, если какие-либо собственные значения имеют величину один или больше.

Мне нужно найти наибольшее абсолютное собственное значение большой разреженной несимметричной матрицы.

Я использовал eigen()функцию R , которая использует алгоритм QR из EISPACK или LAPACK, чтобы найти все собственные значения, а затем я использую, abs()чтобы получить абсолютные значения. Однако мне нужно сделать это быстрее.

Я также попытался использовать интерфейс ARPACK в igraphпакете R. Однако, это дало ошибку для одной из моих матриц.

Окончательная реализация должна быть доступна из R.

Вероятно, будет несколько собственных значений одинаковой величины.

У вас есть какие-нибудь предложения?

РЕДАКТИРОВАТЬ: Точность только должна быть 1e-11. «Типичная» матрица до сих пор была . Я был в состоянии сделать QR-факторизацию на нем. Тем не менее, также возможно иметь гораздо большие. В настоящее время я начинаю читать об алгоритме Арнольди. Я понимаю, что это связано с Lanczsos.386×386

РЕДАКТИРОВАТЬ 2: Если у меня есть несколько матриц, которые я "тестирую", и я знаю, что есть большая подматрица, которая не меняется. Можно ли игнорировать / отменить это?

сила
источник
Смотрите мой ответ здесь: scicomp.stackexchange.com/a/1679/979 . Это актуальная тема исследования, и современные методы могут быть лучше, чем Ланцош. Задача вычисления сингулярных значений эквивалентна проблеме вычисления собственных значений.
dranxo
2
Матрица 400х400! = Большая. Кроме того, что означает наибольшее значение, если «вероятно, будет несколько собственных значений одинаковой величины»? В клочковатой земле: linalg.eig (random.normal (size = (400,400))) занимает около половины секунды. Это слишком медленно?
Meawoppl
@meawoppl да полсекунды слишком медленно. Это потому, что он является частью другого алгоритма, который выполняет этот расчет много раз.
мощность
1
@Power GotCah. У вас есть приближение к собственному вектору. то есть похоже ли это на последнее решение, или вы можете сделать обоснованное предположение о его структуре?
Meawoppl

Ответы:

14

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

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

Если это не относится к вашей проблеме, дайте более конкретную информацию в вашем вопросе. Затем добавьте комментарий к этому ответу, и я обновлю его.

Редактировать: [Это было для старой версии вопроса, asling для наибольшего собственного значения.] Поскольку ваша матрица мала и, по-видимому, плотная, я бы сделал итерацию Арнольди на B = (IA) ^ {- 1}, используя начальный Перестановочная треугольная факторизация IA для дешевого умножения на B. (Или вычисление явного обратного, но это стоит в 3 раза больше факторизации.) Вы хотите проверить, имеет ли B отрицательное собственное значение. Работая с B вместо A, отрицательные собственные значения гораздо лучше разделяются, поэтому, если они есть, вы должны быстро сходиться.

Но мне любопытно, откуда твоя проблема. Несимметричные матрицы обычно имеют сложные собственные значения, поэтому «наибольшее» даже не определено. Таким образом, вы должны знать больше о своей проблеме, что может помочь подсказать, как решить ее еще быстрее и / или более надежно.

Edit2: Это трудно получить с Арнольди в конкретное подмножество интереса. Чтобы надежно получить абсолютно самые большие собственные значения, вы должны выполнить итерацию в подпространстве, используя исходную матрицу, с размером подпространства, соответствующим или превышающим число собственных значений, которое, как ожидается, будет близко к 1 или больше по величине. На маленьких матрицах это будет медленнее, чем алгоритм QR, но на больших матрицах это будет намного быстрее.

Арнольд Ноймайер
источник
Мне нужно проверить, является ли наибольшее собственное значение больше 1. Точность должна быть только до 1e-11. «Типичная» матрица до сих пор была 386 x 386. Мне удалось выполнить QR-факторизацию. Тем не менее, также возможно иметь гораздо большие. В настоящее время я начинаю читать об алгоритме Арнольди. Я понимаю, что это связано с Lanczsos.
сила
Эта информация относится к вашему вопросу - поэтому, пожалуйста, отредактируйте ее, а также добавьте больше информации (почему собственные значения реальны? Или что означает наибольшее значение?) - см. Редактирование моего ответа.
Арнольд Ноймайер
извините, что я не объяснил себе ясно. Я также не объяснил ясно, что собственные значения являются сложными. Я проверяю, имеют ли какие-либо собственные значения величину один или больше.
мощность
1
Это имеет больше смысла, но теперь мой рецепт с работает хорошо, только если плохое собственное значение действительно реально> 1. С другой стороны, новая информация, вероятно, подразумевает, что у вас мало выбора, кроме как вычислить все собственные значения. - Пожалуйста, обновите ваш вопрос, чтобы передать дополнительную информацию! (я-A)-1
Арнольд Ноймайер
1
см. правку 2 в моем ответе
Арнольд Ноймайер,
7

Мощность Итерация (или метод питания), например , то , что описывает Дэн, всегда должен сходиться, хотя и со скоростью ,|λN-1/λN|

Если близко к λ n , оно будет медленным, но вы можете использовать экстраполяцию, чтобы обойти это. Это может показаться сложным, но реализация в псевдокоде приведена в статье.λN-1λN

Pedro
источник
1
что если | λ (n − 1) | = | λ (n) | ?
мощность
@power, тогда обычная итерация Power не будет сходиться. Я не знаю, насколько хорошо методы экстраполяции будут различать различные собственные значения, вам придется прочитать статью для этого.
Педро
2
|λN-1|знак равно|λN|λNλN-1
у вас есть ссылка на академический документ или книгу, которая поддерживает это? Кроме того, что если \ lambda_ {n} является сложным?
мощность
5
Если существует несколько различных собственных значений максимального модуля, степенная итерация сходится только при исключительных обстоятельствах. Как правило, он колеблется несколько непредсказуемым образом.
Арнольд Ноймайер
5

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

Вы можете прочитать больше о новом исследовании здесь:

http://math.berkeley.edu/~strain/273.F10/martinsson.tygert.rokhlin.randomized.decomposition.pdf

http://arxiv.org/abs/0909.4061

Этот код сделает это за вас:

http://cims.nyu.edu/~tygert/software.html

https://bitbucket.org/rcompton/pca_hgdp/raw/be45a1d9a7077b60219f7017af0130c7f43d7b52/pca.m

http://code.google.com/p/redsvd/

https://cwiki.apache.org/MAHOUT/stochastic-singular-value-decomposition.html

Если выбранного вами языка нет, вы можете довольно легко свернуть свой собственный рандомизированный SVD; это требует только умножения вектора матрицы с последующим вызовом готового SVD.

dranxo
источник
4

Здесь вы найдете алгоритмическое введение в алгоритм Якоби-Дэвидсона, который вычисляет максимальное собственное значение.

В этой статье рассматриваются математические аспекты. JD допускает общие (действительные или сложные) матрицы и может использоваться для вычисления диапазонов собственных значений.

Здесь вы можете найти различные реализации библиотек JDQR и JDQZ (включая интерфейс C, на который вы сможете ссылаться из R).

Deathbreath
источник
Я не смог найти никакой литературы, в которой явно говорится, что метод Якоби-Дэвидсона работает для реальной общей матрицы.
мощность
Если в каждой статье явно не указано ограничение и аргумент конвергенции опирается на ограничение, которое не имеет значения.
Смертельное дыхание
Вот еще одно объяснение JD. Рассматриваемые матрицы являются полностью общими. Специальная структура не используется, и результаты, специфичные для эрмитовых матриц, сравниваются и сопоставляются, например, сходимость для общих матриц квадратична, но для эрмитовых матриц кубическая.
Смертельное дыхание
Спасибо за это. Я не нахожу никакого кода C для общей матрицы, поэтому мне придется написать свой собственный. Ссылки на алгоритмы кажутся только для герметических матриц.
сила
1
@power вы также не найдете в литературе результата, в котором говорится, что стандартные реализации QR сходятся для реальной, общей матрицы - это открытая проблема, и, действительно, не так давно был найден контрпример для QR-кода в LAPACK.
Федерико Полони
2

В своем оригинальном сообщении вы говорите:

«Я также пытался использовать интерфейс ARPACK в пакете igraph R. Однако он дал ошибку для одной из моих матриц».

Мне было бы интересно узнать больше об ошибке. Если вы можете сделать эту матрицу общедоступной где-нибудь, мне было бы интересно попробовать ARPACK на ней.

Исходя из того, что я прочитал выше, я ожидаю, что ARPACK очень хорошо справится с извлечением самых больших (или нескольких из самых больших) собственных значений разреженной матрицы. Чтобы быть более конкретным, я ожидал бы, что методы Арнольди будут хорошо работать в этом случае, и это, конечно, то, на чем основан ARPACK.

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

Билл Грин
источник
Я посмотрю, смогу ли я найти свою работу тогда. Я работал над этим год назад.
власть
0

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

Это работает лучше всего, когда наибольшее собственное значение существенно больше, чем любое другое собственное значение. Если другое собственное значение близко по величине к наибольшему, потребуется некоторое время, чтобы сходиться, и может быть трудно определить , сходилось ли оно .

Дэн
источник
1
Спасибо, Дэн, однако: в моих матрицах некоторые другие собственные значения будут иметь такую ​​же (если не такую ​​же) величину, как наибольшее. Ваш метод похож на степенную итерацию и итерацию по Рэлею? Batterson и Smillie (1990) пишут, что для некоторых несимметричных матриц фактор-итерация Рэлея не будет сходиться. Баттерсон, С., Смилли, Дж. (1990) "Фактор Рэлея для несимметричных матриц", Математика вычислений, том 55, номер 191, P 169 - 178
мощность
Если другие собственные значения имеют ту же величину, что и наибольшее ..., не являются ли эти значения тоже "наибольшим"?
Ely
@EMS: они по-прежнему будут "самыми большими собственными значениями", но присутствие более чем одного из них все равно убьет конвергенцию.
Дан
Мне просто интересно, к какому собственному значению вы хотите, чтобы оно сходилось. Такие вещи, как коэффициент Рэлея / метод степеней, подразумеваются, когда существует четкое наибольшее собственное значение. Ваш вопрос просит найти наибольшее собственное значение, но тогда это звучит так, как будто это не совсем определено для вашей проблемы. Я просто введен в заблуждение заголовком поста.
Ely
-1

Пакет RARPACK у меня работает. И это кажется очень быстрым, поскольку это всего лишь интерфейс для ARPACK, стандартного пакета для разреженной линейной алгебры (имеется в виду вычисление нескольких собственных значений и собственных векторов).

HoangDT
источник
Добро пожаловать в SciComp! Как следует из вопроса, ARPACK не работает для OP, поэтому этот ответ не очень полезен.
Кристиан Клэйсон
@HoangDT Этот вопрос предшествовал RARPACK
мощность