sla.varimax {RScaLAPACK}R Documentation

Rotation Methods for Factor Analysis

Description

These functions ‘rotate’ loading matrices in factor analysis using ScaLAPACK.

Usage

sla.varimax(x, normalize = TRUE, eps = 1e-5, NPROWS=0, NPCOLS=0, MB=16)
sla.promax(x, m = 4, NPROWS=0, NPCOLS=0, MB=16)

Arguments

x

A loadings matrix, with p rows and k < p columns

m

The power used the target for promax. Values of 2 to 4 are recommended.

normalize

logical. Should Kaiser normalization be performed? If so the rows of x are re-scaled to unit length before rotation, and scaled back afterwards.

eps

The tolerance for stopping: the relative change in the sum of singular values.

NPROWS

Number of Process Rows in the Process Grid.

NPCOLS

Number of Process Cols in the Process Grid.

MB

Block Size.

Details

These seek a ‘rotation’ of the factors x %*% T that aims to clarify the structure of the loadings matrix. The matrix T is a rotation (possibly with reflection) for varimax, but a general linear transformation for promax, with the variance of the factors being preserved.

Value

A list with components

loadings

The ‘rotated’ loadings matrix, x %*% rotmat.

rotmat

The ‘rotation’ matrix.

Author(s)

Nagiza Samatova (samatovan@ornl.gov), Guruprasad Kora (koragh@ornl.gov), Srikanth Yoginath (yoginathsb@ornl.gov), David Bauer (bauerda@ornl.gov)

References

Hendrickson, A. E. and White, P. O. (1964) Promax: a quick method for rotation to orthogonal oblique structure. British Journal of Statistical Psychology, 17, 65–70.

Horst, P. (1965) Factor Analysis of Data Matrices. Holt, Rinehart and Winston. Chapter 10.

Kaiser, H. F. (1958) The varimax criterion for analytic rotation in factor analysis. Psychometrika 23, 187–200.

Lawley, D. N. and Maxwell, A. E. (1971) Factor Analysis as a Statistical Method. Second edition. Butterworths.

See Also

varimax for the normal varmax, promax for the normal promax, sla.factanal, factanal, Harman74.cor.

Examples

## varimax with normalize = T is the default

v1 <- c(1,1,1,1,1,1,1,1,1,1,3,3,3,3,3,4,5,6)
v2 <- c(1,2,1,1,1,1,2,1,2,1,3,4,3,3,3,4,6,5)
v3 <- c(3,3,3,3,3,1,1,1,1,1,1,1,1,1,1,5,4,6)
v4 <- c(3,3,4,3,3,1,1,2,1,1,1,1,2,1,1,5,6,4)
v5 <- c(1,1,1,1,1,3,3,3,3,3,1,1,1,1,1,6,4,5)
v6 <- c(1,1,1,2,1,3,3,3,4,3,1,1,1,2,1,6,5,4)
m1 <- cbind(v1,v2,v3,v4,v5,v6)
fa <- sla.factanal(m1, factors=3, rotation="sla.promax", NPROWS=2, NPCOLS=2, MB=2 )
sla.varimax(fa$loadings, normalize = FALSE, NPROWS=2, NPCOLS=2, MB=2)
sla.promax(fa$loadings, NPROWS=2, NPCOLS=2, MB=2)

[Package RScaLAPACK version 0.6.1 Index]