nprocess.m

多维数据处理：MATLAB源程序用于处理多维数据

function [Xnew,mX,sX]=nprocess(X,DimX,Cent,Scal,mX,sX,reverse,show);

% $ Version 1.03 $ Date 6. May 1998 $ Drastic error in finding scale parameters corrected $ Not compiled $
% $ Version 1.031 $ Date 25. January 2000 $ Error in scaling part $ Not compiled $
% $ Version 1.032 $ Date 28. January 2000 $ Minor bug$ Not compiled $
%
%
% Copyright, 1998 - 
% This M-file and the code in it belongs to the holder of the
% copyrights and is made public under the following constraints:
% It must not be changed or modified and code cannot be added.
% The file must be regarded as read-only. Furthermore, the
% code can not be made part of anything but the 'N-way Toolbox'.
% In case of doubt, contact the holder of the copyrights.
%
% Rasmus Bro
% Chemometrics Group, Food Technology
% Department of Food and Dairy Science
% Royal Veterinary and Agricultutal University
% Rolighedsvej 30, DK-1958 Frederiksberg, Denmark
% Phone  [phone]
% Fax    [phone]
% E-mail [e-mail]
%
%
% CENTERING AND SCALING OF N-WAY ARRAYS
%
% This m-file works in two ways
%    I.  Calculate center and scale parameters and preprocess data
%    II. Use given center and scale parameters for preprocessing data
% 
% %%% I. Calculate center and scale parameters %%%
% 
%     [Xnew,Means,Scales]=nprocess(X,DimX,Cent,Scal);
% 
%     INPUT
%     X       Data array
%     DimX    Size of X
%     Cent    is binary row vector with as many elements as DimX.
%             If Cent(i)=1 the centering across the i'th mode is performed
%             I.e cnt = [1 0 1] means centering across mode one and three.
%     Scal    is defined likewise. Scal(i)=1, means scaling to standard  
%             deviation one within the i'th mode
% 
%     OUTPUT
%     Xnew    The preprocessed data
%     mX      Sparse vector holding the mean-values 
%     sX      Sparse vector holding the scales
%
% %%% II. Use given center and scale parameters %%%
% 
%     Xnew=nprocess(X,DimX,Cent,Scal,mX,sX);
% 
%     INPUT
%     X       Data array
%     DimX    Size of X
%     Cent    is binary row vector with as many elements as DimX.
%             If Cent(i)=1 the centering across the i'th mode is performed
%             I.e Cent = [1 0 1] means centering across mode one and three.
%     Scal    is defined likewise. Scal(i)=1, means scaling to standard  
%             deviation one within the i'th mode
%     mX      Sparse vector holding the mean-values 
%     sX      Sparse vector holding the scales
%     reverse Optional input
%             if reverse = 1 normal preprocessing is performed (default)
%             if reverse = -1 inverse (post-)processing is performed
%
%     OUTPUT
%     Xnew    The preprocessed data
%
% For convenience this m-file does not use iterative 
% preprocessing, which is necessary for some combinations of scaling
% and centering. Instead the algorithm first standardizes the modes
% successively and afterwards centers. The prior standardization ensures
% that the individual variables are on similar scale (this might be slightly
% disturbed upon centering - unlike for two-way data).
%
%	Copyright
%	Rasmus Bro 1997
%	Denmark
%	E-mail [e-mail]

[I,J]=size(X);

if nargin<4
   error(' Four input arguments must be given')
end

MODE=1;
if nargin>5
  if issparse(mX)&issparse(sX)
    MODE=2;
  end
end

if ~(exist('show')==1)
  show=1;
end

if ~(exist('reverse')==1)
   reverse=1;
end

if ~any([1 -1]==reverse)
  error( 'The input <<reverse>> must be one or minus one')
end

if show~=-1
  if MODE==1
    disp(' Calculating mean and scale and processing data')
  else
    if reverse==1
      disp(' Using given mean and scale values for preprocessing data')
    elseif reverse==-1
      disp(' Using given mean and scale values for postprocessing data')
    end
  end
end

ord=chkpfdim(X,DimX,NaN);
DimX2=DimX;

if MODE==1
     out=0;
     m_idx=[];
     for i=1:ord
        out2=prod(DimX([1:i-1 i+1:ord]));
        out=out+out2;
        m_idx=[m_idx out2];
     end
     mX=sparse(out,1);
     sX=sparse(sum(DimX),1);      

     %Standardize
     for j=ord:-1:2
        X=nshape(X,DimX2,ord);
        if Scal(j)
          if show~=-1
            disp([' Scaling mode ',num2str(j)])
          end
          s=(stdnan(X')').^(-1);
          %s=(stdnan(X')');
          X=X.*(s*ones(1,size(X,2)));
          sX( sum(DimX(1:j-1))+1 : sum(DimX(1:j)) ) = s;
        end
        DimX2=DimX2([ord 1:ord-1]);
     end
     X=nshape(X,DimX2,ord);
     if Scal(1)==1
        if show~=-1
           disp([' Scaling mode ',num2str(1)])
        end
        s=(stdnan(X')').^(-1);
        %s=(stdnan(X')');
        X=X.*(s*ones(1,size(X,2)));
        sX(1:DimX(1))=s;
     end
   
     DimX2=DimX;
     %Center
     for j=ord:-1:2
        X=nshape(X,DimX2,ord);
        if Cent(j)==1
           if show~=-1
              disp([' Centering mode ',num2str(j)])
           end
           m=meannan(X);
           X=X-ones(DimX2(ord),1)*m;
           mX( sum(m_idx(1:j-1))+1 : sum(m_idx(1:j))) = m;
        end
        DimX2=DimX2([ord 1:ord-1]);
     end
    X=nshape(X,DimX2,ord);
     if Cent(1)==1
       if show~=-1
         disp([' Centering mode ',num2str(1)])
       end
       m=meannan(X);
       X=X-ones(DimX(1),1)*m;
       mX(1 : m_idx(1)) = m;
     end
 
  
else  % if MODE==2

     out=0;
     m_idx=[];
     for i=1:ord
        out2=prod(DimX([1:i-1 i+1:ord]));
        out=out+out2;
        m_idx=[m_idx out2];
     end


   if reverse==1
     %Standardize
     for j=ord:-1:2
        X=nshape(X,DimX2,ord);
        if Scal(j)
          if show~=-1
            disp([' Scaling in mode ',num2str(j)])
          end
          s=sX( sum(DimX(1:j-1))+1 : sum(DimX(1:j)) );
          X=X.*(s*ones(1,size(X,2)));
        end
        DimX2=DimX2([ord 1:ord-1]);
     end
     X=nshape(X,DimX2,ord);
     if Scal(1)==1
        if show~=-1
           disp([' Scaling in mode ',num2str(1)])
        end
        s=sX(1:DimX(1));
        X=X.*(s*ones(1,size(X,2)));
     end
   
     DimX2=DimX;
     %Center
     for j=ord:-1:2
        X=nshape(X,DimX2,ord);
        if Cent(j)==1
           if show~=-1
              disp([' Subtracting off-sets in mode ',num2str(j)])
           end
           m=mX( sum(m_idx(1:j-1))+1 : sum(m_idx(1:j)));
           X=X-ones(DimX2(ord),1)*m';
        end
        DimX2=DimX2([ord 1:ord-1]);
     end
     X=nshape(X,DimX2,ord);
     if Cent(1)==1
       if show~=-1
         disp([' Subtracting off-sets in mode ',num2str(1)])
       end
       m=mX(1 : m_idx(1));
       X=X-ones(DimX(1),1)*m';
     end

   elseif reverse==-1

     %Center
     % THIS IS VERY SUB-OPTIMAL WRT. WORKLOAD BUT THINGS WILL BE MUCH EASIER TO IMPLEMENT IN Ver. 5.2
     for j=1:ord
       if Cent(j)==1
         if show~=-1
           disp([' Adding off-sets in mode ',num2str(j)])
         end
         m=reverse*mX( sum(m_idx(1:j-1))+1 : sum(m_idx(1:j)));
         if j==1
           X=X-ones(DimX(j),1)*m';
         else
           for jj=ord:-1:1
             X=nshape(X,DimX([jj+1:ord 1:jj]),ord);
             if jj==j
               X=X-ones(DimX(j),1)*m';
             end
           end
         end
       end
     end

     %Scale
     % THIS IS VERY SUB-OPTIMAL WRT. WORKLOAD BUT THINGS WILL BE MUCH EASIER TO IMPLEMENT IN Ver. 5.2
     for j=1:ord
       if Scal(j)==1
         if show~=-1
           disp([' Rescaling back to original domain in mode ',num2str(j)])
         end
         s=sX( sum(DimX(1:j-1))+1 : sum(DimX(1:j)) ).^(-1);
         if j==1
           X=X.*(s*ones(1,size(X,2)));
         else
           for jj=ord:-1:1
             X=nshape(X,DimX([jj+1:ord 1:jj]),ord);
             if jj==j
                X=X.*(s*ones(1,size(X,2)));
             end
           end
         end
       end
     end


   end
end

Xnew=X;