Скалярное произведение в вещественном пространстве $ (m \times n) $-матриц будем считать заданным равенством $$ \langle A, B \rangle := \sum_{j,k} a_{jk} b_{jk} \ \mbox{при} \ \{A,B\} \subset \mathbb R^{m\times n}. $$ как если бы эти матрицы были предварительно векторизованы, и затем перемножены скалярно как векторы. Эквивалентное определение — через след матрицы: $$ \langle A, B \rangle = \operatorname{Sp} (A \cdot B^{\top}) = \operatorname{Sp} (A^{\top} \cdot B) \, . $$ Как следствие такого способа определения евклидова пространства, вводится понятие евклидовой нормы (или нормы Фробениуса) $$ \| B \|:=\sqrt{\sum_{j,k} b_{jk}^2} = \sqrt{\operatorname{Sp}(B^{\top} B)} = \sqrt{\operatorname{Sp}(B\cdot B^{\top})} \, . $$ Естественно, что расстояние между матрицами $ A $ и $ B $ это обобщение расстояния между точками в стандартном его определении: $$ \|A-B \|:= \sqrt{\sum_{j,k} (a_{jk}-b_{jk})^2} = \sqrt{\operatorname{Sp}((A-B)^{\top}(A- B))} \, . $$
Задача. Для заданной матрицы $ A \in \mathbb R^{m\times n} $ и натурального числа $ \mathfrak r < \operatorname{rank} (A) $ построить матрицу $ \widetilde A \in \mathbb R^{m\times n} $ ранга $ \mathfrak r $ наилучшим образом приближающую матрицу $ A $ в евклидовой норме (норме Фробениуса): $$ \min_{\widetilde A\in \mathbb R^{m\times n}, \operatorname{rank} (\widetilde A)= \mathfrak r} || A- \widetilde A|| \, . $$ Иными словами, среди всех матриц ранга $ \mathfrak r $ найти матрицу, ближайшую к матрице $ A $ (иными словами: найти ортогональную проекцию матрицы $ A $ на множество всех матриц ранга $ \mathfrak r $).
В чем важность этой задачи? — Дело в том, что если произвольная матрица $ A \in \mathbb R^{m\times n} $ задается $ m\cdot n $ элементами, то одноранговая матрица того же порядка задается $ (m+n-1) $ элементами, двухранговая — $ (m+n-1)+ (m+n-3) $ элементами, и т.д. Чем меньше $ \mathfrak r $ — тем меньше требуется ресурсов для хранения.
В чем сложности? —
1. При фиксированном $ \mathfrak r $ найти матрицу $ \widetilde A $.
2. Определить оптимальную величину $ \mathfrak r $, т.е., с одной стороны, чтобы эта величина была поменьше, но, с другой стороны, чтобы соответствующая матрица $ \widetilde A $ давала хорошее приближение матрицы $ A $. Например, чтобы по построенным приближениям $ \widetilde A $ и $ \widetilde B $ матриц $ A $ и $ B $ можно было бы с большой вероятностью правильно идентифицировать эти матрицы.
Представление матрицы $ A $ в виде суммы линейно независимых одноранговых матриц. Такое представление всегда существует в количестве слагаемых равном $ \operatorname{rank} (A) $. К примеру, см. ЗДЕСЬ. И подобных представлений можно привести бесконечно много. Наша цель: найти такое представление, где слагаемые будут ранжированы по «объему вклада» в матрицу $ A $: первая матрица ранга $1$ должна быть максимально близка к матрице $ A $, вторая матрица — максимально близка к разности $ A $ и первой матрицы и т.д.
Требуемое представление получается посредством сингулярного разложения матрицы $ A $. Формальный подход к изложению теории см. ☞ ЗДЕСЬ. Это разложение используется не только для решения поставленной в настоящем разделе задачи.
Если бегло посмотреть на основной результат теории, то увидеть представление матрицы $ A $ в виде суммы каких-то матриц затруднительно: там стоит чистое произведение, т.е. разложение матрицы $ A $ на три сомножителя, из которых один (средний) — простой структуры (диагональная матрица). На самом деле, можно показать (см. ☟ НИЖЕ )эквивалентность этого произведения именно сумме одноранговых матриц.
Обратим внимание на похожесть упоминаемого результата на уже пройденный в курсе алгебры результат о представлении квадратной матрицы $ A $ также в виде произведения трех матриц $$ A=C^{-1} A_{diag} C \, , $$ где матрица $ A_{diag} $ — диагональная, с собственными числами матрицы $ A $ на главной диагонали. Не акцентируем пока внимания на то, что такое представление не для всех матриц возможно; считаем, что случайным образом выбранная матрица $ A $ диагонализуема. Нас, в первую очередь, интересует другое: является ли этот результат частью общего результата о сингулярном разложении?
Для некоторого класса матриц ответ оказывается положительным. Вот именно на основе этих матриц мы и сделаем сначала «заготовку» для общей теории.
Собственные числа, собственные векторы матрицы. Собственные векторы бывают правые и левые.
Вычисление функции от матрицы через интерполяционную формулу Лагранжа.
Симметричные матрицы, свойства их собственных чисел и векторов (в том числе экстремальное). Матрица Грама и ее свойства (она положительно полуопределённа).
Частичная проблема собственных чисел: произвольная матрица и симметричная матрица.
Как произведении $ A=C^{-1} A_{diag} C $ «перегнать» в сумму одноранговых матриц?
Теорема 1. Если все собственные числа $ \lambda_{1},\dots,\lambda_n $ матрицы $ A \in \mathbb C^{n\times n} $ различны, то для произвольного полинома $ g_{}(x) $ справедлива формула
$$ g(A)= \sum_{k=1}^n \frac{g(\lambda_k)}{f_k(\lambda_k)}f_k(A) = \sum_{k=1}^n \frac{g(\lambda_k)}{f'(\lambda_k)}f_k(A) $$ где $ f_k(x) $ означает частное от деления характеристического полинома $ f_{}(x):=\det (A-xE_{n\times n}) $ на $ x- \lambda_k $.
Поскольку этот результат справедлив для произвольного полинома $ g(x) $, то он выполняется и для частного случая $ g(x) \equiv x $: $$ A = \sum_{k=1}^n \lambda_k \frac{f_k(A)}{f_k(\lambda_k)} = \sum_{k=1}^n \lambda_k \frac{f_k(A)}{f^{\prime}(\lambda_k)} \, . $$ То есть мы получили представление матрицы $ A $ в виде линейной комбинации матриц $$ \left\{ M_k(A):=\frac{f_k(A)}{f_k(\lambda_k)} \right\}_{k=1}^n \, . $$ Теперь исследуем свойства этих матриц — они называются ковариантами Фробениуса матрицы $ A $.
Теорема 2. Если алгебраическая кратность собственного числа $ \lambda_k $ равна $ 1 $, то
$$ \operatorname{rank} (M_k(A)) = 1 \, . $$
Любая одноранговая $m\times n $-матрица представима в виде произведения столбца $ U \in \mathbb C^{m\times 1} $ на строчку $ V\in \mathbb C^{1\times n} $: $$ \left(\begin{array}{c} u_1\\ u_2\\ \vdots \\ u_m \end{array} \right) (v_1,\dots,v_n)=\left[u_jv_k\right]_{j=1,\dots,m;k=1,\dots,n} $$ при хотя бы одной паре индексов $ j_{} $ и $ k_{} $ таких, что $ u_j\ne 0, v_k \ne 0 $. Для матрицы $ f_k(A) $ такое представление обеспечивается собственными векторами матрицы $ A $.
Теорема 3. Матрица $ M_k(A) $ представима в виде произведения
$$ M_k(A) = X_k \cdot Y_k \, , $$ где столбец $ X_k $ — правый собственный вектор, а строка $ Y_k $ — левый собственный вектор матрицы $A $, принадлежащие $ \lambda_k $.
Доказательство следует из теоремы Гамильтона-Кэли: $$ f(A)=\mathbb O_{n\times n} \quad \Rightarrow \ (A-\lambda_k E) f_k(A)= \mathbb O_{n\times n} \, . $$ То есть любой ненулевой столбец матрицы $ f_k(A) $ должен быть правым собственным вектором матрицы $ A $, принадлежащим $ \lambda_k $. Аналогично — для любой строки матрицы $ f_k(A) $.
Теорема 4. Если $ \lambda_{j} \ne \lambda_k $, то
$$ M_j(A)M_k(A) = M_k(A)M_j(A) = \mathbb O_{n\times n} \, . $$
Доказательство. $$ M_j(A) = X_j Y_j \, , M_k(A) = X_k Y_k \, . $$ Результат теоремы следует из равенства $ Y_j X_k = 0 $. Оно справедливо поскольку из определений собственных векторов $$ Y_jA=\lambda_j Y_j,\ AX_k=\lambda_k X_k \, , $$ вытекает $$ Y_jAX_k = \left\{\begin{array}{ll} \lambda_j Y_jX_k ;\\ \lambda_k Y_jX_k \end{array} \right. $$ При условии $ \lambda_{j} \ne \lambda_k $ должно выполняться заявленное равенство.
Теорема 5. Матрицы $ M_k(A) $ идемпотентны, т.е.
$$ [M_k(A)]^2=M_k(A) \, . $$ Кроме того, система матриц $ \{M_k(A)\}_{k=1}^n $ линейно независима и $$ \sum_{k=1}^n M_k(A) =E_{n\times n} \, . $$
Таким образом, представление $$ A=\sum_{k=1}^n \lambda_k M_k(A) $$ можно интерпретировать как разложение матрицы $ A $ по базису $ \{M_k(A)\}_{k=1}^n $, где собственные числа играют роль координат.
Теперь проанализируем полученное разложение на предмет ответа на вопрос: какие именно слагаемые вносят в него наибольший вклад? Какие из них наилучшим образом аппроксимируют матрицу $ A $?
Пример 1.
$$ A= \left[ \begin {array}{rrrr} -19873&40167&14897&28232 \\ -11466&23174&8594&16288\\ 14522 &-29351&-10886&-20632\\ -5349&10812&4011&7601 \end {array} \right] $$
Решение. $ f(\lambda):=\det (A-\lambda E) = {\lambda}^{4}-16\,{\lambda}^{3}-27\,{\lambda}^{2}+178\,\lambda-136 $. Спектр матрицы: $$ \lambda_1 =1,\ \lambda_2=2,\ \lambda_3=17,\ \lambda_4=-4 \, . $$ Не вдаваясь пока в вычислительные детали, приведем выражения для матриц $ M_k(A) $ $$ M_1(A) \approx \left[ \begin {array}{rrrr} 356.3&- 757.4&- 370.3& - 705.6\\ 203.6&- 432.8&- 211.6&- 403.2\\ - 254.5& 541& 264.5& 504\\ 95.4375&- 202.875&- 99.1875&- 189 \end {array} \right] , $$ $$ \ M_2(A) \approx \left[ \begin {array}{rrrr} - 751.66& 1469.16& 587.66& 1239.11\\ - 432.66& 845.66& 338.26& 713.24\\ 539.0&- 1053.5&- 421.4&- 888.53\\ - 198.0& 387.0& 154.8& 326.4 \end {array} \right] , \ \dots $$ Теперь оценим разности $ A-\lambda_k M_k $ для каждого $ k \in \{1,2,3,4\} $ вычислив евклидовы нормы (нормы Фробениуса): $$ \| A-\lambda_1 M_1 \| \approx 75041.7, \ \| A-\lambda_2 M_2 \| \approx 82542.5,\ \| A-\lambda_3 M_3 \| \approx 6632.4,\ \| A-\lambda_4 M_4 \| \approx 65668.2 $$ Видим, что наилучшее приближение дает слагаемое, соответствующее максимальному собственному числу. А следующим за ним, в плане качества приближения, оказывается слагаемое, соответствующее собственному числу, абсолютная величина (модуль) которого является максимальной из модулей оставшихся собственных чисел. Возникает гипотеза, что за качество приближения отвечают именно модули собственных чисел. Эта гипотеза тут же опровергается следующим по убыванию модуля собственным числом, но… Делаем скидку на то, что $ \lambda_1 $ и $ \lambda_2 $ близки.
А если взять пару наилучших слагаемых, то ошибка аппроксимации уменьшается еще значительнeй $$ \| A-\lambda_3 M_3 -\lambda_4 M_4 \| \approx 4348.4 \, . $$ ♦
А сколько таких слагаемых можно оставить в сумме $ \sum_{k=1}^n \lambda_k M_k(A) $ чтобы матрица $ A $ стала распознаваемой?
Именно эту задачу мы будем решать в дальнейшем. Поиск решения будет производиться по уже намеченному выше пути. Но придется принять во внимание несколько нетривиальных аспектов.
1. Как правило, возникающие в приложениях матрицы вещественны. И аппроксимировать их требуется также вещественными матрицами. А собственные числа вещественных квадратных матриц не обязательно вещественны (как случилось в предыдущем примере).
2. Собственные числа квадратных матриц не обязательно различны. А все предшествующие рассуждения были справедливы именно в предположении простоты всех собственных чисел.
3. Возникающие в приложениях матрицы не обязательно квадратные. А для неквадратных матриц отсутствует понятие собственного числа.
4. Задача для таких матриц ставится в постановке: аппроксимировать их наилучшим образом матрицами меньшего ранга. Можно ли как-то предсказать ошибку этого наилучшего приближения (без явного построения самого этого приближения)?
$$ A= A^{\top} $$
Среди общего случая квадратных матриц, симметричные матрицы выделяются рядом хороших свойств для собственных чисел и собственных векторов.
1. Все собственные числа $ \lambda_1,\dots, \lambda_n $ вещественны ( ☞ доказательство).
2. Если $ X_k $ — правый собственный вектор матрицы, принадлежащий собственному числу $ \lambda_k $, то $ X_k^{\top} $ — левый собственный вектор, принадлежащий тому же собственному числу: $$ A X_k = \lambda_k X_k \quad \Rightarrow \quad (A X_k)^{\top} = (\lambda_k X_k)^{\top} \quad \Leftrightarrow \quad X_k^{\top} A^{\top} = \lambda_k X_k^{\top} \quad \Leftrightarrow \quad X_k^{\top} A = \lambda_k X_k^{\top} \, . $$
3. Если $ \lambda_j \ne \lambda_k $, то $ X_j^{\top} X_k = 0 $; говорят, что эти собственные векторы ортогональны ( ☞ доказательство).
4. Симметричная матрица всегда диагонализуема, т.е. ее ЖНФ диагональна даже в случае наличия кратных собственных чисел ( ☞ доказательство). Более того, матрица, приводящая к диагональному виду может быть выбрана во множестве ортогональных матриц: $$ \exists P \in \mathbb R^{n\times n} \quad \mbox{такая, что} \ P^{\top}P=E,\ P^{\top}A P= \left( \begin{array}{cccc} \lambda_1 & & & \mathbb O \\ & \lambda_2 & & \\ && \ddots & \\ \mathbb O&& & \lambda_n \end{array} \right) \, . $$
Теорема 6. Симметричную матрицу $ A \ne \mathbb O_{n\times n} $ всегда можно представить в виде линейной комбинации:
$$ A=\sum_{k=1}^{\operatorname{rank}(A)} \lambda_k M_k(A) $$ одноранговых симметричных матриц $$ \{M_k(A):=X_k \cdot X_k^{\top}\}_{k=1}^{\operatorname{rank}(A)} $$ Здесь $ X_k $ — нормированный собственный вектор: $$ AX_k=\lambda_k X_k,\ X_k^{\top}X_k=1 \, . $$ Матрицы $ \{M_k(A)\}_{k=1}^{\operatorname{rank}(A)} $ удовлетворяют условиям $$ \left[ M_k(A) \right]^2= M_k(A),\ \operatorname{Sp} ( M_k(A) ) =1, \ \operatorname{Sp} ( M_j(A) M_k(A) ) = 0 \quad \mbox{при} \ \lambda_j \ne \lambda_k \, . $$
Пример 2. Представить матрицу
$$ A= \left[ \begin {array}{rrr} 0&-1&2\\ -1&3&-1 \\ 2&-1&0\end {array} \right] $$ в виде комбинации матриц из теоремы.
Решение. Собственные числа и нормированные собственные векторы матрицы: $$ \lambda_1=4, X_{1} = \left( \begin{array}{r} 1/\sqrt{6} \\ -2/\sqrt{6} \\ 1/\sqrt{6} \end{array} \right) \ , \ \lambda_2=-2, X_{2} = \left( \begin{array}{r} -1/\sqrt{2} \\ 0 \\ 1/\sqrt{2} \end{array} \right) \ , \ \lambda_3=1, X_{3}= \left( \begin{array}{r} 1/\sqrt{3} \\ 1/\sqrt{3} \\ 1/\sqrt{3} \end{array} \right) \, . $$ Следовательно $$ A= \lambda_1 X_{1} X_{1}^{\top}+ \lambda_2 X_{2} X_{2}^{\top}+ \lambda_3 X_{3} X_{3}^{\top}= $$ $$ = 4 \, \left[ \begin {array}{rrr} 1/6&-1/3&1/6\\ -1/3&2/3& -1/3\\ 1/6&-1/3&1/6 \end {array} \right] - 2\, \left[ \begin {array}{rrr} 1/2&0&-1/2\\ 0&0&0 \\ -1/2&0&1/2\end {array} \right] + 1\, \left[ \begin {array}{rrr} 1/3&1/3&1/3\\ 1/3&1/3&1/3 \\ 1/3&1/3&1/3 \end {array} \right] \, . $$ ♦
Условия на матрицы $ \{M_k(A)\}_{k=1}^n $ из теоремы 6 означают, что эти матрицы образуют ортонормированный «базис». Почему ортонормированный — понятно. А вот почему слово базис приведено в кавычках?
Можно интерпретировать результат теоремы как разложение симметричной матрицы по ортогональному базису подпространства симметричных матриц? Формально это неверно: количество одноранговых матриц в теореме меньше размерности подпространства симметричных матриц ($=n(n-1)/2$). Тем не менее, для конкретной матрицы $ A $ система матриц $ \{M_k(A)\}_{k=1}^n $ может считаться базисной, а собственные числа $\{\lambda_k\}_{k=1}^n $ — координатами матрицы в этом ортонормированном базисе. И в этом базисе мы будем искать решение для следующей задачи.
Задача. Для натурального числа $ \mathfrak r<n $ построить симметричную матрицу $ \widetilde A $ ранга $ \mathfrak r $ наилучшим образом приближающую симметричную матрицу $ A $ в евклидовой норме (норме Фробениуса): $$ \min_{\widetilde A\in \mathbb R^{n\times n}, \operatorname{rank} (\widetilde A)= \mathfrak r} || A- \widetilde A|| \, , $$
Чем хорош ортонормированный базис? — Тем, что длина вектора, координаты которого удалось найти в этом базисе, выражается только через эти координаты.
Теорема 7. Еcли собственные числа $ \lambda_1,\dots, \lambda_{\mathfrak r} $ симметричной матрицы $ A $ отличны от нуля, то матрица
$$ A_{\mathfrak r}:=\sum_{k=1}^{\mathfrak r} \lambda_k M_k(A) $$ имеет ранг $ \mathfrak r $ и $$ || A- A_{\mathfrak r}||= \sqrt{\lambda_{\mathfrak r+1}^2+\dots+ \lambda_n^2} \, . $$
Доказательство. Ввиду теоремы 6, матрицы $ \{ M_k(A) \}_{k=1}^{\mathfrak r} $ попарно ортогональны и отличны от нулевых. Следовательно (см. теорему ЗДЕСЬ), они линейно независимы. Далее, $$ A- A_{\mathfrak r} = \sum_{\mathfrak r+1}^n \lambda_k M_k(A) $$ и $$ \operatorname{Sp} ((A- A_{\mathfrak r})^{\top} (A- A_{\mathfrak r})) = \operatorname{Sp} ( (A- A_{\mathfrak r})^2) = $$ (матрица $ A- A_{\mathfrak r} $ симметрична) $$ =\operatorname{Sp} \left( \sum_{j,k=\mathfrak r+1}^n \lambda_j \lambda_k M_j(A) M_k(A) \right) = $$ (ввиду теоремы 6) $$ = \operatorname{Sp} \left( \sum_{j=\mathfrak r+1}^n \lambda_j^2 M_j^2(A) \right) = \sum_{j=\mathfrak r+1}^n \lambda_j^2 \operatorname{Sp} \left( M_j^2(A) \right) = \sum_{j=\mathfrak r+1}^n \lambda_j^2 \, . $$ ♦
Количество всевозможных комбинаций из $ n $ собственных чисел по $ \mathfrak r $ равно $ C_n^{\mathfrak r} $. Среди них наименьшее значение для $ || A- A_{\mathfrak r}|| $ достигается когда собственные числа пронумерованы в порядке невозрастания: $$ |\lambda_1| \ge |\lambda_2| \ge \dots \ge |\lambda_{\mathfrak r}| \ge |\lambda_{{\mathfrak r}+1}| \ge \dots \ge |\lambda_n | \, . $$ Оказывается, что соответствующая матрица $ A_{\mathfrak r} $ обеспечивает глобальный минимум в поставленной задаче аппроксимации матрицы $ A$ матрицей ранга $ {\mathfrak r} $. Фактически формула $$ A=\sum_{k=1}^n \lambda_k M_k(A) $$ представляет собой разбиение матрицы $ A $ на слагаемые, расположенные по убыванию их «вклада» в матрицу. Если матрица отражает своими элементами какую-то информационную сущность,
то можно попробовать заменить ее хранение (пересылку) на матрицу $ A_{\mathfrak r} $ при подходящем выборе $ {\mathfrak r} < \operatorname{rank} (A) $.
Пример 3. Для матрицы примера 2
$$ \|A-\lambda_1M_1(A) \| =\sqrt{5},\ \|A-\lambda_2M_2(A) \|=\sqrt{17},\ \|A-\lambda_3M_3(A) \|= \sqrt{20} \ ; $$ $$ \|A-\lambda_1M_1(A) -\lambda_2M_2(A) \| =1, \ \|A-\lambda_1M_1(A) -\lambda_3M_3(A) \| = 2, \ \|A-\lambda_2M_2(A) -\lambda_3M_3(A) \| = 4 \, . $$
Каков геометрический смысл этой теоремы?
Задача. Найти координаты точки на плоскости, сумма квадратов расстояний до которой от точек $ P_1=(x_{1},y_1),\dots, P_n=(x_n,y_n) $ минимальна.
Теорема 1. Координаты искомой точки определяются средними значениями координат точек $ \{ P_j\}_{j=1}^n $:
$$ \overline{x}=\frac{x_1+x_2+ \dots+x_n}{n},\quad \overline{y}=\frac{y_1+y_2+ \dots+y_n}{n}\ . $$
Координаты точки, для которой величина
$$ \sum_{j=1}^n m_j \left[(x-x_j)^2+(y-y_j)^2 \right] \ , $$ ($ m_{1},\dots,m_n $ — положительные константы) будет минимальной: $$ \hat x=\frac{m_1x_1+m_2x_2+ \dots+m_nx_n}{m_1+m_2+ \dots+m_n},\quad \hat y=\frac{m_1y_1+m_2y_2+ \dots+m_ny_n}{m_1+m_2+ \dots+m_n}\ . $$
Иными словами: искомая точка совпадает с центром масс (барицентром) системы масс $ \{m_j\}_{j=1}^n $, расположенных в точках $ \{P_j\}_{j=1}^n $.
Задача. Найти прямую на плоскости, сумма квадратов расстояний до которой от точек $ \{P_j\}_{j=1}^n $ минимальна.
Если обозначить $ \delta_j $ расстояние от $ P_j $ до прямой $$ ax+by+c=0 \, , $$ то речь о поиске $$ \min_{(a,b,c)\in \mathbb R^3} (\delta_1^2+\delta_2^2+\dots + \delta_n^2) \, . $$
Теорема 2 [Пирсон]. Искомый минимум равен минимальному положительному корню $ \lambda_{min} $ полинома
$$ f(\lambda):=\left| \begin{array}{cc} \displaystyle \sum_{j=1}^n (x_j-\overline{x})^2 - \lambda & \displaystyle \sum_{j=1}^n (x_j-\overline{x})(y_j-\overline{y}) \\ \displaystyle \sum_{j=1}^n (x_j-\overline{x})(y_j-\overline{y}) & \displaystyle \sum_{j=1}^n (y_j-\overline{y})^2 - \lambda \end{array} \right| \, . $$
Доказательство ☞ ЗДЕСЬ.
Корни квадратного полинома $ f(\lambda) $ из теоремы Пирсона вещественны и неотрицательны. Можно доказать это применением неравенства Коши-Буняковского. С другой стороны, для доказательства можно применить и более общий теоретический результат, который позволит распространить результат теоремы в $\mathbb R^3 $ (и даже в пространства бóльших размерностей). Дело в том, что полином $ f(\lambda) $ является характеристическим полиномом матрицы $$ G=\left(\begin{array}{cc} \displaystyle \sum_{j=1}^n (x_j-\overline{x})^2 & \displaystyle \sum_{j=1}^n (x_j-\overline{x})(y_j-\overline{y}) \\ \displaystyle \sum_{j=1}^n (x_j-\overline{x})(y_j-\overline{y}) & \displaystyle \sum_{j=1}^n (y_j-\overline{y})^2 \end{array} \right) \, , $$ которую можно интерпретировать как матрицу Грама системы векторов $$ \left\{ (x_1-\overline{x},\dots, x_n-\overline{x}), (y_1-\overline{y},\dots, y_n-\overline{y}) \right\} \, . $$ при стандартном способе задания скалярного произведения в $\mathbb R^n $. Доказывается (см. ☞ ЗДЕСЬ), что корни характеристического полинома матрицы Грама (т.е. ее собственные числа) всегда вещественны и неотрицательны.
Если корни полинома $ f(\lambda) $ различны, то прямая, минимизирующая сумму квадратов расстояний до нее от точек $ \{P_j\} $ единственна, она проходит через точку $ (\overline{x},\overline{y}) $, и ее направляющим вектором является собственный вектор матрицы $ G $, принадлежащий максимальному собственному числу этой матрицы.
Пример [Пирсон]. Построить прямую из теоремы $ 2 $ для точек
$$ \begin{array}{c|cccccccccc} x & 0.0 & 0.9 & 1.8 & 2.6 & 3.3 & 4.4 & 5.2 & 6.1 & 6.5 & 7.4 \\ \hline y & 5.9 & 5.4 & 4.4 & 4.6 & 3.5 & 3.7 & 2.8 & 2.8 & 2.4 & 1.5 \end{array} $$
Решение. Здесь $ n=10 $ и $$ \sum_{j=1}^{10} x_j= 38.2,\ \sum_{j=1}^{10} y_j= 37.0 \quad \Rightarrow \quad \overline x =3.82, \overline y = 3.70 \, . $$ $$ G=\left[ \begin {array}{rr} 56.396&- 30.430 \\ - 30.430& 17.220 \end {array} \right] $$ Характеристический полином $$ f(\lambda)=\lambda^2-73.616\,\lambda+45.15422 $$ имеет корни $$ \lambda_1 \approx 0.618573,\ \lambda_2 \approx 72.997427 \, . $$ Собственные векторы, принадлежащие этим собственным числам: $$ \approx \left[ 0.478924, 0.877856\right]^{\top} \quad \mbox{и} \quad \approx \left[ 0.877856, - 0.478924\right]^{\top} $$ соответственно; оба нормированы. Сумма квадратов расстояний от заданных точек до ближайшей к ним прямой равна $ \lambda_1 $. Из следствия к теореме Пирсона получаем уравнение этой прямой $$ 0.478924\, x+0.877856\,y-5.077556=0 \quad \Leftrightarrow \quad y = -0.545561\, x + 5.784082 \, . $$
Эта прямая совпадает с решением Пирсона, но вот значение для минимума у него приведено другое: $ 0.2484 $. Непосредственная проверка подтверждает значение $ \lambda_1 $. ♦
В примере Пирсона исходные данные взяты абстрактными, не имеющими прототипа в реальном мире.
Для симметричных матриц больших порядков проблема вычисления полного разложения в сумму одноранговых симметричных матриц $$ A=\sum_{k=1}^n \lambda_k M_k(A) \quad \mbox{при} \ M_k(A)=X_k X_k^{\top}, \ X_k^{\top} X_k =1 $$ требует нахождения ее спектра, т.е. набора всех ее собственных чисел, а также соответствующих собственных векторов. Однако уже проблема вычисления характеристического полинома матриц порядка $ 100 $ является неприятной. Теоретически алгоритмы существуют, но на практике их не используют. Применяют подход, решающий частичную проблему собственных чисел, т.е. проблему ограничивают вычислением максимальных по модулю собственных чисел и соответствующих собственных векторов.
Пусть собственные числа симметричной матрицы $ A \in \mathbb R^{n\times n} $ занумерованы в порядке убывания модулей: $$ |\lambda_1|>|\lambda_2|> \dots > |\lambda_n | $$ и среди этих модулей нет одинаковых. Тогда $ \lambda_1 $ может быть вычислен методом степенных итераций. Находим его с необходимой точностью (пока не известно какой), а также вычисляем соответствующий ему (правый) собственный вектор $ X_1 \in \mathbb R^{n} $ (столбец), нормируем его $ \| X_1 \|=1 $.
И выбрасываем из разложения матрицы $ A $ первое слагаемое: $$ A_2:=A- \lambda_1 X_1 \cdot X_1^{\top} = \sum_{k=2}^n \lambda_k M_k(A) $$ Если все вычисления производятся точно, то $ \operatorname{rank} (A_2)= \operatorname{rank} (A)-1 $, спектром матрицы $ A_2 $ является набор $ \{\lambda_2,\lambda_3, \dots, \lambda_n, 0 \} $, а ее доминирующим (по модулю) собственным числом становится $ \lambda_2 $. Для его поиска снова применяем метод степенных итераций.
о матрицах $ A\cdot A^{\top} $ и $ A^{\top} A $ для матрицы $ A\in \mathbb R^{m\times n} $.
1. Обе матрицы — симметричны (хотя и не обязательно одинаковых порядков).
2. Матрица $ A\cdot A^{\top} $ — матрица Грама системы строк, а $ A^{\top} A $ — матрица Грама системы столбцов матрицы $ A $. Так если представить матрицу $ A $ как результат конкатенации ее столбцов $ A=[A_{[1]}|A_{[2]}|\dots |A_{[n]} ] $, то $$ A^{\top} A = G(A_{[1]},A_{[2]},\dots ,A_{[n]}) \, ; $$ здесь скалярное произведение задано стандартным способом.
Отсюда, в частности, следует: $$ \operatorname{rank} ( A\cdot A^{\top}) = \operatorname{rank} (A^{\top} A) = \operatorname{rank}(A) \, . $$
3. Характеристические полиномы этих матриц могут различаться только лишь на степень переменной: $$ (-\lambda)^n \det (A \cdot A^{\top} - \lambda E_{m\times m})\equiv (-\lambda)^m \det ( A^{\top} A - \lambda E_{n\times n}) \ . $$
4. Любое ненулевое собственное число одной матрицы вещественно, положительно и является собственным числом другой матрицы. Количество ненулевых собственных чисел совпадает с $ \operatorname{rank} (A) $.
5. Если обозначить $ \{ \sigma_1^2,\dots, \sigma_{\operatorname{rank} (A)}^2 \} $ ненулевые собственные числа матрицы $ A\cdot A^{\top} $ (или $ A^{\top} A $), а соответствующие им ортонормированные собственные векторы $ \{ U_{[1]},\dots U_{[\operatorname{rank} (A)]} \} \subset \mathbb R^m $ (соответственно, $ \{ V_{[1]},\dots V_{[\operatorname{rank} (A)]} \} \subset \mathbb R^n $), то справедливы равенства $$ A\cdot A^{\top} = \sum_{k=1}^{\operatorname{rank} (A)} \sigma_k^2 U_{[k]}U_{[k]}^{\top} \quad \mbox{и} \quad A^{\top} A = \sum_{k=1}^{\operatorname{rank} (A)} \sigma_k^2 V_{[k]}V_{[k]}^{\top} \, . $$
6. Существуют положительно полуопределенные квадратные корни из матриц $ A\cdot A^{\top} $ и $ A^{\top} A $, т.е. вещественные симметричные матрицы $ P \in \mathbb R^{m\times m} $ и $ Q \in \mathbb R^{n\times n}$ такие, что $$ P^2= A\cdot A^{\top},\ Q^2=A^{\top} A \, . $$
Теорема ??. Матрицу $ A \in \mathbb R^{m\times n}, A\ne \mathbb O_{m\times n} $ всегда можно представить в виде линейной комбинации:
$$ A=\sum_{k=1}^{\operatorname{rank}(A)} \sigma_k M_k(A) \qquad \qquad \qquad {\color{Red}{ (1)} } $$ одноранговых матриц $$ \left\{M_k(A):=U_{[k]} V_{[k]}^{\top}\right\}_{k=1}^{\operatorname{rank}(A)} \, . $$ Здесь $ \{\sigma_k>0\} $ — арифметические корни квадратные из ненулевых собственных чисел матрицы $ A \cdot A^{\top} $ . Система столбцов $ \{U_{[k]}\} \subset \mathbb R^{m} $ состоит из ортонормированных собственных векторов, принадлежащих $ \{\sigma_k^2\} $, а система столбцов $ \{V_{[k]}\} \subset \mathbb R^{n} $ состоит из ортонормированных собственных векторов матрицы $ A^{\top} A $, принадлежащих $ \{\sigma_k^2\} $. При этом столбцы $ U_{[k]} $ и $ V_{[k]} $ связаны дополнительными соотношениями $$ A V_{[k]} = \sigma_k U_{[k]}, \ A^{\top} U_{[k]} = \sigma_k V_{[k]} \, \qquad \qquad \qquad {\color{Red}{ (2)} } . $$
Числа $ \{ \sigma_k \} $ называются сингулярными числами матрицы $ A $; столбец $ U_{[k]} $ называют левым сингулярным вектором, а столбец $ V_{[k]} $ — правым сингулярным вектором матрицы $ A $, принадлежащими сингулярному числу $ \sigma_k $.
Представление матрицы $ A $ по формуле $ {\color{Red}{ (1)} } $ называется сингулярным разложением1) матрицы $ A $.
Имеет смысл начинать вычисления с тех сингулярных векторов, которые являются собственными для матрицы наименьшего порядка из $ m\times m $ или $ n \times n $ (да и собственные числа в таком варианте вычислять проще!).
Пример 1. Построить сингулярное разложение матрицы
$$ A=\left( \begin{array}{rrr} 1 & -1 & -2 \\ -7/3 & 1/3 & 2/3\\ 1/3 & -7/3 & -2/3 \\ -5/3 & 5/3 & -2/3 \end{array} \right) \, . $$
Решение. Из двух вариантов $ A \cdot A^{\top} $ и $ A^{\top} A $ выбираем второй: $$ A^{\top} A = \left( \begin{array}{rrr} 28/3 & -16/3 & -8/3 \\ -16/3 & 28/3 & 8/3\\ -8/3 & 8/3 & 16/3 \end{array} \right) \, . $$ Характеристический полином $$ \det (A^{\top} A - \lambda E)= -(\lambda-16)(\lambda-4)^2 \, . $$ Следовательно сингулярными числами матрицы $ A $ являются $ \sigma_1=4,\sigma_2=2,\sigma_3=2 $. Ищем ортонормированные собственные векторы матрицы $ A^{\top} A $. Если бы все собственные числа матрицы $ A^{\top}A $ были простыми, то можно было бы воспользоваться универсальным представлением собственного вектора посредством матрицы, взаимной матрице $ A \cdot A^{\top} - \lambda E $ (см. ЗДЕСЬ). Достаточно найти ее первый столбец: $$ \left( \begin{array}{crr} (\lambda-32/3)(\lambda-4) & * & * \\ 16/3(4-\lambda) & * & * \\ 8/3(4-\lambda) & * & * \end{array} \right) $$ Подстановка $ \lambda=16 $ в него дает собственный вектор в виде $ [64,-64,-32]^{\top} $; после нормирования получим $ V_{[1]}=[2/3,-2/3,-1/3]^{\top} $ . Однако при значении кратного корня $ \lambda = 4 $ столбец вырождается в нулевой, т.е. поиск двух линейно независимых собственных векторов посредством взаимной матрицы принципиально неосуществим. На помощь приходит минимальный аннулирующий полином матрицы $ A^{\top} A $. Таким является $ (\lambda-16)(\lambda-4) $. Матрица $$ A^{\top} A - 16 E= \left( \begin{array}{rrr} -20/3 & -16/3 & -8/3 \\ -16/3 & -20/3 & 8/3 \\ -8/3 & 8/3 & -32/3 \end{array} \right) $$ имеет ранг равный $ 2 $, и любая пара ее столбцов будет линейно независимой системой собственных векторов матрицы, принадлежащих собственному числу $ \lambda=4 $. Возьмем, например, второй и третий столбцы и проведем с ними процедуру ортогонализации Грама-Шмидта. Получим столбцы $ [-16/3,-20/3, 8/3]^{\top} $, $ [-12,-0,-24]^{\top} $, которые после нормирования дают $$ V_{[2]}=\left[ -4/(3\sqrt{5}),-\sqrt{5}/3, 2/\sqrt{5} \right]^{\top}, V_{[3]}=\left[1/\sqrt{5},0, 2/\sqrt{5} \right]^{\top} \, . $$ Нахождение левых сингулярных векторов производится с помощью формулы $ {\color{Red}{ (2)} } $ как $ U_{[j]}= 1/\sigma_j A V_{[j]} $: $$ U_{[1]}=\left[\frac{1}{2},-\frac{1}{2},\frac{1}{2},-\frac{1}{2}\right]^{\top},\ U_{[2]}=\left[-\frac{1}{2\sqrt{5}},\frac{3}{2\sqrt{5}},\frac{3}{2\sqrt{5}},-\frac{1}{2\sqrt{5}}\right]^{\top},\ $$ $$ U_{[3]}=\left[-\frac{3}{2\sqrt{5}},-\frac{1}{2\sqrt{5}},-\frac{1}{2\sqrt{5}},-\frac{3}{2\sqrt{5}}\right]^{\top} \, . $$ Формула $ {\color{Red}{ (1)} } $ дает представление матрицы $ A $ в виде $$ 4 \left( \begin{array}{rrr} 1/3 & -1/3 & -1/6 \\ -1/3 & 1/3 & 1/6 \\ 1/3 & -1/3 & -1/6 \\ -1/3 & 1/3 & 1/6 \end{array} \right)+2 \left( \begin{array}{rrr} -1/3 & -1/6 & -1/3 \\ -1/3 & -1/6 & -1/3 \\ -1/3 & -1/6 & -1/3 \\ -1/3 & -1/6 & -1/3 \end{array} \right) + 2 \left( \begin{array}{rrr} 1/6 & 1/3 & -1/3 \\ -1/6 & -1/3 & 1/3 \\ -1/6 & -1/3 & 1/3 \\ 1/6 & 1/3 & -1/3 \end{array} \right) \, . $$ Впрочем, на практике экономнее хранить не матрицы с фактически одинаковыми строками (столбцами), а именно их представления через $ U_{[j]} $ и $ V_{[j]} $.
Система матриц $\{M_k(A)\} $ ортонормирована.
$$ \operatorname{Sp} \left( (U_{[j]} V_{[j]}^{\top})^{\top} U_{[k]} V_{[k]}^{\top} \right)= \operatorname{Sp} \left( V_{[j]}U_{[j]}^{\top} U_{[k]} V_{[k]}^{\top} \right) = 0 $$ поскольку $ U_{[j]}^{\top} U_{[k]} =0 $.
Часто формула сингулярного разложения матрицы записывается в форме матричного произведения.
Теорема. Составим матрицы из сингулярных чисел и сингулярных векторов матрицы $ A\in \mathbb R^{m\times n} $:
$$ {\color{Red}{\mathbf{\Sigma}}}:= \operatorname{diag}(\sigma_1,\dots, \sigma_{\operatorname{rank} (A)}) = \left( \begin{array}{cccc} \sigma_1 & 0 & \dots & 0 \\ 0 & \sigma_2 & \dots & 0 \\ \vdots && \ddots & \vdots \\ 0 & 0 & \dots & \sigma_{\operatorname{rank} (A)} \end{array} \right) $$ $$ U:=\left[U_{[1]}\mid \dots \mid U_{[\operatorname{rank} (A)]} \right], \ V:=\left[V_{[1]}\mid \dots \mid V_{[\operatorname{rank} (A)]} \right] \ , $$ Справедливо равенство $$ A=U \cdot {\color{Red}{\mathbf{\Sigma}}} \cdot V^{\top} \qquad \qquad {\color{Red}{ (1)^{\prime}} }\, . $$
Эквивалентность этого представления формуле $ {\color{Red}{ (1)} } $ вытекает из следующих соображений: $$ {\color{Red}{\mathbf{\Sigma}}}= \operatorname{diag}(\sigma_1, 0,\dots, 0) + \operatorname{diag}(0, \sigma_2,\dots, 0)+\dots $$ $$ U {\color{Red}{\mathbf{\Sigma}}} V^{\top}=U(\operatorname{diag}(\sigma_1, 0,\dots, 0) + \operatorname{diag}(0, \sigma_2,\dots, 0)+\dots) V^{\top} = $$ $$ =U \operatorname{diag}(\sigma_1, 0,\dots, 0) V^{\top} + \dots = \sigma_1 \left[U_{[1]}\mid \dots \mid U_{[\operatorname{rank} (A)]} \right] \operatorname{diag}(1, 0,\dots, 0) V^{\top} + \dots = $$ $$ = \sigma_1 \left[U_{[1]}\mid \mathbb O \mid \dots \mid \mathbb O \right] V^{\top} + \dots = \sigma_1 \left[U_{[1]}\mid \mathbb O \mid \dots \mid \mathbb O \right] \left[ \begin{array}{c} V_{[1]}^{\top} \\ V_{[2]}^{\top} \\ \vdots \end{array} \right] + \dots = \sigma_1 U_{[1]} V_{[1]}^{\top} + \dots $$ То есть мы получили первое слагаемое формулы $ {\color{Red}{ (1)} } $. Для остальных слагаемых рассуждения аналогичны.
Для симметричной положительно полуопределенной матрицы $ A\in \mathbb R^{n\times n} $ формула сингулярного разложения эквивалентна формуле приведения этой матрицы к диагональному виду с помощью ортогонального преобразования.
Действительно, в качестве сингулярных чисел выступают ненулевые собственные числа матрицы $ A $. Матрицы $ U $ и $ V $ одинаковы, и состоят из ортонормированной системы собственных векторов, принадлежащих этим собственным числам. В частности, при $ \operatorname{rank}(A)=n $ матрица $ U=V $ ортогональна: $ U^{\top} U = E $.
Теорема ??. Упорядочим сингулярные числа матрицы $ A\in \mathbb R^{m\times n} $ по невозрастанию:
$$ \sigma_1 \ge \sigma_2 \ge \dots \ge \sigma_{\operatorname{rank}(A)} >0 \, . $$ Матрицей $ \widetilde A \in \mathbb R^{m\times n} $ ранга $ \mathfrak r , 0<\mathfrak r \le \operatorname{rank}(A) $, наилучшим образом приближающей $ A $ в евклидовой норме (норме Фробениуса): $$ \min_{\widetilde A\in \mathbb R^{m\times n}, \operatorname{rank} (\widetilde A)= \mathfrak r} \| A- \widetilde A\| \, , $$ является матрица $$ A_{\mathfrak r}=\sum_{j=1}^{\mathfrak r} \sigma_j U_{[j]} V_{[j]}^{\top} \, . $$ При такой матрице $ \widetilde A $ величина указанного минимума равна $$ \sqrt{\sum_{j=\mathfrak r+1}^{\operatorname{rank}(A)} \sigma_j^2} \, . $$
Доказательство ЗДЕСЬ
Обозначим
$$ U_{\mathfrak r}:= [U_{[1]}|\dots | U_{[\mathfrak r]}] \in \mathbb R^{n\times \mathfrak r } \, .$$ Тогда матрица $ A_{\mathfrak r} $ из теоремы представима в виде
$$ A_{\mathfrak r} = U_{\mathfrak r} \cdot U_{\mathfrak r}^{\top} A \, . $$ Аналогично, если $$ V_{\mathfrak r}:= [V_{[1]}|\dots | V_{[\mathfrak r]}] \in \mathbb R^{\mathfrak r \times m } \, .$$ то $$ A_{\mathfrak r} = A V_{\mathfrak r} \cdot V_{\mathfrak r}^{\top} \, . $$