A reproducible R-based pipeline for end-to-end analysis of DIA proteomics and bulk RNA-seq transcriptomics data from mouse or human brain tissue. Originally developed for Alzheimer's disease mouse models but generalized for any case-control multi-omics study with multiple genotypes, sexes, and timepoints.
| Step | Analysis |
|---|---|
| QC | Protein abundance matrix loading, boxplot QC, summarization to protein level |
| PCA | Probabilistic PCA (handles missing values) with metadata overlay |
| Differential Abundance | ANOVA + Tukey HSD post-hoc, stratified by sex and age |
| Visualization | Volcano plots, DEP barplots, UpSet overlap plots, expression boxplots |
| Functional Enrichment | KEGG pathway ORA via clusterProfiler |
| Step | Analysis |
|---|---|
| Linear Modeling | abundance ~ Genotype + Age + Sex per protein with covariate control |
| Module Correlation (logFC) | Pearson correlation of ANOVA fold-changes with human AD co-expression modules |
| Module Correlation (LM) | Pearson correlation of LM estimates with human AD co-expression modules |
| BioDomain ORA | GO over-representation analysis annotated with AMP-AD biological domains |
| BioDomain GSEA | Full ranked GSEA annotated with AMP-AD biological domains |
| Sub-BioDomain | Portrait and landscape plots across AD biological sub-domains |
| Step | Analysis |
|---|---|
| Data Inspection | Metadata loading, count matrix QC, transgene merging (optional) |
| Isoform Expression | Visualization of specific gene isoforms e.g. human MAPT (optional) |
| Sex Validation | Xist / Ddx3y expression to confirm sex assignments |
| PCA | DESeq2 VST-based PCA + PCAtools biplot and eigencorplot |
| Differential Expression | DESeq2 pairwise comparisons, stratified by sex and age |
| Visualization | Volcano plots, DEG barplots |
| Functional Enrichment | KEGG pathway ORA via clusterProfiler |
| Step | Analysis |
|---|---|
| Module Correlation (logFC) | Pearson correlation of DESeq2 logFC with sex-stratified AMP-AD modules |
| Linear Modeling | Per-gene LM on TPM data with genotype, sex, and age covariates |
| Module Correlation (LM) | Pearson correlation of LM estimates with AMP-AD modules |
| AD Subtype Correlation | Correlation with Nikhil et al. and Neff et al. AD subtype signatures |
| BioDomain Correlation (LM) | LM estimates vs AMP-AD BioDomain reference fold-changes |
| BioDomain ORA | GO over-representation annotated with AMP-AD biological domains |
| BioDomain GSEA | Full ranked GSEA annotated with AMP-AD biological domains |
brain-multiomics-pipeline/
β
βββ README.md # This file
βββ LICENSE # MIT license
βββ .gitignore # Excludes data, results, and large files
β
βββ notebooks/
β βββ 01_DIA_Proteomics_Pipeline.Rmd # Proteomics QC, PCA, ANOVA, KEGG
β βββ 02_Human_Translation_Analysis.Rmd # Proteomics LM, modules, BioDomain
β βββ 03_RNAseq_Pipeline.Rmd # RNA-seq QC, PCA, DESeq2, KEGG
β βββ 04_RNAseq_Human_Translation.Rmd # RNA-seq LM, modules, subtypes, BioDomain
β
βββ R/
β βββ Helper_Functions_Prot.R # Proteomics utility functions
β βββ Helper_Functions.R # Transcriptomics utility functions
β
βββ data/
β βββ example/
β βββ example_abundance.csv # Example protein abundance matrix
β βββ example_traits.csv # Example sample metadata
β
βββ figures/ # Example output figures
β
βββ results/ # gitignored β outputs saved here locally
β
βββ docs/
βββ input_format.md # Detailed input file format guide
Run order: Notebooks must be run in sequence within each omics type. Notebook 02 requires outputs from 01. Notebook 04 requires outputs from 03. The proteomics and transcriptomics pipelines are independent of each other.
R β₯ 4.2.0
install.packages(c(
"tidyverse", "gt", "UpSetR", "data.table", "broom",
"ggpubr", "ggh4x", "corrplot", "cowplot", "ggrepel",
"matrixStats", "VennDiagram", "readxl", "readr"
))if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install(c(
# Core
"AnnotationDbi", # Gene ID mapping
"clusterProfiler", # Pathway enrichment (ORA and GSEA)
"enrichplot", # Enrichment visualization
"ComplexHeatmap", # Heatmap visualization
# Proteomics specific
"pcaMethods", # Probabilistic PCA (handles missing values)
# Transcriptomics specific
"DESeq2", # Differential expression
"PCAtools", # Advanced PCA and eigencorplot
"DEGreport", # DEG reporting utilities
"ReactomePA", # Reactome pathway analysis
"limma", # Linear models for omics
# Organism annotation β install only what you need
"org.Mm.eg.db", # Mouse
"org.Hs.eg.db", # Human
"org.Rn.eg.db" # Rat
))Species note: Only one organism package is needed per study. Set
organism_dbin the USER CONFIGURATION block of each notebook.
External reference datasets required for human comparison analyses:
| Reference | Source | Used in | Used for |
|---|---|---|---|
| AMP-AD co-expression modules | AMP-AD Knowledge Portal | 02, 04 | Module logFC and LM correlation |
| AMP-AD sex-stratified modules | AMP-AD Knowledge Portal | 04 | Sex-specific module correlation |
| BioDomain annotations | Agora / AMP-AD | 02, 04 | GO term β BioDomain mapping |
| Domain color table | AMP-AD | 02, 04 | BioDomain plot styling |
| Milind et al. AD subtypes | Milind et al. 2021 | 04 | AD subtype correlation |
| Neff et al. AD subtypes | Neff et al. 2021 | 04 | AD subtype correlation |
Step 1: Prepare your input files (see docs/input_format.md):
- Normalized protein abundance matrix (CSV)
- Sample metadata (CSV)
Step 2: Configure and run notebook 01:
# Key settings in the USER CONFIGURATION block:
data_dir <- "../data/"
abundance_file <- "your_abundance_file.csv"
traits_file <- "your_metadata_file.csv"
result_dir <- "../results/"
control_genotype <- "Your_Control_Genotype"
case_genotypes <- c("Case_1", "Case_2", "Case_3")
age_groups <- c(4, 12)
organism_db <- "org.Mm.eg.db" # swap for human
kegg_organism <- "mmu" # "hsa" for human
rmarkdown::render("notebooks/01_DIA_Proteomics_Pipeline.Rmd")Step 3: Configure and run notebook 02 (optional):
# Enable only the analyses you have reference data for
run_module_corr <- TRUE
run_biodomain_ora <- TRUE
run_biodomain_gsea <- TRUE
run_lm_analysis <- TRUE
rmarkdown::render("notebooks/02_Human_Translation_Analysis.Rmd")Step 1: Prepare your input files:
- RSEM gene count matrix (TSV)
- RSEM gene TPM matrix (TSV)
- Sample metadata (CSV)
- Optionally: isoform TPM matrix (TSV)
Step 2: Configure and run notebook 03:
# Key settings in the USER CONFIGURATION block:
data_dir <- "../data/"
gene_counts_file <- "rsem.merged.gene_counts.tsv"
gene_tpm_file <- "rsem.merged.gene_tpm.tsv"
metadata_file <- "metadata_validated.csv"
result_dir <- "../results/"
control_genotype <- "Your_Control_Genotype"
case_genotypes <- c("Case_1", "Case_2", "Case_3")
age_groups <- c(4, 12)
organism_db <- "org.Mm.eg.db"
kegg_organism <- "mmu"
# Transgene handling β set TRUE only if your model has human transgenes
has_transgenes <- FALSE
rmarkdown::render("notebooks/03_RNAseq_Pipeline.Rmd")Step 3: Configure and run notebook 04 (optional):
# Enable only the analyses you have reference data for
run_module_corr_fc <- TRUE
run_module_corr_lm <- TRUE
run_subtype_corr <- TRUE
run_lm_analysis <- TRUE
run_biodomain_ora <- TRUE
run_biodomain_gsea <- TRUE
rmarkdown::render("notebooks/04_RNAseq_Human_Translation.Rmd")ProteinID,Sample1,Sample2,Sample3,...
MAPT|peptide1,12.3,11.8,12.1,...
APP|peptide2,10.5,10.2,10.8,...
- Rows = proteins; columns = samples
- Values should be pre-normalized (e.g., TAMPOR-corrected, log2-scale)
- ProteinID can be
ProteinorProtein|Peptideformat
X,Sex,Genotype,Age
Sample1,Female,Control,4
Sample2,Male,Control,4
Sample3,Female,CaseA,4
X: Sample ID β must exactly match column names in abundance/count matrixSex: Male / FemaleGenotype: group labels matching your config variablesAge: numeric timepoint in months
Standard RSEM TSV output with gene_id as first column and one column per sample. See docs/input_format.md for full format details and edge cases.
# For proteomics notebook 01, set in USER CONFIGURATION:
data_dir <- "../data/example/"
abundance_file <- "example_abundance.csv"
traits_file <- "example_traits.csv"Example data contains 20 proteins and 12 samples across 2 genotypes, 2 sexes, and 2 ages. Designed to verify end-to-end pipeline execution β enrichment results are not biologically meaningful at this scale.
| Output | Location | Generated by | Description |
|---|---|---|---|
| HTML Reports | notebooks/ |
All notebooks | Full rendered analysis with all plots |
| Processed Proteomics Data | data/Processed_DIA_Proteomics.RData |
Notebook 01 | Cleaned abundance matrix + traits |
| ANOVA Results | results/ANOVA_Results_DIA_Proteomics.Rdata |
Notebook 01 | All DEP results by comparison |
| Proteomics Human Translation | results/Human_Translation_Results.Rdata |
Notebook 02 | LM, module, BioDomain results |
| Processed RNA-seq Data | data/ProcessedData_RNAseq.RData |
Notebook 03 | Filtered count + TPM matrices |
| DESeq2 Results | results/DESeq2_Results_RNAseq.RData |
Notebook 03 | All DEG results by comparison |
| LM Results (RNA-seq) | results/LM_Results_RNAseq.RData |
Notebook 04 | LM estimates by gene and variant |
| ORA BioDomain (RNA-seq) | results/ORA_BioDomain_RNAseq.RData |
Notebook 04 | GO-ORA annotated with BioDomains |
| GSEA BioDomain (RNA-seq) | results/GSEA_BioDomain_RNAseq.RData |
Notebook 04 | GO-GSEA annotated with BioDomains |
| RNA-seq Human Translation | results/Human_Translation_RNAseq.RData |
Notebook 04 | All human correlation results |
Proteomics QC: Sample distributions assessed by boxplot after normalization. Proteins with pipe-delimited IDs summarized to protein level by averaging duplicate entries. Probabilistic PCA (pcaMethods::pca()) accommodates missing values common in DIA data.
Differential Abundance (ANOVA): One-way ANOVA per protein with Tukey HSD post-hoc correction, stratified by sex and age group. Only group pairs present in the data are tested β invalid comparisons are automatically skipped.
Transcriptomics QC: Optional transgene merging for models with human gene inserts (configurable via lookup table). Sex validated using Xist and Ddx3y expression. Low-count genes pre-filtered by minimum TPM threshold.
Differential Expression (DESeq2): Pairwise DESeq2 comparisons per sex and age stratum using raw integer counts. Genes annotated with symbols and Entrez IDs via AnnotationDbi.
Linear Modeling: Per-gene/protein linear models estimate the independent contribution of each genotype while controlling for covariates (Age, Sex). Model formula and coefficients are fully user-configurable. Stratified by age group where relevant.
Pathway Enrichment: KEGG over-representation analysis using clusterProfiler::compareCluster(). Gene set enrichment analysis (GSEA) using clusterProfiler::gseGO() on full ranked gene lists.
Human Translation: Mouse model fold-changes and LM estimates correlated (Pearson) with human AMP-AD co-expression module signatures (sex-stratified and combined). GO enrichment results mapped to the AMP-AD BioDomain framework across 19 AD biological domains and sub-domains. AD subtype correlations against Nikhil et al. (ROSMAP, Mayo, MSBB) and Neff et al. (MSBB) reference datasets.
If you use this pipeline, please cite:
Pandey RS. Brain Multi-Omics Analysis Pipeline. GitHub: https://github.com/pandeyravi15/brain-multiomics-pipeline
And the underlying tools:
Issues and pull requests are welcome. Please open an issue first to discuss major changes.
MIT β see LICENSE for details.
Ravi S Pandey
GitHub: @pandeyravi15