ИСПОЛЬЗОВАНИЕ МЕТОДА ГЛАВНЫХ КОМПОНЕНТ В MATLAB: ЭМПИРИЧЕСКОЕ МОДЕЛИРОВАНИЕ, СЖАТИЕ ДАННЫХ, ШУМОПОДАВЛЕНИЕ
Автор: Will Dwinnell
Источник: http://matlabdatamining.blogspot.com/2010/02/putting-pca-to-work.html
В последнем выпуске веб журнала, посвященному методу главных компонент (февраль-26-2010), был сделан обзор метода главных компонент (PCA), и каким образом его можно использовать в MATLAB. Эта статья посвящена использованию метода в трех направлениях:
1. Предварительная обработки временного ряда для эмпирического моделирования;
2. Сжатие данных;
3. Подавления шума.
Ссылаясь на более широкую аудиторию, в данной статье будем использовать для метода PCA только базовые функций MATLAB, и статистический инструмента Toolbox, о котором шла речь в последние публикации, в особенности возможность использования таких инструментов, как princomp и zscore.
Снова будем использовать небольшое количество данных, так же как и в предыдущей статье:
>>A = [269.8 38.9 50.5
272.4 39.5 50.0
270.0 38.9 50.5
270.0 38.9 50.5
272.0 39.3 50.2
269.8 38.9 50.5
269.8 38.9 50.5
268.2 38.6 50.2
268.2 38.6 50.8
267.0 38.2 51.1
267.8 38.4 51.0
273.6 39.6 50.0
271.2 39.1 50.4
269.8 38.9 50.5
270.0 38.9 50.5
270.0 38.9 50.5
];
Рассчитываем параметры примера, и стандартизируем данные в таблице:
>> [n m] = size(A)
n =
15
m =
3
>> AMean = mean(A)
AMean =
269.9733 38.9067 50.4800
>> AStd = std(A)
AStd =
1.7854 0.3751 0.3144
>> B = (A - repmat(AMean,[n 1])) ./ repmat(AStd,[n 1])
B =
-0.0971 -0.0178 0.0636
1.3591 1.5820 -1.5266
0.0149 -0.0178 0.0636
1.1351 1.0487 -0.8905
-0.0971 -0.0178 0.0636
-0.0971 -0.0178 0.0636
-0.9932 -0.8177 -0.8905
-0.9932 -0.8177 1.0178
-1.6653 -1.8842 1.9719
-1.2173 -1.3509 1.6539
2.0312 1.8486 -1.5266
0.6870 0.5155 -0.2544
-0.0971 -0.0178 0.0636
0.0149 -0.0178 0.0636
0.0149 -0.0178 0.0636
Теперь, когда данные в отцентрированы со средним значением 0,0 и стандартным отклонение 1,0, осуществляем eigenanalysis примера ковариаций для определения коэффициенты матрицы, которая формирует основные компоненты:
>> [V D] = eig(cov(B))
V =
0.6505 0.4874 -0.5825
-0.7507 0.2963 -0.5904
-0.1152 0.8213 0.5587
D =
0.0066 0 0
0 0.1809 0
0 0 2.8125
Напомним, что MATLAB EIG функция задает информацию для основных компонент от последнего до первого, когда идет чтение столбца слева направо. Матрица V содержит линейные коэффициенты для основных компонентов. Диагональная матрица D содержит отклонений для основных компонентов. До этого времени, мы проводили анализ основных компонентов. Для того чтобы использовать метод главных компонент, нам необходимо знать, какая доля каждого основного компонента представляет в суммарной дисперсии. Мы можем сделать это путем извлечения и нормализации диагонали матрицы D (мы используем flipud потому что основные компоненты находятся в "обратном" порядке):
>> cumsum(flipud(diag(D))) / sum(diag(D))
ans =
0.9375
0.9978
1.0000
Столбце выше означает, что первый основной компонент содержит 93,75% от общей дисперсии исходных данных, первые две основные компоненты, содержат 99,78% и, конечно, все основные компоненты вместе имею все дисперсии (ровно столько, сколько в оригинале стандартизированных данных).
И наконец, для расчета основных компонентов, просто умножим стандартизированные данные на коэффициент матрицы:
>> PC = B * V
PC =
-0.0571 -0.0003 0.1026
-0.1277 -0.1226 -2.5786
0.0157 0.0543 0.0373
0.0536 0.1326 -1.7779
-0.0571 -0.0003 0.1026
-0.0571 -0.0003 0.1026
0.0704 -1.4579 0.5637
-0.1495 0.1095 1.6299
0.1041 0.2496 3.1841
0.0319 0.3647 2.4306
0.1093 0.2840 -3.1275
0.0892 0.2787 -0.8467
-0.0571 -0.0003 0.1026
0.0157 0.0543 0.0373
0.0157 0.0543 0.0373
Чтобы проверить конденсации дисперсий, вычислим отклонения примера:
>> var(PC)
ans =
0.0066 0.1809 2.8125
Опять же, отметим, что первая главная компонента появляется в последней колонке при использовании EIG функции MATLAB, а столбцы слева имеют все меньше и меньше дисперсию до последнего главного компонента, который хранится в первом столбце.
Применение: Предварительная обработка данных для эмпирического моделирования
Приложение PCA в этом случае есть весьма простым: рассчет основных компонент и выбор из них, а не из исходных данных, зависемостей для построения эмпирических моделей (регрессия, нейронные сети и т.д.). Преимущество этого в том, что PCA сжимает информацию в подмножество новых переменных, меньшие из их будут необходимы при построении модели. На самом деле, это не было бы неразумно просто шагнуть через первые из многих основных компонент для построения модели: Во-первых, используем только первый основной компонент, а затем пробуем первый и второй, затем первый, второй и третий и т.д. Хорошее преимущество этого подхода в том, что все основные компоненты некоррелированы друг с другом.
Как уже говорилось в предыдущей статье, это может или не может хорошо работать, по нескольким причинам: PCA не может быть способным выжать много отклонений, если исходные данные достаточно высоко некоррелированны друг с другом. Кроме того, статистическое отклонение не может быть подобным, как "информацию" для целей построения модели. В последнем случае, даже если этот процесс работает, остается один факт, что PCA нуждается в исходных переменных для расчета основных компонентов, даже если часть из них используется. В любом случае, это способ обработки данных, может принести пользу, поэтому этот случай следует попробовать.
Применение: Сжатие данных
PCA предлагает механизм для выполнения сжатия с потерями данных. Когда данные сжимаются "с потерями", это не может возвращать точные исходные данные. Компромисс в том, что гораздо большего сжатия можно достигнуть, чем "без потерь" сжатия данных (сжатие, в которых исходные данные возвращается точно). Во многих случаях, таких как аудио (MP3) и изображения (JPEG), некоторые потери в точности является приемлемой и большое сжатие есть очень желаемое.
Все схемы сжатия рассчитывают на обнаружение закономерностей в данных. В случае с PCA, регулярность есть линейная зависимость между переменными. С этой стороны, PCA считает, что эти зависимости(отношения) между данными могут быть сжаты. Идея заключается в том, чтобы отказаться от последних основных компонент (тех экспонированных наименьших дисперсии).
В MATLAB, это означает просто отбрасывание столбцов, которые представляют нежелательные основные компонентов. В этом случае, мы будем сохранять только первый основные компоненты:
>> VReduced = V(:,3)
VReduced =
-0.5825
-0.5904
0.5587
>> PCReduced = B * VReduced
PCReduced =
0.1026
-2.5786
0.0373
-1.7779
0.1026
0.1026
0.5637
1.6299
3.1841
2.4306
-3.1275
-0.8467
0.1026
0.0373
0.0373
Декомпрессия осуществляется путем обратного процесса, который мы можно сделать путем перестановки вектор коэффициентов и перемножения:
>> PCReduced * VReduced'
ans =
-0.0598 -0.0606 0.0573
1.5020 1.5224 -1.4406
-0.0217 -0.0220 0.0209
1.0356 1.0497 -0.9933
-0.0598 -0.0606 0.0573
-0.0598 -0.0606 0.0573
-0.3284 -0.3328 0.3150
-0.9494 -0.9623 0.9106
-1.8547 -1.8799 1.7789
-1.4158 -1.4351 1.3580
1.8217 1.8465 -1.7473
0.4932 0.4999 -0.4730
-0.0598 -0.0606 0.0573
-0.0217 -0.0220 0.0209
-0.0217 -0.0220 0.0209
В результате получаем не совсем то же, что и оригинальные стандартизированные данных, но очень близкие. Мы "не стандартизацируем", перевернув оригинальный шаг стандартизации:
>> Z = ((PCReduced * VReduced') .* repmat(AStd,[n 1])) + repmat(AMean,[n 1])
Z =
269.8667 38.8840 50.4980
272.6550 39.4777 50.0270
269.9345 38.8984 50.4866
271.8223 39.3004 50.1677
269.8667 38.8840 50.4980
269.8667 38.8840 50.4980
269.3870 38.7818 50.5790
268.2783 38.5457 50.7663
266.6619 38.2016 51.0393
267.4455 38.3684 50.9070
273.2259 39.5992 49.9306
270.8539 39.0942 50.3313
269.8667 38.8840 50.4980
269.9345 38.8984 50.4866
269.9345 38.8984 50.4866
Опять же, результат очень похож на оригинал, но примерно около 94% дисперсии была сохранена, и мы сжали данные до 33% своего первоначального размера.
Компромисс здесь есть между от сжатием (кол-во основных компонентов сохранились) и сжатия точности (разница сохранилась). В об щем применении, переменных будет больше и разница сжатия, как правило, будет не столь существенна, как в нашем примере. Это означает, что там будет больше урівень сжатия данных, в этом случае число основных компонентов сохраняется.
Применение: шумоподавление
Расширение применения сжатия данных, мы можем использовать метод главных компонент для подавления шума. Основная идея заключается в том, что дисперсия захваченная последними главными компонентами и есть шум который и может бать отвергнут. Предполагая, что переменные несут линейные зависимости, они будут лежать на линии (плоскости, гиперплоскости), а уровни шума будут поднимать их выше над линией. Удаление последних основных компонентов есть средством выравнивания данных в геометрическом смысле и (надеюсь) исключит сам шум.
Этот процесс очень похож на процесс сжатия данных, описанных в предыдущем разделе, за исключением: 1. отбрасываются компоненты имеющие свои коэффициенты равными нулю, вместо того чтобы удалить сразу, и 2. PCA матрицы коэффициентов и обратной перемножения, чтобы ни на шаг обработки, которые (опять же, будем надеяться) снижает уровень шума в данных.
Как и прежде, вычисляем коэффициенты PCA:
>> [V D] = eig(cov(B))
V =
0.6505 0.4874 -0.5825
-0.7507 0.2963 -0.5904
-0.1152 0.8213 0.5587
D =
0.0066 0 0
0 0.1809 0
0 0 2.8125
Решение о ликвидации последних основных компонент, мы устанавливаем его коэффициенты равные нулю:
>> VDenoise = V; VDenoise(:,1) = 0
VDenoise =
0 0.4874 -0.5825
0 0.2963 -0.5904
0 0.8213 0.5587
Эта матрица будет проецировать стандартизированные данные на плоскую поверхность, в данном случае плоскости, так как мы сохранили 2 два вектора. Не желая возится с двумя шагами, умножим эту матрицу на обратную, что в этом случае легко получить, взяв транспонирование:
>> VDenoise = VDenoise * VDenoise'
VDenoise =
0.5769 0.4883 0.0749
0.4883 0.4364 -0.0865
0.0749 -0.0865 0.9867
Эта волшебная матрица, в единстве матрица умножения, исключающая шумы стандартизированных данных:
>> B * VDenoise
ans =
-0.0599 -0.0607 0.0570
1.4422 1.4861 -1.5414
0.0047 -0.0060 0.0654
1.1002 1.0890 -0.8844
-0.0599 -0.0607 0.0570
-0.0599 -0.0607 0.0570
-1.0390 -0.7648 -0.8824
-0.8960 -0.9299 1.0005
-1.7330 -1.8060 1.9839
-1.2380 -1.3270 1.6575
1.9601 1.9307 -1.5141
0.6290 0.5825 -0.2442
-0.0599 -0.0607 0.0570
0.0047 -0.0060 0.0654
0.0047 -0.0060 0.0654
Естественно, что нам все еще нужно умножить назад стандартное отклонение и добавить обратно среднее, для того чтобы вернуть данные в оригинальной шкале:
>> Z = ((B * VDenoise) .* repmat(AStd,[n 1])) + repmat(AMean,[n 1])
Z =
269.8664 38.8839 50.4979
272.5483 39.4640 49.9954
269.9817 38.9044 50.5006
271.9377 39.3151 50.2019
269.8664 38.8839 50.4979
269.8664 38.8839 50.4979
268.1183 38.6198 50.2025
268.3736 38.5579 50.7946
266.8791 38.2293 51.1038
267.7630 38.4090 51.0012
273.4731 39.6308 50.0040
271.0964 39.1251 50.4032
269.8664 38.8839 50.4979
269.9817 38.9044 50.5006
269.9817 38.9044 50.5006
Уровень шума находится в зависимости от числа главных компонент: чем меньше основных компонент сохранилось, тем больше шума. Очевидно, как и все подобные схемы, этот процесс имеет свои ограничения и основное предположение в том, что первоначальные переменные линейно связаны так, что шум, выделяется как отклонение от этой линейности.
Итоги
PCA является мощным инструментом, и быстро вычисляется на современных компьютерах, даже для достаточно больших объемов данных. Хотя существуют пределы тому, что он может сделать, это удобный инструмент, который является недорогим с точки зрения вычислительно времени.
Дополнительная литература
Главная ссылка на PCA см.:
Multivariate Statistical Methods: A Primer, by Manly (ISBN: 0-412-28620-3)
Примечание: первое издание является достаточно простым для понимания и использования PCA, и в настоящее время гораздо дешевле, чем второе или третье издания.
The noise suppression application is described in the article, Vectors help make sense of multiple signals, by Sullivan, Personal Engineering and Instrumentation News (Dec-1997), in which it is referred to as subspace projection.