Skip to contents

1 Overview

This vignette replicates Simulation 2 from the manuscript, in which we examine how between-cluster separation (sep\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 sep∈{0.1,0.3,0.5,0.7,0.9,1.0}\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=50N=50 N_vec = 50
True clusters Ktrue=4K_{true}=4 K_true_vec = 4
Rank per cluster rtrue=4r_{true}=4 r_true_vec = 4
Muscles (channels) M=8M=8 M_vec = 8
Time points T=100T=100 T_each = 100
Between-cluster separation sep=0.1…1.0\text{sep}=0.1\dots1.0 sep_vec
Noise SD σ=2.0\sigma=2.0 noise_vec = 2.0
Candidate ranks r∈[3,5]r\in[3,5] r_candidates
Candidate clusters K∈[1,5]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 (sep≥0.5\text{sep}\ge 0.5) makes cluster boundaries more distinct, leading to ARI close to 1 for MFA and MPCA. When clusters overlap (sep=0.1\text{sep}=0.1), both mixture models struggle, yet still outperform single-subject baselines.

6 Session information

#> R version 4.5.2 (2025-10-31)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.4 LTS
#> 
#> Matrix products: default
#> BLAS:   /opt/R/4.5.2/lib/R/lib/libRblas.so 
#> LAPACK: /opt/R/4.5.2/lib/R/lib/libRlapack.so;  LAPACK version 3.12.1
#> 
#> locale:
#>  [1] LC_CTYPE=ja_JP.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=ja_JP.UTF-8        LC_COLLATE=ja_JP.UTF-8    
#>  [5] LC_MONETARY=ja_JP.UTF-8    LC_MESSAGES=ja_JP.UTF-8   
#>  [7] LC_PAPER=ja_JP.UTF-8       LC_NAME=C                 
#>  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
#> [11] LC_MEASUREMENT=ja_JP.UTF-8 LC_IDENTIFICATION=C       
#> 
#> time zone: Asia/Tokyo
#> tzcode source: system (glibc)
#> 
#> attached base packages:
#> [1] parallel  stats     graphics  grDevices utils     datasets  methods  
#> [8] base     
#> 
#> other attached packages:
#> [1] ggplot2_4.0.2     tidyr_1.3.2       dplyr_1.2.0       synergyMixR_1.0.0
#> 
#> loaded via a namespace (and not attached):
#>  [1] gtable_0.3.6       jsonlite_2.0.0     compiler_4.5.2     tidyselect_1.2.1  
#>  [5] Rcpp_1.1.1         jquerylib_0.1.4    systemfonts_1.3.1  scales_1.4.0      
#>  [9] textshaping_1.0.4  yaml_2.3.12        fastmap_1.2.0      R6_2.6.1          
#> [13] generics_0.1.4     knitr_1.51         htmlwidgets_1.6.4  tibble_3.3.1      
#> [17] desc_1.4.3         RColorBrewer_1.1-3 bslib_0.10.0       pillar_1.11.1     
#> [21] rlang_1.1.7        cachem_1.1.0       xfun_0.56          S7_0.2.1          
#> [25] fs_1.6.6           sass_0.4.10        otel_0.2.0         cli_3.6.5         
#> [29] withr_3.0.2        pkgdown_2.2.0      magrittr_2.0.4     digest_0.6.39     
#> [33] grid_4.5.2         lifecycle_1.0.5    vctrs_0.7.1        evaluate_1.0.5    
#> [37] glue_1.8.0         farver_2.1.2       ragg_1.5.0         rmarkdown_2.30    
#> [41] purrr_1.2.1        tools_4.5.2        pkgconfig_2.0.3    htmltools_0.5.9