by C. A. Andersson
Contents  Data used: howto3.mat,
howto4.mat and howto5.mat Purpose: Learning the advanced aspects of Tucker modeling (e.g. core rotations) Information: R. Henrion Nway principal component analysis theory, algorithms and applications. Chemom.Intell.Lab.Syst. 25:123, 1994 Prerequisites: Knowledge of MATLAB, multiway arrays and the basis of Tucker models 
Using the introduction in the prior chapter as a basis, the features and characteristics of the Tucker3 model will now be discussed. Some of the features, e.g. rotational ambiguity, are known from twoway data analysis. However, many analysts need not regard, or forget, the rotational possibilities of the PCA factors. The nonuniqueness (i.e., rotation, scaling and permutation) of the factors resolved from twoway data results in ambiguous solutions, albeit offers possibilities for exploring the model in many interesting ways. Postprocessing solutions from decomposition/factorization belongs to the discipline of factor analysis (FA). The FA approach will be discussed in relation to the solutions from the Tucker3 model in order to equip the analyst with more ways to explore the estimated model  eventually get a better understanding of his or her data.
For orthogonal models, it can be shown that the size of the squared
core entries are proportional to the variation explained of the model by
the factor combination. In other words, to find the most significant factor
combination, one must find the largest squared entry in the core.
How to do in the Nway Toolbox #3Using Eq. (4) in the previous section we wish, by experiment, to verify that the orthogonality poses no change in the fit. Using the same data set we first estimate an unconstrained model and then a constrained model. You will do this by using the commands below. Use this exciting lesson to get acquainted with the way data are handled and the Nway Toolbox.>load howto3 %Load X and W >%Make the ordinary orthogonality constrained Tucker3 model %The [0 0 0] parameters define an orthogonal model (which is also the default setting) >[Factorso,Go,SSEo]=tucker(X,W,[],[0 0 0]); >[Ao Bo Co]=fac2let(Factorso); %Convert to component matrices > %If the these three matrices exhibit identity the solutions are orthogonal >Ao'*Ao, Bo'*Bo, Co'*Co >%Make an unconstrained model %The [2 2 2] parameters define an unconstrained model >[Factorsu,Gu,SSEu]=tucker(X,W,[],[2 2 2]); >[Au Bu Cu]=fac2let(Factorsu); %Convert to component matrices >%Since the following three matrices are not identity matrices >%the components are oblique (unconstrained)
>Au'*Au, Bu'*Bu, Cu'*Cu
>%The errors of the two model should be equal (to rounding errors)
>[SSEo SSEu] %The errors of the two model should be equal (to rounding errors)

Example #Question:The core from the orthogonal model above has been slightly modified and is listed as an unfolded matrix G below. Remember that w = (3,3,3).
Try to identify the three most significant core elements? From the positions of the significant elements, find out what the corresponding factors are? Is there one particularly important factor in the model? Answer 
From X (r_{1},r_{2},r_{3}), a solution with orthonormal A (r_{1},w_{1}), B (r_{2},w_{2}), C (r_{3},w_{3}) and G (w_{1},w_{2},w_{3}) has been calculated. We may write the model of X as in Eq. (4).
If the component matrices are rotated by arbitrary columnwise orthonormal matrices O (w_{1},w_{1}), P (w_{2},w_{2}) and Q (w_{3},w_{3}) the least squares fit of the model is maintained. Using the rotated factors, the new corresponding, or counterrotated, core is calculated as shown in Eq. (9).
Eq. (9) is the basis for rotation of Tucker3 solution as the rotated core is easily calculated without involving direct rotation of the factors themselves  we will discuss this in the next section. To conclude, Eq. (9) is used in Eq. (10) to formulate the rotated solution.
From Eq. (10) we see that rotation of the factors results in a counterrotation
of the core when applied to X. Hence, any initial solution
to a Tucker3 model may be rotated using arbitrary orthonormal rotation
matrices without affecting the model of X. Thus, the fit
also remains unaffected as postulated above.
How to do in the Nway Toolbox #4In order to demonstrate the rotational abilities of the Tucker3 model we will calculate a model and arbitrarily rotate the factors and the core accordingly. We use r = (8,6,7) and w = (3,2,2). Investigate/plot the differences between the rotated and the unrotated factors. Take a particular look to the differences between the cores. Don't clear the screen before you have answered Example #5. Later in the chapter we will use more intelligent schemes for optimizing the simplicity of cores.>load howto4 %Load X and W >[Factors,G] = tucker(X,W); %Make a Tucker3 model >[A B C] = fac2let(Factors); %Make component matrices >%Lets make three arbitrary rotation matrices, one for each mode >O = orth(rand(W(1),W(1)));P=orth(rand(W(2),W(2)));Q=orth(rand(W(3),W(3))); >Ar = A*O; %Rotate mode 1 >Br = B*P; %Rotate mode 2 >Cr = C*Q; %Rotate mode 3 >Xunf = reshape(X,size(X,1),size(X,2)*size(X,3)); %Unfold X >Gr = Ar'*Xunf*ckron(Cr,Br); %Calculate the new unfolded core >Gr = reshape(Gr,W); %reshape to array >[sum(G(:).^2) sum(Gr(:).^2)] %Must be equal >SSEr = sum(sum(sum( (Xnmodel({Ar,Br,Cr},Gr)).^2 ))); %Calculate the error of the rotated model >SSEo = sum(sum(sum( (Xnmodel({A,B,C},G)).^2 ))); %Calculate the error of the initial model >[SSEo SSEr] %The errors should also be equal 
Example #5Question:Inspect the initial core G and the rotated core Gr from Howto #4 above. Are you free to choose any of these two models for interpretation? If so, what features would you look for and why? Answer 
The analyst may have a desire that is formulated in terms of the components matrices. E.g., it may be desirable to optimize the variance measures of the factors to explore hidden differences or commonalities in selected modes. The perhaps most common measures belong to the family of Orthomax, e.g. Varimax and Quartimax. In addition Quartimin and Covarimin may be useful. Any (orthogonal) twoway rotation can be applied in a straightforward manner to obtain the rotation matrices for the components. These matrices are then applied to the core to complete the rotation of the model.
In chapter 5 of his textbook Kroonenberg (Threemode Principal Component Analysis. Theory and Applications, Leiden:DSWO Press, 1983) suggest that one should seek to diagonalize the planes of the core. To obtain an even simpler core Kiers (TUCKALS core rotations and constrained TUCKALS modelling. Statistica Applicata 4:659667, 1992) propose to diagonalize the core in whole. By the latter method the solutions are rotated in a controlled way until the sum of the squared diagonal elements is at maximum. Note that the core is required to be quadratic, i.e., is required to have w_{1} = w_{2} = w_{3}, otherwise the diagonal is not defined. By this clever proposal the Tucker3 solution will be rotated to have optimal resemblance with a PARAFAC solution with orthogonal factors. Depending on the data and the dimensionality of the model, this rotation may provide a strong simplification of the core. However, in cases where the Tucker3model is needed over a PARAFAC model one often requires the presence of significant offdiagonal elements.
How to do in the Nway Toolbox #5First we will calculate a Tucker3 model with orthogonal component matrices. Then we will rotate the initial core to optimal diagonality and varianceofsquares. Continue with Example #6 before you blank the screen.>load howto5 %Load X and W >[Factors,Go]=tucker(X,W); %Make the Tucker model >[A B C]=fac2let(Factors); %Convert to component matrices >[Gd,Od1,Od2,Od3]=maxdia3(Go); %Rotate to optimum diagonality >[Gv,Ov1,Ov2,Ov3]=maxvar3(Go); %Rotate to optimum varianceofsquares >explcore(Go,7); %Inspect the unrotated solution >explcore(Gd,7); %Inspect the diagonalized solution >explcore(Gv,7); %Inspect the varianceofsquares optimized solution > %Reshape X and G to unfolded matrices
> Xunf = reshape(X,size(X,1),size(X,2)*size(X,3));
> Gounf = reshape(Go,size(Go,1),size(Go,2)*size(Go,3));
> Gdunf = reshape(Gd,size(Go,1),size(Go,2)*size(Go,3));
> Gvunf = reshape(Gv,size(Go,1),size(Go,2)*size(Go,3));
>sum(sum( (Xunf  A*Gounf*kron(C',B')).^2 )) %Find the error >sum(sum( (Xunf  (A*Od1)*Gdunf*kron((C*Od3)',(B*Od2)')).^2 )) %Find the error >sum(sum( (Xunf  (A*Ov1)*Gvunf*kron((C*Ov3)',(B*Ov2)')).^2 )) %Find the error 
Example #6Question:Use some time to discuss the pros and cons of the three models obtained above  which model would you prefer to interprete? Why? What makes you select one model over another? Is there anything in the analysis that suggest to use the PARAFAC model on these data?We have, for some reason become suspicious with regards to the derived transformation matrices derived above. What would you do to check the orthogonality of the transformation matrices and the orthogonality of the rotated factors? Answer 
The Nway tutorial
Copyright © 1998
Changed Jan2001
R. Bro