%\VignetteIndexEntry{Assessing Local Minima in Factor Rotation} %\VignetteEngine{utils::Sweave} %\VignetteEncoding{UTF-8} %\VignettePackage{GPArotation} %\VignetteDepends{GPArotation} %\VignetteKeyword{factor rotation} %\VignetteKeyword{local minima} %\VignetteKeyword{random starts} %\VignetteKeyword{simplimax} \documentclass[english, 10pt]{article} \usepackage{hyperref} \bibliographystyle{apa} \usepackage{natbib} \usepackage[margin=0.7in, letterpaper]{geometry} \begin{document} \SweaveOpts{eval=TRUE, echo=TRUE, results=hide, fig=FALSE} <>= options(continue=" ") pdf.options(pointsize = 8) @ \begin{center} \section*{Assessing Local Minima in Factor Rotation \\ ~~\\ The \texttt{GPFallMinima} Function} \end{center} \begin{center} Author: Coen A. Bernaards \end{center} Many rotation criteria have multiple local minima. The standard \texttt{GPArotation} functions return only the best solution found among all random starts --- the one with the lowest criterion value. While this is the right default behavior, it conceals potentially important information: how many distinct solutions exist, how often each is found, and how different they are from each other. This vignette introduces \texttt{GPFallMinima}, a companion function that runs multiple random starts and retains \emph{all} distinct local minima, not just the global minimum. This allows the analyst to assess the stability of the rotation solution and to examine alternative solutions that may be substantively meaningful. For background on local minima in factor rotation, see \cite{nguwall} and \cite{mansreise}. %------------------------------------------------------------------- \section*{The GPFallMinima Function} %------------------------------------------------------------------- \texttt{GPFallMinima} is not part of the standard \texttt{GPArotation} package. It is a companion utility intended for diagnostic use. Load it alongside \texttt{GPArotation}: <<>>= library("GPArotation") GPFallMinima <- function(A, method = "quartimin", orthogonal = FALSE, randomStarts = 100, eps = 1e-5, maxit = 2000, normalize = FALSE, methodArgs = NULL, minimumInclusion = 2) { engine <- if (orthogonal) GPForth else GPFoblq Qvalues <- numeric(randomStarts) Qconverged <- logical(randomStarts) all_res <- vector("list", randomStarts) for (i in 1:randomStarts) { res <- engine(A, Tmat = Random.Start(ncol(A)), normalize = normalize, eps = eps, maxit = maxit, method = method, methodArgs = methodArgs) Qvalues[i] <- res$Table[nrow(res$Table), 2] Qconverged[i] <- res$convergence all_res[[i]] <- res } converged_idx <- which(Qconverged) nConverged <- length(converged_idx) if (nConverged == 0) stop("No starts converged. Consider increasing maxit or relaxing eps.") Qvalues_conv <- Qvalues[converged_idx] all_res_conv <- all_res[converged_idx] Q_round <- round(Qvalues_conv / eps) * eps Q_unique <- unique(Q_round) minima <- vector("list", length(Q_unique)) for (j in seq_along(Q_unique)) { idx <- which(Q_round == Q_unique[j]) count <- length(idx) proportion <- count / nConverged minima[[j]] <- list( result = all_res_conv[[idx[1]]], f = Q_unique[j], count = count, proportion = proportion ) } ord <- order(sapply(minima, `[[`, "count"), decreasing = TRUE) minima <- minima[ord] keep <- sapply(minima, `[[`, "count") >= minimumInclusion minima <- minima[keep] if (length(minima) == 0) stop("No minima found with count >= minimumInclusion (", minimumInclusion, "). Consider reducing minimumInclusion or increasing randomStarts.") f_values <- sapply(minima, `[[`, "f") f_global <- min(f_values) summary_df <- data.frame( minimum = seq_along(minima), f = round(f_values, 6), deltaF = round(f_values - f_global, 6), count = sapply(minima, `[[`, "count"), proportion = round(sapply(minima, `[[`, "proportion"), 3), isGlobal = f_values == f_global ) result <- list( minima = minima, summary = summary_df, Qvalues = Qvalues_conv, nConverged = nConverged, nStarts = randomStarts, method = method, orthogonal = orthogonal, minimumInclusion = minimumInclusion ) class(result) <- "GPFallMinima" result } print.GPFallMinima <- function(x, ...) { cat("Random start analysis:", x$nConverged, "of", x$nStarts, "starts converged\n") cat("Distinct minima found:", nrow(x$summary), "(minimumInclusion =", x$minimumInclusion, ")\n\n") print(x$summary, row.names = FALSE) cat("\nGlobal minimum: f =", min(x$summary$f), "\n") cat("Access full solutions via $minima[[i]]$result\n") invisible(x) } @ %------------------------------------------------------------------- \section*{Key Arguments} %------------------------------------------------------------------- \begin{itemize} \item \texttt{method} --- the rotation criterion. Criteria with many local minima, such as \texttt{simplimax}, \texttt{geomin}, and \texttt{infomax}, are the most interesting to diagnose. \item \texttt{randomStarts} --- number of random starts to attempt. More starts give a more complete picture of the solution space. For a thorough analysis, 500 or more starts are recommended. \item \texttt{minimumInclusion} --- minimum number of starts required to retain a minimum. The default of 2 filters out singletons that are likely numerical artifacts. With 500 starts, a value of 5 or 10 is more appropriate. \item \texttt{orthogonal} --- \texttt{TRUE} for orthogonal rotation (uses \texttt{GPForth}), \texttt{FALSE} for oblique (uses \texttt{GPFoblq}). Default is \texttt{FALSE}. \end{itemize} Non-converged starts are always discarded before analysis. The \texttt{proportion} column in the summary is calculated over converged starts only, so proportions sum to 1. %------------------------------------------------------------------- \section*{Example: Simplimax on CCAI Climate-Friendly Purchasing Choices Data} %------------------------------------------------------------------- The CCAI Climate-Friendly Purchasing Choices domain (Bi and Barchard, 2024) provides a useful dataset for illustrating local minima behavior. The data have 14 items and 3 factors with strong inter-correlations (0.53--0.59). Bi and Barchard used oblimin rotation in their published analysis. We use simplimax here purely to illustrate how rotation criteria can produce highly complex landscapes with multiple local minima --- not as a recommended analysis of these data. The data illustrate how dramatically rotation criteria can differ in their stability. Oblimin rotation is highly stable on these data --- all 200 random starts converge to the same solution. Simplimax, by contrast, reveals a highly complex landscape: <>= library("GPArotation") data("CCAI", package = "GPArotation") fa_unrotated <- factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "none") A <- loadings(fa_unrotated) # Oblimin: highly stable res_oblimin <- oblimin(A, normalize = TRUE, randomStarts = 200) cat("Oblimin random start diagnostics:\n") res_oblimin$randStartChar # Simplimax: complex landscape set.seed(42) res_ccai <- GPFallMinima(A, method = "simplimax", randomStarts = 200, normalize = TRUE, minimumInclusion = 2) res_ccai @ The contrast is striking. Oblimin converges reliably to a single solution on these data. Simplimax reveals 16 distinct local minima, with only 81 of 200 starts converging (40.5\%) and the global minimum found in just 2.5\% of converged starts. This underscores the importance of using multiple random starts and verifying convergence, particularly when using criteria known to be prone to local minima such as simplimax. The global minimum (f = 0.01080) is clearly separated from the other local minima, with the next best solution having a criterion value three times larger. The sorted loadings plot shows that the local minima produce substantially different factor structures, not merely permutations or sign flips of the same solution. %------------------------------------------------------------------- \section*{Examining Individual Solutions} %------------------------------------------------------------------- Each element of \texttt{res\_ccai\$minima} contains a full \texttt{GPArotation} object in \texttt{result}, so the standard \texttt{print} and \texttt{summary} methods work directly. <>= # Print the most common solution print(res_ccai$minima[[1]]$result) @ <>= # Print solution in row 3 of the summary print(res_ccai$minima[[3]]$result, digits = 2) @ <>= i <- 3 print(res_ccai$minima[[i]]$result, digits = 2) @ <>= # Print the global minimum global <- which(res_ccai$summary$isGlobal) print(res_ccai$minima[[global]]$result, digits = 2) @ <>= # Summary with structure matrix for oblique solution summary(res_ccai$minima[[global]]$result, Structure = TRUE, digits = 2) @ Solutions with small \texttt{deltaF} values may produce loading matrices that are very similar to the global minimum after sorting. Solutions with large \texttt{deltaF} values are likely to produce meaningfully different loading patterns and warrant careful examination. %------------------------------------------------------------------- \section*{Visualizing All Minima} %------------------------------------------------------------------- The sorted absolute loadings plot is a powerful tool for comparing all retained minima at once. Each line represents one distinct solution. Lines that overlap closely indicate practically equivalent solutions; lines that diverge indicate qualitatively different solutions. <<>>= plotSortedLoadings <- function(..., labels = NULL, col = NULL, main = "Sorted Absolute Loadings", ylab = "Absolute loading", xlab = "Rank") { solutions <- list(...) for (i in seq_along(solutions)) { if (!inherits(solutions[[i]], "GPArotation")) stop("Argument ", i, " is not a GPArotation object.") } n <- length(solutions) if (is.null(labels)) labels <- paste("Solution", seq_len(n)) if (is.null(col)) col <- palette.colors(n, palette = "Okabe-Ito") sorted_loadings <- lapply(solutions, function(x) sort(abs(as.vector(x$loadings)), decreasing = FALSE)) all_values <- unlist(sorted_loadings) max_len <- max(sapply(sorted_loadings, length)) plot(NULL, xlim = c(1, max_len), ylim = c(0, max(all_values)), main = main, xlab = xlab, ylab = ylab, las = 1) abline(h = seq(0, 1, by = 0.1), col = "grey90", lty = 1) for (i in seq_len(n)) { lines(seq_along(sorted_loadings[[i]]), sorted_loadings[[i]], col = col[i], lwd = 2) points(seq_along(sorted_loadings[[i]]), sorted_loadings[[i]], col = col[i], pch = 19, cex = 0.6) } legend("topleft", legend = labels, col = col, lwd = 2, pch = 19, bty = "n") invisible(sorted_loadings) } @ <>= do.call(plotSortedLoadings, c(lapply(res_ccai$minima, function(x) x$result), list(labels = paste0("Min ", res_ccai$summary$minimum, " (f=", round(res_ccai$summary$f, 4), ", n=", res_ccai$summary$count, ")")))) @ Lines that are visually indistinguishable correspond to solutions that are practically equivalent regardless of their \texttt{deltaF}. Lines that diverge clearly represent genuinely different factor structures that merit substantive interpretation. %------------------------------------------------------------------- \section*{Practical Guidance} %------------------------------------------------------------------- \begin{itemize} \item \textbf{How many random starts?} For criteria that tend to have many local minima, 100 or more starts are recommended. The goal is for the proportions in the summary to stabilize. \item \textbf{Setting minimumInclusion.} With 100 starts, a value of 2 is reasonable. With 500 starts, use 5 or 10. The goal is to distinguish genuine local minima from numerical noise. \item \textbf{Interpreting deltaF.} Solutions with \texttt{deltaF} $< 0.001$ are typically negligible. Solutions with \texttt{deltaF} $> 0.01$ may produce visibly different loading patterns. \item \textbf{When the global minimum has low proportion.} If the global minimum is found by fewer than 20\% of converged starts, the criterion landscape is complex and the solution may not be stable. Consider trying a different rotation criterion or increasing the number of random starts. \item \textbf{Comparing criteria.} Running \texttt{GPFallMinima} with different rotation criteria on the same dataset can reveal which criteria produce stable solutions and which are prone to local minima for a given data structure. \end{itemize} %------------------------------------------------------------------- \section*{Further Resources} %------------------------------------------------------------------- For detailed discussion of local minima in factor rotation and their implications for substantive interpretation, see \cite{nguwall}. For investigation of local minima using the \textit{fungible} package, see the \texttt{faMain} function. The \textit{psych} package provides \texttt{faRotations} for similar diagnostics. The \texttt{randStartChar} element returned by standard \texttt{GPArotation} functions provides a quick summary of random start results without retaining individual solutions. For routine use, \texttt{randStartChar} is sufficient; \texttt{GPFallMinima} is intended for deeper diagnostic investigation. \begin{thebibliography}{} \bibitem[\protect\citeauthoryear{Bi \& Barchard}{Bi \& Barchard}{2024}]{bibarchard} Bi, Y. and Barchard, K.A. (2024). \newblock Purchasing choices that reduce climate change: An exploratory factor analysis. \newblock \textit{Spectra Undergraduate Research Journal}, \textbf{3}(2), 8--14. \newblock \href{https://doi.org/10.9741/2766-7227.1028} {doi: 10.9741/2766-7227.1028} \bibitem[\protect\citeauthoryear{Mansolf \& Reise}{Mansolf \& Reise}{2016}]{mansreise} Mansolf, M., \& Reise, S. P. (2016). \newblock Exploratory bifactor analysis: The Schmid-Leiman orthogonalization and Jennrich-Bentler analytic rotations. \newblock \textit{Multivariate Behavioral Research}, 51(5), 698--717. \newblock \href{https://doi.org/10.1080/00273171.2016.1215898} {https://doi.org/10.1080/00273171.2016.1215898} \bibitem[\protect\citeauthoryear{Nguyen \& Waller}{Nguyen \& Waller}{2022}]{nguwall} Nguyen, H. V., \& Waller, N. G. (2022). \newblock Local minima and factor rotations in exploratory factor analysis. \newblock \textit{Psychological Methods}. Advance online publication. \newblock \href{https://doi.org/10.1037/met0000467} {https://doi.org/10.1037/met0000467} \end{thebibliography} \end{document}