1 Overview

This vignette replicates Simulation 2 from the manuscript, in which we examine how between-cluster separation (\(\text{sep}\)) influences the ability of

Mixture Factor Analysis (MFA) and Mixture PCA (MPCA)
SingleFA and SinglePCA

to recover latent sub-groups and reconstruct EMG signals.

The procedure varies the separation hyper-parameter \(\text{sep}\in\{0.1,0.3,0.5,0.7,0.9,1.0\}\) while keeping all other generative settings fixed (Table 1). After fitting each candidate model, we record

Bayesian Information Criterion (BIC),
Adjusted Rand Index (ARI) against the true cluster labels, and
Variance Accounted For (VAF).

2 Simulation design

Parameter	Value	R argument
Number of subjects	\(N=50\)	`N_vec = 50`
True clusters	\(K_{true}=4\)	`K_true_vec = 4`
Rank per cluster	\(r_{true}=4\)	`r_true_vec = 4`
Muscles (channels)	\(M=8\)	`M_vec = 8`
Time points	\(T=100\)	`T_each = 100`
Between-cluster separation	\(\text{sep}=0.1\dots1.0\)	`sep_vec`
Noise SD	\(\sigma=2.0\)	`noise_vec = 2.0`
Candidate ranks	\(r\in[3,5]\)	`r_candidates`
Candidate clusters	\(K\in[1,5]\)	`K_candidates`
EM iterations	20	`max_iter`
Random restarts	5	`n_rep_init`

Tip. Use fewer values or a smaller max_iter when prototyping to speed up the vignette build.

3 Running the simulation


set.seed(123)  # reproducibility of meta-parameters

result_dir <- "simulation_sep"

# This runs all combinations in parallel and saves .rds files
simulate_and_run_parallel_save(
  K_true_vec    = 4,

  r_true_vec    = 4,
  N_vec         = 50,
  M_vec         = 8,
  T_each        = 100,
  sep_vec       = seq(0.1, 1, 0.2),  # 0.1 to 1.0
  noise_vec     = 2.0,
  K_candidates  = 1:5,
  r_candidates  = 3:5,
  max_iter      = 20,
  n_rep_init    = 5,
  seed_data_base = 1000,
  mc.cores      = 10,

  output_dir    = result_dir
)

The function automatically uses parallel::mclapply() internally, storing one .rds file per scenario.

4 Aggregating results

After running the simulation, you can load and aggregate the results:

result_dir <- "simulation_sep"  # adjust to your path

# Load all .rds files from the result directory
rds_files <- list.files(result_dir, pattern = "\\.rds$", full.names = TRUE)
df_raw <- do.call(rbind, lapply(rds_files, readRDS))

# Compute summary table
df_result <- compute_ari_vaf_table(df_raw) %>%
  mutate(across(c(sep, noise), as.factor))

head(df_result)

We reshape the data to long format for ggplot convenience:

df_long_ari <- df_result %>%
  select(K_true, r_true, N, M, sep, noise, ARI_MFA, ARI_PCA) %>%
  pivot_longer(c(ARI_MFA, ARI_PCA),
               names_to  = "Method",
               values_to = "ARI")

5 Visualisation

5.1 Clustering accuracy (ARI)

ggplot(df_long_ari,
       aes(x = factor(K_true), y = ARI, fill = Method)) +
  geom_boxplot(position = position_dodge(width = 0.8), outlier.shape = NA) +
  facet_wrap(~ paste("sep =", sep, " noise =", noise), nrow = 2) +
  theme_bw() +
  labs(title = "Adjusted Rand Index (ARI)",
       x = "True K", y = "ARI", fill = "Model")

5.2 Model fit (BIC & VAF)

Observation. Higher separation (\(\text{sep}\ge 0.5\)) makes cluster boundaries more distinct, leading to ARI close to 1 for MFA and MPCA. When clusters overlap (\(\text{sep}=0.1\)), both mixture models struggle, yet still outperform single-subject baselines.

6 Session information

#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04.1 LTS
#> 
#> Matrix products: default
#> BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so;  LAPACK version 3.12.0
#> 
#> locale:
#>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
#>  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
#>  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
#>  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
#> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
#> 
#> time zone: Etc/UTC
#> tzcode source: system (glibc)
#> 
#> attached base packages:
#> [1] parallel  stats     graphics  grDevices utils     datasets  methods  
#> [8] base     
#> 
#> other attached packages:
#> [1] dplyr_1.2.1       tidyr_1.3.2       ggplot2_4.0.3     synergyMixR_1.0.0
#> [5] rmarkdown_2.32   
#> 
#> loaded via a namespace (and not attached):
#>  [1] gtable_0.3.6       jsonlite_2.0.0     compiler_4.6.1     tidyselect_1.2.1  
#>  [5] Rcpp_1.1.2         stringr_1.6.0      jquerylib_0.1.4    scales_1.4.0      
#>  [9] yaml_2.3.12        fastmap_1.2.0      plyr_1.8.9         R6_2.6.1          
#> [13] labeling_0.4.3     generics_0.1.4     knitr_1.52         tibble_3.3.1      
#> [17] maketools_1.3.2    bslib_0.12.0       pillar_1.11.1      RColorBrewer_1.1-3
#> [21] rlang_1.3.0        stringi_1.8.9      cachem_1.1.0       xfun_0.61         
#> [25] sass_0.4.10        sys_3.4.3          S7_0.2.2           otel_0.2.0        
#> [29] cli_3.6.6          withr_3.0.3        magrittr_2.0.5     digest_0.6.39     
#> [33] grid_4.6.1         mvtnorm_1.4-2      mclust_6.1.3       lifecycle_1.0.5   
#> [37] vctrs_0.7.3        evaluate_1.0.5     glue_1.8.1         farver_2.1.2      
#> [41] buildtools_1.0.0   reshape2_1.4.5     purrr_1.2.2        tools_4.6.1       
#> [45] pkgconfig_2.0.3    htmltools_0.5.9

Simulation Study: Effect of Cluster Separation on Model Performance