И другие программы этой серии
SVD/PCA
Наиболее популярным способом сжатия данных в многомерном анализе является метод главных компонент (PCA). С математической точки зрения PCA — это декомпозиция исходной матрицы X, т.е. представление ее в виде произведения двух матриц T и P
X = TPt + E
Матрица T называется матрицей счетов (scores) , матрица P — матрицей нагрузок (loadings), а E — матрицей остатков.
Простейший способ найти матрицы T и P — использовать SVD разложение через стандартную функцию MatLab, называемую svd. function [T, P] = pcasvd(X)
% pcasvd: calculates PCA components.
% The output matrices are T and P.
% T contains scores
% P contains loadings
[U,D,V] = svd(X);
T = U * D;
P = V;
%end of pcasvd
Содержание
5.3 PCA/NIPALS
Для построения PCA счетов и нагрузок, используется рекуррентный алгоритм NIPALS, который на каждом шагу вычисляет одну компоненту. Сначала исходная матрица X преобразуется (как минимум – центрируется; см. раздел 4.3) и превращается в матрицу E0, a=0. Далее применяют следующий алгоритм.
1. Выбрать начальный вектор t
2. pt = tt Ea / ttt
3. p = p / (ptp)½
4. t = Ea p / ptp
5. Проверить сходимость, если нет, то идти на 2
После вычисления очередной (a-ой) компоненты, полагаем ta=t и pa=p. Для получения следующей компоненты надо вычислить остатки Ea+1 = Ea – t pt и применить к ним тот же алгоритм, заменив индекс a на a+1.
Код алгоритма NIPALS может быть написан и самими читателями, в данном же пособии авторы приводят свой вариант. При расчете PCA, можно вводить число главных компонент (переменная numberPC). Если же не известно, сколько необходимо компонент, следует написать в командной строке [P,T] = pcanipals (X) и тогда программа задаст число компонент равным наименьшему из показателей размерности исходной матрицы X. function [T, P] = pcanipals(X, numberPC)
% pcanipals: calculates PCA components.
% The output matrices are T and P.
% T contains scores
% P contains loadings
% calculation of number of components
[X_r, X_c] = size(X); P=[]; T=[];
if lenfth(numberPC) > 0
pc = numberPC{1};
elseif (length(numberPC) == 0) & X_r < X_c
pc = X_r;
else
pc = X_c;
end;
% calculation of scores and loadings for each component
for k = 1:pc
P1 = rand(X_c, 1); T1 = X * P1; d0 = T1\'*T1;
P1 = (T1\' * X/(T1\' * T1))\'; P1 = P1/norm(P1); T1 = X * P1; d = T1\' * T1;
while d - d0 > 0.0001;
P1 = (T1\' * X/(T1\' * T1)); P1 = P1/norm(P1); T1 = X * P1; d0 = T1\'*T1;
P1 = (T1\' * X/(T1\' * T1)); P1 = P1/norm(P1); T1 = X * P1; d = T1\'*T1;
end
X = X - T1 * P1; P = cat(1, P, P1\'); T = [T,T1];
end
О вычислении PCA с помощью надстройки Chemometrics рассказано в пособии Проекционные методы в системе Excel.
Содержание
5.4 PLS1
Самым популярным способом для многомерной калибровки является метод проекции на латентные структуры (PLS).