Qr разложение матрицы

QR-разложение матрицы

Наука о данных и разложение матриц

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

QR—разложение — одно из самых полезных. Ему находится множество важных применений в науке о данных, статистике и анализе данных. Одно из подобных применений — вычисление решения задачи наименьших квадратов.

Содержание статьи

Задача наименьших квадратов

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

Вспомним задачу наименьших квадратов. Нужно решить уравнение ниже:

Проблема состоит в том, что не существует решения для β, потому что обычно, если у нас больше наблюдений, чем переменных, X не имеет обратного значения, следовательно, вычисление ниже невозможно:

Вместо этого попробуем найти некоторое β̂, решающее уравнение неидеально, но с минимально возможной ошибкой. Один из способов— минимизировать следующую целевую функцию, являющуюся функцией от β̂.

Минимизация этой суммы квадратов отклонений и дает имя задаче наименьших квадратов. Взятие производных по β̂ и приравнивание к нулю приведут к нормальным уравнениям и предоставят решение в замкнутой форме.

Это один из способов. Но можно использовать линейную алгебру. И вот здесь QR-разложение подойдет как нельзя лучше.

QR-разложение

Для начала давайте опишем это разложение. QR-разложение позволяет отобразить матрицу как произведение двух отдельных матриц Q и R.

Это означает, что:

Так как R квадратная, до тех пор, пока диагональные элементы не равны нулю, она также обратима. Если столбцы X линейно независимы, это всегда будет верным. Хотя, если в данных есть коллинеарность, все же будут возникать проблемы. Тем не менее — и в этом суть QR-разложения — прямоугольная и необратимая X может быть выражена как две обратимые матрицы! И вот это уже имеет смысл.

Решение задачи наименьших квадратов с помощью QR-разложения.

Теперь, зная, что из себя представляет QR-разложение, решим задачу наименьших квадратов следующим образом:

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

Реализация QR-разложения

Вычисление коэффициентов

Начнем с маленького примера, в котором смоделируем y и X, а затем решим уравнение, используя QR-разложение. Также мы сможем провести двойную проверку QR-разложения — работает ли оно и возвращает ли смоделированный X. Вот смоделированные переменные ответа:

Вот данные, которые мы используем для определения коэффициентов наименьших квадратов. В нашем распоряжении 3 переменные:

2. Убеждаемся в ортогональности Q.

3. И в том, что QR действительно возвращает исходную матрицу X.

Теперь вычисляем актуальные коэффициенты.

Мы получили в точности то же решение, что и для рассчитанных коэффициентов.

Реализация QR-разложения

Процесс Грама-Шмидта

Процесс Грама-Шмидта — это метод вычисления ортогональной матрицы Q, которая состоит из ортогональных или независимых единичных векторов и занимает такое же пространство, что и матрица X.

Выразим это следующим образом:

Получив полный набор ортогональных векторов, просто делим каждый вектор на его нормаль и помещаем их в матрицу:

Зная Q, легко вычисляем R:

Реализация в R и C++

Сравнение реализаций на R и C++

Заключение

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

QR-разложение матрицы

Наука о данных и разложение матриц

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

Содержание статьи

Задача наименьших квадратов

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

Вспомним задачу наименьших квадратов. Нужно решить уравнение ниже:

Проблема состоит в том, что не существует решения для β, потому что обычно, если у нас больше наблюдений, чем переменных, X не имеет обратного значения, следовательно, вычисление ниже невозможно:

Вместо этого попробуем найти некоторое β̂, решающее уравнение неидеально, но с минимально возможной ошибкой. Один из способов— минимизировать следующую целевую функцию, являющуюся функцией от β̂.

Минимизация этой суммы квадратов отклонений и дает имя задаче наименьших квадратов. Взятие производных по β̂ и приравнивание к нулю приведут к нормальным уравнениям и предоставят решение в замкнутой форме.

Это один из способов. Но можно использовать линейную алгебру. И вот здесь QR-разложение подойдет как нельзя лучше.

QR-разложение

Для начала давайте опишем это разложение. QR-разложение позволяет отобразить матрицу как произведение двух отдельных матриц Q и R.

Это означает, что:

Так как R квадратная, до тех пор, пока диагональные элементы не равны нулю, она также обратима. Если столбцы X линейно независимы, это всегда будет верным. Хотя, если в данных есть коллинеарность, все же будут возникать проблемы. Тем не менее — и в этом суть QR-разложения — прямоугольная и необратимая X может быть выражена как две обратимые матрицы! И вот это уже имеет смысл.

Решение задачи наименьших квадратов с помощью QR-разложения.

Теперь, зная, что из себя представляет QR-разложение, решим задачу наименьших квадратов следующим образом:

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

Реализация QR-разложения

Вычисление коэффициентов

Начнем с маленького примера, в котором смоделируем y и X, а затем решим уравнение, используя QR-разложение. Также мы сможем провести двойную проверку QR-разложения — работает ли оно и возвращает ли смоделированный X. Вот смоделированные переменные ответа:

Вот данные, которые мы используем для определения коэффициентов наименьших квадратов. В нашем распоряжении 3 переменные:

2. Убеждаемся в ортогональности Q.

3. И в том, что QR действительно возвращает исходную матрицу X.

Теперь вычисляем актуальные коэффициенты.

Мы получили в точности то же решение, что и для рассчитанных коэффициентов.

Реализация QR-разложения

Процесс Грама-Шмидта

Процесс Грама-Шмидта — это метод вычисления ортогональной матрицы Q, которая состоит из ортогональных или независимых единичных векторов и занимает такое же пространство, что и матрица X.

Выразим это следующим образом:

Получив полный набор ортогональных векторов, просто делим каждый вектор на его нормаль и помещаем их в матрицу:

Реализация в R и C++

Сравнение реализаций на R и C++

Заключение

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

Метод Хаусхолдера (отражений) QR-разложения квадратной матрицы, вещественный точечный вариант

Основные авторы описания: А.В.Фролов

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

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

1.2 Математическое описание алгоритма

1.3 Вычислительное ядро алгоритма

1.4 Макроструктура алгоритма

1.5 Схема реализации последовательного алгоритма

1.6 Последовательная сложность алгоритма

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

1.7 Информационный граф

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

1.8 Ресурс параллелизма алгоритма

Для понимания ресурса параллелизма в разложении матрицы порядка [math]n[/math] методом Хаусхолдера нужно рассмотреть критический путь графа.

Поэтому по грубой (без членов низших порядков) оценке критический путь метода Хаусхолдера будет идти через [math]\frac<2>[/math] умножений и [math]\frac<2>[/math] сложений/вычитаний.

При классификации по высоте ЯПФ, таким образом, метод Хаусхолдера относится к алгоритмам с квадратичной сложностью. При классификации по ширине ЯПФ его сложность будет также квадратичной (без расширения ряда ярусов, связанных с векторными операциями сложения, пришлось бы увеличить вдвое длину критического пути; при таком расширении сложность по ширине ЯПФ станет линейной).

Надо сказать, что здесь в оценках речь идёт именно о классическом способе реализации метода Хаусхолдера. Даже использование схем сдваивания или последовательно-параллельных для вычисления скалярных произведений уменьшает критический путь с квадратичного до степени 3/2 или линейно-логарифмического. Однако все эти широко распространённые методы пока не дали возможности снизить критический путь метода Хаусхолдера до линейного (как, скажем, у метода Гивенса).

1.9 Входные и выходные данные алгоритма

Входные данные: плотная квадратная матрица [math]A[/math] (элементы [math]a_[/math] ).

1.10 Свойства алгоритма

Соотношение последовательной и параллельной сложности, как хорошо видно, является линейным, что даёт определённый стимул для распараллеливания. Однако у наискорейшей ЯПФ ширина квадратична, что указывает на дисбаланс между загруженностями устройств при попытке её реально запрограммировать. Поэтому более практично даже при хорошей (быстрой) вычислительной сети оставить количество устройств (например, узлов кластера) линейным по размеру матрицы, что удвоит критический путь реализуемой ЯПФ.

При этом вычислительная мощность алгоритма, как отношение числа операций к суммарному объему входных и выходных данных, линейна.

Алгоритм в рамках выбранной версии полностью детерминирован.

Вычислительная погрешность в методе отражений (Хаусхолдера) растет линейно, как и в методе Гивенса (вращений).

2 Программная реализация алгоритма

2.1 Особенности реализации последовательного алгоритма

В варианте с кратчайшим критическим путём графа алгоритма (с использованием зависимости между обнуляемым вектором и направляющим вектором отражения) метод Хаусхолдера (отражений) QR-разложения квадратной вещественной матрицы на Фортране 77 можно записать так:

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

2.2 Возможные способы и особенности параллельной реализации алгоритма

Большинство пакетов от LINPACKа и LAPACKa до SCALAPACKa используют для QR-разложения матриц именно метод Хаусхолдера, правда, в различных модификациях (обычно с использованием BLAS). Существует большая подборка исследовательских работ по блочным версиям.

2.3 Результаты прогонов

2.4 Выводы для классов архитектур

Участник:Смирнова Александра/Нахождение собственных чисел квадратной матрицы методом QR разложения (3)

Эта работа успешно выполнена Преподавателю: в основное пространство, в подстраницу Данное задание было проверено и зачтено. Проверено Frolov и VadimVV.

Основные авторы описания:

Содержательно все пункты обсуждались авторами вместе, и их вклад был равным.

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

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

Суть базового QR-алгоритма заключается в итерационном приведении матрицы [math]A[/math] к некоторой подобной ей матрице [math]A_N[/math] при помощи QR-разложения. Матрица [math]A_N[/math] является правой верхней треугольной матрицей, а значит ее диагональ содержит собственные значения. В силу подобия матриц [math]A[/math] и [math]A_N[/math] их наборы собственных значений совпадают. Таким образом задача поиска собственных значений матрицы [math]A[/math] сводится к задаче выведения матрицы [math]A_N[/math] и поиска собственных значений для нее, что является тривиальной задачей.

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

В дальнейшем, в данной статье под модифицированным алгоритмом будет пониматься алгоритм, использующий сдвиги и матрицу с формой Хессенберга. Под базовым алгоритмом будет пониматься алгоритм, не использующий данные приемы.

1.2 Математическое описание алгоритма

1.2.1 QR-разложение матрицы

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

Данное разложение называется QR-разложением.

1.2.1.1 Метод Хаусхолдера (отражений) QR-разложения квадратной матрицы

1.2.1.2 Метод Гивенса (вращений) QR-разложения квадратной матрицы

1.2.2 QR-алгоритм нахождения собственных чисел

1.2.3 Сходимость алгоритма

Предположим, что для [math]\forall m[/math] выполнены следующие условия:

3. Ведущая подматрица порядка [math]m[/math] в [math]X^<-1>[/math] невырождена.

1.2.4 Вещественный вариант QR-алгоритма

В дальнейшем, в данной статье, матрицы, имеющие вышеописанную форму, будут называться квазитреугольными.

1.2.5 Приведение матрицы к форме Хессенберга

Матрицей, имеющей форму Хессенберга, называется такая матрица, у которой все элементы, находящиеся ниже первой поддиагонали, равны нулю ( [math]a_=0[/math] при [math]i\lt j+1[/math] ). Пример такой матрицы приведен ниже:

[math]\begin 1& 2& 3& 4\\ 2& 5& 6& 7\\ 0& 3& 8& 9\\ 0& 0& 4& 1 \end[/math]

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

1.2.6 QR-алгоритм со сдвигами

При этом сохраняется свойство подобия матриц [math]A_k[/math] и [math]A_[/math] :

1.3 Вычислительное ядро алгоритма

1.3.1 Базовый алгоритм

QR-алгоритм обладает двумя вычислительными ядрами, которые повторяются на каждой итерации:

1.3.2 Модифицированный алгоритм

Модифицированный QR-алгоритм обладает тремя вычислительными ядрами, первое вычисляется единожды, второе и третье повторяются на каждой итерации:

1.4 Макроструктура алгоритма

Как уже было описано ранее, базовый QR-алгоритм содержит в себе две макрооперации:

В случае модифицированного QR-алгоритма появляется еще макрооперация приведения матрицы к форме Хессенберга.

1.5 Схема реализации последовательного алгоритма

1.5.1 Базовый алгоритм

Опишем необходимый для реализации цикл при помощи псевдокода:

1.5.2 Модифицированный алгоритм

Опишем необходимый для реализации цикл при помощи псевдокода:

1.6 Последовательная сложность алгоритма

1.6.1 Базовый алгоритм

Подсчитаем сложность одной итерации QR-алгоритма (расчет сложностей для QR-разложения и перемножения матриц представлен в статьях, посвященных данным алгоритмам).

1.6.2 Модифицированный алгоритм

Отличия от базового алгоритма заключаются в следующем:

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

1.7 Информационный граф

1.7.1 Базовый алгоритм

На рисунках 1 и 2 изображен информационный граф QR-алгоритма. Вершины данного графа обозначают следующее:

Вершины V соответствуют операции сравнения с 0. Вершины F соответствуют проверке поддиагональных элементов на соответствие квазитреугольной форме, которая была описана в предыдущем разделе. Вершины & соответствуют логической операции «И».

Вершины -V соответствуют операции вычитания и сравнения с 0. Вершины & соответствуют логической операции «И».

1.7.2 Модифицированный алгоритм

На рисунках 5 и 6 изображен информационный граф модифицированного QR-алгоритма. Дополнительные вершины для этого графа обозначают следующее:

1.8 Ресурс параллелизма алгоритма

1.8.1 Базовый алгоритм

На макроуровне (который можно увидеть на информационном графе QR-алгоритма) алгоритм не обладает ресурсами параллелизма. Все макрооперации внутри итерации, а также сами итерации выполняются последовательно (за исключением операций Triag и Dif, которые могут выполняться параллельно). Основной ресурс параллелизма заложен отдельно в каждой из макроопераций. На каждой итерации алгоритм имеет следующие параллельные характеристики:

1.8.2 Модифицированный алгоритм

Отличия от базового алгоритма заключаются в следующем:

1.9 Входные и выходные данные алгоритма

Объем входных данных:

[math]n^2[/math] (необходимо хранить все элементы матрицы).

Объем выходных данных:

[math]n[/math] (квадратная матрица размера [math]n \times n[/math] имеет ровно [math]n[/math] собственных значений при этом некоторые из них могут быть комплексными).

1.10 Свойства алгоритма

2 Программная реализация алгоритма

2.1 Особенности реализации последовательного алгоритма

2.2 Локальность данных и вычислений

2.3 Возможные способы и особенности параллельной реализации алгоритма

2.4 Масштабируемость алгоритма и его реализации

Масштабируемость алгоритма (рис.7) исследовалась на основе реализации алгоритма в библиотеке ScaLAPACK v2.0.0: функция pdhseqr. Время работы было замерено для различных размеров матрицы, а также различного числа процессоров.

При увеличении числа процессоров время работы программы увеличивается при малом размере матрицы, распараллеливание алгоритма неэффективно. Однако при больших размерах матрицы (более [math] 10^3\times10^3[/math] ) время работы параллельной реализации значительно меньше последовательной.

Метод Гивенса (вращений) QR-разложения квадратной матрицы (вещественный точечный вариант)

Основные авторы описания: А.В.Фролов.

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

номера столбцов: [math]\begin \ _\quad _i\quad _ & \ & _\ \ _j\quad _\end[/math]

Именно [math]t[/math] обычно и хранится в соответствующем элементе массива.

1.2 Математическое описание алгоритма

Для завершения математического описания остаётся записать вычисление [2] матрицы вращения [math]T_[/math] и формулы её умножения на текущую промежуточную матрицу.

1.3 Вычислительное ядро алгоритма

1.4 Макроструктура алгоритма

1.5 Схема реализации последовательного алгоритма

Каждое «обнуление» элемента в позиции [math](j, i)[/math] состоит из двух шагов: а) вычисление параметров матрицы вращения [math]T_[/math] таких, чтобы обнулился элемент в позиции [math](j, i)[/math] ; б) умножение слева матрицы вращения [math]T_[/math] на текущую версию матрицы.

1.6 Последовательная сложность алгоритма

В последовательной версии основная сложность алгоритма определяется прежде всего массовыми операциями вращения. Они, если не учитывать возможную разреженность, составляют (в главном члене) [math]n^3/3[/math] операций комплексного умножения или, если не прибегать к ухищрениям, [math]4n^3/3[/math] вещественных умножений и [math]2n^3/3[/math] вещественных сложений/вычитаний.

При классификации по последовательной сложности, таким образом, метод Гивенса относится к алгоритмам с кубической сложностью.

1.7 Информационный граф

СОДЕРЖАНИЕ

Случаи и определения

Квадратная матрица

Любая вещественная квадратная матрица A может быть разложена как

Прямоугольная матрица

QL, RQ и LQ разложения

Аналогичным образом мы можем определить разложения QL, RQ и LQ, где L является нижней треугольной матрицей.

Вычисление QR-разложения

Использование процесса Грама – Шмидта

Пример

Таким образом, мы имеем

Связь с разложением RQ

Преимущества и недостатки

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

Использование отражений Хаусхолдера

Q может использоваться для отражения вектора таким образом, что все координаты, кроме одной, исчезают.

— верхнетреугольная матрица. Итак, с

Этот метод имеет большую числовую устойчивость, чем метод Грама – Шмидта, описанный выше.

Пример

Теперь обратите внимание:

так что у нас уже есть почти треугольная матрица. Нам нужно только обнулить запись (3, 2).

Тем же способом, что и выше, получаем матрицу преобразования Хаусхолдера

после выполнения прямого суммирования с 1, чтобы убедиться, что следующий шаг в процессе работает правильно.

Теперь мы находим

Или, до четырех десятичных цифр,

Матрица Q ортогональна, а R верхнетреугольная, поэтому A = QR является требуемым QR-разложением.

Преимущества и недостатки

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

Использование вращений Гивенса

Пример

Преимущества и недостатки

QR-разложение с помощью вращений Гивенса является наиболее сложным для реализации, поскольку порядок строк, необходимых для полного использования алгоритма, определить нетривиально. Однако он имеет значительное преимущество в том, что каждый новый нулевой элемент влияет только на строку с обнуляемым элементом ( i ) и строку выше ( j ). Это делает алгоритм вращения Гивенса более эффективным и распараллеливаемым, чем метод отражения Хаусхолдера. а я j <\ displaystyle a_ >

Связь с определителем или произведением собственных значений

Кроме того, поскольку определитель равен произведению собственных значений, мы имеем

Предположим, что QR-разложение для неквадратной матрицы A :

Из свойств SVD и определителя матрицы имеем

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

Вращение колонны

А п знак равно Q р ⟺ А знак равно Q р п Т <\ displaystyle AP = QR \ quad \ iff \ quad A = QRP ^ <\textf>>

Использование для решения линейных обратных задач

По сравнению с прямой обратной матрицей, обратные решения, использующие QR-разложение, более численно стабильны, о чем свидетельствуют их уменьшенные числа обусловленности [Parker, Geophysical Inverse Theory, Ch1.13].

Обобщения

Разложение Ивасавы обобщает QR-разложение на полупростые группы Ли.

Участник:F-morozov/Нахождение собственных чисел квадратной матрицы методом QR разложения (4)

Эта работа успешно выполнена Преподавателю: в основное пространство, в подстраницу Данное задание было проверено и зачтено. Проверено Frolov и VadimVV.

Описание составлено совместно, точно выделить вклад отдельных авторов невозможно.

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

1.2 Математическое описание алгоритма

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

1.2.1 QR-алгоритм

Пусть [math]A_0 = A[/math] — исходная матрица. Для [math]k = 1, 2, \ldots[/math] :

1.2.2 Сходимость QR-алгоритма

Предположим, что для матрицы [math]A \in \mathbb^[/math] выполнены следующие условия:

И пусть QR-алгоритм порождает последовательность матриц

1.2.3 Сокращение затрат на одну итерацию

Одна QR-итерация для матрицы общего вида требует [math]O(n^3)[/math] арифметических операций. Это очень большие затраты, даже если итераций не очень много (обычно их требуется не больше 5 на каждое собственное значение). Для решения данной проблемы используют тот факт, что с помощью отражений или вращений матрицу [math]A[/math] можно привести к унитарно подобной верхней хессенберговой (почти треугольной) матрице:

1.2.4 Уменьшение числа итераций

Такой подход называется QR-алгоритмом со сдвигами и имеет следуюший вид:

[math]A_0 = A[/math] — исходная матрица. Для [math]k = 1, 2, \ldots[/math] :

QR-алгоритм со сдвигами является частным случаем обобщенного QR-алгоритма (QR-алгоритма с мультисдвигами), в котором на каждой итерации используется полином [math]f_k[/math] :

1.3 Вычислительное ядро алгоритма

Основное время работы алгоритма приходится на вычисление QR-разложения на итерациях алгоритма.

Для получения QR-разложения могут использоваться метод Хаусхолдера, метод Гивенса или процесс ортогонализации Грама-Шмидта.

При использовании хессенберговой формы приведение к ней имеет сложность, сравнимую с последующими итерациями QR-алгоритма.

1.4 Макроструктура алгоритма

Основные элементы алгоритма:

1.5 Схема реализации последовательного алгоритма

Пример реализации на языке Python:

1.6 Последовательная сложность алгоритма

Рассмотрим сложность отдельных элементов алгоритма:

Таким образом, каждая итерация имеет кубическую сложность и, если алгоритм сходится через [math]N[/math] итераций, общая сложность составит [math]N \times O(n^3)[/math] арифметических операций.

1.7 Информационный граф

1.8 Ресурс параллелизма алгоритма

Несмотря на то что QR-алгоритм является последовательным, действия, выполняемые на каждой итерации, могут быть эффективно выполнены параллельно:

Таким образом, при классификации по высоте ярусно-параллельной формы сложность равна [math]N \times O(n)[/math] или [math]N \times O(n^2)[/math] в зависимости от используемого метода QR-разложения.

1.9 Входные и выходные данные алгоритма

1.10 Свойства алгоритма

2 Программная реализация алгоритма

2.1 Особенности реализации последовательного алгоритма

2.2 Локальность данных и вычислений

2.3 Возможные способы и особенности параллельной реализации алгоритма

2.4 Масштабируемость алгоритма и его реализации

Исследование проводилось на узлах раздела test суперкомпьютера «Ломоносов» [3] Суперкомпьютерного комплекса Московского университета.

Исследовалась реализация QR-алгоритма из библиотеки Intel MKL версии 11.2.0 (функции LAPACKE_sgehrd и LAPACKE_shseqr). Была написана программа на языке C, использовался компилятор Intel C++. Пример компиляции программы и запуска для матрицы размера [math]8000 \times 8000[/math] :

Оценивалось время вычисления собственных чисел произвольной квадратной матрицы для различных значений порядка матрицы (Matrix size) и числа потоков (Threads). Ввиду того, что количество операций зависит от числа итераций алгоритма, оценить производительность не удается.

Набор и границы значений изменяемых параметров запуска реализации алгоритма:

На рисунке приведен график времени работы выбранной реализации QR-алгоритма в зависимости от изменяемых параметров запуска.

Из рисунка видно, что алгоритм хорошо масштабируется. Время выполнения в зависимости от размера задачи растет значительно быстрее при использовании одного потока по сравнению с использованиемм нескольких потоков. Увеличение числа потоков при фиксированном размере задачи приводит к уменьшению времени выполнения. Для матрицы размера [math]8000 \times 8000[/math] время работы однопоточной реализации более чем в два раза превосходит время работы параллельной.

Исходный код программы на языке C и скрипт для обработки результатов доступны на github.

Нежное введение в матричную факторизацию для машинного обучения

Дата публикации 2018-02-16

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

В этом уроке вы узнаете о разложении матриц и о том, как их вычислять в Python.

После завершения этого урока вы узнаете:

Обзор учебника

Этот урок разделен на 4 части; они есть:

Что такое матричная декомпозиция?

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

Распространенной аналогией для разложения матриц является разложение чисел, такое как разложение 10 на 2 x 5. По этой причине разложение матриц также называется разложением матриц. Подобно факторизации реальных значений, существует много способов разложения матрицы, поэтому существует целый ряд различных методов разложения матрицы.

Далее мы подробнее рассмотрим каждый из этих методов.

LU Matrix Разложение

Разложение LU предназначено для квадратных матриц и разбивает матрицу на компоненты L и U.

Или без точечной записи.

Факторы L и U являются треугольными матрицами. Факторизация, возникающая в результате исключения, равна A = LU.

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

Вариант этого разложения, который численно более устойчив для решения на практике, называется разложением LUP, или разложением LU с частичным поворотом.

Строки родительской матрицы переупорядочены, чтобы упростить процесс декомпозиции, а дополнительная P-матрица указывает способ перестановки результата или возврата результата в исходный порядок. Есть и другие варианты LU.

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

Декомпозиция LU может быть реализована в Python с помощью функции lu (). Более конкретно, эта функция вычисляет разложение LPU.

В приведенном ниже примере сначала определяется квадратная матрица 3 × 3. Рассчитывается разложение LU, а затем исходная матрица восстанавливается по компонентам.

При выполнении примера сначала печатается определенная матрица 3 × 3, затем компоненты разложения P, L и U, а затем, наконец, восстанавливается исходная матрица.

Разложение QR-матрицы

QR-разложение предназначено для m x n матриц (не ограничено квадратными матрицами) и разбивает матрицу на Q и R компоненты.

Или без точечной записи.

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

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

QR-разложение может быть реализовано в NumPy с помощью функции qr (). По умолчанию функция возвращает матрицы Q и R с меньшими или «уменьшенными» размерами, что является более экономичным. Мы можем изменить это, чтобы вернуть ожидаемые размеры m x m для Q и m x n для R, указав аргумент mode как «complete», хотя это не требуется для большинства приложений.

В приведенном ниже примере определяется матрица 3 × 2, вычисляется QR-разложение, а затем восстанавливается исходная матрица из разложенных элементов.

При выполнении примера сначала печатается определенная матрица 3 × 2, затем элементы Q и R, а затем, наконец, восстановленная матрица, которая соответствует тому, с чего мы начали.

Холесский Разложение

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

Что касается наших интересов в машинном обучении, мы сосредоточимся на разложении Холецкого для вещественных матриц и проигнорируем случаи при работе с комплексными числами.

Разложение определяется следующим образом:

Или без точечной записи:

Разложение также может быть записано как произведение верхней треугольной матрицы, например:

Разложение Холецкого используется для решения линейных наименьших квадратов для линейной регрессии, а также для методов моделирования и оптимизации.

При разложении симметричных матриц разложение Холецкого почти в два раза эффективнее разложения LU и должно быть предпочтительным в этих случаях.

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

Разложение Cholesky может быть реализовано в NumPy путем вызова функции cholesky (). Функция возвращает только L, так как мы можем легко получить доступ к транспонированию L по мере необходимости.

Приведенный ниже пример определяет симметричную и положительно определенную матрицу 3 × 3 и вычисляет разложение Холецкого, после чего восстанавливается исходная матрица.

При выполнении примера сначала печатается симметричная матрица, затем нижняя треугольная матрица из разложения, за которым следует восстановленная матрица.

расширения

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

Если вы исследуете какое-либо из этих расширений, я хотел бы знать.

Дальнейшее чтение

Этот раздел предоставляет больше ресурсов по теме, если вы хотите углубиться.

книги

статьи

Резюме

В этом уроке вы обнаружили матричные декомпозиции и способы их вычисления в Python.

В частности, вы узнали:

У вас есть вопросы?
Задайте свои вопросы в комментариях ниже, и я сделаю все возможное, чтобы ответить.

Иллюстрированный самоучитель по MathCAD 12

QR-разложение

Среди матричных разложений особую роль играют ортогональные преобразования, обладающие свойством сохранения нормы вектора. Напомним из курса линейной алгебры, что матрица Q называется ортогональной, если Q T Q=I, где I – единичная матрица. Свойство сохранения нормы вектора при ортогональных преобразованиях, т. е. |Qx|=|x|, дает рецепт поиска псевдорешения вырожденных СЛАУ, а именно, замену исходной задачи минимизации невязки с «плохой» матрицей |Аx-b| задачей |Q T (Ax-b)|, в которой матрица уже будет «хорошей» благодаря специальному построению матрицы Q. Таким образом, ортогональные разложения используются при решении любых систем (в том числе с прямоугольной матрицей А, причем как переопределенных, так и недоопределенных).

Одним из важнейших вариантов ортогональных разложений некоторой матрицы А является QR-разложение вида A=QR, где Q – ортогональная матрица, а R – верхняя треугольная матрица:

Результатом действия функции qr (А) является матрица, составленная из матриц Q и R соответственно (подобно рассмотренному в предыдущем разделе LU-разложению). Выделить матрицы QR-разложения несложно при помощи встроенной функции выделения подматрицы submatrix (листинг 8.22).

Листинг 8.22. QR-разложение сингулярной матрицы:

Если бы исходная СЛАУ Aх=b не была вырожденной, то можно было сразу записать: Q T Q-Rx=Q T b, откуда следует (благодаря ортогональности матрицы Q): Rx=Q T b. Так как матрица к – треугольная, решение данной системы получается по формулам прямого хода. Если использовать нашу пользовательскую функцию решения треугольной системы (см. разд. 8.3. Т), то результат исходной СЛАУ запишется в виде одной строки кода: trg(R,Q T b) (соответствующий пример решения невырожденной СЛАУ вы найдете на компакт-диске).

Обратимся теперь к проблеме решения СЛАУ с сингулярной квадратной NxN или прямоугольной NxM матрицей А. Если ее ранг равен k (он может быть меньше N и м, как в листинге 8.22, где N=M=S, a k=2), то получающаяся треугольная матрица R имеет следующую структуру:

Где R1 – верхняя треугольная матрица, a R2 – просто некоторая матрица, а нули обозначают, в общем случае, нулевые матрицы соответствующих размеров.

Если система вырожденная, то она имеет бесконечное множество псевдорешений (векторов, минимизирующих норму невязки). При помощи QR-разложения можно сразу выписать одно из них (правда, не обладающее минимальной нормой). В нашем примере последняя строка матрицы R содержит одни нули, поэтому последняя компонента вектора псевдорешения х может быть абсолютно любой. Если положить х2=0, то остальные компоненты х определятся из треугольной СЛАУ R1x=Q T b, как это проиллюстрировано листингом 8.23.

Участник:Elijah/Нахождение собственных чисел квадратной матрицы методом QR разложения

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

QR-алгоритм — это численный метод в линейной алгебре, предназначенный для решения полной проблемы собственных значений, то есть отыскания всех собственных чисел и собственных векторов матрицы. Был разработан в конце 1950-х годов независимо В. Н. Кублановской и Дж. Фрэнсисом. [1]

Данный алгоритм хорошо подходит для поиска собственных значений матриц общего вида, порядок которых составляет порядке тысячи.

Суть алгоритма состоит в итеративном приведении матрицы A к верхнетреугольном виду с помощью последовательных QR разложений и матричный умножений. Данные матрицы на каждой итерации подобны между собой, так как данные умножения равносильны ортогональным преобразованиям. В следствии подобия, собственные значения этих матриц равны между собой.

1.2 Математическое описание алгоритма

1.2.1 QR-алгоритм нахождения собственных чисел

[math] A_ = R_k Q_k = Q_k^ <-1>Q_k R_k Q_k = Q_k^ <-1>A_k Q_k = Q_k^A_k Q_k, [/math]

то есть все матрицы [math]A_k[/math] являются подобными, то есть их собственные значения равны.

Пусть все диагональные миноры матрицы [math]A[/math] вырождены. Тогда последовательность матриц [math]A[/math] при [math] k \rightarrow \infty[/math] сходится по форме к клеточному правому треугольному виду, соответствующему клеткам с одинаковыми по модулю собственными значениями. [1]

1.2.2 Сходимость QR-алгоритма нахождения собственных чисел

Сходимость QR-алгоритма определяется как сходимость поддиагональной подматрицы матрицы [math]A_k[/math] ( [math]A_k = R_Q_[/math] ) к нулевой матрице.

Предположим, что для матрицы [math]A \in \mathbb^[/math] выполнены следующие условия:

Пусть QR-алгоритм порождает последовательность матриц вида:

1.2.3 Приведение матрицы к матрице Хессенберга

Для каждого [math] k = 1,2. n-1 [/math]

При этом матрица [math]P[/math] не вычисляется и используется напрямую. Вместо этого используются следующие вычисления

Таким образом перемножение матриц [math]P^A^<(k)>[/math] и [math]B^P^<(k)>[/math] имеет сложность [math]O(n^2) [/math] (умножения матрицы на вектор) вместо [math]O(n^3)[/math] (перемножение матриц)

1.2.4 Использование сдвигов

Такой подход называется QR-алгоритмом со сдвигами. Пусть [math]A_0 = A[/math] — исходная матрица. Тогда для [math]k = 1, 2, \ldots[/math] :

QR-алгоритм со сдвигами является частным случаем обобщенного QR-алгоритма, в котором на каждой итерации используется некоторый полином [math]f_k[/math] :

1.3 Вычислительное ядро алгоритма

1.3.1 Базовый алгоритм

У базового QR алгоритма есть два вычислительных ядра:

Матричное умножение может быть выполнено любым классическим способом: для плотных матриц [math]A[/math] (элементы [math]a_[/math] ) и [math]B[/math] (элементы [math]b_[/math] ) вычисляется плотная матрица [math]C[/math] с элементами [math]c_[/math] по формуле:

[math] \begin c_ = \sum_^ a_ b_, \quad i \in [1, m], \quad i \in [1, l]. \end [/math]

1.3.2 Алгоритм с преобразованием к матрице Хессенберга

Для QR алгоритма с использования матрицы Хессенберга:

Для преобразования исходной матрицы к матрице Хессенберга используются преобразования Хаусхолдера (см. раздел 1.2.3)

Более подробное описание вычислительных ядер можно найти по соответствующим ссылкам в разделе 1.4.

1.4 Макроструктура алгоритма

Как уже описано в описании ядра алгоритма, базовая версия QR-алгоритма на каждой итерации использует следующие алгоритмы:

В случае использования матрицы Хессенберга:

1.5 Схема реализации последовательного алгоритма

1.5.1 Базовый алгоритм

Последовательная реализация простейшего алгоритма состоит из некоторого числа итераций

Каждая итерация состоит из следующих действий (с начальным условием [math] A_0 = A[/math] ):

1.5.2 Алгоритм с преобразованием к матрице Хессенберга

Затем производится некоторое число итераций, состоящих из следующих действий (с начальным условием [math] A_0 = H[/math] ):

1.6 Последовательная сложность алгоритма

1.6.1 Последовательная сложность базового алгоритма

1.6.2 Последовательная сложность алгоритма с приведением матрицы к матрице Хессенбрга

Данный алгоритм состоит из двух этапов:

1.6.2.1 Сложность алгоритма приведения матрицы к матрице Хессенбрга

На каждой итерации алгоритма приведения матрицы к матрице Хессенбрга требуется вычислить:

1.6.2.2 Сложность алгоритма QR-разложения методом Гивенса

1.6.2.3 Итоговая сложность

Таким образом, общая вычислительная сложность алгоритма вычисления собственных значений с приведением матрицы к матрице Хессенбрга равна [math]O(n^3) + N * O(n^2)[/math]

1.7 Информационный граф

Макрограф алгоритма изображён на рисунке 1. Он представляет из себя верхний уровень алгоритма, на котором в цикле итеративно вызываются следующие операции:

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

В случае оптимизации алгоритма с приведением матрицы к форме Хессенберга, инициализация матрицы (int matrices (1)), будет включать в себя преобразование хаусхолдера, имеющее следующий граф алгоритма:

1.8 Ресурс параллелизма алгоритма

Все итерации базового QR-алгоритма производятся последовательно, поэтому на верхнем уровне алгоритм чисто последователен.

Основной ресурс параллелизма представлен на нижнем уровне, при реализации различных операций, используемых в алгоритме, таких как QR разложение методом вращений и матричное умножение.

1.8.1 Параллельная сложность базового QR-алгоритма

Посчитаем параллельную сложность каждой операции по отдельности, а затем, по полученным данным, и всего алгоритма целиком:

1.8.2 Параллельная сложность алгоритма с приведением матрицы к форме Хессенберга:

Снова вычислим параллельную сложность каждой операции по отдельности, а затем, по полученным данным, и всего алгоритма целиком:

1.9 Входные и выходные данные алгоритма

1.10 Свойства алгоритма

Соотношение последовательной и параллельной сложности для базового является квадратичным, что даёт хороший стимул для распараллеливания. Для оптимизированного алгоритма данное соотношение линейно.

Вычислительная погрешность растет линейно, из-за использования метода вращений для QR разложения.

2 Программная реализация

2.1 Особенности реализации последовательного алгоритма

Основным типом данных для QR алгоритма является матрица. В данной статье представлена реализация на языке C++, что накладывает определенные ограничения на формат хранени матриц: все матрицы хранятся в одномерном массиве в построчном (ROW_MAJOR) формате.

В простейшем случае QR-алгоритм выглядит следующим образом:

В последовательном случае функции dgeqrf, dorgqr и dgemm представляют собой последовательные функции стандарта BLAS, вычисляющие QR разложение(geqrf), матрицу Q после разложения (dorgqr) и матричное умножение (dgemm). Так же здесь используются простейшие функции set_Q_and_R и check_if_upper_triangular. Функция set_Q_and_R получает после QR разложения из матрицы A данные о матрицах Q и R. Функция check_if_upper_triangular проверяет, приведена ли матрица A к треугольном виду, и если да, происходит выход из цикла.

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

2.2 Локальность данных и вычислений

2.3 Возможные способы и особенности параллельной реализации алгоритма

Как уже говорилось в разделе 2.1, при использовании различных библиотек можно реализовать данный алгоритм практически на всех архитектурах. Далее более подробно будут рассмотрены примеры для некоторых из них. Особенности параллельных реализаций для базовых макроопераций (QR разложения и матричного умножения), уже описаны в соответствующих статьях (см. ссылки в разделе 1.4)

Реализация на многоядерных CPU (одном узле)

Для распарарллеливания на многоядерных CPU используются функции LAPACKE_sgeqrf, LAPACKE_sorgqr и cblas_sgemm. Распараллеливание оставшихся функций может быть произведено с использованием openMP, так как внутренние циклы в них независимы по данным.

Исходный код предложенной реализации представлен далее:

Реализация на GPU (и multi-GPU)

Для реализации алгоритма на GPU и multi-GPU можно использовать аналогичные функции dgemm и geqrf из библиотеки magma.

2.4 Масштабируемость алгоритма и его реализации

2.4.1 Масштабируемость реализации алгоритма

Для данного алгоритма была проверена масштабируемость оптимизированной версии с приведением матрицы к форме Хессенберга. Тесты проводились для случайно сгенерированной матрицы размера 8000×8000 на 1-14 параллельных вычислительных потоках.

Код, для которого проверлась масштабируемость, реализован с помощью функций библиотеки intel MKL, выглядит следующим образом:

Из данных графика видно, что алгоритм хорошо масштабируется при небольшом числе процессов (1-4), однако дальнейшее прибавление процессов не ведет к линейному ускорению.

2.5 Динамические характеристики и эффективность реализации алгоритма

2.6 Выводы для классов архитектур

2.7 Существующие реализации алгоритма

Последовательный QR алгоритм, а так же различные базовые функции для его реализации, можно найти в библиотеке LAPACK и её аналогах и расширениях.

Базовая версия QR алгоритма может быть реализована описанным в разделе 2.1 способом при помощи функций:

Для реализации оптимизированного алгоритма с приведением матрицы к матрице Хессенберга можно использовать следующие функции:

Подробная документация по упомянутым функциям может быть найдена здесь.

Параллельный QR алгоритм может быть реализован с помощью библиотек ScaLAPACK и PlaLAPACK.

Участник:Elijah/Нахождение собственных чисел квадратной матрицы методом QR разложения

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

QR-алгоритм — это численный метод в линейной алгебре, предназначенный для решения полной проблемы собственных значений, то есть отыскания всех собственных чисел и собственных векторов матрицы. Был разработан в конце 1950-х годов независимо В. Н. Кублановской и Дж. Фрэнсисом. [1]

Данный алгоритм хорошо подходит для поиска собственных значений матриц общего вида, порядок которых составляет порядке тысячи.

Суть алгоритма состоит в итеративном приведении матрицы A к верхнетреугольном виду с помощью последовательных QR разложений и матричный умножений. Данные матрицы на каждой итерации подобны между собой, так как данные умножения равносильны ортогональным преобразованиям. В следствии подобия, собственные значения этих матриц равны между собой.

1.2 Математическое описание алгоритма

1.2.1 QR-алгоритм нахождения собственных чисел

[math] A_ = R_k Q_k = Q_k^ <-1>Q_k R_k Q_k = Q_k^ <-1>A_k Q_k = Q_k^A_k Q_k, [/math]

то есть все матрицы [math]A_k[/math] являются подобными, то есть их собственные значения равны.

Пусть все диагональные миноры матрицы [math]A[/math] вырождены. Тогда последовательность матриц [math]A[/math] при [math] k \rightarrow \infty[/math] сходится по форме к клеточному правому треугольному виду, соответствующему клеткам с одинаковыми по модулю собственными значениями. [1]

1.2.2 Сходимость QR-алгоритма нахождения собственных чисел

Сходимость QR-алгоритма определяется как сходимость поддиагональной подматрицы матрицы [math]A_k[/math] ( [math]A_k = R_Q_[/math] ) к нулевой матрице.

Предположим, что для матрицы [math]A \in \mathbb^[/math] выполнены следующие условия:

Пусть QR-алгоритм порождает последовательность матриц вида:

1.2.3 Приведение матрицы к матрице Хессенберга

Для каждого [math] k = 1,2. n-1 [/math]

При этом матрица [math]P[/math] не вычисляется и используется напрямую. Вместо этого используются следующие вычисления

Таким образом перемножение матриц [math]P^A^<(k)>[/math] и [math]B^P^<(k)>[/math] имеет сложность [math]O(n^2) [/math] (умножения матрицы на вектор) вместо [math]O(n^3)[/math] (перемножение матриц)

1.2.4 Использование сдвигов

Такой подход называется QR-алгоритмом со сдвигами. Пусть [math]A_0 = A[/math] — исходная матрица. Тогда для [math]k = 1, 2, \ldots[/math] :

QR-алгоритм со сдвигами является частным случаем обобщенного QR-алгоритма, в котором на каждой итерации используется некоторый полином [math]f_k[/math] :

1.3 Вычислительное ядро алгоритма

1.3.1 Базовый алгоритм

У базового QR алгоритма есть два вычислительных ядра:

Матричное умножение может быть выполнено любым классическим способом: для плотных матриц [math]A[/math] (элементы [math]a_[/math] ) и [math]B[/math] (элементы [math]b_[/math] ) вычисляется плотная матрица [math]C[/math] с элементами [math]c_[/math] по формуле:

[math] \begin c_ = \sum_^ a_ b_, \quad i \in [1, m], \quad i \in [1, l]. \end [/math]

1.3.2 Алгоритм с преобразованием к матрице Хессенберга

Для QR алгоритма с использования матрицы Хессенберга:

Для преобразования исходной матрицы к матрице Хессенберга используются преобразования Хаусхолдера (см. раздел 1.2.3)

Более подробное описание вычислительных ядер можно найти по соответствующим ссылкам в разделе 1.4.

1.4 Макроструктура алгоритма

Как уже описано в описании ядра алгоритма, базовая версия QR-алгоритма на каждой итерации использует следующие алгоритмы:

В случае использования матрицы Хессенберга:

1.5 Схема реализации последовательного алгоритма

1.5.1 Базовый алгоритм

Последовательная реализация простейшего алгоритма состоит из некоторого числа итераций

Каждая итерация состоит из следующих действий (с начальным условием [math] A_0 = A[/math] ):

1.5.2 Алгоритм с преобразованием к матрице Хессенберга

Затем производится некоторое число итераций, состоящих из следующих действий (с начальным условием [math] A_0 = H[/math] ):

1.6 Последовательная сложность алгоритма

1.6.1 Последовательная сложность базового алгоритма

1.6.2 Последовательная сложность алгоритма с приведением матрицы к матрице Хессенбрга

Данный алгоритм состоит из двух этапов:

1.6.2.1 Сложность алгоритма приведения матрицы к матрице Хессенбрга

На каждой итерации алгоритма приведения матрицы к матрице Хессенбрга требуется вычислить:

1.6.2.2 Сложность алгоритма QR-разложения методом Гивенса

1.6.2.3 Итоговая сложность

Таким образом, общая вычислительная сложность алгоритма вычисления собственных значений с приведением матрицы к матрице Хессенбрга равна [math]O(n^3) + N * O(n^2)[/math]

1.7 Информационный граф

Макрограф алгоритма изображён на рисунке 1. Он представляет из себя верхний уровень алгоритма, на котором в цикле итеративно вызываются следующие операции:

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

В случае оптимизации алгоритма с приведением матрицы к форме Хессенберга, инициализация матрицы (int matrices (1)), будет включать в себя преобразование хаусхолдера, имеющее следующий граф алгоритма:

1.8 Ресурс параллелизма алгоритма

Все итерации базового QR-алгоритма производятся последовательно, поэтому на верхнем уровне алгоритм чисто последователен.

Основной ресурс параллелизма представлен на нижнем уровне, при реализации различных операций, используемых в алгоритме, таких как QR разложение методом вращений и матричное умножение.

1.8.1 Параллельная сложность базового QR-алгоритма

Посчитаем параллельную сложность каждой операции по отдельности, а затем, по полученным данным, и всего алгоритма целиком:

1.8.2 Параллельная сложность алгоритма с приведением матрицы к форме Хессенберга:

Снова вычислим параллельную сложность каждой операции по отдельности, а затем, по полученным данным, и всего алгоритма целиком:

1.9 Входные и выходные данные алгоритма

1.10 Свойства алгоритма

Соотношение последовательной и параллельной сложности для базового является квадратичным, что даёт хороший стимул для распараллеливания. Для оптимизированного алгоритма данное соотношение линейно.

Вычислительная погрешность растет линейно, из-за использования метода вращений для QR разложения.

2 Программная реализация

2.1 Особенности реализации последовательного алгоритма

Основным типом данных для QR алгоритма является матрица. В данной статье представлена реализация на языке C++, что накладывает определенные ограничения на формат хранени матриц: все матрицы хранятся в одномерном массиве в построчном (ROW_MAJOR) формате.

В простейшем случае QR-алгоритм выглядит следующим образом:

В последовательном случае функции dgeqrf, dorgqr и dgemm представляют собой последовательные функции стандарта BLAS, вычисляющие QR разложение(geqrf), матрицу Q после разложения (dorgqr) и матричное умножение (dgemm). Так же здесь используются простейшие функции set_Q_and_R и check_if_upper_triangular. Функция set_Q_and_R получает после QR разложения из матрицы A данные о матрицах Q и R. Функция check_if_upper_triangular проверяет, приведена ли матрица A к треугольном виду, и если да, происходит выход из цикла.

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

2.2 Локальность данных и вычислений

2.3 Возможные способы и особенности параллельной реализации алгоритма

Как уже говорилось в разделе 2.1, при использовании различных библиотек можно реализовать данный алгоритм практически на всех архитектурах. Далее более подробно будут рассмотрены примеры для некоторых из них. Особенности параллельных реализаций для базовых макроопераций (QR разложения и матричного умножения), уже описаны в соответствующих статьях (см. ссылки в разделе 1.4)

Реализация на многоядерных CPU (одном узле)

Для распарарллеливания на многоядерных CPU используются функции LAPACKE_sgeqrf, LAPACKE_sorgqr и cblas_sgemm. Распараллеливание оставшихся функций может быть произведено с использованием openMP, так как внутренние циклы в них независимы по данным.

Исходный код предложенной реализации представлен далее:

Реализация на GPU (и multi-GPU)

Для реализации алгоритма на GPU и multi-GPU можно использовать аналогичные функции dgemm и geqrf из библиотеки magma.

2.4 Масштабируемость алгоритма и его реализации

2.4.1 Масштабируемость реализации алгоритма

Для данного алгоритма была проверена масштабируемость оптимизированной версии с приведением матрицы к форме Хессенберга. Тесты проводились для случайно сгенерированной матрицы размера 8000×8000 на 1-14 параллельных вычислительных потоках.

Код, для которого проверлась масштабируемость, реализован с помощью функций библиотеки intel MKL, выглядит следующим образом:

Из данных графика видно, что алгоритм хорошо масштабируется при небольшом числе процессов (1-4), однако дальнейшее прибавление процессов не ведет к линейному ускорению.

2.5 Динамические характеристики и эффективность реализации алгоритма

2.6 Выводы для классов архитектур

2.7 Существующие реализации алгоритма

Последовательный QR алгоритм, а так же различные базовые функции для его реализации, можно найти в библиотеке LAPACK и её аналогах и расширениях.

Базовая версия QR алгоритма может быть реализована описанным в разделе 2.1 способом при помощи функций:

Для реализации оптимизированного алгоритма с приведением матрицы к матрице Хессенберга можно использовать следующие функции:

Подробная документация по упомянутым функциям может быть найдена здесь.

Параллельный QR алгоритм может быть реализован с помощью библиотек ScaLAPACK и PlaLAPACK.

Участник:Смирнова Александра/Нахождение собственных чисел квадратной матрицы методом QR разложения (3)

Эта работа успешно выполнена Преподавателю: в основное пространство, в подстраницу Данное задание было проверено и зачтено. Проверено Frolov и VadimVV.

Основные авторы описания:

Содержательно все пункты обсуждались авторами вместе, и их вклад был равным.

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

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

Суть базового QR-алгоритма заключается в итерационном приведении матрицы [math]A[/math] к некоторой подобной ей матрице [math]A_N[/math] при помощи QR-разложения. Матрица [math]A_N[/math] является правой верхней треугольной матрицей, а значит ее диагональ содержит собственные значения. В силу подобия матриц [math]A[/math] и [math]A_N[/math] их наборы собственных значений совпадают. Таким образом задача поиска собственных значений матрицы [math]A[/math] сводится к задаче выведения матрицы [math]A_N[/math] и поиска собственных значений для нее, что является тривиальной задачей.

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

В дальнейшем, в данной статье под модифицированным алгоритмом будет пониматься алгоритм, использующий сдвиги и матрицу с формой Хессенберга. Под базовым алгоритмом будет пониматься алгоритм, не использующий данные приемы.

1.2 Математическое описание алгоритма

1.2.1 QR-разложение матрицы

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

Данное разложение называется QR-разложением.

1.2.1.1 Метод Хаусхолдера (отражений) QR-разложения квадратной матрицы

1.2.1.2 Метод Гивенса (вращений) QR-разложения квадратной матрицы

1.2.2 QR-алгоритм нахождения собственных чисел

1.2.3 Сходимость алгоритма

Предположим, что для [math]\forall m[/math] выполнены следующие условия:

3. Ведущая подматрица порядка [math]m[/math] в [math]X^<-1>[/math] невырождена.

1.2.4 Вещественный вариант QR-алгоритма

В дальнейшем, в данной статье, матрицы, имеющие вышеописанную форму, будут называться квазитреугольными.

1.2.5 Приведение матрицы к форме Хессенберга

Матрицей, имеющей форму Хессенберга, называется такая матрица, у которой все элементы, находящиеся ниже первой поддиагонали, равны нулю ( [math]a_=0[/math] при [math]i\lt j+1[/math] ). Пример такой матрицы приведен ниже:

[math]\begin 1& 2& 3& 4\\ 2& 5& 6& 7\\ 0& 3& 8& 9\\ 0& 0& 4& 1 \end[/math]

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

1.2.6 QR-алгоритм со сдвигами

При этом сохраняется свойство подобия матриц [math]A_k[/math] и [math]A_[/math] :

1.3 Вычислительное ядро алгоритма

1.3.1 Базовый алгоритм

QR-алгоритм обладает двумя вычислительными ядрами, которые повторяются на каждой итерации:

1.3.2 Модифицированный алгоритм

Модифицированный QR-алгоритм обладает тремя вычислительными ядрами, первое вычисляется единожды, второе и третье повторяются на каждой итерации:

1.4 Макроструктура алгоритма

Как уже было описано ранее, базовый QR-алгоритм содержит в себе две макрооперации:

В случае модифицированного QR-алгоритма появляется еще макрооперация приведения матрицы к форме Хессенберга.

1.5 Схема реализации последовательного алгоритма

1.5.1 Базовый алгоритм

Опишем необходимый для реализации цикл при помощи псевдокода:

1.5.2 Модифицированный алгоритм

Опишем необходимый для реализации цикл при помощи псевдокода:

1.6 Последовательная сложность алгоритма

1.6.1 Базовый алгоритм

Подсчитаем сложность одной итерации QR-алгоритма (расчет сложностей для QR-разложения и перемножения матриц представлен в статьях, посвященных данным алгоритмам).

1.6.2 Модифицированный алгоритм

Отличия от базового алгоритма заключаются в следующем:

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

1.7 Информационный граф

1.7.1 Базовый алгоритм

На рисунках 1 и 2 изображен информационный граф QR-алгоритма. Вершины данного графа обозначают следующее:

Вершины V соответствуют операции сравнения с 0. Вершины F соответствуют проверке поддиагональных элементов на соответствие квазитреугольной форме, которая была описана в предыдущем разделе. Вершины & соответствуют логической операции «И».

Вершины -V соответствуют операции вычитания и сравнения с 0. Вершины & соответствуют логической операции «И».

1.7.2 Модифицированный алгоритм

На рисунках 5 и 6 изображен информационный граф модифицированного QR-алгоритма. Дополнительные вершины для этого графа обозначают следующее:

1.8 Ресурс параллелизма алгоритма

1.8.1 Базовый алгоритм

На макроуровне (который можно увидеть на информационном графе QR-алгоритма) алгоритм не обладает ресурсами параллелизма. Все макрооперации внутри итерации, а также сами итерации выполняются последовательно (за исключением операций Triag и Dif, которые могут выполняться параллельно). Основной ресурс параллелизма заложен отдельно в каждой из макроопераций. На каждой итерации алгоритм имеет следующие параллельные характеристики:

1.8.2 Модифицированный алгоритм

Отличия от базового алгоритма заключаются в следующем:

1.9 Входные и выходные данные алгоритма

Объем входных данных:

[math]n^2[/math] (необходимо хранить все элементы матрицы).

Объем выходных данных:

[math]n[/math] (квадратная матрица размера [math]n \times n[/math] имеет ровно [math]n[/math] собственных значений при этом некоторые из них могут быть комплексными).

1.10 Свойства алгоритма

2 Программная реализация алгоритма

2.1 Особенности реализации последовательного алгоритма

2.2 Локальность данных и вычислений

2.3 Возможные способы и особенности параллельной реализации алгоритма

2.4 Масштабируемость алгоритма и его реализации

Масштабируемость алгоритма (рис.7) исследовалась на основе реализации алгоритма в библиотеке ScaLAPACK v2.0.0: функция pdhseqr. Время работы было замерено для различных размеров матрицы, а также различного числа процессоров.

При увеличении числа процессоров время работы программы увеличивается при малом размере матрицы, распараллеливание алгоритма неэффективно. Однако при больших размерах матрицы (более [math] 10^3\times10^3[/math] ) время работы параллельной реализации значительно меньше последовательной.

Тензорные разложения и их применения. Лекция в Яндексе

Предыдущая лекция с Data Fest была посвящена алгоритмам, необходимым для построения нового вида поиска. Сегодняшний доклад тоже в некотором смысле про разные алгоритмы, а точнее про математику, лежащую в основе множества из них. О матричных разложениях зрителям рассказал доктор наук и руководитель группы вычислительных методов «Сколтеха» Иван Оселедец.

Под катом — расшифровка и большинство слайдов.

Есть небольшое — относительно machine learning — сообщество, которое занимается линейной алгеброй. Есть линейная алгебра, а в ней есть люди, которые занимаются мультилинейной алгеброй, тензорными разложениями. Я к ним отношусь. Там получено достаточно большое количество интересных результатов, но в силу способности сообщества самоорганизовываться и не выходить наружу, некоторые из этих результатов неизвестны. Кстати, это относится и к линейной алгебре тоже. Если посмотреть на статьи ведущих конференций NIPS, ICML и т. д., то можно увидеть достаточно много статей, связанных с тривиальными фактами линейной алгебры. Так происходит, но тем не менее.

Будем говорить о тензорных разложениях, но прежде надо сказать про матричные разложения.

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

Надо сказать, что есть LU-разложение, разложение Гаусса, QR-разложение, но я буду в основном говорить о сингулярном разложении, которое почему-то, например, в курсах линейной алгебры в наших университетах незаслуженно забывается. Читают такую ерунду с точки зрения вычислительной науки, как жорданова форма, которая нигде не используется. Сингулярное разложение обычно где-то в конце, опционально и т. д. Это неправильно. Если взять книжку Голуба и ван Лоуна, «Матричные вычисления», то она начинается ни с LU-, ни с QR-, а с сингулярного разложения. Это представление матрицы в виде произведения унитарной, диагональной и унитарной.

Про матричные разложения хорошо то, что про них действительно много чего известно. Для базовых разложений существуют эффективные алгоритмы, они устойчивые, существует ПО, которое вы можете вызвать в MATLAB, в Python, написаны они и на Фортране. Корректно расценивать это как единичную операция, которую можно посчитать, и вы уверены, что получите точный результат — что вам не нужно делать никаких градиентных спусков. Фактически это один квант операции. Если вы свели свою задачу к вычислению матричных разложений, то вы её, фактически, решили. Это, например, лежит в основе различных спектральных методов. Например, если есть ЕМ-алгоритм, он итерационный, у него медленная сходимость. А если вам удалось получить спектральный метод, то вы молодцы и у вас будет решение в один шаг.

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

Есть посмотреть на то, что такое сингулярное разложение и записать его в индексном виде, мы получим следующую формулу. Это приводит к классической старой идее — к разделению переменных, если вернуться к первым курсам университета. Раскладываем x от t, а здесь отделяем I от j, один дискретный индекс отделяем от другого, то есть представляем функцию от двух переменных в виде произведения функций от одной переменной. Ясно, что практически никакая функция в таком виде не приближается, но если вы напишете сумму таких функций, вы получите достаточно широкий класс.

Фактически, основное применение сингулярного разложения — это построение таких малоранговых аппроксимаций. Это матрица ранга r, и, естественно, нас интересует случай, когда ранг существенно меньше, чем общий размер матрицы.

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

Тем не менее, она не закрыта, но хочется сделать больше, чем два индекса.

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

А что будет, если вместо А(i, j) мы напишем А (i, jk), и тоже попробуем разделить переменные. Отсюда совершенно естественным образом возникает задача тензорных разложений. У нас есть тензор — под тензором я понимаю многомерный массив, простую двумерную таблицу, — и мы хотим эти многомерные массивы, грубо говоря, сжать, построить малопараметрическое представление и, обычно красивые слова говорят, — преодолеть проклятье размерности. На самом деле этот термин — он, как всегда, появился совершенно в другой связи с совершенно другой задачей, в работе Беллмана в 1961 году. Но сейчас он используется в том смысле, что n в степени d — это очень плохо. Если вы хотите увеличить d, то у вас будет экспоненциальный рост и т. д.

Факторизация в самом общем смысле: вы пытаетесь выбрать некий класс тензоров, который возникает в вашей прикладной задаче, и построить такое малопараметрическое представление, что по нему вы достаточно быстро сможете восстановить сам тензор. А параметров мало. В принципе, можно сказать, что нейронные сети в каком-то смысле решают задачу аппроксимации тензоров, потому что пытаются построить отображение из пространства параметров куда-то еще.

Если в задачах факторизации, то есть в представления в виде произведения простых объектов, есть матрицы и произведения матриц с тензорами, то определить тензоры — задача уже более сложная. Но если мы посмотрим на индексную запись неполного сингулярного разложения малоранговой аппроксимации, у нас совершенно естественным образом возникнет аналогичная формула. Вот она. Мы пытаемся приблизить тензор объектами, которые являются произведениями функции от одной дискретной переменной, функция t1 умножить на функцию t2. Идея совершенно естественная. Если слагаемое только одно, ничего вы так не приблизите. А если вы разрешаете брать достаточно небольшое количество параметров, возникает очень интересный класс тензоров, который явно или неявно используется во многих задачах методов. Можно сказать, что примерно такой формат используется в квантовой химии для представления волновой функции. Там еще антисимметричность, вместо произведения возникает определитель — но целая огромная область науки основана, фактически, на такой сепарабельной аппроксимации.

Поэтому, естественно, интересно изучить свойства, понять, что здесь происходит, и т. д.

Эта сущность — представление в виде суммы произведения функций от одной переменной — называется каноническим разложением. Впервые оно было предложено где-то в 1927 году. Потом использовалось далеко не в математике. Были работы в статистике, в журналах — таких как Chemometrics, Psychomatrix, — где люди получали кубы данных. Есть признаки — что-то получилось. Таким образом они строили многофакторные модели. Эти матрицы факторов можно интерпретировать как некие факторы. Из них, в свою очередь, можно делать какие-то выводы.

Если посчитать число параметров, мы получим d-n-r параметров. Если r маленький, то все хорошо. Но проблема, которую я пытаюсь решать — это что в принципе посчитать такую аппроксимацию, даже если известно, что она существует, — задача сложная. Можно даже сказать, что она NP-сложная в общем случае — как и задача вычисления ранга. Для матрицы это не так. Ранг матрицы вы можете посчитать — например, с помощью Гауссова исключения это можно сделать за полиномиальное число операций. Здесь всё не так.

Есть пример тензора размера 9х9х9, для которого до сих пор неизвестно точное значение минимального числа слагаемых. Известно, что оно не больше 23 и не меньше 21. Суперкомпьютеры, простые полиномиальные: это вопрос разрешимости системы полиномиальных уравнений. Но такая задача может быть очень сложной. А тензор имеет практическое значение. Ранг этого тензора связан с показателем сложности алгоритма умножения матриц. Штрассен-логарифм двоичный от семи. Соответственно, на самом деле задача вычисления экспоненты сводится к вычислению канонического ранга некоторого тензора.

Во многих приложениях это представление для какого-то конкретного тензора работает хорошо. Но в общем случае эта задача плохая.

Недавно появилась работа достаточно уважаемых людей, где они это представление с разделяющими переменными пытаются связать с глубокими нейронными сетями и доказывают некие результаты. На самом деле они, если их расшифровать, являются достаточно интересными алгебраическими результатами по тензорным разложениям. Если перевести результаты на язык тензорных разложений, получится, что одно тензорное разложение, каноническое, лучше другого. То, которое лучше, я еще не показал, а каноническое показал. Каноническое — плохое в этом смысле. Такой вывод обосновывается тем, что каноническая архитектура — shallow, неглубокая, а должна быть глубокая, более обширная и т. д.

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

У такой модели есть одно хорошее свойство: если решение все-таки найдено, оно единственное с точностью до тривиальных вещей. Например, не возникает проблемы матричных разложений, когда вы можете вставить s, s-1. Если каноническое разложение посчитано, оно единственное, и это хорошо. Но к сожалению, посчитать его иногда достаточно сложно.

Есть другое разложение Таккера. Это был какой-то психометрик. Идея в том, чтобы ввести некий здесь склеивающий множитель — чтобы не было связи каждого с каждым. Все хорошо, это разложение становится устойчивым, наилучшее приближение всегда существует. Но если вы попытаетесь его использовать, скажем, для d=10, то возникнет вспомогательный двумерный массив, который нужно будет хранить, и проклятье размерности, хотя и в ослабленном виде, останется.

Да, например, для трехмерных-четырехмерных задач разложение Таккера прекрасно работает. Существует достаточно большое количество работ, где оно применяется. Но тем не менее, наша конечная цель состоит не в этом.

Мы хотим получить нечто, где по умолчанию нет экспоненциального числа параметров, но где мы всё можем посчитать. Исходная идея была достаточно простая: если для матриц всё хорошо, давайте лепить матрицы из тензоров, многомерных массивов. Как это делать? Достаточно просто. У нас d индексов — давайте разбивать их на группы, делать матрицы, считать разложения.

Разбиваем. Одну часть индексов объявляем строчными индексами, другую часть — столбцовыми. Может, как-то переставляем. В MATLAB и Python это все делается командой reshape и transpose, считаем малоранговую аппроксимацию. Другой вопрос — как мы это делаем.

Если сделать это аккуратно, то — я только одну фразу скажу — возникнет один дополнительный индекс. Его надо считать новым индексом. Пусть у нас был девятимерный массив и мы разбили индексы на 5+4, разложили, получили один дополнительный индекс, получили один шестимерный и один пятимерный тензор и продолжаем их разлагать. В таком виде никакой экспоненциальной сложности не возникает, и получается представление, которое имеет вот такое число параметров. А если нам кто-то сказал, чтобы все получающиеся по дороге матрицы имели малый ранг? Получается такая рекурсивная структура. Она достаточно противная, программировать ее достаточно противно, мне было лень это делать.

В какой-то момент я понял, что можно выбрать самый простой вид данной структуры — который, тем не менее, достаточно мощный. Это то, что стало тензорным поездом, tensor-train. Мы такое предложили, имя прижилось. Потом оказалось, что мы были, естественно, не первыми, кто это придумал. В физике твердого тела оно было известно под названием matrix product states. Но я на каждом докладе повторяю, что хорошую вещь надо открыть как минимум дважды — иначе не очень понятно, хорошая ли она. По крайней мере, это представление открывалось минимум два, а то и два с половиной раза или даже три.

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

Мы представляем тензор в виде произведения объектов, которые зависят только от одного индекса. Фактически, речь идет о разделении переменных, но за одним маленьким исключением: теперь эти маленькие штучки — матрицы. Мы перемножаем, чтобы получить значения тензоров точки, эти матрицы. Первый элемент — строчка, потом матрица, матрица, столбец — мы получаем число. А указанные матрицы зависят от параметров, их можно хранить как трехмерный тензор. У нас есть коллекция матриц, вначале — коллекция строчек. Мы говорим: нам нужна нулевая строчка отсюда, пятая матрица отсюда, третья отсюда — мы их просто перемножаем. Ясно, что вычисление элемента занимает dr 2 операций. В индексном виде тоже. Компактный такой вид.

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

Здесь — некие свойства. У нас теперь не один ранг. В матрице только один ранг, столбцовый или строчный, но они одинаковые для матрицы, так повезло. Для тензоров рангов много. Число r вначале — это канонический ранг. Здесь у нас ранг d-1. Для каждого d свой ранг, у нас много рангов. Тем не менее, все эти определяющие сложность ранги и память для хранения такого представления можно посчитать, по крайней мере, теоретически — как ранги некоторых матриц, так называемых разверток. Берем многомерный массив, превращаем в матрицы, считаем ранги — ага, есть разложение с такими рангами.

Почему это называется форматом? Если у нас объекты хранятся в таком формате, алгебраические операции сложения, умножения, вычисления нормы можно сделать напрямую — не выходя из формата. При этом все равно что-то происходит. Например, если мы сложили два тензора, то результат будет иметь ранг, который не превосходит сумму рангов. Если проделывать такое много раз в каком-нибудь итерационном процессе, то ранг станет очень большим — даже несмотря на то, что по d зависимость линейная, а по r — квадратичная. Спокойно можно работать с рангом 100-200. Физики на кластерах умеют работать с рангами 4000-5000, но больше нельзя, потому что надо хранить достаточно большие плотные матрицы.

Но есть очень простой и красивый алгебраический алгоритм, позволяющий уменьшать ранги. Я имею некую допустимую точность в норме Фрабениуса и хочу найти минимальные ранги, которые дают точность. Но я одновременно уменьшаю ранги. Округление. Когда вы работаете с числами с плавающей точкой, вы же не храните все 100 цифр. Вы оставляете только 16 цифр. Мы находим те параметры, такой размер этих матриц, которые дают допустимую точность. Есть надежный робастный алгоритм, который это делает, — он занимает 20-30 строк в MATLAB или на Python.

Есть и еще одна красивая вещь. Физики ее точно не знали. Этим результатом я горжусь, он мой с моим тогдашним шефом, Евгением Евгеньевичем Тартышниковым. Если мы знаем, что данный нам тензор имеет такую структуру, то мы его точно можем восстановить по dnr 2 элементов. То есть мы, не вычисляя всех элементов, можем посмотреть на небольшое их количество и восстановить весь тензор, не используя ни базисов, ничего. Это интересный факт, он для матриц не насколько хорошо известный, каким мог бы быть. Есть d = 2, матрица ранга r. Ее можно восстановить по r столбцам и r строкам.

Когда у нас есть такое многообразие, мы можем минимизировать функционал с ограничением на ранги. Стандартная ситуация: доступа к элементам у нас нет, но есть некий функционал качества, который мы хотели бы минимизировать. Для этого возникает целое отдельное направление исследований вопроса, как оптимизировать с указанными ограничениями. Речь идет не о выпуклом множестве, но оно имеет много интересных топологических свойств, требующих отдельного изучения.

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

Пробежимся по тому, как считать ТТ-ранги. Берем матрицу, первые k индексов считаем строчными индексами, последние d-k — столбцовыми. Делаем огромную матрицу, считаем ранги. Тогда существует разложение вот с такими рангами.

Теорема конструктивная и дает алгоритм вычисления. Фактически, речь идет просто о последовательном вычислении сингулярных разложений небольших матриц. И можно доказать оценку устойчивости такого алгоритма. Там возникает множитель — корень из d-1 относительно сингулярных чисел, которые мы по дороге отбрасываем. Но в целом, если мы работаем с точностью 10-6, это ни на что не влияет.

Операции имеют линейную по d сложность. Вот моя любимая операция — ее можно написать короче всего. Если мы берем два тензора и хотим их поэлементно перемножить, то ядра просто перемножаются кронекоровым образом. Размер ядер был 10х10, значит будет 100х100. Если сделать так пару раз, уже будет 10000х10000, поэтому нужно округлять.

Округление можно сделать за dnr 3 операций. Можно точно построить это разложение, если кто-то дал нам уже неоптимальное представление.

То же самое для матриц. Если вам кто-то дал матрицу ранга 100 и эту факторизацию, построить факторизацию ранга 10 можно дешево, не считая всю матрицу.

Проекцию на множество тензоров с ограниченными рангами тоже можно вычислить устойчивым образом. Это очень важно в различных оптимизационных алгоритмах.

Крестовая аппроксимация. Я говорил, что матрицы и тензоры малого ранга можно восстанавливать по подмножеству элементов. Это незаслуженно малоизвестный факт. Сейчас в Америке он приобрел некую популярность, но тем не менее, всё надо переоткрывать и т. д. Однако он давно был известен и в России в группе Тартышникова, и в Германии.

Если у вас есть матрица ранга r, вы можете точно восстановить ее по r столбцам и r строкам, используя вот такую замечательную формуле. Взяли столбцы, взяли подматрицу, которая невырождена, и вот она формула.

Если на нее внимательно посмотреть, выяснится, что это формула гауссова исключения. Можно спокойно сказать, что все это придумал Гаусс и не мучиться, но тем не менее формула замечательная. Если матрица — миллион на миллион и вы знаете, что она ранга 5, то вам в ней достаточно посчитать 5 столбцов и 5 строчек, чтобы восстановить ее целиком.

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

Существуют достаточно эффективные эмпирические алгоритмы, которые, естественно, не будут работать в общем случае. Какой-то нехороший человек возьмет и здесь один элемент испортит. Вы такой элемент никогда не найдете, не просмотрев всю матрицу. Но если матрица обладает неким свойством гладкости и к ней нехорошие люди не допускаются, то для классов таких матриц можно доказать, что у вас все будет хорошо. И на практике эти методы работают совершенно блестяще. Устоявшийся алгоритм: вам дали матрицу, они это делают.

Совершенно то же самое можно сделать для тензоров. Если мы знаем, что тензоры имеют малые ранги в смысле ТТ-формата, то можем потыкать в некие подмножества и восстановить весь тензор. Естественно, тогда эффект может быть огромным. Вы сможете считать тензоры, число элементов в которых формально больше, чем число атомов во Вселенной. Причем делаться это будет на ноутбуке. Целиком никто никогда не посчитает, но тем не менее.

Существуют, как я сказал, эффективные эмпирические алгоритмы.

Самое серьезное направление исследований — решение оптимизационных задач, к которым сводится большинство практических. Когда у вас есть некий функционал, неизвестный объект может быть проиндексирован d-индексами. d-индекс введен либо естественным образом, либо искусственны. Это отдельная тема. Можно брать вектор длины 2 d и превращать его в тензор 2х2х2х2. И у него будут малые ранги, даже если взятый вектор состоит просто из значений одномерной функции. Такова одна из красивых вещей. Можно искусственно делать тензоры, вводить виртуальные размерности в данные, где их не было, и получать малые ранги.

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

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

Есть два пакета: ttpy для Python и TT-Toolbox для MATLAB. Документированы они фигово, но все же работают, и все алгоритмы, которые есть у нас, имеются и там. Базовые, продвинутые, придуманные нами и не только нами.

Есть ощущение, что можно что-то задействовать, но требуется серьезная доработка. Нельзя взять условный TensorFlow и применить. Мы взяли, применили — ничего не работает, надо думать. Когда что-то не работает, это замечательная ситуация. Значит, ты должен понять, почему оно не работает.

На самом деле, число параметров очень маленькое: 2 на d, которое 60, и на ранг, который 50. Формально мы решаем задачу, где очень много неизвестных.

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

Дальше есть пример, когда ТТ был применен напрямую к нейронным сетям. Это TensorNet, работа Саши Новикова и группы Дмитрия Ветрова, работа с NIPS прошлого года. Просто вместо полносвязного слоя взят слой, который имеет вот такую структуру. У нас получается нелинейный структурированный слой, но при этом число параметров падает в огромное количество раз, а точность меняется несильно.

Можно сжимать и сверточные слои. Была у нас работа с группой Виктора Лемпицкого. Он сегодня будет читать блестящую лекцию, всем рекомендую.

Там мы ТТ не применяли, применили каноническое разложение. Сверточный слой можно представить как свертку с четырехмерным тензором, который мы просто засунули в ПО, считающее малоранговые аппроксимации. Оказалось, что можно восемь раз сжать. Там есть детали, напрямую оно не сработало, но результат следующий: можно действительно сжать восемь раз, и без особой потери точности.

Вот пример из биологии. Что мне нравится в таких вещах: можно применять практически один и тот же код в биологии, в машинном обучении, в химии, да еще и не зная предметного языка, а переводя всё на язык многомерных массивов и запуская черный ящик. Конечно, этому ящику требуется доработка, но всё же примеров применения в совершенно разных областях у нас уже достаточно много. Речь идет о неком биологическом цикле, где задача шестимерная. Он возникает, когда медведь спит — ему нужно быстро выработать энергию, чтобы прогреться. Это футильный бесполезный цикл, работа группы Кристофа Шваба из Цюриха. Владимир Козеев, который был студентом нашей кафедры ЭВМ, когда я там работал. Есть метод Монте-Карло, решающий некое стохастическое дифференциальное уравнение. 1500 ядер — 10 5 секунд. В MATLAB это все считалось тысячу секунд. Хорошее ускорение.

Считали колебательные состояния молекул. Было 12 степеней свободы, 12-мерная задача, примерно 30 12 формальных неизвестных. По сравнению с некоторыми хитрыми и достаточно сложными химическими методами удалось добиться ускорения в десятки раз.

Последний пример связан с рекомендательными системами. Мой аспирант Женя Фролов, большой эксперт в рекомендательных системах, может часами про них говорить. Часто поговорить ему о них не с кем, поэтому если кто-то интересуется рекомендательными системами — напишите ему, пожалуйста. Недавно нашу статью приняли на RecSys. Он создал некий фреймворк для рекомендательных систем под названием Polara. Там мы тензоры применили к классической задаче построения рекомендаций, где у вас есть пользователи, фильмы и рейтинг. В этой области есть некая проблема: все современные методы рекомендаций не учитывают негативный фидбек. Здесь пример: человек говорит, что ему не нравится фильм «Лицо со шрамом», а ему SVD рекомендует посмотреть еще «Крестного отца». «История игрушек» ему бы больше понравилась. Статья написана очень хорошо. Идея такая: вместо рассмотрения матрицы «пользователь — фильм» мы рассматриваем матрицу «пользователь — фильм — рейтинг» и пишем единицу, если пользователь поставил фильму оценку. В каком-то смысле рейтинги становятся выровненными. Мы к указанному трехмерному тензору уже строим факторизацию, делаем практически то, что делается стандартно, и получаем что-то вроде парадигмы: если вам не нравится это, вам, возможно, понравится нечто другое. Мы можем давать осмысленные рекомендации на основе одного негативного фидбека, что достаточно сложно.

Есть классы задач, когда у вас имеется экспоненциальная точность — то есть когда ошибка по числу параметров бывает экспоненциальной. Есть надежный набор алгоритмов. Можно вложить и сказать, что это тоже deep learning. В данной работе товарища из Израиля как раз говорится, что речь идет про arithmetic circuit. Но для него есть надежный набор алгоритмов, можно считать.

Участник:Odbaev/Нахождение собственных чисел квадратной матрицы методом QR разложения (2)

Эта работа успешно выполнена Преподавателю: в основное пространство, в подстраницу Данное задание было проверено и зачтено. Проверено Frolov и VadimVV.

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

QR-алгоритм — это численный метод в линейной алгебре, предназначенный для решения полной проблемы собственных значений, то есть отыскания всех собственных чисел и собственных векторов матрицы. Был разработан в конце 1950-х годов независимо В.Н. Кублановской и Дж. Фрэнсисом.

1.2 Математическое описание алгоритма

1.2.1 Базовый QR-алгоритм

[math] A_ = R_k Q_k = Q_k^ <-1>Q_k R_k Q_k = Q_k^ <-1>A_k Q_k = Q_k^A_k Q_k, [/math]

то есть все матрицы [math]A_k[/math] являются подобными и их собственные значения равны.

Пусть все диагональные миноры матрицы [math]A[/math] не вырождены. Тогда последовательность матриц [math]A_k[/math] при [math]k \rightarrow \infty[/math] сходится по форме к клеточному правому треугольному виду, соответствующему клеткам с одинаковыми по модулю собственными значениями. [1]

1.2.2 Приведение матрицы к форме Хессенберга

Любая матрица [math]A[/math] может быть приведена к подобной матрице, имеющей форму Хессенберга. Подобная манипуляция позволяет значительно ускорить дальнейший процесс QR-разложения. Одним из способов приведения матрицы к форме Хессенберга является метод Хаусхолдера.

1.2.3 QR-алгоритм со сдвигами

При этом матрицы [math]A_k[/math] и [math]A_[/math] подобны, поэтому набор собственных чисел не изменяется.

1.3 Вычислительное ядро алгоритма

1.3.1 Базовый QR-алгоритм

Основными вычислительными ядрами алгоритма являются:

Существуют различные алгоритмы вычисления QR-разложения. Можно выделить такие, как метод Гивенса и метод Хаусхолдера.

Подробное описание вычислительных ядер содержится в соответствующих статьях, указанных в пункте 1.4 (Макроструктура алгоритма).

1.3.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

Основными вычислительными ядрами данного варианта алгоритма являются:

1.4 Макроструктура алгоритма

1.4.1 Базовый QR-алгоритм

QR-алгоритм на каждой итерации использует следующие макрооперации:

1.4.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

В данном варианте макроструктура алгоритма выглядит следующим образом:

1.5 Схема реализации последовательного алгоритма

1.5.1 Базовый QR-алгоритм

QR-алгоритм является итерационным.

На каждой итерации [math] k [/math] :

Алгоритм выполняется до сходимости матрицы [math] A_k[/math] к треугольному виду.

1.5.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

Как уже было описано ранее, данный алгоритм состоит из следуюших этапов:

1.6 Последовательная сложность алгоритма

1.6.1 Базовый QR-алгоритм

Распишем последовательную сложность каждого действия:

Таким образом, на каждой итерации последовательный алгоритм выполняет [math]O(n^3)[/math] операций. Для всего алгоритма потребуется выполнить [math]N * O(n^3)[/math] операций.

1.6.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

Распишем последовательную сложность каждой операции:

1.7 Информационный граф

На Рисунке 1 изображен макрограф описываемого QR-алгоритма, представляющего собой последовательность итераций алгоритма. На рисунке 2 представлен граф, на котором более подробно отображено содержимое каждой итерации. Информационные графы QR-разложения и матричного умножения представлены в соответствующих статьях.

QR-разложение

-разложение матрицы — представление матрицы в виде произведения унитарной (или ортогональной матрицы) и верхнетреугольной матрицы.

Определение

Матрица размера с комплексными элементами может быть представлена в виде:

где — унитарная матрица размера , а — верхнетреугольная матрица размера .

В частном случае, когда матрица состоит из вещественных чисел, является ортогональной матрицей (то есть , где — единичная матрица).

По аналогии, можно определить варианты этого разложения: -, -, и -разложения, где — нижнетреугольная матрица.

Свойства

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

Алгоритмы

-разложение может быть получено различными методами. Проще всего оно может быть вычислено, как побочный продукт в процессе Грама — Шмидта. На практике следует использовать модифицированный алгоритм Грама ― Шмидта, поскольку классический алгоритм обладает плохой численной устойчивостью.

Альтернативные алгоритмы для вычисления -разложения основаны на отражениях Хаусхолдера и вращениях Гивенса.

Полезное

Смотреть что такое «QR-разложение» в других словарях:

РАЗЛОЖЕНИЕ — РАСПАД И РАЗЛОЖЕНИЕ В словарь общерусского литературного языка впиталось много научных и специальных терминов. Выйдя за пределы профессиональной речи, эти термины расширяют свои значения и вовлекаются в новые фразеологические контексты.… … История слов

РАЗЛОЖЕНИЕ — РАЗЛОЖЕНИЕ, разложения, мн. нет, ср. (книжн.). 1. Действие по гл. разложить в 6, 7 и 8 знач. разлагать. Разложение воды на составные части. Разложение множителя. Разложение окиси на ртуть и кислород. Разложение армии врага. 2. Состояние по гл.… … Толковый словарь Ушакова

разложение — См. испытание … Словарь синонимов

РАЗЛОЖЕНИЕ — • РАЗЛОЖЕНИЕ, в биологии естественная деградация органического вещества с образованием более простых веществ, например, углекислого газа и воды. Разложение обычно является результатом деятельности таких организмов, как бактерии и грибки. За счет… … Научно-технический энциклопедический словарь

разложение в ряд Фурье — преобразование в ряд Фурье гармоническое разложение — [Л.Г.Суменко. Англо русский словарь по информационным технологиям. М.: ГП ЦНИИС, 2003.] Тематики информационные технологии в целом Синонимы преобразование в ряд Фурьегармоническое… … Справочник технического переводчика

разложение по степеням — разложение в ряд — [Л.Г.Суменко. Англо русский словарь по информационным технологиям. М.: ГП ЦНИИС, 2003.] Тематики информационные технологии в целом Синонимы разложение в ряд EN power expansion … Справочник технического переводчика

разложение Карунена-Лоэва — Представление случайного процесса второго порядка в виде суперпозиции ряда детерминированных функций, причем коэффициенты такого разложения являются взаимно некоррелированными случайными переменными. Данное разложение широко используется в теории … Справочник технического переводчика

разложение на множители — Разложение целого числа на его наибольшие сомножители (главные факторы). [http://www.rfcmd.ru/glossword/1.8/index.php?a=index d=23] Тематики защита информации EN factoring … Справочник технического переводчика

разложение — 1. Разложение органического вещества на его более простые составляющие. 2. Ослабление и разрушение горных пород вследствие химического выветривания. Syn.: распад … Словарь по географии

РАЗЛОЖЕНИЕ НА МНОЖИТЕЛИ — многочлена представление его в виде произведения двух или большего числа многочленов низших степеней. Напр.: х2 1 = (х 1)(х + 1) … Большой Энциклопедический словарь

РАЗЛОЖЕНИЕ СИЛЫ — замена одной силы, приложенной к телу, системой сил, производящей такое же механическое воздействие на тело, как и данная сила … Большой Энциклопедический словарь

В линейная алгебра, а QR-разложение, также известный как QR-факторизация или же QU факторизация это разложение из матрица А в продукт А = QR из ортогональная матрица Q и верхнетреугольная матрица р. QR-разложение часто используется для решения линейный метод наименьших квадратов проблема и является основой для конкретной алгоритм собственных значений, то QR-алгоритм.

Содержание

Случаи и определения

Квадратная матрица

Любая настоящая квадратная матрица А может быть разложен как

куда Q является ортогональная матрица (его столбцы ортогональный единичные векторы смысл Q Т Q = Q Q Т = я < Displaystyle Q ^ < extf> Q = QQ ^ < extf> = I> ) и р это верхний треугольная матрица (также называется прямоугольной матрицей, отсюда и название). Если А является обратимый, то факторизация единственна, если нам потребуются диагональные элементы р быть позитивным.

Если А имеет п линейно независимый столбцы, затем первый п столбцы Q для мужчин ортонормированный базис для пространство столбца из А. В более общем плане первый k столбцы Q образуют ортонормированный базис для охватывать из первых k столбцы А для любого 1 ≤k ≤ п. [1] Дело в том, что любой столбец k из А зависит только от первого k столбцы Q отвечает за треугольную формур. [1]

Прямоугольная матрица

В более общем плане мы можем разложить на множители сложные м×п матрица А, с м ≥ п, как продукт м×м унитарная матрица Q и м×п верхнетреугольная матрица р. Как внизу (м−п) ряды м×п верхнетреугольная матрица полностью состоит из нулей, часто бывает полезно разбить р, или оба р и Q:

Голуб и Ван Лоан (1996), §5.2) вызов Q₁р₁ в тонкая QR-факторизация из А; Трефетен и Бау называют это уменьшенная QR-факторизация. [1] Если А полный классифицировать п и потребуем, чтобы диагональные элементы р₁ положительны тогда р₁ и Q₁ уникальны, но в целом Q₂ не является. р₁ тогда равна верхнему треугольному множителю Разложение Холецкого из А * А (= А Т А если А реально).

QL, RQ и LQ разложения

Аналогичным образом мы можем определить разложения QL, RQ и LQ с L быть ниже треугольная матрица.

Вычисление QR-разложения

Существует несколько методов реального вычисления QR-разложения, например, с помощью Процесс Грама – Шмидта, Преобразования домовладельцев, или же Гивенса вращения. У каждого есть ряд преимуществ и недостатков.

Использование процесса Грама – Шмидта

Пример

Таким образом, мы имеем

Связь с разложением RQ

Преимущества и недостатки

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

Использование отражений Хаусхолдера

Q может использоваться для отражения вектора таким образом, что все координаты, кроме одной, исчезают.

(Stoer & Bulirsch, 2002 г., п. 225) и заменить транспонирование сопряженным транспонированием при построении Q ниже.

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

Это можно повторить для А‘ (получен из Q₁А удалив первую строку и первый столбец), в результате получится матрица Хаусхолдера Q′₂. Обратите внимание, что Q′₂ меньше чем Q₁. Поскольку мы хотим, чтобы он действительно работал на Q₁А вместо А′ Нам нужно развернуть его в верхний левый угол, заполнив 1, или в общем:

— верхнетреугольная матрица. Итак, с

У этого метода больше числовая стабильность чем вышеописанный метод Грама – Шмидта.

В следующей таблице указано количество операций в k-й шаг QR-разложения с помощью преобразования Хаусхолдера, предполагающий квадратную матрицу размером п.

Пример

Теперь обратите внимание:

так что у нас уже есть почти треугольная матрица. Нам нужно только обнулить запись (3, 2).

Возьмите (1, 1) незначительный, а затем снова примените процесс к

Тем же способом, что и выше, получаем матрицу преобразования Хаусхолдера

после выполнения прямого суммирования с 1, чтобы убедиться, что следующий шаг в процессе работает правильно.

Теперь мы находим

Или, до четырех десятичных цифр,

Преимущества и недостатки

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

Использование вращений Гивенса

QR разложения также могут быть вычислены с помощью серии Гивенса вращения. Каждый поворот обнуляет элемент в поддиагонали матрицы, образуя р матрица. Конкатенация всех вращений Гивенса образует ортогональную Q матрица.

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

Пример

Преимущества и недостатки

QR-разложение с помощью вращений Гивенса является наиболее сложным для реализации, поскольку порядок строк, необходимых для полного использования алгоритма, определить нетривиально. Однако у него есть существенное преимущество в том, что каждый новый нулевой элемент а я j < displaystyle a_ > влияет только на строку с обнуляемым элементом (i) и строку выше (j). Это делает алгоритм вращения Гивенса более эффективным и распараллеливаемым, чем метод отражения Хаусхолдера.

Связь с определителем или произведением собственных значений

Кроме того, поскольку определитель равен произведению собственных значений, мы имеем

Предположим, что QR-разложение для неквадратной матрицы А:

Из свойств СВД и определитель матрицы, имеем

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

Вращение колонны

А п = Q р ⟺ А = Q р п Т < displaystyle AP = QR quad iff quad A = QRP ^ < extf>>

Использование для решения линейных обратных задач

По сравнению с прямой обратной матрицей, обратные решения, использующие QR-разложение, более численно устойчивы, о чем свидетельствует их уменьшенная числа условий [Паркер, Геофизическая обратная теория, глава 1.13].

Обобщения

Разложение Ивасавы обобщает QR-разложение на полупростые группы Ли.

В линейная алгебра, а QR-разложение, также известный как QR-факторизация или же QU факторизация это разложение из матрица А в продукт А = QR из ортогональная матрица Q и верхнетреугольная матрица р. QR-разложение часто используется для решения линейный метод наименьших квадратов проблема и является основой для конкретной алгоритм собственных значений, то QR-алгоритм.

Содержание

Случаи и определения

Квадратная матрица

Любая настоящая квадратная матрица А может быть разложен как

куда Q является ортогональная матрица (его столбцы ортогональный единичные векторы смысл Q Т Q = Q Q Т = я < Displaystyle Q ^ < extf> Q = QQ ^ < extf> = I> ) и р это верхний треугольная матрица (также называется прямоугольной матрицей, отсюда и название). Если А является обратимый, то факторизация единственна, если нам потребуются диагональные элементы р быть позитивным.

Если А имеет п линейно независимый столбцы, затем первый п столбцы Q для мужчин ортонормированный базис для пространство столбца из А. В более общем плане первый k столбцы Q образуют ортонормированный базис для охватывать из первых k столбцы А для любого 1 ≤k ≤ п. [1] Дело в том, что любой столбец k из А зависит только от первого k столбцы Q отвечает за треугольную формур. [1]

Прямоугольная матрица

В более общем плане мы можем разложить на множители сложные м×п матрица А, с м ≥ п, как продукт м×м унитарная матрица Q и м×п верхнетреугольная матрица р. Как внизу (м−п) ряды м×п верхнетреугольная матрица полностью состоит из нулей, часто бывает полезно разбить р, или оба р и Q:

Голуб и Ван Лоан (1996), §5.2) вызов Q₁р₁ в тонкая QR-факторизация из А; Трефетен и Бау называют это уменьшенная QR-факторизация. [1] Если А полный классифицировать п и потребуем, чтобы диагональные элементы р₁ положительны тогда р₁ и Q₁ уникальны, но в целом Q₂ не является. р₁ тогда равна верхнему треугольному множителю Разложение Холецкого из А * А (= А Т А если А реально).

QL, RQ и LQ разложения

Аналогичным образом мы можем определить разложения QL, RQ и LQ с L быть ниже треугольная матрица.

Вычисление QR-разложения

Существует несколько методов реального вычисления QR-разложения, например, с помощью Процесс Грама – Шмидта, Преобразования домовладельцев, или же Гивенса вращения. У каждого есть ряд преимуществ и недостатков.

Использование процесса Грама – Шмидта

Пример

Таким образом, мы имеем

Связь с разложением RQ

Преимущества и недостатки

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

Использование отражений Хаусхолдера

Q может использоваться для отражения вектора таким образом, что все координаты, кроме одной, исчезают.

(Stoer & Bulirsch, 2002 г., п. 225) и заменить транспонирование сопряженным транспонированием при построении Q ниже.

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

Это можно повторить для А‘ (получен из Q₁А удалив первую строку и первый столбец), в результате получится матрица Хаусхолдера Q′₂. Обратите внимание, что Q′₂ меньше чем Q₁. Поскольку мы хотим, чтобы он действительно работал на Q₁А вместо А′ Нам нужно развернуть его в верхний левый угол, заполнив 1, или в общем:

— верхнетреугольная матрица. Итак, с

У этого метода больше числовая стабильность чем вышеописанный метод Грама – Шмидта.

В следующей таблице указано количество операций в k-й шаг QR-разложения с помощью преобразования Хаусхолдера, предполагающий квадратную матрицу размером п.

Пример

Теперь обратите внимание:

так что у нас уже есть почти треугольная матрица. Нам нужно только обнулить запись (3, 2).

Возьмите (1, 1) незначительный, а затем снова примените процесс к

Тем же способом, что и выше, получаем матрицу преобразования Хаусхолдера

после выполнения прямого суммирования с 1, чтобы убедиться, что следующий шаг в процессе работает правильно.

Теперь мы находим

Или, до четырех десятичных цифр,

Преимущества и недостатки

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

Использование вращений Гивенса

QR разложения также могут быть вычислены с помощью серии Гивенса вращения. Каждый поворот обнуляет элемент в поддиагонали матрицы, образуя р матрица. Конкатенация всех вращений Гивенса образует ортогональную Q матрица.

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

Пример

Преимущества и недостатки

QR-разложение с помощью вращений Гивенса является наиболее сложным для реализации, поскольку порядок строк, необходимых для полного использования алгоритма, определить нетривиально. Однако у него есть существенное преимущество в том, что каждый новый нулевой элемент а я j < displaystyle a_ > влияет только на строку с обнуляемым элементом (i) и строку выше (j). Это делает алгоритм вращения Гивенса более эффективным и распараллеливаемым, чем метод отражения Хаусхолдера.

Связь с определителем или произведением собственных значений

Кроме того, поскольку определитель равен произведению собственных значений, мы имеем

Предположим, что QR-разложение для неквадратной матрицы А:

Из свойств СВД и определитель матрицы, имеем

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

Вращение колонны

А п = Q р ⟺ А = Q р п Т < displaystyle AP = QR quad iff quad A = QRP ^ < extf>>

Использование для решения линейных обратных задач

По сравнению с прямой обратной матрицей, обратные решения, использующие QR-разложение, более численно устойчивы, о чем свидетельствует их уменьшенная числа условий [Паркер, Геофизическая обратная теория, глава 1.13].

Обобщения

Разложение Ивасавы обобщает QR-разложение на полупростые группы Ли.

QR-разложение

>

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

Резюме

Расширения

Можно вычислить разложение матрицы RQ или даже разложение QL и LQ, где матрица L имеет нижнюю треугольную форму.

Методы

Есть несколько методов для достижения этой декомпозиции:

У каждого из них есть свои достоинства и недостатки. Поскольку QR-разложение не уникально, разные методы дадут разные результаты.

Метод домохозяина

Пример

Мы сейчас только на диагональных нулях в 1- м столбце.

Чтобы повторить процесс, мы берем основную подматрицу

Тем же методом получаем

В итоге получаем

Стоимость и преимущества

Стоимость этого метода для матрицы размера n × n составляет: 4 / 3 n 3 Эта стоимость относительно высока (метод Холецкого для положительно определенных симметричных матриц находится в 1 / 3 п 3 ). Однако у метода Хаусхолдера есть значительное преимущество в том, что он намного более устойчив в численном отношении, ограничивая деление небольшим числом. Метод Гивенса, несмотря на еще более высокую стоимость, будет обеспечивать еще большую стабильность.

Метод Шмидта

Затем мы переставляем уравнения так, чтобы a _i были слева, используя тот факт, что e _i являются единичными векторами:

Пример

Возьмем матрицу примера

Напомним, что ортогональная матрица Q удовлетворяет

Затем мы можем вычислить Q с помощью средств Грама-Шмидта следующим образом:

В этом случае мы имеем:

Связь с разложением RQ

Метод Гивенса

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

Пример

Возьмем тот же пример

Произведение G ₁ A отменяет коэффициент a ₃₁ :

Участник:Odbaev/Нахождение собственных чисел квадратной матрицы методом QR разложения (2)

Эта работа успешно выполнена Преподавателю: в основное пространство, в подстраницу Данное задание было проверено и зачтено. Проверено Frolov и VadimVV.

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

QR-алгоритм — это численный метод в линейной алгебре, предназначенный для решения полной проблемы собственных значений, то есть отыскания всех собственных чисел и собственных векторов матрицы. Был разработан в конце 1950-х годов независимо В.Н. Кублановской и Дж. Фрэнсисом.

1.2 Математическое описание алгоритма

1.2.1 Базовый QR-алгоритм

[math] A_ = R_k Q_k = Q_k^ <-1>Q_k R_k Q_k = Q_k^ <-1>A_k Q_k = Q_k^A_k Q_k, [/math]

то есть все матрицы [math]A_k[/math] являются подобными и их собственные значения равны.

Пусть все диагональные миноры матрицы [math]A[/math] не вырождены. Тогда последовательность матриц [math]A_k[/math] при [math]k \rightarrow \infty[/math] сходится по форме к клеточному правому треугольному виду, соответствующему клеткам с одинаковыми по модулю собственными значениями. [1]

1.2.2 Приведение матрицы к форме Хессенберга

Любая матрица [math]A[/math] может быть приведена к подобной матрице, имеющей форму Хессенберга. Подобная манипуляция позволяет значительно ускорить дальнейший процесс QR-разложения. Одним из способов приведения матрицы к форме Хессенберга является метод Хаусхолдера.

1.2.3 QR-алгоритм со сдвигами

При этом матрицы [math]A_k[/math] и [math]A_[/math] подобны, поэтому набор собственных чисел не изменяется.

1.3 Вычислительное ядро алгоритма

1.3.1 Базовый QR-алгоритм

Основными вычислительными ядрами алгоритма являются:

Существуют различные алгоритмы вычисления QR-разложения. Можно выделить такие, как метод Гивенса и метод Хаусхолдера.

Подробное описание вычислительных ядер содержится в соответствующих статьях, указанных в пункте 1.4 (Макроструктура алгоритма).

1.3.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

Основными вычислительными ядрами данного варианта алгоритма являются:

1.4 Макроструктура алгоритма

1.4.1 Базовый QR-алгоритм

QR-алгоритм на каждой итерации использует следующие макрооперации:

1.4.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

В данном варианте макроструктура алгоритма выглядит следующим образом:

1.5 Схема реализации последовательного алгоритма

1.5.1 Базовый QR-алгоритм

QR-алгоритм является итерационным.

На каждой итерации [math] k [/math] :

Алгоритм выполняется до сходимости матрицы [math] A_k[/math] к треугольному виду.

1.5.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

Как уже было описано ранее, данный алгоритм состоит из следуюших этапов:

1.6 Последовательная сложность алгоритма

1.6.1 Базовый QR-алгоритм

Распишем последовательную сложность каждого действия:

Таким образом, на каждой итерации последовательный алгоритм выполняет [math]O(n^3)[/math] операций. Для всего алгоритма потребуется выполнить [math]N * O(n^3)[/math] операций.

1.6.2 QR-алгоритм с приведением к форме Хессенберга и сдвигами

Распишем последовательную сложность каждой операции:

1.7 Информационный граф

На Рисунке 1 изображен макрограф описываемого QR-алгоритма, представляющего собой последовательность итераций алгоритма. На рисунке 2 представлен граф, на котором более подробно отображено содержимое каждой итерации. Информационные графы QR-разложения и матричного умножения представлены в соответствующих статьях.

QR-разложение матрицы

Наука о данных и разложение матриц

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

QR—разложение — одно из самых полезных. Ему находится множество важных применений в науке о данных, статистике и анализе данных. Одно из подобных применений — вычисление решения задачи наименьших квадратов.

Содержание статьи

Задача наименьших квадратов

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

Вспомним задачу наименьших квадратов. Нужно решить уравнение ниже:

Проблема состоит в том, что не существует решения для β, потому что обычно, если у нас больше наблюдений, чем переменных, X не имеет обратного значения, следовательно, вычисление ниже невозможно:

Вместо этого попробуем найти некоторое β̂, решающее уравнение неидеально, но с минимально возможной ошибкой. Один из способов— минимизировать следующую целевую функцию, являющуюся функцией от β̂.

Минимизация этой суммы квадратов отклонений и дает имя задаче наименьших квадратов. Взятие производных по β̂ и приравнивание к нулю приведут к нормальным уравнениям и предоставят решение в замкнутой форме.

Это один из способов. Но можно использовать линейную алгебру. И вот здесь QR-разложение подойдет как нельзя лучше.

QR-разложение

Для начала давайте опишем это разложение. QR-разложение позволяет отобразить матрицу как произведение двух отдельных матриц Q и R.

Это означает, что:

Так как R квадратная, до тех пор, пока диагональные элементы не равны нулю, она также обратима. Если столбцы X линейно независимы, это всегда будет верным. Хотя, если в данных есть коллинеарность, все же будут возникать проблемы. Тем не менее — и в этом суть QR-разложения — прямоугольная и необратимая X может быть выражена как две обратимые матрицы! И вот это уже имеет смысл.

Решение задачи наименьших квадратов с помощью QR-разложения.

Теперь, зная, что из себя представляет QR-разложение, решим задачу наименьших квадратов следующим образом:

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

Реализация QR-разложения

Вычисление коэффициентов

Начнем с маленького примера, в котором смоделируем y и X, а затем решим уравнение, используя QR-разложение. Также мы сможем провести двойную проверку QR-разложения — работает ли оно и возвращает ли смоделированный X. Вот смоделированные переменные ответа:

Вот данные, которые мы используем для определения коэффициентов наименьших квадратов. В нашем распоряжении 3 переменные:

2. Убеждаемся в ортогональности Q.

3. И в том, что QR действительно возвращает исходную матрицу X.

Теперь вычисляем актуальные коэффициенты.

Мы получили в точности то же решение, что и для рассчитанных коэффициентов.

Реализация QR-разложения

Процесс Грама-Шмидта

Процесс Грама-Шмидта — это метод вычисления ортогональной матрицы Q, которая состоит из ортогональных или независимых единичных векторов и занимает такое же пространство, что и матрица X.

Выразим это следующим образом:

Получив полный набор ортогональных векторов, просто делим каждый вектор на его нормаль и помещаем их в матрицу:

Зная Q, легко вычисляем R:

Реализация в R и C++

Сравнение реализаций на R и C++

Заключение

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

В линейная алгебра, а QR-разложение, также известный как QR-факторизация или же QU факторизация это разложение из матрица А в продукт А = QR из ортогональная матрица Q и верхнетреугольная матрица р. QR-разложение часто используется для решения линейный метод наименьших квадратов проблема и является основой для конкретной алгоритм собственных значений, то QR-алгоритм.

Содержание

Случаи и определения

Квадратная матрица

Любая настоящая квадратная матрица А может быть разложен как

куда Q является ортогональная матрица (его столбцы ортогональный единичные векторы смысл Q Т Q = Q Q Т = я < Displaystyle Q ^ < extf> Q = QQ ^ < extf> = I> ) и р это верхний треугольная матрица (также называется прямоугольной матрицей, отсюда и название). Если А является обратимый, то факторизация единственна, если нам потребуются диагональные элементы р быть позитивным.

Если А имеет п линейно независимый столбцы, затем первый п столбцы Q для мужчин ортонормированный базис для пространство столбца из А. В более общем плане первый k столбцы Q образуют ортонормированный базис для охватывать из первых k столбцы А для любого 1 ≤k ≤ п. [1] Дело в том, что любой столбец k из А зависит только от первого k столбцы Q отвечает за треугольную формур. [1]

Прямоугольная матрица

В более общем плане мы можем разложить на множители сложные м×п матрица А, с м ≥ п, как продукт м×м унитарная матрица Q и м×п верхнетреугольная матрица р. Как внизу (м−п) ряды м×п верхнетреугольная матрица полностью состоит из нулей, часто бывает полезно разбить р, или оба р и Q:

Голуб и Ван Лоан (1996), §5.2) вызов Q₁р₁ в тонкая QR-факторизация из А; Трефетен и Бау называют это уменьшенная QR-факторизация. [1] Если А полный классифицировать п и потребуем, чтобы диагональные элементы р₁ положительны тогда р₁ и Q₁ уникальны, но в целом Q₂ не является. р₁ тогда равна верхнему треугольному множителю Разложение Холецкого из А * А (= А Т А если А реально).

QL, RQ и LQ разложения

Аналогичным образом мы можем определить разложения QL, RQ и LQ с L быть ниже треугольная матрица.

Вычисление QR-разложения

Существует несколько методов реального вычисления QR-разложения, например, с помощью Процесс Грама – Шмидта, Преобразования домовладельцев, или же Гивенса вращения. У каждого есть ряд преимуществ и недостатков.

Использование процесса Грама – Шмидта

Пример

Таким образом, мы имеем

Связь с разложением RQ

Преимущества и недостатки

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

Использование отражений Хаусхолдера

Q может использоваться для отражения вектора таким образом, что все координаты, кроме одной, исчезают.

(Stoer & Bulirsch, 2002 г., п. 225) и заменить транспонирование сопряженным транспонированием при построении Q ниже.

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

Это можно повторить для А‘ (получен из Q₁А удалив первую строку и первый столбец), в результате получится матрица Хаусхолдера Q′₂. Обратите внимание, что Q′₂ меньше чем Q₁. Поскольку мы хотим, чтобы он действительно работал на Q₁А вместо А′ Нам нужно развернуть его в верхний левый угол, заполнив 1, или в общем:

— верхнетреугольная матрица. Итак, с

У этого метода больше числовая стабильность чем вышеописанный метод Грама – Шмидта.

В следующей таблице указано количество операций в k-й шаг QR-разложения с помощью преобразования Хаусхолдера, предполагающий квадратную матрицу размером п.

Пример

Теперь обратите внимание:

так что у нас уже есть почти треугольная матрица. Нам нужно только обнулить запись (3, 2).

Возьмите (1, 1) незначительный, а затем снова примените процесс к

Тем же способом, что и выше, получаем матрицу преобразования Хаусхолдера

после выполнения прямого суммирования с 1, чтобы убедиться, что следующий шаг в процессе работает правильно.

Теперь мы находим

Или, до четырех десятичных цифр,

Преимущества и недостатки

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

Использование вращений Гивенса

QR разложения также могут быть вычислены с помощью серии Гивенса вращения. Каждый поворот обнуляет элемент в поддиагонали матрицы, образуя р матрица. Конкатенация всех вращений Гивенса образует ортогональную Q матрица.

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

Пример

Преимущества и недостатки

QR-разложение с помощью вращений Гивенса является наиболее сложным для реализации, поскольку порядок строк, необходимых для полного использования алгоритма, определить нетривиально. Однако у него есть существенное преимущество в том, что каждый новый нулевой элемент а я j < displaystyle a_ > влияет только на строку с обнуляемым элементом (i) и строку выше (j). Это делает алгоритм вращения Гивенса более эффективным и распараллеливаемым, чем метод отражения Хаусхолдера.

Связь с определителем или произведением собственных значений

Кроме того, поскольку определитель равен произведению собственных значений, мы имеем

Предположим, что QR-разложение для неквадратной матрицы А:

Из свойств СВД и определитель матрицы, имеем

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

Вращение колонны

А п = Q р ⟺ А = Q р п Т < displaystyle AP = QR quad iff quad A = QRP ^ < extf>>

Использование для решения линейных обратных задач

По сравнению с прямой обратной матрицей, обратные решения, использующие QR-разложение, более численно устойчивы, о чем свидетельствует их уменьшенная числа условий [Паркер, Геофизическая обратная теория, глава 1.13].

Обобщения

Разложение Ивасавы обобщает QR-разложение на полупростые группы Ли.

Метод Хаусхолдера (отражений) QR-разложения квадратной матрицы, вещественный точечный вариант

Основные авторы описания: А.В.Фролов

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

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

1.2 Математическое описание алгоритма

1.3 Вычислительное ядро алгоритма

1.4 Макроструктура алгоритма

1.5 Схема реализации последовательного алгоритма

1.6 Последовательная сложность алгоритма

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

1.7 Информационный граф

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

1.8 Ресурс параллелизма алгоритма

Для понимания ресурса параллелизма в разложении матрицы порядка [math]n[/math] методом Хаусхолдера нужно рассмотреть критический путь графа.

Поэтому по грубой (без членов низших порядков) оценке критический путь метода Хаусхолдера будет идти через [math]\frac<2>[/math] умножений и [math]\frac<2>[/math] сложений/вычитаний.

При классификации по высоте ЯПФ, таким образом, метод Хаусхолдера относится к алгоритмам с квадратичной сложностью. При классификации по ширине ЯПФ его сложность будет также квадратичной (без расширения ряда ярусов, связанных с векторными операциями сложения, пришлось бы увеличить вдвое длину критического пути; при таком расширении сложность по ширине ЯПФ станет линейной).

Надо сказать, что здесь в оценках речь идёт именно о классическом способе реализации метода Хаусхолдера. Даже использование схем сдваивания или последовательно-параллельных для вычисления скалярных произведений уменьшает критический путь с квадратичного до степени 3/2 или линейно-логарифмического. Однако все эти широко распространённые методы пока не дали возможности снизить критический путь метода Хаусхолдера до линейного (как, скажем, у метода Гивенса).

1.9 Входные и выходные данные алгоритма

Входные данные: плотная квадратная матрица [math]A[/math] (элементы [math]a_[/math] ).

1.10 Свойства алгоритма

Соотношение последовательной и параллельной сложности, как хорошо видно, является линейным, что даёт определённый стимул для распараллеливания. Однако у наискорейшей ЯПФ ширина квадратична, что указывает на дисбаланс между загруженностями устройств при попытке её реально запрограммировать. Поэтому более практично даже при хорошей (быстрой) вычислительной сети оставить количество устройств (например, узлов кластера) линейным по размеру матрицы, что удвоит критический путь реализуемой ЯПФ.

При этом вычислительная мощность алгоритма, как отношение числа операций к суммарному объему входных и выходных данных, линейна.

Алгоритм в рамках выбранной версии полностью детерминирован.

Вычислительная погрешность в методе отражений (Хаусхолдера) растет линейно, как и в методе Гивенса (вращений).

2 Программная реализация алгоритма

2.1 Особенности реализации последовательного алгоритма

В варианте с кратчайшим критическим путём графа алгоритма (с использованием зависимости между обнуляемым вектором и направляющим вектором отражения) метод Хаусхолдера (отражений) QR-разложения квадратной вещественной матрицы на Фортране 77 можно записать так:

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

2.2 Возможные способы и особенности параллельной реализации алгоритма

Большинство пакетов от LINPACKа и LAPACKa до SCALAPACKa используют для QR-разложения матриц именно метод Хаусхолдера, правда, в различных модификациях (обычно с использованием BLAS). Существует большая подборка исследовательских работ по блочным версиям.

2.3 Результаты прогонов

2.4 Выводы для классов архитектур

Метод QR-разложения

Лабораторная работа №3

· Знакомство с алгоритмом QR-разложения матрицы коэффициентов;

· Применение метода QR-разложения к решению систем линейных уравнений;

· Использование возможностей системы MATHCAD для выполнения QR-разложения.

Найти решение системы линейных уравнений с помощью метода QR-разложения, где

– матрица коэффициентов размера ,

, — столбец неизвестных и столбец свободных членов соответственно.

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

Пусть к шагу с номером найдены матрицы и , такие, что матрица имеет вид ,

— ортогональная матрица и . На первом шаге . Составим квадратную матрицу

, где .

, , .

Возьмем число таким, что ( ), а знак числа противоположен знаку числа . Если , то возьмем положительным. Пусть . Очевидно, что – единичный вектор. Положим , — ортогональная матрица (это легко проверить). Преобразование с матрицей будет преобразованием отражения относительно плоскости с нормальным вектором и будет вектор переводить в вектор .

В произведении первый столбец получается умножением матрицы на столбец и поэтому станет равным

.

Образуем матрицу следующего вида

,

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

.

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

Преимуществом QR-разложения является то, что элементы матрицы R не могут сильно превышать по модулю элементы матрицы A. Действительно, , т.е. каждый столбец матрицы R получается умножением ортогональной матрицы на соответствующей столбец матрицы A. Так как при умножении ортогональной матрицы на столбец вторая норма столбца не меняется, то нормы столбцов матрицы R совпадают с нормами соответствующих столбцов матрицы A. Норма каждого столбца ортогональной матрицы равна 1. Поэтому все элементы матрицы Q по модулю не больше 1.

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

С помощью QR-разложения можно найти разложение прямоугольной матрицы коэффициентов. Если матрица А размера m×n, где , то матрица будет иметь размер , а матрица .

Недостатком метода служит то, что его реализация требует в два раза больше операций, чем LU-разложение. Кроме того QR-разложение требует дополнительную память для хранения матрицы Q, в то время, как в LU-разложении матрицы L и U могут формироваться в памяти компьютера на месте, занимаемом матрицей A. Впрочем, недостатки, как и преимущества, сказываются только при больших значениях n.

Ход лабораторной работы:

1. Ввести матрицу коэффициентов A (n×n) и столбец свободных членов b.

2. На первом шаге .

3. Пусть .

4. Создать матрицу .

5. Создать векторы .

6. Найти число .

7. Найти единичный вектор .

8. Найти матрицы и . Увеличить j на 1: . Если перейти к шагу 9, иначе к шагу 4.

9. Ввести обозначение , создать матрицу .

10. Найти решение системы .

11. Найти решение системы .

12. Выполнить проверку.

1) Для проверки разложения:

Вычислить произведение матриц QR, сравнить с матрицей A.

2) Для проверки решения:

Посмотреть выполняется ли равенство ?

Найти решение системы линейных уравнений , где

Получим QR-разложение матрицы коэффициентов:

На первом шаге :

Создаем векторы , :

Находим число . Функция определяет знак элемента :

Найдем единичный вектор :

, — ортогональная матрица:

Составляем матрицу следующего вида . На первом шаге :

Найдем матрицу . Заметим, что матрица единичная, следовательно :

Увеличиваем j на 1. Составляем матрицу . Воспользуемся функцией Mathcad submatrix:

Найдем матрицу : Матрица (см. описание метода):

Аналогично находим матрицы и :

Составим матрицу :

Решим систему уравнений :

Решив вторую систему , получим:

Проверка:

Требования к отчету:

1. Отчет должен быть представлен в электронном или бумажном виде;

2. Отчет должен содержать:

· Расчеты и проверку.

3. Ответы на вопросы:

· Какова точность найденного решения;

· Преимущество метода QR-разложения по сравнению с

Иллюстрированный самоучитель по MathCAD 12

QR-разложение

Листинг 8.23. Поиск одного из псевдорешений вырожденной СЛАУ посредством QR-разложения (продолжение листинга 8.22):

Для того чтобы выбрать из всего множества псевдорешений (минимизирующих невязку исходной СЛАУ) нормальное псевдорешение (т. е. обладающее минимальной нормой), необходимо решить соответствующую задачу минимизации (см. разд. 8.2.3). Если построено QR-разложение, сделать это намного легче. Если произвольную компоненту решения обозначить переменной у: х₂=у, то она определится из соответствующей задачи оптимизации (листинг 8.24, 1-3 строки), а остальные составляющие самого решения х – из треугольной СЛАУ R1x=Q T b-R2y. В последней строке листинга выводится полученное нормальное псевдорешение, а также его норма и соответствующая норма невязки (которые полезно сравнить с результатом прошлого листинга). Вспомогательный рис. 8.13 помогает понять структуру минимизируемой функции из листинга 8.24.

Листинг 8.24. Нормальное псевдорешение вырожденной СЛАУ (продолжение листингов 8.22 и 8.23):

Рис. 8.13. Норма псевдорешения в зависимости от у (продолжение листинга 8.24)

Резюмируя практические аспекты применения QR-разложения, надо отметить, что алгоритмы решения СЛАУ на его основе практически одинаковы как для хорошо обусловленных, так и для сингулярных систем.

Примечание
На прилагаемом к книге компакт-диске вы найдете дополнительные примеры построения QR-разложения как для квадратной, так и прямоугольной матрицы А.

QR-разложение

-разложение матрицы — представление матрицы в виде произведения унитарной (или ортогональной матрицы) и верхнетреугольной матрицы.

Определение

Матрица размера с комплексными элементами может быть представлена в виде:

где — унитарная матрица размера , а — верхнетреугольная матрица размера .

В частном случае, когда матрица состоит из вещественных чисел, является ортогональной матрицей (то есть , где — единичная матрица).

По аналогии, можно определить варианты этого разложения: -, -, и -разложения, где — нижнетреугольная матрица.

Свойства

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

Алгоритмы

-разложение может быть получено различными методами. Проще всего оно может быть вычислено, как побочный продукт в процессе Грама — Шмидта. На практике следует использовать модифицированный алгоритм Грама ― Шмидта, поскольку классический алгоритм обладает плохой численной устойчивостью.

Альтернативные алгоритмы для вычисления -разложения основаны на отражениях Хаусхолдера и вращениях Гивенса.

Полезное

Смотреть что такое «QR-разложение» в других словарях:

РАЗЛОЖЕНИЕ — РАСПАД И РАЗЛОЖЕНИЕ В словарь общерусского литературного языка впиталось много научных и специальных терминов. Выйдя за пределы профессиональной речи, эти термины расширяют свои значения и вовлекаются в новые фразеологические контексты.… … История слов

РАЗЛОЖЕНИЕ — РАЗЛОЖЕНИЕ, разложения, мн. нет, ср. (книжн.). 1. Действие по гл. разложить в 6, 7 и 8 знач. разлагать. Разложение воды на составные части. Разложение множителя. Разложение окиси на ртуть и кислород. Разложение армии врага. 2. Состояние по гл.… … Толковый словарь Ушакова

разложение — См. испытание … Словарь синонимов

РАЗЛОЖЕНИЕ — • РАЗЛОЖЕНИЕ, в биологии естественная деградация органического вещества с образованием более простых веществ, например, углекислого газа и воды. Разложение обычно является результатом деятельности таких организмов, как бактерии и грибки. За счет… … Научно-технический энциклопедический словарь

разложение в ряд Фурье — преобразование в ряд Фурье гармоническое разложение — [Л.Г.Суменко. Англо русский словарь по информационным технологиям. М.: ГП ЦНИИС, 2003.] Тематики информационные технологии в целом Синонимы преобразование в ряд Фурьегармоническое… … Справочник технического переводчика

разложение по степеням — разложение в ряд — [Л.Г.Суменко. Англо русский словарь по информационным технологиям. М.: ГП ЦНИИС, 2003.] Тематики информационные технологии в целом Синонимы разложение в ряд EN power expansion … Справочник технического переводчика

разложение Карунена-Лоэва — Представление случайного процесса второго порядка в виде суперпозиции ряда детерминированных функций, причем коэффициенты такого разложения являются взаимно некоррелированными случайными переменными. Данное разложение широко используется в теории … Справочник технического переводчика

разложение на множители — Разложение целого числа на его наибольшие сомножители (главные факторы). [http://www.rfcmd.ru/glossword/1.8/index.php?a=index d=23] Тематики защита информации EN factoring … Справочник технического переводчика

разложение — 1. Разложение органического вещества на его более простые составляющие. 2. Ослабление и разрушение горных пород вследствие химического выветривания. Syn.: распад … Словарь по географии

РАЗЛОЖЕНИЕ НА МНОЖИТЕЛИ — многочлена представление его в виде произведения двух или большего числа многочленов низших степеней. Напр.: х2 1 = (х 1)(х + 1) … Большой Энциклопедический словарь

РАЗЛОЖЕНИЕ СИЛЫ — замена одной силы, приложенной к телу, системой сил, производящей такое же механическое воздействие на тело, как и данная сила … Большой Энциклопедический словарь

Инструменты сайта

Основное

Навигация

Информация

Действия

Содержание

Ортогональная матрица

Пример. Матрицы

Одно из подмножеств таких матриц имеют специальное название.

Произведение

Теорема. Произведение ортогональных матриц является ортогональной матрицей.

Определитель

Множество ортогональных матриц одинакового порядка образует группу относительно операции умножения.

Будет ли эта группа коммутативной, т.е. абелевой?

Норма

Теорема Кэли

Выражение через экспоненциал кососимметричной матрицы

Ортогональная матрица третьего порядка

Теорема [Эйлер]. Произвольное вращение твердого тела вокруг неподвижной точки, закрепленной в начале координат, задается ортогональной матрицей

Ортогональная матрица четвертого порядка: кватернионы

Характеристический полином, собственные числа

Произведение возратных полиномов будет возвратным полиномом. ♦

Пример. Найти характеристический полином и спектр матрицы

Вещественная жорданова нормальная форма

Применения

QR-разложение матрицы

Источники

[1]. Rodrigues O. Des lois géométriques qui régissent les déplacements d’un système solide dans l’espace, et de la variation des coordonnées provenant de ces déplacements considérés indépendamment des causes qui peuvent les produire. Journal de Mathématiques Pures et Appliquées. 1840. V. 5, 380–440

[2]. Проскуряков И.В. Сборник задач по линейной алгебре. М.Наука. 1974; задача N 896, задача 1570.

Алгоритм QR-разложения. Ортогональные матрицы и матрицы плоского вращения

Задача наименьших квадратов. Решение методом QR-разложения.

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

Эмпирические данные представляют собой значения неизвестной функции, полученные в результате эксперимента на сетке узлов, в общем случае, неравномерной.

Узлы сетки могут представлять собой моменты времени или пространственные координаты линейной электрической или механической системы.

Координаты узлов и значения функции в этих узлах объединяются в набор точек , где – координаты узлов, – значения функции.

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

Пример: . Такая функция часто используется при отслеживании дрейфа временных рядов в экономике.

В общем случае такая функция может не обеспечить требуемую точность восстановления зависимости. Поэтому в качестве аппроксимирующих функций используются общие многочлены по системе линейно-независимых функций:

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

В случае степенных функций .

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

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

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

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

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

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

Евклидова норма: или квадрат этой нормы: .

Задача наименьших квадратов возникает из задачи минимизации квадрата евклидовой нормы

) Ортогональное преобразование любого вектора не изменяет его длины (евклидова норма):

2) Ортогональное преобразование не изменяет углов между векторами в n-мерном пространстве.

Угол между векторами в n-мерном пространстве:

Для того, чтобы выражение имело смысл необходимо, чтобы его правая часть не превышала по модулю 1.

Это следует из неравенства Коши-Шварца:

Сохранение углов между векторами следует из равенства:

QR-разложение может быть осуществлено методами вращения и отражения.

Рассмотрим вращение вектора на плоскости.

Матрица вращения задается в виде: , – угол вращения.

Свойство ортогональной матрицы – сохранение угла между векторами.

Видно, что матрица вращения – ортогональная матрица:

Если принять, что или , то .

Рассмотрим систему линейных алгебраических уравнений второго порядка:

Найдем матрицу Q такую, что

, где

Рассмотрим систему уравнений с матрицей .

Плоской матрицей вращения называется матрица, имеющая вид:

Можно подтвердить, что матица Q является ортогональной матрицей с определителем, равным 1.

Применение указанной матрицы к i-му столбцу матрицы A: , дает вектор , имеющий в j-ой позиции 0. (Верхний индекс обозначает номер вектора).

Применяя к исходной матрице указанные плоские матрицы вращения получим матрицу:

С помощью указанных матриц вращения все элементы матрицы R ниже главной диагонали становятся равными нулю.

Для исключения соответствующих элементов, коэффициенты c и s определяются выражениями:

Произведение ортогональных матриц является ортогональной матрицей.

Чтобы решить задачу, матрица A дополняется матрицей , матрица A является произвольной.

Учитывая, что ортогональное преобразование вектора невязки:

второе слагаемое не зависит от коэффициентов многочлена, линейное значение первого слагаемого сводится к решению системы уравнений: , где R – верхняя треугольная матрица.

Решение задачи наименьших квадратов при , сводится к задаче решения системы алгебраических уравнений с верхней треугольной матрицей:

Чтобы применить метод QR-разложения к решению задачи наименьших квадратов, нужно привести матрицу A к квадратной форме:

матрица B – произвольная.

Исходное уравнение:

Матрица является квадратной . К этой системе можно применить метод QR-разложения.

Применяя метод вращения, уравнение записывается: ,

размерность вектора , размерность вектора ,

– ортогональная матрица, – верхняя треугольная матрица.

Разобьем матрицу R на блоки:

Умножая матрицу R справа на можем записать:

От неизвестных параметров зависит только первое слагаемое нормы невязки.

Минимальное значение этого слагаемого, если матрица A максимальный размер, определяется из уравнения: .

Таким образом, задача наименьших квадратов решается в два этапа.

На первом этапе осуществляется QR-разложение расширенной матрицы и определяются ее подматрицы и .

На втором этапе решается задача решения системы линейных алгебраических уравнений, матрица которой представлена в QR форме.

Матрица Q является ортогональной матрицей, т.е. матрицей, транспонирование которой совпадает с обратной матрицей.

Матрица R – верхняя треугольная матрица, решение которой осуществляется методом обратной подстановки.

Это следует из неравенства Коши-Шварца:

Сохранение углов между векторами следует из равенства:

QR-разложение может быть осуществлено методами вращения и отражения.

Рассмотрим вращение вектора на плоскости.

Матрица вращения задается в виде: , – угол вращения.

Свойство ортогональной матрицы – сохранение угла между векторами.

Видно, что матрица вращения – ортогональная матрица:

Если принять, что или , то .

Рассмотрим систему линейных алгебраических уравнений второго порядка:

Найдем матрицу Q такую, что

, где

Рассмотрим систему уравнений с матрицей .

Плоской матрицей вращения называется матрица, имеющая вид:

Можно подтвердить, что матица Q является ортогональной матрицей с определителем, равным 1.

Применение указанной матрицы к i-му столбцу матрицы A: , дает вектор , имеющий в j-ой позиции 0. (Верхний индекс обозначает номер вектора).

Применяя к исходной матрице указанные плоские матрицы вращения получим матрицу:

С помощью указанных матриц вращения все элементы матрицы R ниже главной диагонали становятся равными нулю.

Для исключения соответствующих элементов, коэффициенты c и s определяются выражениями:

Произведение ортогональных матриц является ортогональной матрицей.

Чтобы решить задачу, матрица A дополняется матрицей , матрица A является произвольной.

Учитывая, что ортогональное преобразование вектора невязки:

второе слагаемое не зависит от коэффициентов многочлена, линейное значение первого слагаемого сводится к решению системы уравнений: , где R – верхняя треугольная матрица.

Решение задачи наименьших квадратов при , сводится к задаче решения системы алгебраических уравнений с верхней треугольной матрицей:

Чтобы применить метод QR-разложения к решению задачи наименьших квадратов, нужно привести матрицу A к квадратной форме:

матрица B – произвольная.

Исходное уравнение:

Матрица является квадратной . К этой системе можно применить метод QR-разложения.

Применяя метод вращения, уравнение записывается: ,

размерность вектора , размерность вектора ,

– ортогональная матрица, – верхняя треугольная матрица.

Разобьем матрицу R на блоки:

Умножая матрицу R справа на можем записать:

От неизвестных параметров зависит только первое слагаемое нормы невязки.

Минимальное значение этого слагаемого, если матрица A максимальный размер, определяется из уравнения: .

Таким образом, задача наименьших квадратов решается в два этапа.

На первом этапе осуществляется QR-разложение расширенной матрицы и определяются ее подматрицы и .

На втором этапе решается задача решения системы линейных алгебраических уравнений, матрица которой представлена в QR форме.

Матрица Q является ортогональной матрицей, т.е. матрицей, транспонирование которой совпадает с обратной матрицей.

Матрица R – верхняя треугольная матрица, решение которой осуществляется методом обратной подстановки.

Дата добавления: 2015-04-24 ; Просмотров: 4686 ; Нарушение авторских прав? ;

Нам важно ваше мнение! Был ли полезен опубликованный материал? Да | Нет

Документация

Матричное разложение для решения линейных систем

Описание

Можно использовать decomposition объект dA со многими из тех же операторов вы можете использовать на исходной матрице коэффициентов A :

Комплексное сопряженное транспонирование dA’

Создание

Синтаксис

Описание

dA = decomposition( A ) возвращает разложение матричного A то, что можно использовать, чтобы решить линейные системы более эффективно. Тип разложения автоматически выбран на основе свойств входной матрицы.

Входные параметры

A — Матрица коэффициентов
матрица

Матрица коэффициентов. Матрица коэффициентов появляется в системе линейных уравнений слева как Ax = b или справа как xA = b.

Типы данных: single | double
Поддержка комплексного числа: Да

type — Тип разложения
‘auto’ (значение по умолчанию) | ‘qr’ | ‘cod’ | ‘lu’ | ‘ldl’ | ‘chol’ | ‘triangular’ | ‘permutedTriangular’ | ‘banded’ | ‘hessenberg’ | ‘diagonal’

Тип разложения, заданный как одна из опций в этих таблицах.

Эти опции работают на любую матрицу коэффициентов.

Матричное разложение A

‘auto’ (значение по умолчанию)

Q унитарно, R верхний треугольный, и P матрица перестановок.

Разложения QR дают решение методом наименьших квадратов.

R верхний треугольный, и оба Q и Z имейте ортонормированные столбцы.

Завершитесь ортогональные разложения дают минимальное решение методом наименьших квадратов нормы.

Для квадратных содействующих матриц также можно использовать эти опции.

Матричное разложение A

L является нижним треугольным, U верхний треугольный, и P матрица перестановок.

P ( R \ A ) Q = L U

P и Q матрицы сочетания и R диагональная матрица масштабирования.

L нижняя треугольная матрица с 1 с на диагонали, D диагональная матрица и P матрица перестановок.

P * S A S P = L D L *

S масштабирующаяся матрица.

A должно быть симметричным.

L является нижним треугольным.

P матрица перестановок.

A должен быть симметричен положительный определенный.

T является треугольным.

A должно быть треугольным.

T является треугольным, и оба P и Q матрицы сочетания.

A должно быть сочетание треугольной матрицы.

P матрица перестановок и оба L и U соединены.

Самый эффективный для матриц с низкой пропускной способностью. Смотрите bandwidth для получения дополнительной информации.

P матрица перестановок, L соединен, и U верхний треугольный.

A должен иметь нули ниже первой поддиагонали.

D является диагональным.

A должно быть диагональным.

triangularFlag — Отметьте, чтобы использовать только верхний или нижний треугольный фрагмент матрицы коэффициентов
‘upper’ | ‘lower’

‘triangular’ — Если и верхняя и нижняя треугольная матрица хранится в той же матрице, то используйте triangularFlag задавать только одну из треугольных матриц.

‘chol’ и ‘ldl’ — Используйте triangularFlag постараться не симметрировать почти симметричную матрицу коэффициентов.

Аргументы в виде пар имя-значение

Пример: dA = decomposition(A,’qr’,’CheckCondition’,false) выполняет разложение QR A и выключает предупреждения об условии матрицы коэффициентов, когда это используется, чтобы решить линейную систему.

Общие параметры

‘CheckCondition’ — Переключитесь, чтобы проверять условие матрицы коэффициентов
true (значение по умолчанию) | false

Переключитесь, чтобы проверять условие матрицы коэффициентов, заданной как разделенная запятой пара, состоящая из ‘CheckCondition’ и любой логический 1 TRUE ) или логический 0 ложь ). Если CheckCondition true и матрица коэффициентов плохо обусловливается или низкого ранга, затем решая линейные системы с помощью mldivide (\) или mrdivide (/) производит предупреждения.

Типы данных: логический

‘RankTolerance’ ‘RankTolerance’
неотрицательный скаляр

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

Примечание

Параметры разреженной матрицы

‘BandDensity’ — Порог плотности полосы
0.5 (значение по умолчанию) | скаляр

Плотность полосы задана как: (# ненули в полосе) / (# элементы в полосе). Значение 1.0 указывает, чтобы никогда не использовать полосный решатель.

‘LDLPivotTolerance’ — Допуск центра к LDL-разложению
0.01 (значение по умолчанию) | скаляр

Этим допуском центра является тот же тот ldl использование для действительных разреженных матриц.

‘LUPivotTolerance’ — Допуск центра к LU-факторизации
[0.1 0.001] (значение по умолчанию) | скаляр | вектор

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

Этим допуском центра является тот же тот lu использование для разреженных матриц.

Свойства

MatrixSize — Размер матрицы коэффициентов
вектор

Это свойство доступно только для чтения.

Размер матрицы коэффициентов, возвращенной как двухэлементный вектор-строка.

Типы данных: double

Type — Тип разложения
‘qr’ | ‘cod’ | ‘lu’ | ‘ldl’ | ‘chol’ | ‘triangular’ | ‘permutedTriangular’ | ‘banded’ | ‘hessenberg’ | ‘diagonal’

Это свойство доступно только для чтения.

Типы данных: char

CheckCondition — Переключитесь, чтобы проверять условие матрицы коэффициентов
true (значение по умолчанию) | false

Переключитесь, чтобы проверять условие матрицы коэффициентов, заданной как любой логический 1 TRUE ) или логический 0 ложь ). Если CheckCondition true и матрица коэффициентов плохо обусловливается или низкого ранга, затем решая линейные системы с помощью mldivide (\) или mrdivide (/) производит предупреждения.

Типы данных: логический

Datatype — Тип данных матрицы коэффициентов
‘double’ | ‘single’

Это свойство доступно только для чтения.

Типы данных: char

IsConjugateTransposed — Индикатор, что матрица коэффициентов является транспонированным сопряженным комплексным числом
false (значение по умолчанию) | true

Это свойство доступно только для чтения.

Типы данных: логический

IsReal — Индикатор, что матрица коэффициентов действительна
true | false

Это свойство доступно только для чтения.

Индикатор, что матрица коэффициентов действительна, возвратился как любой логический 1 TRUE ) или логический 0 ложь ). Значение false указывает, что матрица коэффициентов содержит комплексные числа.

Типы данных: логический

IsSparse — Индикатор, что матрица коэффициентов разреженна
true | false

Это свойство доступно только для чтения.

Индикатор, что матрица коэффициентов разреженна, возвратился как любой логический 1 TRUE ) или логический 0 ложь ).

Типы данных: логический

ScaleFactor — Мультипликативный масштабный коэффициент для матрицы коэффициентов
1 (значение по умолчанию) | скаляр

Это свойство доступно только для чтения.

Типы данных: double
Поддержка комплексного числа: Да

Функции объекта

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

Также можно проверять число обусловленности или ранг базовой матрицы decomposition объекты. Поскольку различные алгоритмы используются, результаты использования этих функций на decomposition объект может отличаться по сравнению с использованием тех же функций непосредственно на матрице коэффициентов.

Только форма с 1 входом rank(dA) поддерживается.

Осуществляет ту же проверку условия та обратная косая черта \ использование, чтобы определить, выдать ли предупреждение.

Примеры

Решение линейной системы с несколькими правыми сторонами

Покажите как использование decomposition объекты могут повысить эффективность решения Ax = b со многими правыми сторонами.

Теперь используйте decomposition объект решить ту же задачу.

Выбор Decomposition Type

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

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

Решите линейную систему с помощью вектора из единиц для правой стороны.

Задайте тип разложения, чтобы использовать ‘qr’ метод вместо ‘ldl’ по умолчанию метод. Это обеспечивает обратную косую черту (\), чтобы найти решение методом наименьших квадратов к проблеме вместо того, чтобы возвратить вектор NaN s.

Использование треугольного фрагмента матрицы

Задайте ‘upper’ использовать только верхний треугольный фрагмент входной матрицы в разложении.

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

Выключение матричных предупреждений условия

Создайте матрицу коэффициентов, которая плохо обусловливается. В этой матрице, составляя в среднем вместе первые два столбца производит третий столбец.

Решите линейную систему Ax = b использование вектора 1 с для правой стороны. mldivide производит предупреждение о создании условий матрицы коэффициентов.

Теперь создайте decomposition объект для матрицы и решает ту же линейную систему. Задайте ‘CheckCondition’ как false так, чтобы mldivide не проверяет условие матрицы коэффициентов. Даже при том, что то же решение возвращено, mldivide не отображает предупреждающее сообщение.

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

Расширенные возможности

Распределенные массивы
Большие массивы раздела через объединенную память о вашем кластере с помощью Parallel Computing Toolbox™.

Указания и ограничения по применению:

Если A разреженно, затем вводит ‘chol’ и ‘qr’ не поддержаны.

Если A разрежен прямоугольный, затем decomposition не поддержан.

Для получения дополнительной информации смотрите функции MATLAB Запуска с Распределенными Массивами (Parallel Computing Toolbox).

Смотрите также

Введенный в R2017b

Открытый пример

У вас есть модифицированная версия этого примера. Вы хотите открыть этот пример со своими редактированиями?

Документация MATLAB

Поддержка

© 1994-2020 The MathWorks, Inc.

1. Если смысл перевода понятен, то лучше оставьте как есть и не придирайтесь к словам, синонимам и тому подобному. О вкусах не спорим.

2. Не дополняйте перевод комментариями “от себя”. В исправлении не должно появляться дополнительных смыслов и комментариев, отсутствующих в оригинале. Такие правки не получится интегрировать в алгоритме автоматического перевода.

4. Не имеет смысла однотипное исправление перевода какого-то термина во всех предложениях. Исправляйте только в одном месте. Когда Вашу правку одобрят, это исправление будет алгоритмически распространено и на другие части документации.

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

В линейная алгебра, а QR-разложение, также известный как QR-факторизация или же QU факторизация это разложение из матрица А в продукт А = QR из ортогональная матрица Q и верхнетреугольная матрица р. QR-разложение часто используется для решения линейный метод наименьших квадратов проблема и является основой для конкретной алгоритм собственных значений, то QR-алгоритм.

Содержание

Случаи и определения

Квадратная матрица

Любая настоящая квадратная матрица А может быть разложен как

куда Q является ортогональная матрица (его столбцы ортогональный единичные векторы смысл Q Т Q = Q Q Т = я < Displaystyle Q ^ < extf> Q = QQ ^ < extf> = I> ) и р это верхний треугольная матрица (также называется прямоугольной матрицей, отсюда и название). Если А является обратимый, то факторизация единственна, если нам потребуются диагональные элементы р быть позитивным.

Если А имеет п линейно независимый столбцы, затем первый п столбцы Q для мужчин ортонормированный базис для пространство столбца из А. В более общем плане первый k столбцы Q образуют ортонормированный базис для охватывать из первых k столбцы А для любого 1 ≤k ≤ п. [1] Дело в том, что любой столбец k из А зависит только от первого k столбцы Q отвечает за треугольную формур. [1]

Прямоугольная матрица

В более общем плане мы можем разложить на множители сложные м×п матрица А, с м ≥ п, как продукт м×м унитарная матрица Q и м×п верхнетреугольная матрица р. Как внизу (м−п) ряды м×п верхнетреугольная матрица полностью состоит из нулей, часто бывает полезно разбить р, или оба р и Q:

Голуб и Ван Лоан (1996), §5.2) вызов Q₁р₁ в тонкая QR-факторизация из А; Трефетен и Бау называют это уменьшенная QR-факторизация. [1] Если А полный классифицировать п и потребуем, чтобы диагональные элементы р₁ положительны тогда р₁ и Q₁ уникальны, но в целом Q₂ не является. р₁ тогда равна верхнему треугольному множителю Разложение Холецкого из А * А (= А Т А если А реально).

QL, RQ и LQ разложения

Аналогичным образом мы можем определить разложения QL, RQ и LQ с L быть ниже треугольная матрица.

Вычисление QR-разложения

Существует несколько методов реального вычисления QR-разложения, например, с помощью Процесс Грама – Шмидта, Преобразования домовладельцев, или же Гивенса вращения. У каждого есть ряд преимуществ и недостатков.

Использование процесса Грама – Шмидта

Пример

Таким образом, мы имеем

Связь с разложением RQ

Преимущества и недостатки

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

Использование отражений Хаусхолдера

Q может использоваться для отражения вектора таким образом, что все координаты, кроме одной, исчезают.

(Stoer & Bulirsch, 2002 г., п. 225) и заменить транспонирование сопряженным транспонированием при построении Q ниже.

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

Это можно повторить для А‘ (получен из Q₁А удалив первую строку и первый столбец), в результате получится матрица Хаусхолдера Q′₂. Обратите внимание, что Q′₂ меньше чем Q₁. Поскольку мы хотим, чтобы он действительно работал на Q₁А вместо А′ Нам нужно развернуть его в верхний левый угол, заполнив 1, или в общем:

— верхнетреугольная матрица. Итак, с

У этого метода больше числовая стабильность чем вышеописанный метод Грама – Шмидта.

В следующей таблице указано количество операций в k-й шаг QR-разложения с помощью преобразования Хаусхолдера, предполагающий квадратную матрицу размером п.

Пример

Теперь обратите внимание:

так что у нас уже есть почти треугольная матрица. Нам нужно только обнулить запись (3, 2).

Возьмите (1, 1) незначительный, а затем снова примените процесс к

Тем же способом, что и выше, получаем матрицу преобразования Хаусхолдера

после выполнения прямого суммирования с 1, чтобы убедиться, что следующий шаг в процессе работает правильно.

Теперь мы находим

Или, до четырех десятичных цифр,

Преимущества и недостатки

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

Использование вращений Гивенса

QR разложения также могут быть вычислены с помощью серии Гивенса вращения. Каждый поворот обнуляет элемент в поддиагонали матрицы, образуя р матрица. Конкатенация всех вращений Гивенса образует ортогональную Q матрица.

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

Пример

Преимущества и недостатки

QR-разложение с помощью вращений Гивенса является наиболее сложным для реализации, поскольку порядок строк, необходимых для полного использования алгоритма, определить нетривиально. Однако у него есть существенное преимущество в том, что каждый новый нулевой элемент а я j < displaystyle a_ > влияет только на строку с обнуляемым элементом (i) и строку выше (j). Это делает алгоритм вращения Гивенса более эффективным и распараллеливаемым, чем метод отражения Хаусхолдера.

Связь с определителем или произведением собственных значений

Кроме того, поскольку определитель равен произведению собственных значений, мы имеем

Предположим, что QR-разложение для неквадратной матрицы А:

Из свойств СВД и определитель матрицы, имеем

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

Вращение колонны

А п = Q р ⟺ А = Q р п Т < displaystyle AP = QR quad iff quad A = QRP ^ < extf>>

Использование для решения линейных обратных задач

По сравнению с прямой обратной матрицей, обратные решения, использующие QR-разложение, более численно устойчивы, о чем свидетельствует их уменьшенная числа условий [Паркер, Геофизическая обратная теория, глава 1.13].

Обобщения

Разложение Ивасавы обобщает QR-разложение на полупростые группы Ли.

Иллюстрированный самоучитель по MatLab

LU- и QR-разложения

Так называемые LU- и QR-разложения реализуются следующими матричными функциями:

Функция lu выражает любую квадратную матрицу X как произведение двух треугольных матриц, одна из которых (возможно, с перестановками) – нижняя треугольная матрица, а другая – верхняя треугольная матрица [В MATLAB 6 аргументом (входным аргументом) функции lu может быть и полная прямоугольная матрица. – Примеч. ред.]. Иногда эту операцию называют LR-разложением. Для выполнения этой операции служит следующая функция:

Функция qr выполняет QR-разложение матрицы. Эта операция полезна для квадратных и треугольных матриц. Она выполняет QR-разложение, вычисляя произведение унитарной [Квадратная матрица с комплексными элементами, обладающая тем свойством, что обратная матрица ее комплексно сопряженной матрицы равна транспонированной, т. е. (А*)»-А’. – Примеч. ред. ] матрицы и верхней треугольной матрицы. Функция используется в следующих формах: [ Квадратная матрица с комплексными элементами, обладающая тем свойством, что обратная матрица ее комплексно сопряженной матрицы равна транспонированной, т. е. (А*)»-А’. – Примеч. ред.]

Метод Гивенса (вращений) QR-разложения квадратной матрицы (вещественный точечный вариант)

Основные авторы описания: А.В.Фролов.

Содержание

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

номера столбцов: [math]\begin \ _\quad _i\quad _ & \ & _\ \ _j\quad _\end[/math]

Именно [math]t[/math] обычно и хранится в соответствующем элементе массива.

1.2 Математическое описание алгоритма

Для завершения математического описания остаётся записать вычисление [2] матрицы вращения [math]T_[/math] и формулы её умножения на текущую промежуточную матрицу.

1.3 Вычислительное ядро алгоритма

1.4 Макроструктура алгоритма

1.5 Схема реализации последовательного алгоритма

Каждое «обнуление» элемента в позиции [math](j, i)[/math] состоит из двух шагов: а) вычисление параметров матрицы вращения [math]T_[/math] таких, чтобы обнулился элемент в позиции [math](j, i)[/math] ; б) умножение слева матрицы вращения [math]T_[/math] на текущую версию матрицы.

1.6 Последовательная сложность алгоритма

В последовательной версии основная сложность алгоритма определяется прежде всего массовыми операциями вращения. Они, если не учитывать возможную разреженность, составляют (в главном члене) [math]n^3/3[/math] операций комплексного умножения или, если не прибегать к ухищрениям, [math]4n^3/3[/math] вещественных умножений и [math]2n^3/3[/math] вещественных сложений/вычитаний.

При классификации по последовательной сложности, таким образом, метод Гивенса относится к алгоритмам с кубической сложностью.

1.7 Информационный граф

QR-алгоритм

Авторы описания: Смирнова А.С. (Раздел 2), Фролов А.В. (общая редактура и правка)

Содержание

1 Общее описание метода

Суть базового QR-алгоритма заключается в итерационном приведении матрицы [math]A[/math] к некоторой унитарно подобной ей матрице [math]A_N[/math] при помощи QR-разложения. Матрица [math]A_N[/math] является правой верхней треугольной матрицей, а значит ее диагональ содержит собственные значения. В силу подобия матриц [math]A[/math] и [math]A_N[/math] их наборы собственных значений совпадают. Таким образом задача поиска собственных значений матрицы [math]A[/math] сводится к задаче выведения матрицы [math]A_N[/math] и поиска собственных значений для нее, что является тривиальной задачей.

2 Математическое описание базового метода

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

Пусть [math]A_0 = A[/math] — исходная матрица. Для [math]k = 0, 1, 2, \ldots[/math] выполняется QR-разложение:

В специальной литературе по численным методам приводится доказательство [2] сходимости получающихся матриц к клеточной правой треугольной матрице с диагональными клетками, размеры которых зависят от размеров канонических ящиков Жордана исходной матрицы. Таким образом, проблема собственных значений матрицы сводится к проблемам собственных значений матриц меньшего размера. QR-алгоритм в узком смысле и есть эта процедура. Однако обычно к нему относят не только её, но и весь комплекс приёмов ускорения этого итерационного метода.

2.1 Особенности вещественного несимметричного варианта QR-алгоритма

3 Ускорение сходимости QR-алгоритма

3.1 Использование хессенберговой формы матрицы

Поэтому естественным приёмом ускорения итераций QR-алгоритма является начальное унитарно подобное приведение матрицы к хессенбергову виду. После этого все процедуры QR-алгоритма следует проводить уже с матрицами хессенберговой формы.

3.2 Сдвиги QR-алгоритма

При этом сохраняется свойство подобия матриц [math]A_k[/math] и [math]A_[/math] :

3.3 Неявное проведение итераций

Q-теорема [6] помогает не вычислять в явном виде требуемые матрицы с их последующим перемножением, а выполнять эти операции неявно. При этом в соответствии с особенностями вещественного случая несимметричных матриц удаётся выполнять т.н. неявный двойной сдвиг с использованием двух комплексно сопряжённых сдвигов (например, из правила Фрэнсиса), оставаясь в рамках вещественной арифметрики, с помощью выполнения двусторонними преобразованиями Хаусхолдера т.н. вытеснения выступа. В вещественном случае аналогично (но с преимуществами симметрии) выполняется неявный одинарный сдвиг, вытеснением выступа двусторонними преобразованиями вращения (Гивенса).

3.4 Обнуление поддиагональных элементов (дефляция)

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

4 QR-алгоритм в различных вариантах

В классическом (без приёмов ускорения) виде QR-алгоритм не применяется из-за медленности. Варианты метода, в основном, используются для матриц разного вида:

4.1 Метод для матриц общего вида

QR-алгоритм, используемый в SCALAPACK для решения проблемы собственных значений у матриц общего вида в настоящее время использует тот подход к выбору сдвигов, который не подтверждён теоретически, но имел лучшие практические показатели на момент создания библиотеки LAPACK. Это выполнение процедуры «двойного сдвига» с величинами, равными собственным значениям последней диагональной клетки размера 2х2 у всех диагональных блоков, больших размером. Поскольку существуют матрицы, для которых этот выбор сдвигов не даёт сходимости, то после 10й итерации при отсутствии уменьшения поддиагональных элементов выполняется некий «особый» сдвиг.

4.2 Метод для симметричных и эрмитовых матриц

Qr разложение матрицы

Совет
Никогда не применяйте для обозначения скалярного произведения символ х,который является общеупотребительным символом векторного произведения (см. разд. 9.1.7).

Рис. 9.3. Оператор векторизации
Листинг 9.16. Использование векторизациидля перемножения элементов вектора

Оператор векторизации можно использовать только с векторами и матрицами одинакового размера.
Большинство неспецифических функций MathCAD не требуют векторизации для проведения одной и той же операции над всеми элементами вектора. Например, аргументом тригонометрических функций по определениюявляется скаляр. Если попытаться вычислить синус векторной величины,MathCAD осуществит векторизацию по умолчанию, вычислив синус каждого элемента и выдав в качестве результата соответствующий вектор. Пример показан в листинге 9.17.
Листинг 9.17. Векторизация необязательнадля большинства функций MathCAD

9.1.12. Символьные операции с матрицами
Все матричные и векторные операторы, о которых шла речь выше, допустимо использовать в символьных вычислениях. Мощь символьных операцийзаключается в возможности проводить их не только над конкретными числами, но и над переменными. Несколько примеров приведены в листинге 9.18.
Листинг 9.18. Примеры символьных операций над векторами и матрицами

Cовет
Смело используйте символьный процессор в качестве мощного математического справочника. Например, когда вы хотите вспомнить какое-либо определение из области линейной алгебры (так, правила перемножения и обращенияматриц показаны в первых строках листинга 9.18).

Примечание
О вложенных массивах читайте в разд. «Создание тензора» гл. 4.

Пример использования функции CreateSpace показан на рис. 9.4. Заметьте, дляпостроения графика спирали не потребовалось никакого дополнительногокода, кроме определения параметрической зависимости в вектор-функции F!

Примечание
Размер NXM матрицы А для функции geninv должен быть таким, чтобы N>M.

Примечание
Выделить из матрицы один столбец или строку можно и с помощью функции submatrix.

Примечание
Число элементов вектора и индекс его последнего элемента совпадают, еслииндексы нумеруются с 1, т. е. системная константа ORIGIN равна 1 (см. гл. 4).

Примечание
Если элементы матриц или векторов комплексные, то сортировка ведется подействительной части, а мнимая часть игнорируется.

Совет
В большинстве задач неважно, какую норму использовать. Как видно, в обычных случаях разные нормы дают примерно одинаковые значения, хорошо отражая порядок величины матричных элементов. Определение остальных нормзаинтересованный читатель отыщет в справочниках по линейной алгебре илив Центре Ресурсов MathCAD.

Примечание
Как нетрудно понять, матрицы А и в из предыдущего листинга 9.30 обладаютодинаковыми числами обусловленности, т.к. в=100-А, и, следовательно, обематрицы определяют одну и ту же систему уравнений.

Примечание
К системам линейных уравнений сводится множество, если не сказать большинство, задач вычислительной математики. Один из таких примеров приведен в разд. «Разностные схемы»гл. 12.

Примечание
Соответствующая матрице А и вектору B система уравнений выписана явнов листинге 9.35.

Листинг 9.33. Решение СЛАУ

Листинг 9.34. Символьное решение СЛАУ (продолжение листинга 9.33)

В некоторых случаях, для большей наглядности представления СЛАУ, егоможно решить точно так же, как систему нелинейных уравнений (см. гл. 8). Пример численного решения СЛАУ из предыдущих листингов показан влистинге 9.35. Не забывайте, что при численном решении всем неизвестнымтребуется присвоить начальные значения (это сделано в первой строке листинга 9.35). Они могут быть произвольными, т. к. решение СЛАУ с невырожденной матрицей единственно.

Примечание
Фактически, треугольное разложение матрицы системы линейных уравненийпроизводится при ее решении численным методом Гаусса.

В линейная алгебра, а QR-разложение, также известный как QR-факторизация или же QU факторизация это разложение из матрица А в продукт А = QR из ортогональная матрица Q и верхнетреугольная матрица р. QR-разложение часто используется для решения линейный метод наименьших квадратов проблема и является основой для конкретной алгоритм собственных значений, то QR-алгоритм.

Содержание

Случаи и определения

Квадратная матрица

Любая настоящая квадратная матрица А может быть разложен как

куда Q является ортогональная матрица (его столбцы ортогональный единичные векторы смысл Q Т Q = Q Q Т = я < Displaystyle Q ^ < extf> Q = QQ ^ < extf> = I> ) и р это верхний треугольная матрица (также называется прямоугольной матрицей, отсюда и название). Если А является обратимый, то факторизация единственна, если нам потребуются диагональные элементы р быть позитивным.

Если А имеет п линейно независимый столбцы, затем первый п столбцы Q для мужчин ортонормированный базис для пространство столбца из А. В более общем плане первый k столбцы Q образуют ортонормированный базис для охватывать из первых k столбцы А для любого 1 ≤k ≤ п. [1] Дело в том, что любой столбец k из А зависит только от первого k столбцы Q отвечает за треугольную формур. [1]

Прямоугольная матрица

В более общем плане мы можем разложить на множители сложные м×п матрица А, с м ≥ п, как продукт м×м унитарная матрица Q и м×п верхнетреугольная матрица р. Как внизу (м−п) ряды м×п верхнетреугольная матрица полностью состоит из нулей, часто бывает полезно разбить р, или оба р и Q:

Голуб и Ван Лоан (1996), §5.2) вызов Q₁р₁ в тонкая QR-факторизация из А; Трефетен и Бау называют это уменьшенная QR-факторизация. [1] Если А полный классифицировать п и потребуем, чтобы диагональные элементы р₁ положительны тогда р₁ и Q₁ уникальны, но в целом Q₂ не является. р₁ тогда равна верхнему треугольному множителю Разложение Холецкого из А * А (= А Т А если А реально).

QL, RQ и LQ разложения

Аналогичным образом мы можем определить разложения QL, RQ и LQ с L быть ниже треугольная матрица.

Вычисление QR-разложения

Существует несколько методов реального вычисления QR-разложения, например, с помощью Процесс Грама – Шмидта, Преобразования домовладельцев, или же Гивенса вращения. У каждого есть ряд преимуществ и недостатков.

Использование процесса Грама – Шмидта

Пример

Таким образом, мы имеем

Связь с разложением RQ

Преимущества и недостатки

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

Использование отражений Хаусхолдера

Q может использоваться для отражения вектора таким образом, что все координаты, кроме одной, исчезают.

(Stoer & Bulirsch, 2002 г., п. 225) и заменить транспонирование сопряженным транспонированием при построении Q ниже.

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

Это можно повторить для А‘ (получен из Q₁А удалив первую строку и первый столбец), в результате получится матрица Хаусхолдера Q′₂. Обратите внимание, что Q′₂ меньше чем Q₁. Поскольку мы хотим, чтобы он действительно работал на Q₁А вместо А′ Нам нужно развернуть его в верхний левый угол, заполнив 1, или в общем:

— верхнетреугольная матрица. Итак, с

У этого метода больше числовая стабильность чем вышеописанный метод Грама – Шмидта.

В следующей таблице указано количество операций в k-й шаг QR-разложения с помощью преобразования Хаусхолдера, предполагающий квадратную матрицу размером п.

Пример

Теперь обратите внимание:

так что у нас уже есть почти треугольная матрица. Нам нужно только обнулить запись (3, 2).

Возьмите (1, 1) незначительный, а затем снова примените процесс к

Тем же способом, что и выше, получаем матрицу преобразования Хаусхолдера

после выполнения прямого суммирования с 1, чтобы убедиться, что следующий шаг в процессе работает правильно.

Теперь мы находим

Или, до четырех десятичных цифр,

Преимущества и недостатки

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

Использование вращений Гивенса

QR разложения также могут быть вычислены с помощью серии Гивенса вращения. Каждый поворот обнуляет элемент в поддиагонали матрицы, образуя р матрица. Конкатенация всех вращений Гивенса образует ортогональную Q матрица.

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

Пример

Преимущества и недостатки

QR-разложение с помощью вращений Гивенса является наиболее сложным для реализации, поскольку порядок строк, необходимых для полного использования алгоритма, определить нетривиально. Однако у него есть существенное преимущество в том, что каждый новый нулевой элемент а я j < displaystyle a_ > влияет только на строку с обнуляемым элементом (i) и строку выше (j). Это делает алгоритм вращения Гивенса более эффективным и распараллеливаемым, чем метод отражения Хаусхолдера.

Связь с определителем или произведением собственных значений

Кроме того, поскольку определитель равен произведению собственных значений, мы имеем

Предположим, что QR-разложение для неквадратной матрицы А:

Из свойств СВД и определитель матрицы, имеем

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

Вращение колонны

А п = Q р ⟺ А = Q р п Т < displaystyle AP = QR quad iff quad A = QRP ^ < extf>>

Использование для решения линейных обратных задач

По сравнению с прямой обратной матрицей, обратные решения, использующие QR-разложение, более численно устойчивы, о чем свидетельствует их уменьшенная числа условий [Паркер, Геофизическая обратная теория, глава 1.13].

Обобщения

Разложение Ивасавы обобщает QR-разложение на полупростые группы Ли.

Материал из Википедии — свободной энциклопедии

Определение

Алгоритмы

Пример QR-разложения

Далее, применяем алгоритм ортогонализации Грама — Шмидта и нормируем полученные вектора, получаем следующий набор:

R = ( 26 3 2 13 + 9 26 11 2 13 0 20 2 1807 + 99 3614 56 2 1807 0 0 31 139 ) <\displaystyle R=<\begin<\sqrt <26))&3<\sqrt <\frac <2><13))>+<\frac <9><\sqrt <26))>&11<\sqrt <\frac <2><13))>\\0&20<\sqrt <\frac <2><1807))>+<\frac <99><\sqrt <3614))>&56<\sqrt <\frac <2><1807))>\\0&0&<\frac <31><\sqrt <139))>\end — искомая верхнетреугольная матрица.

Примечания

Литература

Tell your friends about Wikiwand!

Gmail Facebook Twitter Link

Suggest as cover photo

Would you like to suggest this photo as the cover photo for this article?

Thank you for helping!

Your input will affect cover photo selection, along with input from other users.

Thanks for reporting this video!

This browser is not supported by Wikiwand 🙁
Wikiwand requires a browser with modern capabilities in order to provide you with the best reading experience.
Please download and use one of the following browsers:

An extension you use may be preventing Wikiwand articles from loading properly.

If you’re using HTTPS Everywhere or you’re unable to access any article on Wikiwand, please consider switching to HTTPS (https://www.wikiwand.com).

An extension you use may be preventing Wikiwand articles from loading properly.

If you are using an Ad-Blocker, it might have mistakenly blocked our content. You will need to temporarily disable your Ad-blocker to view this page.

MATLAB Language
Матричные разложения

Синтаксис

Разложение Холецкого

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

Это возвращает верхнюю треугольную матрицу. Нижняя получается транспозицией.

Наконец, мы можем проверить правильность разложения.

QR-разложение

Этот метод разложит матрицу в верхнюю треугольную и ортогональную матрицу.

Это вернет верхнюю треугольную матрицу, а следующее вернет обе матрицы.

На следующем графике будет отображаться время выполнения qr зависящее от квадратного корня элементов матрицы.

Разложение LU

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

Однако довольно часто этот метод не работает или плохо работает, поскольку он нестабилен. Например

Достаточно добавить матрицу перестановок, для которой PA = LU:

Далее мы построим время выполнения `lu ‘, зависящее от квадратного корня элементов матрицы.

Разложение Шура

Мы также показываем время выполнения schur зависящее от квадратного корня из матричных элементов:

Разложение сингулярного значения

Учитывая m раз n матрицу A с n больше m. Разложение сингулярного значения

вычисляет матрицы U, S, V.

Матрица U состоит из левых сингулярных собственных векторов, являющихся собственными векторами A*A.’ в то время как V состоит из правых сингулярных собственных значений, являющихся собственными векторами оператора A.’*A Матрица S имеет квадратные корни из собственных значений A*A.’ и A.’*A на его диагонали.

Если m больше n, можно использовать

для выполнения разложения особых величин с экономичным размером.

СОДЕРЖАНИЕ

Случаи и определения

Квадратная матрица

Любая вещественная квадратная матрица A может быть разложена как

Прямоугольная матрица

QL, RQ и LQ разложения

Аналогичным образом мы можем определить разложения QL, RQ и LQ, где L является нижней треугольной матрицей.

Вычисление QR-разложения

Использование процесса Грама – Шмидта

Пример

Таким образом, мы имеем

Связь с разложением RQ

Преимущества и недостатки

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

Использование отражений Хаусхолдера

Q может использоваться для отражения вектора таким образом, что все координаты, кроме одной, исчезают.

— верхнетреугольная матрица. Итак, с

Этот метод имеет большую числовую устойчивость, чем метод Грама – Шмидта, описанный выше.

Пример

Теперь обратите внимание:

так что у нас уже есть почти треугольная матрица. Нам нужно только обнулить запись (3, 2).

Тем же способом, что и выше, получаем матрицу преобразования Хаусхолдера

после выполнения прямого суммирования с 1, чтобы убедиться, что следующий шаг в процессе работает правильно.

Теперь мы находим

Или, до четырех десятичных цифр,

Матрица Q ортогональна, а R верхнетреугольная, поэтому A = QR является требуемым QR-разложением.

Преимущества и недостатки

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

Использование вращений Гивенса

Пример

Преимущества и недостатки

QR-разложение с помощью вращений Гивенса является наиболее сложным для реализации, поскольку порядок строк, необходимых для полного использования алгоритма, определить нетривиально. Однако он имеет значительное преимущество в том, что каждый новый нулевой элемент влияет только на строку с обнуляемым элементом ( i ) и строку выше ( j ). Это делает алгоритм вращения Гивенса более эффективным и распараллеливаемым, чем метод отражения Хаусхолдера. а я j <\ displaystyle a_ >

Связь с определителем или произведением собственных значений

Кроме того, поскольку определитель равен произведению собственных значений, мы имеем

Предположим, что QR-разложение для неквадратной матрицы A :

Из свойств SVD и определителя матрицы имеем

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

Вращение колонны

А п знак равно Q р ⟺ А знак равно Q р п Т <\ displaystyle AP = QR \ quad \ iff \ quad A = QRP ^ <\textf>>

Использование для решения линейных обратных задач

По сравнению с прямой обратной матрицей, обратные решения, использующие QR-разложение, более численно стабильны, о чем свидетельствуют их уменьшенные числа обусловленности [Parker, Geophysical Inverse Theory, Ch1.13].

Обобщения

Разложение Ивасавы обобщает QR-разложение на полупростые группы Ли.

Как решить линейную регрессию с помощью линейной алгебры

Дата публикации 2018-03-05

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

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

После завершения этого урока вы узнаете:

Обзор учебника

Этот урок разделен на 6 частей; они есть:

Линейная регрессия

Модель предполагает, что y является линейной функцией или взвешенной суммой входной переменной.

Или указано с коэффициентами.

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

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

Матричная формулировка линейной регрессии

Линейная регрессия может быть задана с использованием матричной записи; например:

Или без точечной записи.

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

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

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

Эта формулировка имеет уникальное решение, если входные столбцы независимы (например, некоррелированы).

В матричной записи эта проблема формулируется с использованием так называемого нормального уравнения:

Это может быть переупорядочено, чтобы указать решение для b как:

Это может быть решено напрямую, хотя наличие обратной матрицы может быть численно сложным или нестабильным.

Набор данных линейной регрессии

Чтобы изучить матричную формулировку линейной регрессии, давайте сначала определим набор данных как контекст.

Мы будем использовать простой 2D-набор данных, в котором данные легко визуализировать в виде точечной диаграммы, а модели легко представить в виде линии, которая пытается соответствовать точкам данных.

Приведенный ниже пример определяет матричный набор данных 5 × 2, разбивает его на компоненты X и y и строит набор данных как график рассеяния.

При запуске примера сначала печатается определенный набор данных.

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

Решить напрямую

То есть, учитывая X, каков набор коэффициентов b, который при умножении на X даст y. Как мы видели в предыдущем разделе, нормальные уравнения определяют, как вычислять b напрямую.

Это можно вычислить непосредственно в NumPy, используя функцию inv () для вычисления обратной матрицы.

После того как коэффициенты рассчитаны, мы можем использовать их для прогнозирования результатов с учетом X.

Объединяя это с набором данных, определенным в предыдущем разделе, полный пример приведен ниже.

При выполнении примера выполняется вычисление и выводится вектор коэффициента b.

Затем создается точечная диаграмма набора данных с линейным графиком для модели, показывающим разумное соответствие данным

Решить с помощью QR-разложения

QR-разложение является популярным подходом для решения линейного уравнения наименьших квадратов.

Переходя ко всему выводу, коэффициенты могут быть найдены с использованием элементов Q и R следующим образом:

Подход все еще включает матричную инверсию, но в этом случае только на более простой R-матрице.

QR-разложение можно найти с помощью функции qr () в NumPy. Расчет коэффициентов в NumPy выглядит следующим образом:

Связав это с набором данных, полный пример приведен ниже.

Выполнение примера сначала печатает решение коэффициента и наносит на график данные с моделью.

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

Решить через разложение по сингулярному значению

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

После разложения коэффициенты могут быть найдены путем вычисления псевдообращения входной матрицы X и умножения ее на выходной вектор y.

Где псевдообратный рассчитывается следующим образом:

Инверсия матриц не определена для матриц, которые не являются квадратными. […] Когда A имеет больше столбцов, чем строк, то решение линейного уравнения с использованием псевдообратного представления дает одно из многих возможных решений.

Мы можем получить U и V из операции SVD. D ^ + можно рассчитать, создав диагональную матрицу из Sigma и вычислив обратную величину каждого ненулевого элемента в Sigma.

Мы можем вычислить SVD, затем псевдообратную вручную. Вместо этого NumPy предоставляет функцию pinv (), которую мы можем использовать напрямую.

Полный пример приведен ниже.

При выполнении примера печатается коэффициент и данные отображаются красной линией, показывающей прогнозы из модели.

Фактически, NumPy предоставляет функцию для замены этих двух шагов в функции lstsq (), которую вы можете использовать напрямую.

расширения

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

Если вы исследуете какое-либо из этих расширений, я хотел бы знать.

Дальнейшее чтение

Этот раздел предоставляет больше ресурсов по теме, если вы хотите углубиться.

книги

статьи

Учебники

Резюме

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

В частности, вы узнали:

У вас есть вопросы?
Задайте свои вопросы в комментариях ниже, и я сделаю все возможное, чтобы ответить.

Предисловие

2. Вычисления с плавающей точкой

4. Базовые процедуры линейной алгебры

6. Процедуры IMSL для систем линейных уравнений и вычисления определителей

8. Собственные значения и собственные векторы

10. Процедуры и операции библиотеки IMSL Fortran 90 MP

Приложение 2. Список представленных в пособии процедур математической библиотеки IMSL

Предметный указатель

1.1. Состав библиотеки

1.2. Вызов процедур IMSL

1.3. Выделение памяти

1.4. Соглашения об именах

1.5. Вычислительная сложность алгоритмов

1.6. Оптимизация кода

1.7. Учет особенностей машинной арифметики

1.8. Обработка ошибок

1.9. Соглашения, действующие при описании процедур библиотеки IMSL

1.10. Вспомогательная функция ru_doswin

2. Вычисления с плавающей точкой

2.1. Двоичные представления чисел

2.2. Исключения

2.3. Управление вычислениями с плавающей точкой

2.4. Обработка исключений

3. Вспомогательные процедуры

3.2. Машинные константы

3.3. Вывод результатов

3.4. Обработка ошибок

3.5. Процедуры даты и времени

3.6. Задание и чтение настроек процедур IMSL

3.7. Некоторые специальные процедуры

3.8. Версия библиотеки IMSL

4. Базовые процедуры линейной алгебры

4.1. Введение

4.2. Некоторые сведения о векторах и матрицах

4.3. Представление матриц в программах

4.4. Базовые процедуры линейной алгебры уровня 1

4.5. Базовые процедуры линейной алгебры уровней 2 и 3

4.6. Дополнительные операции над векторами и матрицами

5. Решение систем линейных уравнений

5.1. Постановка задачи

5.2. Метод исключения Гаусса

5.3. Линейные системы с симметрической матрицей

6. Процедуры IMSL для систем линейных уравнений и вычисления определителей

6.1. Введение

6.2. Системы линейных уравнений с вещественной несимметрической матрицей

6.3. Процедуры с вещественной симметрической матрицей

6.4. Системы линейных уравнений с вещественной матрицей Теплица

6.5. Системы линейных уравнений с комплексной матрицей

6.6. Вычисление определителей

7. Матричные разложения и их пересчет

7.1. Введение

7.2. Преобразования Хаусхолдера и Гивенса

7.3. Методы вычисления QR-разложения

7.4. QR-разложение процедурами IMSL

7.5. Пересчет матричных разложений

7.6. Вычисление сингулярного разложения

8. Собственные значения и собственные векторы

8.1. Подходы к решению проблемы собственных значений

8.2. Процедуры IMSL для вычисления собственных значений

8.3. Оценка точности решения проблемы собственных значений

9. Решение прямоугольных систем Ax ≈ b с ограничениями и без них

9.1. Метод наименьших квадратов

9.2. Подпрограммы IMSL для решения переопределенной и недоопределенной проблемы Ax ≈ b

9.3. Проблема Ax ≈ b с линейными ограничениями

10. Процедуры и операции библиотеки IMSL Fortran 90 MP

10.1. Введение

10.2. Операции библиотеки IMSL 90

10.3. Функции библиотеки IMSL 90

10.4. Решение систем линейных уравнений

10.5. Сингулярные и собственные значения

10.6. Решение линейных систем методом наименьших квадратов

Приложение 1. Вывод русского текста в DOS-окно

Приложение 2. Список представленных в пособии процедур математической библиотеки IMSL

1. Характеристики процедур IMSL

1.1. Состав библиотеки

1.2. Вызов процедур IMSL

1.3. Выделение памяти

1.4. Соглашения об именах

1.5. Вычислительная сложность алгоритмов

1.6. Оптимизация кода

1.7. Учет особенностей машинной арифметики

1.7.1. Разрывы между вещественными числами

1.7.2. Ошибки округления

1.7.3. Устранение переполнения и исчезновения порядка

1.7.4. Оценка точности результата

1.8. Обработка ошибок

1.9. Соглашения, действующие при описании процедур библиотеки IMSL

1.10. Вспомогательная функция ru_doswin

2. Вычисления с плавающей точкой

2.1. Двоичные представления чисел

2.1.1. Представление целых чисел

2.1.2. Вещественные числа с плавающей точкой

2.1.3. Виды вещественных чисел с плавающей точкой

2.1.4. характеристики модели представления чисел

2.1.5. Абсолютная и относительная ошибки

2.2. Исключения

2.3. Управление вычислениями с плавающей точкой

2.3.1. Значения статуса и управляющего слова

2.3.2. Чтение статуса

2.3.3. Изменение управляющего слова

2.4. Обработка исключений

2.4.1. Обработка исключений с плавающей точкой

2.4.2. Обработка математических исключений

3. Вспомогательные процедуры

3.2. Машинные константы

3.2.1. Целые машинные константы

3.2.2. Вещественные машинные константы

3.2.3. Проверка на NaN

3.3. Вывод результатов

3.3.1. Список, вызовы и параметры процедур

3.3.2. Вывод прямоугольной матрицы с нумерацией строк и столбцов

3.3.3. Вывод прямоугольной матрицы по заданному формату и с заданными именами строк и столбцов

3.3.3.1. Вывод вещественной матрицы

3.3.3.2. Вывод комплексной матрицы

3.3.3.3. Вывод целочисленной матрицы

3.3.4. Управление выводом

3.3.5. Управление размером страницы

3.3.6. Задание номера устройства В/В

3.4. Обработка ошибок

3.4.1. Виды ошибок IMSL

3.4.2. Примеры ошибок IMSL

3.4.3. Изменение характера реагирования на информационные ошибки

3.4.4. Использование системы обработки ошибок

3.5. Процедуры даты и времени

3.6. Задание и чтение настроек процедур IMSL

3.7. Некоторые специальные процедуры

3.7.1. Разложение числа на простые множители

3.7.2. Дополнительные символьные функции

3.7.3. Вычисление SQRT(a**2 + b**2) без исчезновения порядка и переполнения

3.7.4. Математические и физические константы

3.7.5. Перевод величин из одной системы единиц в другую

3.8. Версия библиотеки IMSL

4. Базовые процедуры линейной алгебры

4.1. Введение

4.2. Некоторые сведения о векторах и матрицах

4.2.1. Обозначение векторов и матриц

4.2.2. Операции над векторами

4.2.3. Операции над матрицами

4.2.4. Элементарные преобразования

4.2.5. Определитель и миноры матрицы

4.2.6. Виды матриц

4.2.7. Умножение матриц

4.2.8. Умножение матрицы на вектор

4.2.9. Блочные матрицы

4.2.10. Двумерные матрицы вращения и отражения

4.2.11. Линейная независимость и базис

4.2.12. Ортогональность векторов

4.2.13. Ранг матрицы

4.2.14. Векторные нормы

4.2.15. Матричные нормы

4.2.16. Нормы и ортогональные преобразования

4.3. Представление матриц в программах

4.3.1. Представление всех элементов матрицы

4.3.2. Ленточные матрицы

4.3.3. Ленточные симметрические матрицы

4.3.4. Ленточные эрмитовы матрицы

4.3.5. Ленточные треугольные матрицы

4.3.6. Кодиагональное представление ленточных симметрических матриц

4.3.7. Кодиагональное представление ленточных эрмитовых матриц

4.3.8. Разреженные матрицы

4.4. Базовые процедуры линейной алгебры уровня 1

4.4.1. Замечания для программиста

4.4.2. Список процедур уровня 1

4.4.3. Присваивание вектору скаляра

4.4.4. Копирование вектора

4.4.5. Масштабирование вектора

4.4.6. Умножение вектора на скаляр

4.4.7. Сумма вектора и скаляра

4.4.8. Вычитание из скаляра вектора

4.4.9. Присваивание вида y = αx + y

4.4.10. Обмен данными векторов

4.4.11. Скалярное произведение векторов

4.4.12. Скалярное произведение векторов с аккумулятором двойной точности

4.4.13. Сумма скаляра и скалярного произведения векторов с аккумулятором двойной точности

4.4.14. Скалярное произведение векторов с аккумулятором двойной точности, предаваемым в качестве параметра

4.4.15. Произведение Адамара

4.4.16. Сумма произведений трех векторов

4.4.17. Сумма элементов вектора

4.4.18. Сумма абсолютных значений элементов вектора

4.4.19. Евклидова или 2-норма вектора

4.4.20. Произведение элементов вектора

4.4.21. Индекс минимального элемента вектора

4.4.22. Индекс максимального элемента вектора

4.4.23. Индекс минимального по модулю элемента вектора

4.4.24. Индекс максимального по модулю элемента вектора

4.4.25. Построение плоского вращения Гивенса

4.4.26. Применение плоского вращения Гивенса

4.4.27. Построение быстрого вращения Гивенса

4.4.28. Применение быстрого вращения Гивенса

4.5. Базовые процедуры линейной алгебры уровней 2 и 3

4.5.1. Замечания для программиста

4.5.2. Список процедур уровней 2 и 3

4.5.3. Умножение матрицы общего вида на вектор

4.5.4. Умножение ленточной матрицы на вектор

4.5.5. Умножение эрмитовой матрицы на вектор

4.5.6. Умножение эрмитовой ленточной матрицы на вектор

4.5.7. Умножение вещественной симметрической матрицы на вектор

4.5.8. Умножение симметрической ленточной матрицы на вектор

4.5.9. Умножение треугольной матрицы на вектор

4.5.10. Умножение треугольной ленточной матрицы на вектор

4.5.11. Умножение обратной треугольной матрицы на вектор

4.5.12. Умножение обратной треугольной ленточной матрицы на вектор

4.5.13. Преобразования ранга 1 матрицы общего вида

4.5.14. Преобразования ранга 1, 2, k и 2k симметрической матрицы

4.5.15. Преобразования ранга 1, 2, k и 2k эрмитовой матрицы

4.5.16. Произведение матриц общего вида

4.5.17. Произведение с симметрической матрицей

4.5.18. Произведение с эрмитовой матрицей

4.5.19. Произведение с треугольной матрицей

4.5.20. Вычисления с треугольной матрицей

4.6. Дополнительные операции над векторами и матрицами

4.6.1. Процедуры, которые могут быть заменены встроенными возможностями Фортрана

4.6.2. Преобразование матриц

4.6.3. Билинейная и квадратичная формы

4.6.4. Оценка матричного степенного ряда

4.6.5. Умножение ленточной матрицы на вектор

4.6.6. Сложение ленточных матриц

4.6.7. Вычисление нормы матрицы

4.6.8. Вычисление расстояния между точками

4.6.9. Вычисления с повышенной точностью

5. Решение систем линейных уравнений

5.1. Постановка задачи

5.2. Метод исключения Гаусса

5.2.1. Идея метода исключения

5.2.2. LU-разложение матрицы

5.2.2.1. Вычисление LU-разложения матрицы

5.2.2.2. Применение LU-разложения для решения линейных систем

5.2.3. Повышение качества разложения

5.2.3.1. Схема частичного выбора

5.2.3.2. Полный выбор ведущего элемента

5.2.4. Устойчивость метода исключений Гаусса

5.2.5. Некоторые приложения метода исключений Гаусса

5.2.5.1. Вычисление обратной матрицы

5.2.5.2. Вычисление определителя

5.2.6. Число обусловленности матрицы

5.2.6.1. Понятие числа обусловленности

5.2.6.2. Оценка числа обусловленности матрицы

5.2.7. Оценка точности и уточнение решения

5.2.7.1. Точность решения

5.2.7.2. Итерационное уточнение решения

5.2.8. Системы с разреженными матрицами

5.3. Линейные системы с симметрической матрицей

5.3.1. Положительно определенные системы

5.3.2. Получение разложения Холецкого

5.3.3. Неотрицательно определенные системы

5.3.3.1. Свойства неотрицательно определенных матриц

5.3.3.2. LLT-разложение симметрических неотрицательно определенных матриц

5.3.3.3. Симметричный выбор диагонального элемента

5.3.4. Симметрические неопределенные системы

6. Процедуры IMSL для систем линейных уравнений и вычисления определителей

6.1. Введение

6.2. Системы линейных уравнений с вещественной несимметрической матрицей

6.2.1. Список процедур с вещественной несимметрической матрицей

6.2.2. Возможные ошибки процедур

6.2.3. Параметры процедур с вещественными несимметрическими и комплексными неэрмитовыми матрицами

6.2.4. Процедуры с вещественной несимметрической общего вида матрицей

6.2.4.1. Задание вещественной несимметрической общего вида матрицы

6.2.4.2. Вычисление LU-разложения и числа обусловленности матрицы

6.2.4.3. Вычисление LU-разложения матрицы

6.2.4.4. Решение системы линейных уравнений с применением ранее найденного LU-разложения

6.2.4.5. Итерационное уточнение решения линейной системы

6.2.4.6. Решение системы линейных уравнений без уточнения корней

6.2.4.7. Решение системы линейных уравнений с итерационным уточнением корней

6.2.4.8. Вычисление обратной матрицы

6.2.5. Пример описания процедуры с вещественной несимметрической ленточной матрицей

6.2.5.1. Задание ленточной матрицы

6.2.5.2. Вычисление LU-разложения и числа обусловленности матрицы

6.2.6. Процедуры с вещественной несимметрической трехдиагональной матрицей

6.2.6.1. Вычисление LDU-разложения и решение линейной системы Ax = b с трехдиагональной матрицей

6.2.6.2. Решение трехдиагональной системы линейных уравнений

6.2.7. Процедуры с вещественной треугольной матрицей

6.2.7.1. Задание треугольной матрицы

6.2.7.2. Оценка числа обусловленности треугольной матрицы

6.2.7.3. Вычисление обратной треугольной матрицы

6.2.7.4. Решение треугольной системы линейных уравнений

6.2.7.5. Решение треугольной (возможно вырожденной) системы линейных уравнений и вычисление обобщенного обращения верхней треугольной матрицы

6.2.8. Процедуры с вещественной несимметрической разреженной матрицей

6.2.8.1. Вычисление LU-разложения разреженной матрицы

6.2.8.2. Решение системы линейных уравнений с применением ранее найденного LU-разложения

6.2.8.3. Решение системы линейных уравнений

6.3. Процедуры с вещественной симметрической матрицей

6.3.1. Список процедур

6.3.2. Возможные ошибки процедур с вещественной симметрической матрицей

6.3.3. Параметры процедур с вещественной симметрической и комплексной эрмитовой матрицами

6.3.4. Процедуры с вещественной симметрической положительно определенной матрицей

6.3.4.1. Об использовании процедур

6.3.4.2. Примеры для процедур с симметрической положительно определенной матрицей

6.3.5. Процедуры с вещественной симметрической положительно определенной ленточной матрицей

6.3.5.1. Об использовании процедур

6.3.5.2. RTDR-разложение Холецкого симметрической положительно определенной ленточной матрицы и решение системы Ax = b

6.3.5.3. Примеры для процедур с симметрической положительно определенной ленточной матрицей

6.3.6. Процедуры с вещественной симметрической положительно определенной разреженной матрицей

6.3.6.1. Метод сопряженных градиентов с аппроксимацией и обратной связью

6.3.6.2. Метод сопряженных градиентов с аппроксимацией Якоби и обратной связью

6.3.6.3. Вычисление LLT-разложения Холецкого симметрической положительно определенной разреженной матрицы

6.3.6.4. Решение системы линейных уравнений с симметрической положительно определенной разреженной матрицей с использованием ранее найденного LLT-разложения Холецкого

6.3.6.5. Решение системы линейных уравнений с симметрической положительно определенной разреженной матрицей методом Холецкого

6.3.7. Разложение Холецкого симметрической неотрицательно определенной матрицы

6.3.8. Разложение Холецкого симметрической неотрицательно определенной матрицы с применением симметричного выбора

6.3.9. Процедуры с вещественной симметрической неопределенной матрицей

6.3.9.1. Об использовании процедур с симметрической неопределенной матрицей

6.3.9.2. Примеры для процедур с симметрической неопределенной матрицей

6.3.10. Разложение симметрической матрицы с добавлением к ней диагональной матрицы

6.4. Системы линейных уравнений с вещественной матрицей Теплица

6.5. Системы линейных уравнений с комплексной матрицей

6.5.1. Процедуры с комплексной неэрмитовой матрицей

6.5.1.1. Список процедур

6.5.1.2. Системы линейных уравнений с циклической матрицей

6.5.1.3. Примеры процедур с комплексной неэрмитовой матрицей

6.5.2. Процедуры с комплексной эрмитовой матрицей

6.5.2.1. Список процедур

6.5.2.2. Возможные ошибки процедур

6.5.2.3. Примеры для процедур с комплексной эрмитовой матрицей

6.6. Вычисление определителей

6.6.1. Список процедур

6.6.2. Параметры процедур

6.6.3. Описание подпрограммы LFDRG

6.6.4. Примеры вычисления определителей

7. Матричные разложения и их пересчет

7.1. Введение

7.2. Преобразования Хаусхолдера и Гивенса

7.2.1. Преобразование Хаусхолдера

7.2.1.1. Матрица и вектор Хаусхолдера

7.2.1.2. Вычисление вектора Хаусхолдера

7.2.1.3. Умножение матрицы на матрицу Хаусхолдера

7.2.1.4. Ошибки округления

7.2.2. Преобразование Гивенса

7.2.2.1. Матрица вращения Гивенса

7.2.2.2. Умножение матрицы на матрицу Гивенса

7.2.2.3. Компактное представление матриц Гивенса

7.3. Методы вычисления QR-разложения

7.3.1. QR-разложение прямоугольной матрицы с применением преобразований Хаусхолдера

7.3.2. QR-разложение симметрической матрицы

7.3.2.1. Метод Хаусхолдера приведения симметрической матрицы к трехдиагональному виду

7.3.2.2. QR-разложение симметрической трехдиагональной матрицы

7.4. QR-разложение процедурами IMSL

7.4.1. QR-разложение прямоугольной матрицы

7.4.2. Восстановление ортогональной матрицы Q из QR-разложения прямоугольной матрицы

7.5. Пересчет матричных разложений

7.5.1. Алгоритм пересчета QR-разложения

7.5.2. Перечень процедур

7.5.3. Пересчет разложения Холецкого

7.5.3.2. Пересчет разложения матрицы A + xx T

7.5.4. Пересчет QR-разложения

7.6. Вычисление сингулярного разложения

7.6.1. Понятие сингулярного разложения

7.6.2. Сингулярное разложение комплексной матрицы

7.6.3. Сингулярное разложение вещественной матрицы

7.6.4. Вычисление обобщенной обратной матрицы

8. Собственные значения и собственные векторы

8.1. Подходы к решению проблемы собственных значений

8.1.1. Понятия собственного значения и собственного вектора

8.1.2. Преобразования подобия

8.1.3. Некоторые свойства собственных значений и собственных векторов

8.1.4. Собственные значения и собственные векторы вещественной симметрической и комплексной эрмитовой матрицы

8.1.4.2. Метод обратных итераций

8.2. Процедуры IMSL для вычисления собственных значений

8.2.1. Организация материала

8.2.2. Обычная проблема собственных значений

8.2.2.1. Список, вызовы, параметры и ошибки процедур

8.2.2.3. Вычисление всех собственных значений и собственных векторов вещественной симметрической матрицы

8.2.2.4. Вычисление нескольких наибольших или наименьших собственных значений и им соответствующих собственных векторов вещественной симметрической матрицы

8.2.2.5. Вычисление собственных значений и собственных векторов вещественной симметрической матрицы из заданного интервала

8.2.2.6. Вычисление всех собственных значений и собственных векторов вещественной несимметрической матрицы

8.2.2.7. Вычисление всех собственных значений и собственных векторов комплексной эрмитовой матрицы

8.2.2.8. Вычисление всех собственных значений и собственных векторов комплексной неэрмитовой матрицы

8.2.4.2. Вычисление всех собственных значений и собственных векторов комплексной матрицы Хессенберга

8.2.2.9. Примеры применения подпрограмм, вычисляющих собственные значения и собственные векторы

8.2.3. Обобщенная проблема собственных значений

8.2.3.1. Список, вызовы, параметры и ошибки процедур

8.2.3.2. Вычисление всех собственных значений и собственных векторов вещественной симметрической обобщенной системы Az = λBz

8.2.3.3. Вычисление всех собственных значений и собственных векторов вещественной обобщенной системы Az = λBz

8.3. Оценка точности решения проблемы собственных значений

8.3.1. Список, вызовы, параметры и ошибки процедур

9. Решение прямоугольных систем Ax ≈ b с ограничениями и без них

9.1. Метод наименьших квадратов

9.1.1. Постановка задачи

9.1.2. Применение QR-разложения в методе наименьших квадратов

9.1.3. Вырожденная задача наименьших квадратов

9.1.4. Решение задачи наименьших квадратов с помощью сингулярного разложения

9.2. Подпрограммы IMSL для решения переопределенной и недоопределенной проблемы Ax ≈ b

9.2.1. Перечень подпрограмм

9.2.2. Проблема Ax ≈ b без ограничений

9.2.2.1. Решение проблемы Ax ≈ b методом наименьших квадратов с итерационным уточнением корней

9.2.2.2. Решение проблемы Ax ≈ b методом наименьших квадратов без итерационного уточнения корней

9.2.2.3. Решение проблемы Ax ≈ b методом наименьших квадратов с использованием блочных преобразований Хаусхолдера

9.2.2.4. вычисление преобразований координат и завершение решения проблемы Ax ≈ b методом наименьших квадратов

9.3. Проблема Ax ≈ b с линейными ограничениями

10. Процедуры и операции библиотеки IMSL Fortran 90 MP

10.1. Введение

10.2. Операции библиотеки IMSL 90

10.3. Функции библиотеки IMSL 90

10.4. Решение систем линейных уравнений

10.4.1. Процедура LIN_SOL_GEN для линейных систем с матрицей общего вида

10.4.2. Процедура LIN_SOL_SELF для линейных систем с самоприсоединенной матрицей

10.4.3. Процедура LIN_SOL_TRI для линейных систем с трехдиагональной матрицей

10.5. Сингулярные и собственные значения

10.5.1. Процедура LIN_SVD для вычисления сингулярного разложения

10.5.2. Процедура LIN_EIG_SELF для вычисления собственных значений самоприсоединенной матрицы

10.5.3. Процедура LIN_EIG_GEN для вычисления собственных значений матрицы общего вида

10.5.4. Процедура LIN_GEIG_GEN для обобщенной проблемы собственных значений

10.6. Решение линейных систем методом наименьших квадратов

10.6.1. Процедура LIN_SOL_LSQ для решения прямоугольных линейных систем методом наименьших квадратов

10.6.2. Процедура LIN_SOL_SVD, использующая сингулярное разложение для решения прямоугольных линейных систем

Приложение 1. Вывод русского текста в DOS-окно

Приложение 2. Список представленных в пособии процедур математической библиотеки IMSL