Перейдем к рассмотрению метода МГК (PCA). Вначале найдем константу \(\mu\), которая наилучшим образом описывает исходные данные \[ \varepsilon \left( \mu \right) = \sum\limits_{i = 1}^n {\left( {x_i - \mu } \right)^2 \to \mathop {\min }\limits_\mu }. \] Для нахождения минимума приравняем производную нулю и найдем значение\(\mu\)доставляющее минимум \[ \frac{d}{d\mu }\varepsilon \left( \mu \right) = - 2\sum\limits_{i = 1}^n {\left( {x_i - \mu } \right) = 0 \Rightarrow } \sum\limits_{i = 1}^n {\mu = } \sum\limits_{i = 1}^n {x_i \Rightarrow } \mu n = \sum\limits_{i = 1}^n {x_i \Rightarrow } \mu = \frac{1}{n}\sum\limits_{i = 1}^n {x_i } . \] Далее проведем центрирование данных, то есть, переопределим исходные значения \(x_{new} = x_{old} - \mu\), вычтя из каждого полученное значение \(\mu \).
Ясно, что новые данные имеют математическое ожидание равное нулю \[ E\left( {X - E\left( X \right)} \right) = E\left( X \right) - E\left( X\right) = 0. \] По сути дела, мы сделали параллельный перенос в существующей системе координат.
В дальнейшем будем считать, что исходные данные центрированы и мы хотим найти самое точное представление данных \(D = \left\{ {x_1 ,...,x_n } \right\}\) в некотором подпространстве \(W\), которое имеет размерность \(k \lt n\).
Чтобы минимизировать ошибку, нужно взять частные производные и учесть ортогональность\(\left\{ {e_1 ,...,e_k } \right\}\). Вначале упростим соотношение (1) \[ \varepsilon \left( {e_1 ,...,e_k ,\alpha _{1,1} ,...,\alpha _{n,k} } \right) = \sum\limits_{j = 1}^n {\left\| {x_j - \sum\limits_{i = 1}^k {\alpha _{j,i} } e_i } \right\|_2^2 = } \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - 2\sum\limits_{j = 1}^n {x_j^T \sum\limits_{i = 1}^k {\alpha _{j,i} } e_i + \sum\limits_{j = 1}^n {\sum\limits_{i = 1}^k {\alpha _{j,i}^2 } = } } } \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - 2\sum\limits_{j = 1}^n {\sum\limits_{i = 1}^k {\alpha _{j,i} x_j^T } e_i + \sum\limits_{j = 1}^n {\sum\limits_{i = 1}^k {\alpha _{j,i}^2 } .} } } \] Тогда \[ \frac{\partial }{\partial \alpha _{m,\ell} }\varepsilon \left( {e_1 ,...,e_k,\alpha _{1,1} ,...,\alpha _{n,k} } \right) = - 2x_m^T e_\ell + 2\alpha _{m,\ell}. \] Необходимое и достаточное условие экстремума будет иметь вид \[ - 2x_m^T e_\ell + 2\alpha _{m,\ell} = 0 \Rightarrow \alpha _{m,\ell} = x_m^T e_\ell . \] Таким образом, ошибка (1) примет вид \[ \varepsilon \left( {e_1 ,...,e_k } \right) = \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - 2\sum\limits_{j = 1}^n {\sum\limits_{i = 1}^k {\left( {x_j^T e_i } \right)x_j^T } e_i + \sum\limits_{j = 1}^n {\sum\limits_{i = 1}^k {\left( {x_j^T e_i } \right)^2} .} } } \] Упрощая это соотношение, получаем \begin{equation}\label{eq2.2} \varepsilon \left( {e_1 ,...,e_k } \right) = \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - \sum\limits_{j = 1}^n {\sum\limits_{i = 1}^k {\left( {x_j^T e_i } \right)^2} .} } \end{equation} Учитывая, что\(\left\langle {a,b} \right\rangle = a^Tb \) и \(\left\langle {b,a} \right\rangle = \left\langle {a,b} \right\rangle\), получаем \[ \left( {a^Tb} \right)^2 = \left( {a^Tb} \right)\left( {a^Tb} \right) = \left( {b^Ta} \right)\left( {a^Tb} \right) = b^T\left( {aa^T} \right)b, \] следовательно, \[ \varepsilon \left( {e_1 ,...,e_k } \right) = \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - \sum\limits_{i = 1}^k {e_i^T \left( {\sum\limits_{j = 1}^n {\left( {x_j x_j^T } \right)} } \right)e_i } } = \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - \sum\limits_{i = 1}^k {e_i^T Se_i } } , \] где \(S = \sum\limits_{j = 1}^n {\left( {x_j x_j^T } \right)}\) является ковариационной матрицей.
Следующим шагом будет минимизация \[\varepsilon \left( {e_1 ,...,e_k }\right) = \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - \sum\limits_{i = 1}^k {e_i^T Se_i } }\] при условии \(e_i^T e_i = 1\) для всех \(i.\) Используя метод неопределенных множителей Лагранжа (Lagrange), введем множители \(\lambda _1 ,...,\lambda _k\) и, замечая, что \(\sum\limits_{j =1}^n {\left\| {x_j } \right\|_2^2 } \equiv Const\), выпишем функцию цели \[ L\left( {e_1 ,...,e_k } \right) = \sum\limits_{i = 1}^k {e_i^T Se_i } - \sum\limits_{i = 1}^k {\lambda _i \left( {e_i^T e_i - 1} \right)} . \] Замечая, что \(\frac{d}{dX}\left( {X^TX} \right) = \frac{d}{dX}\left\langle {X,X} \right\rangle = 2X\) и если A симметрическая матрица, то \(\frac{d}{dX}\left( {X^TAX} \right) = 2AX\), получаем, \[ \frac{\partial }{\partial e_m }L \left( {e_1 ,...,e_k } \right) = 2Se_m - 2\lambda _m e_m = 0, \] то есть, \(Se_m = \lambda _m e_m \). Таким образом, необходимо найти решение уравнения \(\left( {S - \lambda I} \right)e = 0\) (здесь \(I\)-единичная матрица), что эквивалентно тому, что \(\lambda _m\)и\(e_m\) есть собственные значения и собственные векторы ковариационной матрицы S.
При этом ошибка приобретает вид \begin{equation}\label{eq2.3} \varepsilon \left( {e_1 ,...,e_k } \right) = \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - \sum\limits_{i = 1}^k {\lambda _i\left\| {e_i } \right\|_2^2 } } = \sum\limits_{j = 1}^n {\left\| {x_j } \right\|_2^2 - \sum\limits_{i = 1}^k {\lambda _i} } . \end{equation} Минимизация выписанной ошибки состоит в выборе базиса \(W\) из \(k\) собственных векторов матрицы \(S\), которые соответствуют \(k\) наибольшим собственным значениям. Большее собственное значение \(S\) дает большую вариацию в направлении соответственного собственного вектора.
Этот результат можно переформулировать следующим образом - проекция \(X\) на подпространство размерности \(k\), которая обеспечивает наибольшую вариацию. Таким образом МГК может трактоваться следующим образом - берем ортогональный базис и вращаем его пока на одном из направлений не получим максимальную вариацию. Фиксируем это направление и вращаем остальные, пока не найдем второе направление и так далее.
Пусть \(\left\{ {e_1 ,...,e_n } \right\}\) все собственные вектора матрицы \(S\), сортированные в порядке уменьшения соответственного собственного значения, тогда для любого \[ x_i = \sum\limits_{j = 1}^n {\alpha _{i,j} } e_j = \underbrace {\alpha _{i,1} e_1 + ... + \alpha _{i,k} e_k }_{approximation} + \overbrace {\alpha _{i,k + 1} e_{k + 1} + ... + \alpha _{i,n} e_n }^{error} \] коэффициенты \(\alpha _{m,\ell} = x_m^T e_\ell\) являются координатами главных компонент, и чем больше значение \(k\), тем получаем лучшую аппроксимацию. При этом главные компоненты располагаются в порядке значимости, более важные имеют меньший номер.

Описанный метод определения главных компонент является достаточно ресурсоемким и неустойчивым, особенно в случае, если собственные значения матрицы близки к нулю.
Более эффективным является использование итерационного метода определения главных компонент. Для этой цели рассмотрим задачу (1) с другой точки зрения.
Для случая \(i=1\) задача (1) сводится к определению одной компоненты, которая наилучшим образом восстанавливает все исходные данные \(\left\{ {x_1,...,x_n } \right\}\) \begin{equation} \label{eq2.4} \varepsilon \left( {e_1 ,\alpha _{1,1} ,...,\alpha _{n,1} } \right) = \sum\limits_{j = 1}^n {\left\| {x_j - \alpha _{j,1} e_1 } \right\|_2^2 } \to \min \end{equation} по всем \(e_1\) и \(\left\{ {\alpha _{i,1} } \right\}_{i = 1}^n\) при условии \(\|e_1\|=1\).
Если \(\left\{ {\tilde {\alpha }_{i,1} } \right\}_{i = 1}^n\) и \(\tilde {e}_1\) есть решение этой задачи и \(\Delta x_j = x_j - \tilde {\alpha }_{j,1}\tilde {e}_1\)- ошибка восстановления данных одной первой главной компонентой, то решая задачу \[ \sum\limits_{j = 1}^n {\left\| {\Delta x_j - \alpha _{j,2} e_2 } \right\|_2^2 } \to \min \] по всем \(e_2\) и \(\left\{ {\alpha _{i,2} } \right\}_{i = 1}^n\) при условии \(\|e_2\|=1\), получаем вторую главную компоненту \(\left\{ {\tilde {\alpha }_{i,2} } \right\}_{i = 1}^n\) и соответствующий вектор \(\tilde {e}_2\) и т.д.
При фиксированных или \(\left\{ {\alpha _{i,1} } \right\}_{i = 1}^n\) или \(e_1\) задача (1) решается методом наименьших квадратов. В силу того, что функция цели представляет собой квадратичный функционал, необходимое и достаточное условия экстремума совпадают. Таким образом, решение задачи сводится к поиску решения уравнения \[ \frac{\partial }{\partial e_1 }\varepsilon \left( {e_1 ,\alpha _{1,1} ,...,\alpha _{n,1} } \right) = - 2\sum\limits_{j = 1}^n {\left( {x_j - \alpha _{j,1} e_1 } \right)\alpha _{j,1} } = - 2\left( {\sum\limits_{j =1}^n {x_j \alpha _{j,1} } - \sum\limits_{j = 1}^n {\alpha _{j,1}^2 e_1 } } \right). \] Отсюда получаем \[ \bar{e}_1 = \frac{\sum\limits_{j = 1}^n {x_j \alpha _{j,1} } }{\sum\limits_{j =1}^n {\alpha _{j,1}^2 } }, \] учитывая условие нормирования единицей, имеем \[ e_1 = \frac{\bar{e}_1}{\|\bar{e}_1\|}. \] Следующий шаг будем делать исходя из предположения, что в задаче (1) нам известен вектор \(e_1\) и требуется найти экстремум по \(\left\{ {\alpha_{i,1} } \right\}_{i = 1}^n\) \[ \frac{\partial }{\partial \alpha _{\nu ,1} }\varepsilon \left( {e_1 ,\alpha_{1,1} ,...,\alpha _{n,1} } \right) = - 2\left( {x_\nu - \alpha _{\nu ,1} e_1 } \right)e_1 = - 2\left( {\left\langle {x_\nu ,e_1 } \right\rangle - \alpha _{\nu ,1} \left\langle {e_1 ,e_1 } \right\rangle } \right) = 0, \] то есть \[ \alpha _{\nu ,1} = \frac{\left\langle {x_\nu ,e_1 } \right\rangle }{\left\langle {e_1 ,e_1 } \right\rangle }= \left\langle {x_\nu ,e_1 } \right\rangle (\left\langle {e_1 ,e_1 } \right\rangle =1), \] где, как обычно, \(\left\langle {x,y} \right\rangle\)- скалярное произведение векторов \(x\) и \(y\).
Далее, считая найденные \(\left\{ {\alpha _{i,1} } \right\}_{i = 1}^n\) известными, повторяем весь процесс, пока не произойдет стабилизация ошибки.
Полученные \(\left\{ {\tilde{\alpha}_{i,1} } \right\}_{i = 1}^n\) будем считать первой главной компонентой.
Тогда \(\Delta x_j = x_j - \tilde {\alpha }_{j,1} \tilde {e}_1\)- ошибка восстановления данных одной первой главной компонентой.
Применяя этот алгоритм к ошибке восстановления \(\Delta x_j\), находим вторую главную компоненту \(\left\{ {\tilde{\alpha}_{i,2} } \right\}_{i = 1}^n\) вместе с коэффициентами \(e_2\) и т.д.

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

Область

Потребление юридическими потребителями млн.кВт.ч.

Общее потребление бытовыми потребителями млн.кВт.ч.

Потребление городскими бытовыми потребителями млн.кВт.ч.

Потребление сельскими бытовыми потребителями млн.кВт.ч.

Потребление бытовыми потребителями с электроплитами, электрообогревом в городской и сельской местности млн.кВт.ч.

---

Вначале проведем нормирование данных.

Для этого вычислим для каждого столбца данных минимальное \(x^{min}_j\), максимальное \(x^{max}_j\) и среднее \(x^{av}_j\) значения.

Потребление юридическими потребителями млн.кВт.ч.

Общее потребление бытовыми потребителями млн.кВт.ч.

Потребление городскими бытовыми потребителями млн.кВт.ч.

Потребление сельскими бытовыми потребителями млн.кВт.ч.

Потребление бытовыми потребителями с электроплитами, электрообогревом в городской и сельской местности млн.кВт.ч.

Далее нормируем исходные данные, для того, чтобы они имели равноправные значения \[ \tilde{x}_{i,j}=\frac{x_{i,j}-x^{av}_j}{x^{max}_j-x^{min}_j}, i=0,1,...,n-1. \]

Область

Потребление юридическими потребителями млн.кВт.ч.

Общее потребление бытовыми потребителями млн.кВт.ч.

Потребление городскими бытовыми потребителями млн.кВт.ч.

Потребление сельскими бытовыми потребителями млн.кВт.ч.

Потребление бытовыми потребителями с электроплитами, электрообогревом в городской и сельской местности млн.кВт.ч.

Возьмем стартовые значения углов направлений первой главной компоненты \(e^1_j=1/\sqrt{5}, j=1,2,3,4,5.\) и найдем первое приближение главной компоненты \[ \alpha^1_{1,i}=\sum_{j=1}^5e^1_j\tilde{x}_{i,j}, i=0,1,...,n-1. \] Далее найдем \[ \bar{e}_j=\frac{\sum_{i=0}^{n-1}\tilde{x}_{i,j}\alpha^1_{1,i}}{\sum_{i=0}^{n-1}\left(\alpha^1_{1,i}\right)^2} \] и нормируем единицей \(e^2_j=\bar{e}_j/\|\bar{e}_j\|=\) , после чего найдем второе приближение главной компоненты \[ \alpha^2_{1,i}=\sum_{j=1}^5e^2_j\tilde{x}_{i,j}, i=0,1,...,n-1. \] После десяти итераций получим Используя полученные значения, можно поставить каждой области в соответствие точку на оси, получая, тем самым рейтинг области. Для получения более полной информации, найдем вторую главную компоненту. Для этого найдем восстановление \(\tilde{x}_{i,j}\approx e^{10}_j\alpha^{10}_{1,i}\), и вычислим ошибку восстановления \[\Delta\tilde{x}_{i,j}=\tilde{x}_{i,j}- e^{10}_j\alpha^{10}_{1,i}\]

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

В результате получим следующую таблицу.

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

И, наконец, найдем три главные компоненты

Визуализируем полученные данные -

Цвет является результатом восприятия в сетчатке глаза светового потока в спектре, имеющего длину волны от 400 nm до 700 nm. Сетчатка человеческого глаза имеет три типа рецепторов, каждый из которых отвечает за восприятия света определенной части спектра. Четвертый тип рецепторов, присутствующих в сетчатке - палочки, которые эффективны только на чрезвычайно низких светлых уровнях (образно говоря, ночное зрение). Хотя для зрения они достаточно важны, никакой роли для восприятия цвета они не играют. Поскольку для цветового восприятия существуют ровно три типа рецепторов, то именно три цветовых компонента необходимы для описания цвета. Из всех существующих цветовых моделей наиболее популярны следующие представления цветов (цветовые пространства)

• RGB - модель смешивания цветов (базовая компьютерная модель),
• CMYK - модель вычитания цветов (используется при цветной печати),
• YIQ, YUV, YCrCb (используются в видео - системах).

Наиболее распространенной и часто используемой в компьютерной графике моделью является модель смешивания трех базовых цветов - красного, зеленого и синего. Смешивание базовых цветов в разных пропорциях позволяет получить всю цветовую палитру. Три цветовых компоненты можно рассматривать как три линейно независимых вектора, описывающих трехмерное цветовое пространство или цветовой куб. Если размер каждой компоненты не более 255 (то есть лежит в одной байте), то точка RGB (0, 0, 0) представляет собой черный цвет, белому цвету соответствует точка RGB (255, 255, 255), а сама диагональ соответствует значениям серого от черного до белого.
Тот факт, что цветовые компоненты независимы, говорит о том, что внесение помех в одну из компонент никах не скажется на остальных. Оно вроде бы и неплохо. Но. Человеческий глаз воспринимает эти цвета по разному, зеленый получше, синий похуже, поэтому если пожертвовать качеством синего, то на человеческом восприятии это скажется очень мало.
В ряде задач, связанных с хранением и передачей изображений (и видео в том числе), независимость цветовых компонент R,G,B является не преимуществом, а недостатком. В этом случае применяют, так называемые неравноправные цветовые компоненты (подробнее об этом можно посмотреть здесь), среди которых цветовая модель YCrCb, используемвая в JPEG. \[ \left\{ \begin{array}{lrl} Y=&\frac{1}{23}\left(7R+14G+2B\right), \\ Cr=&\frac{4}{69}\left(8R-7G-B\right),\\ Cb=&-\frac{4}{23}\left(R+2G-3B\right),\\ \end{array} \right. \Longleftrightarrow \left\{ \begin{array}{ll} R&=Y+\frac{3}{2}Cr, \\ G&=Y-\frac{1}{4}(3Cr+Cb),\\ B&=Y+\frac{7}{4}Cb. \end{array} \right. \] Y- люминисцентная составляющая (освещенность), Cr- красная хроматическая составляющая (теплые тона), Cb - синяя хроматическая составляющая (холодные тона).
Ценность использования неравноправных компонент состоит в том, что вклад каждой из них в формирование изображений, разный. Этот вклад можно оценить величиной нормы каждого из частотных доменов. Чем больше норма домена, тем больший вклад он делает в формирование изображения.
Можно подходить к проблеме построения неравноправного цветового пространства с другой точки зрения, исходя из максимальной информативности каждой компоненты. Применим метод главных компонент для получения неравноправной трехкомпонентной модели оптимальной с точки зрения минимизации среднеквадратичной ошибки восстановления исходного изображения. То есть, первая их полученных цветовых компонент будет нести наибольшую информацию об изображении среди всех полученных цветовых компонент, а вторая будет содержать наибольшую информацию среди оставшихся.
Таким образом, в нашей терминологии, задача (1) примет вид \[ \left\| {R - \sum\limits_{i = 1}^3 {\alpha _{i} } e_{r,i} } \right\|_2^2 + \left\| {G - \sum\limits_{i = 1}^3 {\alpha _{i} } e_{g,i} } \right\|_2^2 + \left\| {B - \sum\limits_{i = 1}^3 {\alpha _{i} } e_{b,i} } \right\|_2^2 \to \min , \] где минимум берется по всем \(\alpha _{i}\) и \(e_{color,i} \quad i=1,2,3, color=\{r,g,b\}.\)
Применяя метод главных компонент, получаем \[ \left( {{\begin{array}{*{20}c} e_1 \hfill \\ e_2 \hfill \\ e_3 \hfill \\ \end{array} }} \right) = \left( {\begin{array}{rrr} 0.767785 &0.45439 &0.4517034 \\ - 0.6164395 &0.716085 &0.3274513 \\ - 0.174667 &- 0.5298553& 0.8299064 \\ \end{array}} \right) \] Здесь в первом столбце матрицы стоят коэффициенты \(e_{r,1} , e_{g,1} , e_{b,1}\), во втором -\(e_{r,2} , e_{g,2} , e_{b,2}\) и в третьем -\(e_{r,3} ,e_{g,3} ,e_{b,3}\). Тогда восстановление тестового изображения по одной компоненте можно записать в виде

R_i,j=0.767785 × α_i,j, G_i,j=0.45439 × α_i,j, B_i,j=0.4517034 × α_i,j,

где α_i,j- значения первой главной компоненты \(\alpha_1\), соответствующие пикселю с координатами \((i,j)\). О вкладе разных главных компонент в формирование изображения можно судить по нормам частотных доменов.

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

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

Для получения синтезированных псевдоцветных изображений используются методы сокращения размерности сигнала, одним из которых является колоризация синтезированного изображения с учетом естественного разбиения сигнала на инфракрасный, видимый и ультрафиолетовый диапазоны, а также связь между ними. В основе предлагаемой технологии лежит модифицированный метод главных компонент. Пусть \(X_i (i=1,\ldots ,N)\) многозональный сигнал, причем значения \(X_{i} (i=1,\ldots ,N_1)\) соответствуют инфракрасному диапазону, \(X_{i} (i=N_{1},\ldots ,N_{2})\) - видимому и \(X_{i} (i=N_{2},\ldots ,N)\) ультрафиолетовому.
Метод главных компонент (МГК) позволяет выделить часть всего сигнала, наиболее характерно отражающего поведение всего сигнала в целом, что реализуется путем решения экстремальной задачи \[ \min\left\{\left.\sum_{k=0}^{n-1}\left\|X^k-\sum_{\mu=0}^{m-1}\alpha^\mu_kY^\mu\right\|_2 \right| Y^\mu,\alpha^\mu_k: \sum_{k=0}^{n-1}\left(\alpha^\mu_k\right)^2 =1 \right\}. \] Для решения этой задачи используем итерационный алгоритм.
Традиционно для колоризации синтезированного изображения используется пространство цветов YCrCb, где Y - величина освещенности, Cr - составляющая теплых тонов и Cb - составляющая холодных тонов. При этом в качестве Y берется первая главная компонента, в качестве Cr - вторая и вместо Cb используется третья главная компонента. Этот подход кажется достаточно надуманным и не всегда соответствующим реальному восприятию. Наши дальнейшие рассуждения посвящены построению естественной колоризации синтезированного изображения.
Пусть векторы \(\alpha^{i }(i=1,\ldots ,N_{1})\), соответствующие главной компоненте \(Y^R\) кадров \(X_{i} (i=1,\ldots ,N_{1})\), векторы \(\alpha^{i }(i= N_{1},\ldots ,N_{2})\), соответствующие главной компоненте \(Y^G\) кадров \(X_{i} (i= N_{1},\ldots ,N_{2})\) и \(\alpha^{i }(i= N_{2},\ldots ,N)\) векторы, соответствующие \(Y^B\) кадров \(X_{i} (i= N_{2},\ldots ,N).\)
На следующем шаге найдем векторы \(\alpha^R,\alpha^G,\alpha^B\), доставляющие решение задачи \begin{equation} \label{eq_multi2} \sum\limits_{k = 1}^N{\left({X_k - \alpha_k^R Y^R - \alpha_k^G Y^G - \alpha_k^B Y^B} \right)^2 \to \mathop{\min }\limits_{\alpha_k^R ,\alpha _k^G ,\alpha_k^B } }. \end{equation} Эта задача сводится к системе из \(3N\) линейных уравнений с \(3N\) неизвестными
Решая ее, получаем векторы \[\tilde{\alpha }^R = \left({\tilde{\alpha }_1^R ,\tilde{\alpha }_2^R ,...,\tilde{\alpha }_N^R } \right),\tilde {\alpha }^G = \left({\tilde{\alpha }_1^G ,\tilde{\alpha }_2^G ,...,\tilde {\alpha }_N^G } \right)$ и $\tilde{\alpha }^B = \left({\tilde{\alpha }_1^B ,\tilde{\alpha }_2^B ,...,\tilde{\alpha }_N^B } \right).\] Используя найденные векторы, проведем последовательное уточнение главных компонент. Вначале рассмотрим задачу \[ \sum\limits_{k = 1}^N{\left({X_k - \tilde{\alpha }_k^R Y^R - \tilde {\alpha }_k^G Y^G - \tilde{\alpha }_k^B Y^B} \right)^2 \to \mathop{\min }\limits_{Y^R} } , \] решением которой будет \[ \tilde{Y}^R = \frac{\sum\limits_{k = 1}^N{\tilde{\alpha }_k^R \left({X_k - \tilde{\alpha }_k^G Y^G - \tilde{\alpha }_k^B Y^B} \right)} }{\sum\limits_{k = 1}^N{\left({\tilde{\alpha }_k^R } \right)^2} }. \] Следующим шагом будет уточнение \(Y^{G}\), для чего рассмотрим следующую экстремальную задачу \[ \sum\limits_{k = 1}^N{\left({X_k - \tilde{\alpha }_k^R \tilde{Y}^R - \tilde{\alpha }_k^G Y^G - \tilde{\alpha }_k^B Y^B} \right)^2 \to \mathop {\min }\limits_{Y^G} } , \] решением которой будет \[ \tilde{Y}^G = \frac{\sum\limits_{k = 1}^N{\tilde{\alpha }_k^G \left({X_k - \tilde{\alpha }_k^R \tilde{Y}^R - \tilde{\alpha }_k^B Y^B} \right)} }{\sum\limits_{k = 1}^N{\left({\tilde{\alpha }_k^G } \right)^2} } \] и, наконец, уточним последнюю компоненту \[ \tilde{Y}^B = \frac{\sum\limits_{k = 1}^N{\tilde{\alpha }_k^B \left({X_k - \tilde{\alpha }_k^R \tilde{Y}^R - \tilde{\alpha }_k^G \tilde{Y}^G} \right)} }{\sum\limits_{k = 1}^N{\left({\tilde{\alpha }_k^B } \right)^2} }. \] Далее найдем нормированные \(\alpha ^R,\alpha ^G,\alpha ^B\) векторы \[ \alpha_i^R = \frac{\tilde{\alpha }_i^R }{\sum\limits_{k = 1}^N{\left( {\tilde{\alpha }_k^R } \right)^2} }, \quad \alpha_i^G = \frac{\tilde{\alpha }_i^G }{\sum\limits_{k = 1}^N{\left( {\tilde{\alpha }_k^G } \right)^2} }, \quad \alpha_i^B = \frac{\tilde{\alpha }_i^B }{\sum\limits_{k = 1}^N{\left( {\tilde{\alpha }_k^B } \right)^2} } \] и полагаем \(Y^R = \tilde{Y}^R,Y^G = \tilde{Y}^G,Y^B = \tilde{Y}^B.\)
Следующим этапом переходим к задаче (\ref{eq_multi2}), повторяя этот процесс итерационно, пока компоненты \(Y^R,Y^G,Y^B\) не стабилизируются.
Результатом выполнения алгоритма будут три главные компоненты \(Y^R,Y^G,Y^B\), а также три вектора \(\alpha^R,\alpha^G,\alpha^B\).
Для получения колоризованного синтезированного изображения, возьмем первую главную компоненту в качестве носителя красного цвета, вторую - зеленого и третью - синего цвета. Эта последовательность цветов определена тем, что первые \(N_{1}\) были получены в инфракрасном спектре, а последние \(N-N_{2}\) в ультрафиолетовом.
Заметим, что по полученным компонентам легко восстановить (с некоторой точностью) все изображения исходного многоканального сигнала: \[X_k \approx \tilde{X}_k = \alpha_k^R Y^R + \alpha_k^G Y^G + \alpha _k^B Y^B (k=1,2,\ldots ,N).\] Естественным является проблема улучшения качества восстановления многозонального изображения.
Примером решения этой проблемы может служить нелинейная конструкция восстановления многозонального изображения \[ \sum\limits_{k = 1}^N{\left({X_k - \alpha_k^R Y^R - \alpha_k^G Y^G - \alpha_k^B Y^B - \beta_k^R \left({Y^R} \right)^2 - \beta_k^G \left( {Y^G} \right)^2 - \beta_k^B \left({Y^B} \right)^2} \right)^2 \to \mathop {\min }\limits} \] На следующей диаграмме приведены значения ошибки восстановления по каждому кадру многозонального изображения. Использование квадратичной функции для восстановления сигнала усложняется тем фактом, что нужно анализировать на какой ветви параболы лежит искомое решение. Для решения этой проблемы используем метод улучшения восстановления, использующий линейную комбинацию главных компонент и компонент, получаемых применением оператора дифференцирования (дискретного аналога). Для этого найдем свертку каждой найденной главной компоненты с разностным фильтром \[ D = \left\{{d_{i,j} } \right\}_{i,j = - 1}^1 = \left( {{\begin{array}{*{20}r} 0 \hfill & 1 \hfill & 0 \hfill \\ 1 \hfill &{ - 4} \hfill & 1 \hfill \\ 0 \hfill & 1 \hfill & 0 \hfill \\ \end{array} }} \right), \] то есть, если \(T^\chi = \left\{{t_{i,j}^\chi } \right\},\chi \in \left\{ {R,G,B} \right\}\), то \[ DT^\chi = \left\{{dt_{i,j}^\chi } \right\} = \left\{{\sum\limits_{\nu ,\mu = - 1}^1{d_{\nu ,\mu } t_{i + \nu ,j + \mu }^\chi } } \right\}. \] Восстановление кадров многозонального изображения будем искать в виде решения оптимизационной задачи \begin{equation} \label{eq_multi3} \sum\limits_{k = 1}^N{\left({X_k - \alpha_k^R Y^R - \alpha_k^G Y^G - \alpha_k^B Y^B - \beta_k^R DY^R - \beta_k^G DY^G - \beta_k^B DY^B} \right)^2 \to \mathop{\min }\limits_{\beta_k^R ,\beta_k^G ,\beta_k^B } }. \end{equation} Понятно, что не все кадры многозонального изображения восстанавливаются с одним и тем же качеством - в случае наличия на кадрах оригинального сигнала уникальных отклонений, например, теплового шума в кадрах инфракрасного диапазона, ошибка восстановления будет большей.
Этот метод обработки многозонального сигнала дает возможность обрабатывать не все пикселы исходного сигнала, а избирательно, что позволяет фильтровать помехи, порождаемые, например, тепловым шумом. Этот факт позволяет получать более точные спектральные характеристики сигнала.
На следующих рисунках приведены спектрограммы точки многозонального изображения по оригинальным (желтым цветом) и данным, восстановленным по трем компонентам \(Y^R,Y^G, Y^B\)(синим цветом). Во втором случае взята точка с импульсным тепловым шумом, что существенно искажает спектрограмму. Поэтому целесообразно проводить предварительную обработку кадров. Приведем схему предварительной обработки. После получения восстановленных кадров \(\tilde{X}_k= \left\{{\tilde{x}_i } \right\}_{i = 1}^M \) \((k=1,2,\ldots ,N)\), находим среднеквадратичную и интегральную ошибки восстановления оригинальных кадров \(X_k= \left\{{x_i } \right\}_{i = 1}^M \) \[ \varepsilon_k^1 = \frac{1}{M}\sum\limits_{i = 1}^M{\left|{x_i - \tilde {x}_i } \right|} , \quad \varepsilon_k^2 = \sqrt{\frac{1}{M}\sum\limits_{i = 1}^M{\left|{x_i - \tilde{x}_i } \right|^2} } . \] Если для кадра величина \(\varepsilon_k^2/\varepsilon_k^1\) больше заданного порога, то этот кадр содержит пикселы с большими значениями шума. Эти пикселы изымаются из изображения и подлежат специальному анализу. Месторасположения вынутых точек заполняются средними значениями точек окрестности, получаемыми, например, в результате решения стационарной задачи теплопроводности.
К исправленному таким образом многозональному изображению применим весь алгоритм, описанный выше, что позволит получить более корректное восстановление спектрограмм.
Заметим, этот метод восстановления многозонального изображения не является единственным, можно использовать различные комбинации полученных главных компонент.
Использование модификаций метода главных компонент позволяет повысить качество восстанавливаемых кадров многозонального изображения по синтезированному изображению с возможностью фильтрации спектрограмм от импульсного шума.