if (!requireNamespace("BiocManager", quietly = TRUE)) {
install.packages("BiocManager")
}
BiocManager::install(c("OMAWorkflows", "mia"))Filtering input for Joint-RPCA
A reproducible workflow using mia
2026-08-03
Source:vignettes/filterRPCAInput-validation.qmd
Introduction
Robust principal component analysis (RPCA) methods require appropriate filtering of the raw count table before robust centered log-ratio transformation. Low-depth samples and low-abundance or low-prevalence features can otherwise contribute little information and increase computational cost.
The mia::filterRPCAInput() function provides Gemelli-style filters for SummarizedExperiment objects. Samples are retained when their total count is strictly greater than min.sample.count. Features are retained when both their total count and prevalence are strictly greater than the corresponding thresholds.
This vignette demonstrates the filtering step with executable R code. The complete cross-language comparison of mia and Gemelli is maintained as a separate validation report. See the OMAWorkflows package website for additional workflows.
Installation
Install OMAWorkflows and mia from Bioconductor:
Load the packages used in the example:
Example data
The GlobalPatterns dataset contains a raw count assay and is small enough for a package vignette.
Filter the RPCA input
This example retains samples with more than 10,000 reads and features with more than 100 total reads and prevalence greater than 10%.
min_sample_count <- 10000
min_feature_count <- 100
min_feature_frequency <- 0.10
filtered_tse <- filterRPCAInput(
tse,
assay.type = "counts",
min.sample.count = min_sample_count,
min.feature.count = min_feature_count,
min.feature.frequency = min_feature_frequency
)
data.frame(
object = c("Before filtering", "After filtering"),
features = c(nrow(tse), nrow(filtered_tse)),
samples = c(ncol(tse), ncol(filtered_tse))
) object features samples
1 Before filtering 19216 26
2 After filtering 6213 26
Verify the filtering rules
The following calculation independently reconstructs the expected identifiers. The checks also demonstrate that all three thresholds use strict greater-than comparisons.
expected_samples <- colSums(count_matrix) > min_sample_count
expected_features <- rowSums(count_matrix) > min_feature_count
feature_prevalence <-
rowSums(count_matrix > 0) / ncol(count_matrix)
expected_features <- expected_features &
feature_prevalence > min_feature_frequency
stopifnot(
identical(
colnames(filtered_tse),
colnames(tse)[expected_samples]
),
identical(
rownames(filtered_tse),
rownames(tse)[expected_features]
)
)
data.frame(
criterion = c(
"Sample total count",
"Feature total count",
"Feature prevalence"
),
threshold = c(
min_sample_count,
min_feature_count,
min_feature_frequency
),
comparison = rep("strictly greater than", 3L)
) criterion threshold comparison
1 Sample total count 1e+04 strictly greater than
2 Feature total count 1e+02 strictly greater than
3 Feature prevalence 1e-01 strictly greater than
Individual filters can be disabled with NULL. This is useful when a threshold has already been applied upstream.
unfiltered_tse <- filterRPCAInput(
tse,
assay.type = "counts",
min.sample.count = NULL,
min.feature.count = NULL,
min.feature.frequency = NULL
)
stopifnot(identical(dim(unfiltered_tse), dim(tse)))Complete mia-Gemelli validation
The complete cross-language validation includes filtering and rCLR agreement, component alignment, sample and feature loadings, explained variation, sample distances, and cross-validation trajectories.
View the complete Joint-RPCA validation report or inspect the reproducibility files.
References
Martino, C., Shenhav, L., Marotz, C. A., et al. (2020). Context-aware dimensionality reduction deconvolutes gut microbial community dynamics. Nature Biotechnology, 39, 165-168. doi:10.1038/s41587-020-0660-7.
Session information
R version 4.6.1 (2026-06-24)
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 24.04.4 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.26.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] stats4 stats graphics grDevices utils datasets methods
[8] base
other attached packages:
[1] mia_1.21.6 TreeSummarizedExperiment_2.21.0
[3] Biostrings_2.81.6 XVector_0.53.0
[5] SingleCellExperiment_1.35.2 MultiAssayExperiment_1.39.0
[7] SummarizedExperiment_1.43.0 Biobase_2.73.2
[9] GenomicRanges_1.65.1 Seqinfo_1.3.0
[11] IRanges_2.47.2 S4Vectors_0.51.6
[13] BiocGenerics_0.59.10 generics_0.1.4
[15] MatrixGenerics_1.25.0 matrixStats_1.5.0
loaded via a namespace (and not attached):
[1] tidyselect_1.2.1 viridisLite_0.4.3
[3] dplyr_1.2.1 vipor_0.4.7
[5] farver_2.1.2 viridis_0.6.5
[7] S7_0.2.2 fastmap_1.2.0
[9] lazyeval_0.2.3 bluster_1.23.0
[11] digest_0.6.39 rsvd_1.0.5
[13] DirichletMultinomial_1.55.0 lifecycle_1.0.5
[15] cluster_2.1.8.3 tidytree_0.4.8
[17] magrittr_2.0.5 compiler_4.6.1
[19] rlang_1.3.0 tools_4.6.1
[21] igraph_2.3.3 yaml_2.3.12
[23] knitr_1.51 S4Arrays_1.13.0
[25] DelayedArray_0.39.3 plyr_1.8.9
[27] RColorBrewer_1.1-3 abind_1.4-8
[29] BiocParallel_1.47.0 purrr_1.2.2
[31] grid_4.6.1 beachmat_2.29.0
[33] ggplot2_4.0.3 MASS_7.3-66
[35] scales_1.4.0 cli_3.6.6
[37] vegan_2.7-5 rmarkdown_2.31
[39] crayon_1.5.3 treeio_1.37.0
[41] otel_0.2.0 reshape2_1.4.5
[43] DelayedMatrixStats_1.35.0 scuttle_1.23.1
[45] DBI_1.3.0 ggbeeswarm_0.7.3
[47] ape_5.8-1 stringr_1.6.0
[49] splines_4.6.1 parallel_4.6.1
[51] yulab.utils_0.2.4 vctrs_0.7.3
[53] Matrix_1.7-6 jsonlite_2.0.0
[55] BiocSingular_1.29.0 BiocNeighbors_2.7.2
[57] ggrepel_0.9.8 decontam_1.33.0
[59] irlba_2.3.7 beeswarm_0.4.0
[61] scater_1.41.2 tidyr_1.3.2
[63] glue_1.8.1 codetools_0.2-20
[65] stringi_1.8.7 gtable_0.3.6
[67] ScaledMatrix_1.21.0 tibble_3.3.1
[69] pillar_1.11.1 rappdirs_0.3.4
[71] htmltools_0.5.9 R6_2.6.1
[73] sparseMatrixStats_1.25.0 evaluate_1.0.5
[75] lattice_0.22-9 DECIPHER_3.9.2
[77] Rcpp_1.1.2 permute_0.9-10
[79] gridExtra_2.3.1 SparseArray_1.13.2
[81] nlme_3.1-170 mgcv_1.9-4
[83] xfun_0.60 fs_2.1.0
[85] pkgconfig_2.0.3