Documentation of varmaxt

Global Index (short | long) | Local contents | Local Index (short | long)

Function Synopsis

[Frot,AT,Cscor,Vrot,h] = varmaxt(Fm,L,norm,A)

Help text

 USE: [Frot,AT,Cscor,Vrot,h] = varmaxt(Fm,L,norm,A)
 This is a program to perform a varimax rotation of factor 
 loadings.
 Input:
   Fm = MxN (L =< N <= M) input factor loadings
        of which an M by L subset is rotated (if eigenvectors are used
        remember to un-normalize them so that Fm'*Fm is a diagonal
        matrix with the eigenvalues on the diagonal.
   L  = number of PCs to use in rotation (the program uses the L
        largest PCs assuming the input loadings (eigenvectors)
        are arranged in descending order).
   norm = 'Y' if loadings are to be normalized by communalities.
   A  = Mx1 vector of area weights (if not specified it is set to 
	 a vector of ones.

 Output:
   Frot  = MxL rotated factor loadings.  
   AT    = LxL orthogonal rotation matrix.
   Cscor = MxL score coefficient matrix.
   Vrot  = Lx1 vector giving the fractional variance explained by
           each loading vector.

 To calculate scores do S=P(:,1:L)*AT (where P are the first L principal 
 components) or S=D'*diag(A)*Cscor, where D is the M by N data matrix.

Listing of function varmaxt

function [Frot,AT,Cscor,Vrot,h] = varmaxt(Fm,L,norm,A)
suma=sum(A);
fprintf('Performing Factor Rotation \n')
sqrt2 = 1./sqrt(2.);
L1 = L-1; nit = -1; ncm = 0;
%... Copy a subset of F into Fr:
[M,N]=size(Fm);
Frot = Fm(:,1:L);
AT = eye(L);
if (nargin == 3), A=ones(M,1); end
%... Calculate communalities and normalize loadings:
h=(sum((Frot.^2)'))';
fvar = sum(h .* A)
if norm == 'Y'
     Frot=diag(1 ./ sqrt(h))*Frot;
end
tvi = 0.;
% Iterate to achieve the Varimax criterion until conversion:
while ncm < 6 
%... Calculate the variamx criterion:
    for l=1:L
        sf1(l)=sum((Frot(:,l).^2).*A);
        sf2(l)=sum((Frot(:,l).^4).*A);
    end    
    tvf = sum(sf2*suma - sf1.^2)/(suma*suma);
    if nit >= 0 
       if abs(tvf-tvi) <= 0.000001
          ncm = ncm +1;
       end
    end
    nit = nit + 1;
    tvi = tvf;
    fprintf('Iteration No. %f Varimax criterion = %f \n',nit,tvf)
%... Rotate columns i and j:
    for i=1:L1
        L2 = i+1;
        for j=L2:L
            p = Frot(:,i);
            q = Frot(:,j);             
	    pa = AT(:,i);
            qa = AT(:,j);             
            u = (p+q) .* (p-q);
            v = 2*(p .* q);
            s = (u+v) .* (u-v);
            t = 2*(u .* v);
            a = sum(u .* A);
            b = sum(v .* A);
            c = sum(s .* A);
            d = sum(t .* A);
            xn = d - 2*(a*b)/suma;
            xo = c-(a*a-b*b)/suma;
            xr = sqrt(xn*xn + xo*xo);
            if xr > .001
                cos4t = xo/xr;
                cos2t = sqrt((1. + cos4t)/2.);
                cos1t = sqrt((1. + cos2t)/2.);
                sin1t = sqrt( 1. - cos1t^2  );
                if sin1t > .001
                    if xn < 0.
                        sin1t = -sin1t;
                    end
%... Rotate columns i and j:
			Frot(:,i) = cos1t*p + sin1t*q;
                    	Frot(:,j) = cos1t*q - sin1t*p;
%... Rotate AT in the same way
			AT(:,i) = cos1t*pa + sin1t*qa;
			AT(:,j) = cos1t*qa - sin1t*pa;
                end
            end  
        end     % Go to next j
    end         % Go to next i
end             % Go for the next iteration until convergence    
%... Un-normalize rotated loadings:
if norm == 'Y'
     Frot = diag(sqrt(h))*Frot;
end
%... Recalculate communalities:
for m=1:M
    Vrot = Frot(m,:).^2;
    h(m,1) = sum(Vrot');
end
fvar = sum(h .* A)
%... Calculate partial variance explained by each PC
for l=1:L
    p = Frot(:,l) .^2;
    q = p .* A;
    Vrot(l) = sum(q);
end
[pvar(L:-1:1),k] = sort(Vrot); 
Frot(:,L:-1:1) = Frot(:,k);
AT(:,L:-1:1) =AT(:,k);
for l=1:L
    p = Frot(:,l) .^2;
    q = p .* A;
    Vrot(l) = sum(q);
end
Vrot
%... Calculating score-coefficient matrix by least-square method:
Cscor = Frot*inv(Frot'*diag(A)*Frot);
%