Участник:OscarFox/Обобщенный прокрустов анализ
Автор: Фостенко О.А.
Свойства и структура алгоритма
Общее описание алгоритма
Обобщенный прокрустов анализ (Generalized Procrustes Analysis, GPA) - это метод, предназначенный для сравнения и выравнивания нескольких наборов из одного и того же числа точек (форм) в пространстве одной размерности. Цель GPA - найти оптимальные преобразования, включающие вращения, масштабирование и сдвиги, для наилучшего согласования форм друг с другом.
Формулировка похожей задачи встречается в работе B. F. Green 1952 года[1], аналитическое решение задачи выравнивания двух множеств точек друг относительно друга с учетом только вращений (ортогональная задача Прокруста) было предложено P. H. Schönemann в 1966 году[2], в 1970 году им же в соавторстве с R. M. Carroll[3] предложено аналитическое решение для модификации, учитывающей, помимо вращений, смещения и масштабирование (классический прокрустов анализ), а непосредственно алгоритм GPA был разработан и опубликован J. C. Gower в 1975 году[4].
GPA - это итеративный алгоритм: он последовательно согласовывает все формы относительно некоторой «средней формы», обновляя её до тех пор, пока изменения не станут незначительными. К «средней форме» набор из точек приводится с помощью классического прокрустова анализа. Иллюстрация выравнивания, выполняемого с помощью как классического, так и обобщенного прокрустова анализа, приведена на рисунке 1.
Математическое описание алгоритма
Описание данных и формул
Входные данные:
- [math]\displaystyle{ \{X_1, X_2,\dots,X_m\} }[/math] - набор входных форм, то есть матриц одного размера: [math]\displaystyle{ X_i \in \mathbb{R}^{n \times d}, \; i=1,\dots,m }[/math]. Каждая матрица [math]\displaystyle{ X_i }[/math] представляет собой набор из [math]\displaystyle{ n }[/math] точек в [math]\displaystyle{ d }[/math]-мерном пространстве. Точки с одинаковыми индексами (т.е. строки матриц с одинаковыми номерами) предполагаются соответствующими друг другу.
Выходные данные:
- Выровненные формы [math]\displaystyle{ \{ X_1', X_2',\dots,X_m'\} }[/math], где [math]\displaystyle{ X_i'\in \mathbb{R}^{n \times d}, \; i=1,\dots,m }[/math].
Вспомогательные понятия:
- [math]\displaystyle{ CPA(X, Y): \mathbb{R}^{n \times d} \times \mathbb{R}^{n \times d} \rightarrow \mathbb{R}^{n \times d} }[/math] - функция, выполняющая классический прокрустов анализ, вращая, смещая и масштабируя форму [math]\displaystyle{ X }[/math] так, чтобы она наилучшим образом подходила под форму [math]\displaystyle{ Y }[/math]. Возвращает матрицу из [math]\displaystyle{ \mathbb{R}^{n \times d} }[/math], отвечающую преобразованной форме [math]\displaystyle{ X }[/math].
Описание алгоритма (GPA)
Пусть [math]\displaystyle{ X_i' = s_i X_i R_i + \mathbf{1}_n t_i }[/math]. Здесь [math]\displaystyle{ s_i\in \mathbb{R} }[/math] есть константа масштаба, [math]\displaystyle{ R_i\in \mathbb{R}^{d \times d} }[/math] - ортогональная матрица, отвечающая вращению, [math]\displaystyle{ \mathbf{1}_n\in \mathbb{R}^{n \times 1} }[/math] - вектор-столбец из единиц, а [math]\displaystyle{ t_i\in \mathbb{R}^{1 \times n} }[/math] - вектор переноса.
В процессе работы GPA итеративно минимизирует следующую функцию:
[math]\displaystyle{ S = \displaystyle \sum_{i\lt j}^{m} \| X_i' - X_j' \|_F^2 }[/math]
по всевозможным смещениям [math]\displaystyle{ t_i }[/math], вращениям [math]\displaystyle{ R_i }[/math] и константам масштаба [math]\displaystyle{ s_i }[/math] для всех [math]\displaystyle{ i=1,\dots,m }[/math], с условием, что [math]\displaystyle{ \displaystyle \sum_{i=1}^{m} s_i^2 \| X_i \|_F^2 = \sum_{i=1}^{m} \| X_i \|_F^2 }[/math] (чтобы избежать зануления всех [math]\displaystyle{ s_i }[/math]).
Минимизация [math]\displaystyle{ S }[/math] сводится к минимзации [math]\displaystyle{ S_r = \displaystyle \sum_{i=1}^{m} tr(s_i^2 X_i' X_i'^T) - m \cdot tr(Y Y^T) }[/math], где [math]\displaystyle{ Y }[/math] - средняя форма (консенсус).
Шаг 1. Инициализация.
Задаем начальные значения выравниваемых форм, центрируя их:
Вычисляем смещение для центрирования форм: [math]\displaystyle{ \bar{X_i} = \frac{1}{n} \mathbf{1}_n^T X_i }[/math]
Центрируем формы: [math]\displaystyle{ X_i' := X_i - \mathbf{1}_n \bar{X_i}^T, \; i = 1,\dots,m. }[/math]
Далее нормируем масштаб форм (в совокупности):
[math]\displaystyle{ \displaystyle X_i' := \lambda X_i', \; \lambda = \sqrt{\frac{m}{ \sum_{i=1}^{m} \| X_i \|_F^2}}, \; i=1,\dots,m. }[/math]
Задаем начальную среднюю форму, выбирая одну из центрированных форм входного набора (например, первую): [math]\displaystyle{ Y := X_1' }[/math].
Шаг 2. Инициализация функции ошибки.
Выполнить подгонку каждой формы (без учета первой - так как выбрали ее за изначальный консенсус) [math]\displaystyle{ X_i' }[/math] к [math]\displaystyle{ Y }[/math] и оценить функцию ошибки:
[math]\displaystyle{ \forall i=2,\dots,m \; X_i' := CPA(X_i', Y) }[/math];
[math]\displaystyle{ S_r := m(1 - tr(Y Y^T)) }[/math];
Также инициализируем [math]\displaystyle{ s_i := 1, \; i=1,\dots,m }[/math].
Шаг 3. Выравнивание форм относительно средней.
Для всех форм выполняется их выравнивание относительно средней формы (консенсуса):
[math]\displaystyle{ X_i' }[/math] присваивается новое значение: [math]\displaystyle{ X_i' := CPA(X_i', M^{(k)}), \; i=1,\dots,m }[/math].
, \; i=2,\dots,m;</math> [math]\displaystyle{ X_1' := \frac{X_1}{\left\|X_1\right\|_F} }[/math], где [math]\displaystyle{ \left\|X_1\right\|_F = \sqrt{tr(X_1 X_1^T)} }[/math] - норма Фробениуса. Нормировка выполняется, так как масштаб существенно влияет в том числе на расстояния между формами, а следовательно и на число итераций перед достижением критерия останова (как правило, меньше масштаб - меньше итераций, что нежелательно).
Для каждой формы [math]\displaystyle{ X_i }[/math] вычисляется её геометрический центр (среднее по строкам):
[math]\displaystyle{ \bar{X_i} = \frac{1}{n} \mathbf{1}_n^T X_i, \; i=1,\dots,m, }[/math]
[math]\displaystyle{ \sum{j=1}^{n} x_{ij}, }[/math] где [math]\displaystyle{ x_{ij} }[/math] — j-я точка формы [math]\displaystyle{ X_i }[/math] (строка матрицы).
Затем форма центрируется по этому среднему: [math]\displaystyle{ X_i \leftarrow X_i - \mathbf{1}_n \bar{x}_i^T, }[/math] где [math]\displaystyle{ \mathbf{1}_n \in \mathbb{R}^{n \times 1} }[/math] — вектор единиц.
[math]\displaystyle{ \bar{X}_i = \frac{1}{n} \mathbf{1}_n^T X_i, \quad i=1,\dots,m }[/math] [math]\displaystyle{ X_i^{(0)} = X_i - \mathbf{1}_n \bar{X}_i }[/math]
Программная реализация алгоритма
Литература
- ↑ Green B. F. The orthogonal approximation of an oblique structure in factor analysis // Psychometrika. — 1952. — Vol. 17, no. 4. — P. 429–440. [1]
- ↑ Schönemann P. H. A Generalized Solution of the Orthogonal Procrustes Problem // Psychometrika. — 1966. — Vol. 31, no. 1. — P. 1–10. [2]
- ↑ Schönemann P. H., Carroll R. M. Fitting one matrix to another under choice of a central dilation and a rigid motion // Psychometrika. — 1970. — Vol. 35, no. 2. — P. 245–255. [3]
- ↑ Gower J. C. Generalized procrustes analysis // Psychometrika. — 1975. — Vol. 40, no. 1. — P. 33–51. [4]