PCA с использованием princomp в MATLAB (для распознавания лиц)

Я пытаюсь уменьшить размерность с помощью MATLAB princomp, но я не уверен, что делаю это правильно.

Вот мой код только для тестирования, но я не уверен, что правильно делаю проекцию:

A = rand(4,3)
AMean = mean(A)
[n m] = size(A)
Ac = (A - repmat(AMean,[n 1]))
pc = princomp(A)
k = 2; %Number of first principal components
A_pca = Ac * pc(1:k,:)'  %Not sure I'm doing projection right
reconstructedA = A_pca * pc(1:k,:)
error = reconstructedA- Ac

И мой код для распознавания лиц с использованием набора данных ORL:

%load orl_data 400x768 double matrix (400 images 768 features)
%make labels
orl_label = [];
for i = 1:40
    orl_label = [orl_label;ones(10,1)*i];
end

n = size(orl_data,1);
k = randperm(n);
s = round(0.25*n); %Take 25% for train

%Raw pixels
%Split on test and train sets
data_tr = orl_data(k(1:s),:);
label_tr = orl_label(k(1:s),:);
data_te = orl_data(k(s+1:end),:);
label_te = orl_label(k(s+1:end),:);

tic
[nn_ind, estimated_label] = EuclDistClassifier(data_tr,label_tr,data_te);
toc

rate = sum(estimated_label == label_te)/size(label_te,1)

%Using PCA
tic
pc = princomp(data_tr);
toc

mean_face = mean(data_tr);
pc_n = 100;
f_pc = pc(1:pc_n,:)';
data_pca_tr = (data_tr - repmat(mean_face, [s,1])) * f_pc;
data_pca_te = (data_te - repmat(mean_face, [n-s,1])) * f_pc;

tic
[nn_ind, estimated_label] = EuclDistClassifier(data_pca_tr,label_tr,data_pca_te);
toc

rate = sum(estimated_label == label_te)/size(label_te,1)

Если я выберу достаточно основных компонентов, это даст мне равные показатели признания. Если я использую небольшое количество основных компонентов (PCA), то скорость использования PCA будет хуже.

Вот несколько вопросов:

  1. Является ли функция princomp лучшим способом вычислить первые k главных компонентов с помощью MATLAB?
  2. Использование проектируемых функций PCA по сравнению с необработанными функциями не дает дополнительной точности, а дает только меньший размер векторных функций? (быстрее сравнивать векторы признаков).
  3. Как автоматически выбрать min k (количество главных компонентов), которые дают такую ​​же точность по сравнению с необработанным вектором признаков?
  4. Что делать, если у меня очень большой набор образцов, могу ли я использовать только часть из них с сопоставимой точностью? Или я могу вычислить PCA для некоторого набора, а затем «добавить» какой-то другой набор (я не хочу повторно вычислять pca для set1 + set2, но каким-то образом итеративно добавляю информацию из set2 в существующий PCA из set1)?

Я также пробовал версию с графическим процессором, просто используя gpuArray:

%Test using GPU
tic
A_cpu = rand(30000,32*24);
A = gpuArray(A_cpu);
AMean = mean(A);
[n m] = size(A)
pc = princomp(A);
k = 100;
A_pca = (A - repmat(AMean,[n 1])) * pc(1:k,:)';
A_pca_cpu = gather(A_pca);
toc
clear;

tic
A = rand(30000,32*24);
AMean = mean(A);
[n m] = size(A)
pc = princomp(A);
k = 100;
A_pca = (A - repmat(AMean,[n 1])) * pc(1:k,:)';
toc
clear;

Работает быстрее, но для больших матриц не подходит. Может я ошибаюсь?

Если я использую большую матрицу, это дает мне:

Ошибка при использовании gpuArray Недостаточно памяти на устройстве.


person mrgloom    schedule 13.04.2013    source источник


Ответы (1)


"Является ли функция princomp лучшим способом вычислить первые k главных компонентов с помощью MATLAB?"

Он вычисляет полный SVD, поэтому на больших наборах данных он будет медленным. Вы можете значительно ускорить это, указав необходимое количество измерений в начале и вычислив частичный svd. Функции Matlab для частичного svd - svds.

Если svds недостаточно быстр для вас, есть более современная реализация:

http://cims.nyu.edu/~tygert/software.html (matlab версия: http://code.google.com/p/framelet-mri/source/browse/pca.m)

(см. статью с описанием алгоритма http://cims.nyu.edu/~tygert/blanczos.pdf)

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

>> A = rand(40,30); %random rank-30 matrix
>> [U,S,V] = pca(A,2); %compute a rank-2 approximation to A
>> norm(A-U*S*V',2)/norm(A,2) %relative error               

ans =

    0.1636

>> [U,S,V] = pca(A,25); %compute a rank-25 approximation to A
>> norm(A-U*S*V',2)/norm(A,2) %relative error                 

ans =

    0.0410

Когда у вас есть большие данные и разреженная матрица, вычисление полного SVD часто невозможно, поскольку факторы никогда не будут разреженными. В этом случае вы должны вычислить частичный SVD, чтобы он уместился в памяти. Пример:

>> A = sprandn(5000,5000,10000);
>> tic;[U,S,V]=pca(A,2);toc;
no pivots
Elapsed time is 124.282113 seconds.
>> tic;[U,S,V]=svd(A);toc;   
??? Error using ==> svd
Use svds for sparse singular values and vectors.

>> tic;[U,S,V]=princomp(A);toc;
??? Error using ==> svd
Use svds for sparse singular values and vectors.

Error in ==> princomp at 86
    [U,sigma,coeff] = svd(x0,econFlag); % put in 1/sqrt(n-1) later

>> tic;pc=princomp(A);toc;     
??? Error using ==> eig
Use eigs for sparse eigenvalues and vectors.

Error in ==> princomp at 69
        [coeff,~] = eig(x0'*x0);
person dranxo    schedule 07.07.2013
comment
Какое потребление памяти у этих методов? - person mrgloom; 08.07.2013
comment
Для полного SVD на матрице MxN (то есть с использованием princomp или svd) вам нужно будет хранить плотные матрицы U и V, поэтому 2 * M N. Это недопустимо, когда входные данные большие (и поэтому хранятся в разреженной матрице). Использование svds или pca.m требует только сохранения k max (M, N), где k - количество необходимых вам измерений. Если ваши данные действительно большие, вы можете использовать реализацию PCA в Mahout (это просто реализация статьи, указанной в моем ответе) builds.apache.org/job/Mahout-Quality/javadoc/org/apache/mahout/ - person dranxo; 08.07.2013