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

Сингулярное разложение - одна одна из лучших вещей, которые есть в линейной алгебре!

И знаете, с этим трудно не согласиться! Сингулярное разложение матриц (Singular Value Decomposition) используется в алгоритмах сжатия данных, теории автоматического управления, машинного обучения и много где ещё. Поскольку SVD не является главным объектом исследования данной статьи, я предоставляю читателю возможность самостоятельно ознакомиться с примерами ниже. Если же сингулярное разложение не вызывает вопросов, то перейдём к самому интересному!

Немного о сингулярном разложении

Рассмотрим некоторую матрицу M вместе с её сингулярным разложением:

\begin{gather} M = U \Sigma V^T \\ M \in K^{m \times n}, U \in K^{m \times m}, \Sigma \in K^{m \times n}, V \in K^{n \times n} \end{gather}

В данной записи нас интересует лишь матрица \Sigma, поскольку именно она содержит mсингулярных чисел, расположенных на её главной диагонали. Матрицы U, Vсодержат в себе левые и правые сингулярные вектора соответственно. Элементы матриц должны принадл��жать некоторому полю K- вещественному или комплексному.

Геометрический смысл показанного уравнения прекрасно показан ниже:

Геометрический смысл сингулярного разложения
Геометрический смысл сингулярного разложения

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

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

Пример сжатия изображения с помощью SVD
Пример сжатия изображения с помощью SVD

Автор привёл потрясающий пример того, как взяв 5, 10, 15 и 20 наибольших сингулярных чисел матрицы можно восстановить большую часть информации об исходном изображении. Даже при использовании 5 сингулярных чисел изображение становится различимым, а ведь изначально их было аж 360!

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

Первичное голосование

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

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

Пример фотографии звёздного неба
Пример фотографии звёздного неба

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

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

Теперь в системе координат камеры нужно провести луч до каждой из потенциальных звёзд. Любой луч задаётся единичным вектором (обозначим его как p_i), координаты которого можно посчитать, используя не очень сложные математические формулы. Я допускаю, что некоторым читателям они покажутся знакомыми.

Как посчитать вектора

Предположим, что у нас есть ширина \omega и высота h исходного изображения, а также фокусное расстояние камеры f_{pix} в пикселях. Координаты центроида, соответствующего единичному вектору к потенциальной звезде, обозначим как(x_i, y_i). Тогда:

\begin{gather} cx = w / 2 \\ cy = h / 2 \\[4pt] p_i = \begin{bmatrix} (x_i - cx) / f_{pix} \\ (y_i - cy)/f_{pix} \\ 1 \end{bmatrix}\end{gather}

Для использования вектораp_i в дальнейших расчётах его осталось лишь нормализовать.

Для всех центроидов, найденных на изображении, необходимо выбрать тройки и составить из них треугольники. И здесь может возникнуть вопрос: а разве треугольников не будет слишком много? И правда, их количество равно C_N^{3}, где N - количество потенциальных кандидатов. То есть для N = 40 потребуется составить почти 10 тысяч треугольников, что непозволимо дорого, если алгоритм работает, например, на микроконтроллере. По этой причине для составления треугольников лучше использовать N_{max} самых ярких центроидов.

Данная статья начиналась с матриц. Мы вынуждены вернуться к ним. Для каждой выбранной тройки единичных векторов (p_i, p_j, p_k) составим матрицу M таким образом, чтобы каждый из векторов являлся столбцом:

M = \begin{bmatrix} p_{i_x} & p_{j_x} & p_{k_x} \\ p_{i_y} & p_{j_y} & p_{k_y} \\ p_{i_z} & p_{j_z} & p_{k_z} \end{bmatrix} = U \begin{bmatrix} \sigma_1 & 0 & 0 \\ 0 & \sigma_2 & 0 \\ 0 & 0 & \sigma_3\end{bmatrix}V^T, \; \sigma_1 > \sigma_2 > \sigma_3

Два наибольших сингулярных числа \sigma_1 и \sigma_2будут "кодировать" достаточно информации о треугольнике. Но на этом работа с выделением ключевых признаков треугольника не закончена. Логично предположить, что среди огромного числа троек звёзд на небе наверняка найдется несколько (на самом деле, может оказаться очень много) тех, которые образуют очень похожие треугольники. Для решения этой проблемы вводятся дополнительные параметры \Delta_1 и \Delta_2 - шаги дискретизации сингулярных чисел, которые намеренно выбираются очень малыми. После дискретизации числа округляются вниз до ближайшего целого числа. Рассмотрим простой пример, иллюстрирующий вышесказанное:

\begin{gather} \sigma_1 = 1.56155, \sigma_2 = 0.12598 \\ \Delta_1 = 0.0001, \Delta_2 = 0.001 \\ \sigma_{d_1} = \sigma_1 / \Delta_1 = 15615 \\ \sigma_{d_2} = \sigma_2/\Delta_2 = 125 \end{gather}

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

Я не уверен, поэтому лучше промолчать, чем выдать неверно идентифицированные звёзды!

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

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

Финишная прямая (этап валидации)

Осталась лишь одна небольшая проблема, которую очень легко продемонстрировать на примитивном примере.

Пример расположения точек
Пример расположения точек

Пусть у нас есть некоторый набор точек A, B, C, D, образующий два треугольника с одной общей стороной.

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

Я думаю, что вершины нужно обозначить вот так!

Обычно расстояние между звёздами называется угловым и измеряется в градусах. Внутри каталога есть информация о том, какое угловое расстояние (отрезок какой длины) между каждой из вершин (звёзд). Нужно лишь проверить, что в каталоге расстояние между точками A и D является таким же, как на кадре, ну или хотя бы с некоторой погрешностью. Если это так, то звезда скорее была идентифицирована верно, и ей нужно присвоить валидационный голос. Иными словами, угловые расстояния должны быть согласованы между собой. Ведь если окажется, что в каталоге между звёздами A и D угловое расстояние 150 градусов, а мы насчитали 10, то кандидаты точно являются ложными.

Введение понятия "валидационный голос" является ключевым, ведь именно количество таких голосов показывает, насколько алгоритм уверен в правильности идентификации. И правда, если мы имеем N победителей голосования, а каждый из кандидатов на этапе верификации набрал ровно N дополнительных голосов, то можно быть почти полностью уверенным в том, что алгоритм действительно верно идентифицировал каждую из звёзд. Почему я употребил слово "почти"? Потому что всё ещё есть вероятность, что где-то в другой части звёздного неба возможно сделать кадр, на котором есть M очень похожих треугольников. При плохом выборе шага дискретизации и при малом количестве звёзд на кадре такая ситуация вполне возможна.

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

Подмечу несколько свойств описанного алгоритма:

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

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

  • Гибкость настройки - я уже упомянул, что изменяя шаги дискретизации \Delta_1 и \Delta_2 перед нами открывается возможность повышать робастность системы. Например, их увеличение может поспособствовать тому, что алгоритм будет способен идентифицировать большее количество звёзд на кадре. В то же время это снижает количество уникальных пар дискретизированных сингулярных чисел в каталоге, то есть растёт вероятность банально перепутать наборы треугольников.

Заключение

Наверное, если перейти к результатам практического использования данного алгоритма, то данную статью вполне можно будет превратить если не в диплом, то в полноценную научную статью. А ведь среди результатов есть весьма интересные. Например, кажется очевидным, что просто положить в каталог всевозможные комбинации троек звёзд нельзя. Ну разве что если вы не пытаетесь реализовать алгоритм, имея за спиной целый дата-центр, наполненный современным оборудованием. Треугольники нужно генерировать таким образом, чтобы они равномерно покрывали небо. При этом не должно быть слишком вытянутых. Более того, само распределение треугольников по небу может в конечном итоге очень сильно влиять на качество идентификации. При неграмотном подходе к генерации каталога треугольников возможно возникновение вот такой ситуации:

Пример плохого распределения треугольников по небу
Пример плохого распределения треугольников по небу

Если получится так, что в кадре окажутся звезды A, D, F, но не окажется B и C, то алгоритм идентификации скажет, что ничего идентифицировать удалось. Хотя казалось бы, треугольники покрывают площадь достаточно равномерно, при этом нету ни одного слишком узкого или длинного. Звёзды в кадре могут быть очень яркими и даже образовывать известное созвездие, которое вы периодически в ясную погоду видите у себя из окна. Однако в каталоге его не окажется, а значит, и алгоритму эти звёзды в такой конфигурации неизвестны. Если все 5 звезд попадут в один кадр, то с большой долей вероятности они все будут идентифицированы.

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

Использованные источники:

  1. https://www.mdpi.com/1424-8220/21/9/3084#:~:text=The multilayer voting algorithm involves,based on star pairs voting.

  2. https://dmicz.github.io/machine-learning/svd-image-compression/

  3. https://ru.wikipedia.org/wiki/Сингулярное_разложение