Random Starting Matrix for Factor Rotation

Generates a random orthogonal matrix for use as an initial rotation matrix (Tmat) in gradient projection rotation functions.

Usage

Random.Start(k = 2L)

Arguments

k: a positive integer specifying the dimension of the square orthogonal matrix to generate (default 2).

Value

A \(k \times k\) orthogonal matrix drawn uniformly from the Haar measure on the orthogonal group O(k). Columns have unit length and are mutually orthogonal.

Details

The function generates a random orthogonal matrix using QR decomposition of a matrix of standard normal variates, with a sign correction applied to the diagonal of R to ensure uniform sampling from the Haar measure. This follows the approach of Stewart (1980) and Mezzadri (2007).

The naive approach of qr.Q(qr(matrix(rnorm(k*k), k))) does not guarantee uniform sampling from the Haar measure. The sign correction Q %*% diag(sign(diag(R))) is required to achieve this. This was updated in GPArotation 2024.2-1 following a suggestion by Yves Rosseel.

For oblique rotation a random starting transformation matrix can be generated by normalizing the columns of a random matrix: X %*% diag(1/sqrt(diag(crossprod(X)))) where X <- matrix(rnorm(k*k), k).

References

Stewart, G.W. (1980). The efficient generation of random orthogonal matrices with an application to condition estimators. SIAM Journal on Numerical Analysis, 17(3), 403–409. doi: 10.1137/0717034

Mezzadri, F. (2007). How to generate random matrices from the classical compact groups. Notices of the American Mathematical Society, 54(5), 592–604. arXiv:math-ph/0609050

Author

Coen A. Bernaards and Robert I. Jennrich with some R modifications by Paul Gilbert. Updated following a suggestion by Yves Rosseel.

Examples

  # Generate a 5 x 5 random orthogonal matrix
  Random.Start(5)
#>            [,1]        [,2]       [,3]       [,4]        [,5]
#> [1,]  0.5688814  0.70609178  0.2707621  0.2681352 -0.18055389
#> [2,]  0.4221066 -0.68054718  0.3700148  0.4622146 -0.09015679
#> [3,]  0.6377921 -0.11761046 -0.6815294 -0.1508173  0.30358011
#> [4,] -0.2650616  0.15584326 -0.1708915  0.7240478  0.59330096
#> [5,] -0.1454767 -0.01306089 -0.5441422  0.4092310 -0.71770825

  # Verify orthogonality: t(Q) %*% Q should be the identity matrix
  Q <- Random.Start(4)
  round(t(Q) %*% Q, 10)
#>      [,1] [,2] [,3] [,4]
#> [1,]    1    0    0    0
#> [2,]    0    1    0    0
#> [3,]    0    0    1    0
#> [4,]    0    0    0    1

  # Use as starting matrix for rotation
  data("Thurstone", package = "GPArotation")
  simplimax(box26, Tmat = Random.Start(3))
#> Oblique rotation method Simplimax (k=26) converged.
#> Loadings:
#>         [,1]   [,2]   [,3]
#>  [1,]  0.716  0.722  0.488
#>  [2,]  0.723 -0.481  0.432
#>  [3,]  0.751 -0.126 -0.521
#>  [4,]  0.899  0.105  0.615
#>  [5,]  0.917  0.366 -0.072
#>  [6,]  0.874 -0.416 -0.099
#>  [7,]  0.883  0.388  0.583
#>  [8,]  0.864 -0.174  0.542
#>  [9,]  0.876  0.537  0.142
#> [10,]  0.978  0.222 -0.256
#> [11,]  0.842 -0.478  0.083
#> [12,]  0.857 -0.338 -0.258
#> [13,] -0.008  0.992 -0.015
#> [14,]  0.008 -0.992  0.015
#> [15,]  0.000  0.661  0.875
#> [16,]  0.000 -0.661 -0.875
#> [17,] -0.014 -0.397  0.837
#> [18,]  0.014  0.397 -0.837
#> [19,]  0.875  0.005  0.606
#> [20,]  0.902  0.347 -0.128
#> [21,]  0.882 -0.405 -0.080
#> [22,]  0.870  0.001  0.589
#> [23,]  0.883  0.317 -0.133
#> [24,]  0.874 -0.386 -0.051
#> [25,]  1.004  0.017  0.176
#> [26,]  0.970 -0.092  0.132
#> 
#>                  [,1]  [,2]  [,3]
#> SS loadings    14.909 5.459 5.041
#> Proportion Var  0.573 0.210 0.194
#> Cumulative Var  0.573 0.783 0.977
#> 
#> Phi:
#>        [,1]   [,2]   [,3]
#> [1,]  1.000 -0.037 -0.170
#> [2,] -0.037  1.000 -0.199
#> [3,] -0.170 -0.199  1.000