Skip to contents

Estimates the optimal number of principal components (r) using a global Probabilistic PCA (PPCA) model on pooled data. This function is designed for use before fitting mixture MPCA models, as it provides a cluster-independent estimate of r based on the overall data structure.

Usage

select_r_global_mpca(
  list_of_data,
  rvec = 1:6,
  model = c("PPCA", "FA"),
  criterion = c("CV", "ParallelAnalysis"),
  folds = 5,
  cv_metric = c("negLogLik", "frobenius"),
  n_permutations = 100,
  seed = NULL,
  verbose = FALSE
)

Arguments

list_of_data

A list of N matrices, each (T_i x M).

rvec

Vector of candidate r values to evaluate. Default is 1:6.

model

Model type for global r estimation. Options: "PPCA" (Probabilistic PCA, default), "FA" (Factor Analysis). PPCA is recommended as it is more stable and theory-consistent for mixture models.

criterion

Selection criterion. Options: "CV" (Cross-Validation, default), "ParallelAnalysis" (Parallel Analysis). CV is the primary criterion as it is stable and works for both FA and PPCA. Parallel Analysis is optional and lower priority.

folds

Number of folds for cross-validation. Default is 5. CV is performed at the subject level (not time-window level) to respect the generative structure.

cv_metric

Metric for CV evaluation. Options: "negLogLik" (negative log-likelihood, default), "frobenius" (Frobenius norm of reconstruction error). Negative log-likelihood is the principled choice for probabilistic models.

n_permutations

Number of permutations for Parallel Analysis. Default is 100. Only used when criterion = "ParallelAnalysis".

seed

Random seed for reproducibility (CV fold assignment and PA permutations). Default is NULL (no seed set).

verbose

Logical. If TRUE, print detailed diagnostic information during r selection. Default is FALSE.

Value

A list with:

r_global

The selected optimal r value.

summary

A data.frame with columns r and criterion-specific metrics (e.g., cv_error, cv_se).

criterion

The criterion used for selection.

model

The model type used (FA or PPCA).

cv_curve

For CV: ggplot object showing CV error vs r (if ggplot2 available).

vaf_curve

For CV: ggplot object showing VAF vs r (if ggplot2 available).

pa_curve

For PA: ggplot object showing observed vs permuted eigenvalues (if ggplot2 available).

diagnostics

Additional diagnostic information.

Details

This function implements model-based r estimation for use before mixture MPCA fitting. The key design principles are:

Why Global Model? The global r must be estimated using PPCA (or FA) on pooled data, NOT using MPCA(K>1). Estimating r from a mixture model is cluster-dependent and unstable, as different cluster configurations yield different r estimates. A global model provides a cluster-independent baseline.

Cross-Validation Design: CV is performed at the subject level, not time-window level. Time-window CV would leak subject-level covariance structure and provide overly optimistic estimates. Subject-level CV properly reflects the generative structure where subjects are the independent units.

For each fold:

  1. Split subjects into training and test sets

  2. Fit PPCA/FA on training data (pooled across subjects)

  3. Evaluate reconstruction error or log-likelihood on test data

  4. Average across folds to get CV error for each r

The optimal r is selected by minimizing CV error (or using the one-standard-error rule for parsimony).

Parallel Analysis: PA compares observed eigenvalues to those from permuted data. The optimal r is where observed eigenvalues exceed the 95th percentile of permuted eigenvalues. Note: PA is theoretically questionable for FA (as it assumes isotropic noise), so CV is preferred.

Examples

if (FALSE) { # \dontrun{
# Basic usage with default PPCA and CV
result <- select_r_global_mpca(
  list_of_data, rvec = 2:6,
  model = "PPCA", criterion = "CV", folds = 5
)
print(result$r_global)
print(result$summary)

# With FA model
result_fa <- select_r_global_mpca(
  list_of_data, rvec = 2:6,
  model = "FA", criterion = "CV"
)

# Using Parallel Analysis
result_pa <- select_r_global_mpca(
  list_of_data, rvec = 2:6,
  criterion = "ParallelAnalysis", n_permutations = 100
)
} # }