

%%
function [a_tr, a_te]=find_lpc(pos,k_values,tr_all,te_all)

% pos: voxel_number by 3 matrix 
% k_values: the number of neighbor voxels ( p_values in the algorithm )
% tr_all, te_all: train and test data (N by voxel_number matrix, N:sample size
%%

p=k_values;

a_tr=[];
a_te=[];

for i=1:size(tr_all,2) % if size(tr_all)~=size(te_all), change here

    temp=1:size(pos,1);


    [neighbor_index] = find_nn(pos(temp~=i,:), pos(i,:), k_values);

    d=setdiff(temp,neighbor_index{1,1});


    x1=tr_all(:,i);
    x2=tr_all(:,neighbor_index{1,1});

    x=[x1 x2];
    [a{1,i},g] = lpc(x',p);
    a1=a{1,i};
    a_tr=[a_tr a1(:,2:end)];
    clear x;

    d=setdiff(temp,neighbor_index{1,1});

    x1=te_all(:,i);
    x2=te_all(:,neighbor_index{1,1});

    x=[x1 x2];
    [a{1,i},g] = lpc(x',p);

    a1=a{1,i};

    a_te=[a_te a1(:,2:end)];

    clear x;
end
   