lp: Lp Rotation

In GPArotation: Gradient Projection Factor Rotation

L^p Rotation

Description

Performs L^p rotation to obtain sparse loadings.

Usage

     GPForth.lp(A, Tmat=diag(rep(1, ncol(A))), p=1, normalize=FALSE, eps=1e-05, 
            maxit=10000, gpaiter=5) 
     GPFoblq.lp(A, Tmat=diag(rep(1, ncol(A))), p=1, normalize=FALSE, eps=1e-05, 
            maxit=10000, gpaiter=5)  

Arguments

A	Initial factor loadings matrix to be rotated.
Tmat	Initial rotation matrix.
p	Component-wise L^p where 0 < p =< 1.
normalize	Not recommended for L^p rotation.
eps	Convergence is assumed when the norm of the gradient is smaller than eps.
maxit	Maximum number of iterations allowed in the main loop.
gpaiter	Maximum iterations for GPA rotation. The goal is to decrease the objective value, not optimize the inner loop. Warnings may appear, but they can be ignored if the main loop converges.

Details

These functions optimize an L^p rotation objective, where 0 < p =< 1. A smaller p promotes sparsity in the loading matrix but increases computational difficulty. For guidance on choosing p, see the Concluding Remarks in the references.

Since the L^p function is nonsmooth, a different optimization method is required compared to smooth rotation criteria. Two new functions, GPForth.lp and GPFoblq.lp, replace GPForth and GPFoblq for orthogonal and oblique L^p rotations, respectively.

The optimization method follows an iterative reweighted least squares (IRLS) approach. It approximates the nonsmooth objective function with a smooth weighted least squares function in the main loop and optimizes it using GPA in the inner loop.

Normalization is not recommended for L^p rotation. Its use may have unexpected results.

Value

A GPArotation object, which is a list containing:

loadings	Rotated loadings matrix, with one column per factor. If randomStarts were used, this contains the loadings with the lowest criterion value.
Th	Rotation matrix, satisfying loadings %*% t(Th) = A.
Table	Matrix recording iteration details during optimization.
method	String indicating the rotation objective function.
orthogonal	Logical indicating whether the rotation is orthogonal.
convergence	Logical indicating whether convergence was achieved. Convergence is controlled element-wise by tolerance.
Phi	Covariance matrix of rotated factors, t(Th) %*% Th.

Author(s)

Xinyi Liu, with minor modifications for GPArotation by C. Bernaards

References

Liu, X., Wallin, G., Chen, Y., & Moustaki, I. (2023). Rotation to sparse loadings using L^p losses and related inference problems. Psychometrika, 88(2), 527–553.

See Also

lpT, lpQ, vgQ.lp.wls

Examples

  data("WansbeekMeijer", package = "GPArotation")
  fa.unrotated <- factanal(factors = 2, covmat = NetherlandsTV, rotation = "none")
  
  options(warn = -1)
  
  # Orthogonal rotation
  # Single start from random position
  fa.lpT1 <- GPForth.lp(loadings(fa.unrotated), p = 1)
  # 10 random starts
  fa.lpT <- lpT(loadings(fa.unrotated), Tmat=Random.Start(2), p = 1, randomStarts = 10)
  print(fa.lpT, digits = 5, sortLoadings = FALSE, Table = TRUE, rotateMat = TRUE)
  
  p <- 1
  # Oblique rotation
  # Single start
  fa.lpQ1 <- GPFoblq.lp(loadings(fa.unrotated), p = p)
  # 10 random starts
  fa.lpQ <- lpQ(loadings(fa.unrotated), p = p, randomStarts = 10)
  summary(fa.lpQ, Structure = TRUE)


  # this functions ensures consistent ordering of factors of a
  # GPArotation object for cleaner comparison
  # Inspired by fungible::orderFactors and fungible::faSort functions
  sortFac <- function(x){
    # Only works for object of class GPArotation 
    if (!inherits(x, "GPArotation")) {stop("argument not of class GPArotation")}
    cln <- colnames(x$loadings) 
    # ordering for oblique slightly different from orthogonal
    ifelse(x$orthogonal, vx <- colSums(x$loadings^2), 
      vx <- diag(x$Phi %*% t(x$loadings) %*% x$loadings) )
    # sort by squared loadings from high to low
    vxo <- order(vx, decreasing = TRUE)
    vx <- vx[vxo]
    # maintain the right sign
    Dsgn <- diag(sign(colSums(x$loadings^3))) [ , vxo]
    x$Th <- x$Th %*% Dsgn 
    x$loadings <- x$loadings %*% Dsgn
    if (match("Phi", names(x))) { 
            # If factor is negative, reverse corresponding factor correlations
            x$Phi <- t(Dsgn) %*% x$Phi %*% Dsgn
      }
    colnames(x$loadings) <- cln 
    x 
  } 
  
  # seed set to see the results of sorting
  set.seed(1020)
  fa.lpQ1 <- lpQ(loadings(fa.unrotated),p=1,randomStarts=10)
  fa.lpQ0.5 <- lpQ(loadings(fa.unrotated),p=0.5,randomStarts=10)
  fa.geo <- geominQ(loadings(fa.unrotated), randomStarts=10)
  
  # with ordered factor loadings
  res <- round(cbind(sortFac(fa.lpQ1)$loadings, sortFac(fa.lpQ0.5)$loadings, 
    sortFac(fa.geo)$loadings),3)
  print(c("oblique-  Lp 1           Lp 0.5         geomin")); print(res)

  # without ordered factor loadings
  res <- round(cbind(fa.lpQ1$loadings, fa.lpQ0.5$loadings, fa.geo$loadings),3)
  print(c("oblique-  Lp 1           Lp 0.5         geomin")); print(res)

  options(warn = 0)
  

GPArotation documentation built on April 13, 2025, 1:08 a.m.