Bio Causal Genomics Transcriptome Wide Association

Version Compatibility

Reference examples tested with: FUSION (head of gusevlab/fusion_twas, scripts dated 2023+), MetaXcan / S-PrediXcan / S-MultiXcan 0.7.5+ (hakyimlab/MetaXcan), PrediXcan model files from PredictDB (GTEx v8 elastic-net + MASHR), UTMOST (head of Joker-Jerome/UTMOST), pyfocus 0.8+ (bogdanlab/focus), MA-FOCUS (head of mancusolab/ma-focus), TIGAR-V2 (head of yanglab-emory/TIGAR), PLINK 1.9 + PLINK 2.0, R 4.3+, Python 3.9-3.11.

Before using code patterns, verify installed versions match. If versions differ:

• R: Rscript --version; for FUSION scripts inspect --help flags directly in the source
• Python: pip show metaxcan pyfocus then SPrediXcan.py --help, SMulTiXcan.py --help, focus finemap --help
• CLI: plink2 --version; FUSION ships as R scripts not a binary

If a script throws an error about an argument that has moved (e.g. --gwas_file vs --gwas-file) or a model database schema change, introspect the installed script with --help and adapt rather than retrying. PredictDB model file paths change with GTEx version; pin the version explicitly in scripts.

Transcriptome-Wide Association

"Find genes whose predicted tissue expression is associated with my GWAS trait" -> Train SNP -> expression prediction models on a reference eQTL panel, apply the per-gene SNP weights to GWAS summary statistics or genotypes, and produce a gene-level Z-score equivalent to a weighted sum of SNP Z-scores. The output is a gene-by-tissue association, but TWAS is NOT direct evidence of causal mediation: an LD-tagged eQTL signal produces the same statistical association as a truly causal one, and the dominant failure modes are LD-induced false positives at gene-dense loci, tissue mis-specification, and ancestry mismatch between GWAS and prediction weights.

• CLI (sumstat TWAS, R): FUSION.assoc_test.R --sumstats g.sumstats --weights weights.pos --weights_dir wgt/ --ref_ld_chr 1KG/EUR. --chr 22 --out chr22.dat
• CLI (S-PrediXcan, Python): SPrediXcan.py --model_db_path gtex_v8.db --covariance gtex_v8.cov --gwas_file g.txt --output_file out.csv
• CLI (S-MultiXcan joint): SMulTiXcan.py --models_folder mashr_models/ --gwas_folder gwas/ --metaxcan_folder spredixcan_per_tissue/ --output joint.csv
• CLI (UTMOST cross-tissue): joint test across tissues via UTMOST's per-tissue GBJ / GBJ2 step
• CLI (FOCUS fine-mapping): focus finemap chr22.dat 1KG_EUR focus.db --p-threshold 5e-8 --out chr22.focus
• CLI (MA-FOCUS multi-ancestry): focus finemap with colon-separated per-ancestry sumstats / LD / weights and ancestry codes in --locations

TWAS, cis-eQTL MR, and coloc operate on overlapping evidence: TWAS asks "is the gene's predicted expression associated with the trait?"; cis-eQTL MR asks "does the eQTL effect on expression mediate the trait effect under IV assumptions?"; coloc asks "do the GWAS and eQTL share a causal variant?". Strong causal claims require triangulation, not single-method significance.

Algorithmic Taxonomy

| Tool | Model | Input | Output | Strength | Fails when | |------|-------|-------|--------|----------|------------| | FUSION (Gusev 2016 Nat Genet 48:245) | Weighted sum of SNP Z-scores; per-gene multi-model (BLUP, lasso, elnet, top1, bslmm) selected by cross-validation R^2 | GWAS sumstats + pre-computed weights (.pos + per-gene RData) | TWAS Z, p; conditional joint analysis | Mature, ENCODE-style, pre-trained weights for many tissues (GTEx, CMC, METSIM, YFS, NTR) | LD-induced false positives at gene-dense loci; needs ancestry-matched LD reference; weights are heritability-thresholded so low-h2 genes drop | | PrediXcan (Gamazon 2015 Nat Genet 47:1091) | Elastic-net SNP -> expression prediction; individual-level genotypes | PLINK genotypes + GWAS phenotype + GTEx prediction DB | Per-gene association test with full regression machinery | Most flexible (allows covariates, interactions, binary outcomes); same TWAS interpretation | Requires individual-level data; biobank-scale compute | | S-PrediXcan (Barbeira 2018 Nat Commun 9:1825) | Summary-statistic equivalent of PrediXcan; Z-score weighted sum analogous to FUSION | GWAS sumstats + PredictDB model + covariance | Per-gene Z, p | Public PredictDB models (GTEx v8 elastic-net + MASHR-EUR); minimal compute | Pre-trained weights ancestry-specific (EUR primarily); covariance file must match model DB | | S-MultiXcan (Barbeira 2019 PLoS Genet 15:e1007889) | Joint multi-tissue test combining per-tissue S-PrediXcan via PCA-regularised regression | Folder of per-tissue S-PrediXcan outputs + model folder | Single joint p per gene + per-tissue significance | Boosts power when causal tissue is unknown; standard for transcriptome-wide screens | Joint test cannot pinpoint causal tissue; correlated tissues produce ill-conditioned regression | | UTMOST (Hu 2019 Nat Genet 51:568) | Cross-tissue elastic-net (group lasso) for joint tissue prediction + GBJ test | Per-tissue eQTL data + GWAS sumstats | Cross-tissue joint statistic | Often better-powered than S-MultiXcan at cross-tissue genes | Computationally heavier than MetaXcan; tissue weights are less interpretable | | MOSTWAS (Bhattacharya 2021 PLoS Genet 17:e1009398) | Mediator-aware TWAS adding distal trans-mediating SNPs to cis-only models | GWAS sumstats + MOSTWAS weights | Per-gene Z, p (TWAS + distal mediator extension) | Recovers signal at genes with non-cis genetic regulation | Trans-mediation models need large reference panel; weights less broadly available | | kTWAS (Cao 2021 Brief Bioinform 22:bbaa270) | Kernel-based TWAS using SKAT-style aggregation | GWAS sumstats or genotypes + per-gene SNP set | Per-gene p | Robust to non-linear and rare-variant contributions | Loses the eQTL-weighting interpretability; less standardised | | EpiXcan (Zhang 2019 Nat Commun 10:3300) | Adds epigenome-derived per-SNP prior to weight training | eQTL + epigenome + GWAS sumstats | Per-gene Z, p (epigenome-informed) | Higher prediction R^2 in epigenome-rich tissues | Requires matched epigenome data for the prediction tissue | | TIGAR-V2 (Tang 2021 NAR 49:e64) | Dirichlet process regression (non-parametric Bayes) for SNP -> expression | Reference eQTL + GWAS sumstats | Per-gene Z, p | Captures non-elastic-net effect structures; more flexible weight learning | Slower training; benefits depend on locus genetics | | FOCUS (Mancuso 2019 Nat Genet 51:675) | Probabilistic gene-level fine-mapping over TWAS Z-scores using gene-by-gene predicted-expression correlation as analog of LD | FUSION/S-PrediXcan TWAS Z + ancestry-matched LD reference | Per-gene PIP + credible gene set | Resolves co-regulated gene clusters into a probabilistic causal gene; standard add-on after TWAS | Requires the same prediction-weight panel that produced TWAS Z; PIPs depend on prior | | MA-FOCUS (Lu 2022 AJHG 109:1388-1404) | Multi-ancestry FOCUS; joint gene fine-mapping across ancestries with shared causal-gene assumption | Per-ancestry TWAS sumstats + per-ancestry weights + per-ancestry LD | Cross-ancestry PIP | Smaller credible gene sets when AFR/EAS contribute non-EUR LD information | Trans-ethnic gene-effect heterogeneity violated; weights must be ancestry-matched | | JEPEG / JEPEG-Mix (Lee 2015 / 2016) | Gene-based test combining eQTL and functional weights | GWAS sumstats + JEPEG annotation database | Per-gene p | Lightweight gene-burden alternative to TWAS | Less granular than full PrediXcan/FUSION machinery; minimally updated |

Methodology evolves; the FUSION-vs-PrediXcan landscape has been stable but probabilistic fine-mapping (FOCUS, MA-FOCUS) and integrative methods (MOSTWAS, EpiXcan, OmicsXcan) continue to advance. Verify against the current PredictDB release notes (predictdb.org) and the latest FUSION weight panels (gusevlab.org/projects/fusion) before locking on a tissue or model.

S-PrediXcan and FUSION are mathematically near-identical: both compute a weighted sum of GWAS SNP Z-scores using per-gene SNP-expression weights and an LD-aware variance correction. The practical differences reduce to (a) the weight panel (FUSION elastic-net vs PredictDB MASHR), (b) the LD reference, and (c) the per-gene heritability threshold; the supplement of Barbeira 2018 works through the algebra. Method choice should therefore be driven by panel availability and ancestry-match, not by the underlying algorithm.

PredictDB Model Choice and GTEx Versioning

| Release | Cohort | Status (2026) | When to use | |---------|--------|----------------|-------------| | GTEx v8 (838 donors, 49 tissues, 2020) | EUR-dominant (~85%) | Current PredictDB standard | Default; pre-trained MASHR + elastic-net DBs available | | GTEx v9 (2023) | Expanded harmonization | Not migrated into PredictDB | Do not use until PredictDB rebuilds | | GTEx v10 (2024 AnVIL release) | Re-aligned to GRCh38 v44 | Limited harmonization; not PredictDB-default | Wait for community-validated weight panels |

Operational rule: Use GTEx v8 unless there is an explicit biological reason to deviate (tissue not in v8, ancestry-specific panel preferred). In methods, pin exactly: "GTEx v8 MASHR-EUR, PredictDB release 2022-01".

MASHR vs Elastic-Net Models

PredictDB ships two cross-validated model families per tissue (Barbeira 2020 Genome Biol 21:36):

| Model | Construction | Per-gene SNP count | When to use | |-------|--------------|---------------------|-------------| | MASHR | Cross-tissue posterior mean from DAP-G fine-mapped SNPs | ~10x sparser | Primary discovery; higher per-gene R^2 in most genes; standard for S-MultiXcan | | Elastic-net | Per-tissue lasso/ridge mix (alpha = 0.5) | Denser | Tissues where MASHR's cross-tissue prior is mis-specified (ovary, testis, isolated-organ traits) |

Operational rule: Never mix MASHR and elastic-net within a single S-MultiXcan run; the inter-tissue covariance and condition number assumptions break. Pick one family and apply consistently across all tissues.

Decision Tree by Experimental Scenario

| Scenario | Recommended workflow | Why | |----------|---------------------|-----| | GWAS summary stats only, EUR, single hypothesis tissue (e.g. liver for LDL) | S-PrediXcan with GTEx v8 MASHR-EUR weights, or FUSION with GTEx liver | Standard pre-trained pipeline; minimal compute | | GWAS summary stats only, tissue unknown a priori | S-MultiXcan (standard GTEx v8) OR UTMOST (custom panel) -- see S-MultiXcan vs UTMOST table | Joint multi-tissue inflates power; tissue prioritisation requires LDSC-SEG separately | | Multiple TWAS hits at one locus (gene-dense region) | Run TWAS then FOCUS for probabilistic fine-mapping | LD ties co-regulated genes; FOCUS PIP distinguishes likely causal gene | | Multi-ancestry GWAS (EUR + EAS + AFR) | Per-ancestry S-PrediXcan with matched weights, then MA-FOCUS to combine | Single-ancestry weights miscalibrated in other ancestries; joint fine-mapping shrinks credible set | | Individual-level genotypes available (UKB) | PrediXcan (full regression) | Allows covariates, interactions, binary outcomes natively | | Low-N tissue (GTEx N --ref-ld-chr-cts --w-ld-chr ` against 200+ tissue-specific gene expression annotations

1. CELLEX (Timshel 2020 eLife 9:e55851): single-cell tissue / cell-type prioritization on the same GWAS
2. MAGMA gene-property analysis (de Leeuw 2016 PLoS Comput Biol 12:e1004219): cheaper substitute when LDSC unavailable

Operational rule: Primary TWAS tissue = the tissue with FDR-significant enrichment in at least two of the three methods. Run secondary tissues in S-MultiXcan for cross-tissue replication. Bonferroni for tissue selection alone: 0.05 / 200 annotations = 2.5e-4.

S-MultiXcan vs UTMOST

| Method | Use case | Rationale | |--------|----------|-----------| | S-MultiXcan (Barbeira 2019) | Standard GTEx v8 analysis | Pre-computed MASHR weights; lower compute barrier; PCA-regularised inter-tissue regression | | UTMOST (Hu 2019 Nat Genet 51:568) | Custom eQTL panel with cross-tissue retraining | Higher power at genes with shared cross-tissue eQTL architecture; group-lasso enforces sparsity across tissues |

Benchmarks: Hu 2019, Barbeira 2019. Choose by panel availability first; the methods recover overlapping but non-identical gene sets.

Single-Cell and Cell-Type-Resolved TWAS

Bulk-tissue TWAS averages over cell composition; sc-eQTL TWAS recovers cell-type-specific regulation but at lower per-cell-type power.

| Panel | Reference | Cells / tissue | |-------|-----------|----------------| | OneK1K | Yazar 2022 Science 376:eabf3041 | PBMC, ~982 donors, 14 cell types | | HipSci iPSC-eQTL | Kilpinen 2017 Nature 546:370 | iPSC, ~317 donors | | BLUEPRINT | Chen 2016 Cell 167:1398 | Monocytes, neutrophils, T cells |

Tooling state (2026): No fully pre-built scPrediXcan equivalent to MASHR; sc-eQTL weights are panel-specific. Train custom PredictDB or use TIGAR-V2's Bayesian DPR on the sc-eQTL matrix.

Operational rule: Run standard bulk-tissue TWAS first; run sc-eQTL TWAS in the prioritized cell type as a secondary analysis; require concordance between bulk and sc results before nominating a cell-type-specific gene. Upstream sc preprocessing: cross-reference single-cell/preprocessing.

Per-Tool Failure Modes

LD-induced TWAS false positives (most common pitfall)

Trigger: Two or more genes at the same locus have correlated cis-eQTLs (shared causal eQTL SNP or LD-linked eQTL SNPs).

Mechanism: TWAS Z-scores are linear combinations of SNP Z-scores weighted by per-gene SNP effects. When two genes share many high-weight SNPs (e.g. nearby genes regulated by the same enhancer or LD-tagged independent eQTLs), their TWAS Z-scores are positively correlated. A single causal GWAS variant therefore produces significant Z at multiple co-regulated genes (Wainberg 2019 Nat Genet 51:592; Mancuso 2019 Nat Genet 51:675).

Symptom: A GWAS lead locus shows 3-10 genes all passing genome-wide TWAS significance (p = 0.8 as candidate causal; report co-significant genes with PIP = 0.7) for the candidate causal gene. FUSION's FUSION.post_process.R --joint performs conditional analysis as a lighter-weight alternative.

Tissue mis-specification

Trigger: Running TWAS in a tissue that does not host the causal regulatory effect (e.g. whole blood for a psychiatric trait; pancreas for an LDL trait).

Mechanism: A gene's cis-eQTL effect size varies across tissues; a wrong-tissue model has weaker per-gene prediction R^2 and lower power. Conversely, eQTL effects in the wrong tissue can still tag the GWAS signal via LD and produce spurious associations not present in the causal tissue.

Symptom: Strong TWAS signal in a tissue biologically irrelevant to the trait; null in the expected tissue; tissue-prioritisation methods (LDSC-SEG, Finucane 2018 Nat Genet 50:621; CELL-TYPE-SPECIFIC LDSC, Calderon 2017) disagree with the TWAS tissue.

Fix: Run S-MultiXcan to combine tissues if causal tissue is unknown. For prioritisation, use LDSC-SEG / CELLEX / EWCE on the GWAS sumstats independently of TWAS, and report TWAS in the prioritised tissues. Never report a single-tissue TWAS hit as causal without independent tissue evidence (single-cell eQTL, chromatin accessibility in matched cell type).

Ancestry mismatch in prediction weights

Trigger: Running TWAS on a non-EUR GWAS using GTEx (~ 85% EUR) weights, or vice versa.

Mechanism: Cis-eQTL effect sizes and LD structure are ancestry-specific; prediction weights trained in one ancestry transfer with reduced R^2 and biased Z-scores in another. Power is lost preferentially at loci where the causal eQTL is not shared across ancestries (Patel 2022 Genome Med 14:31).

Symptom: Genome-wide TWAS hit count much lower than expected given GWAS power; non-EUR-specific GWAS loci fail to produce TWAS hits; per-gene prediction R^2 substantially reduced.

Fix: Use ancestry-matched prediction panels where available: MESA multi-ethnic eQTL (Mogil 2018 PLoS Genet), eQTLGen-Asian,

…

Source & license

This open-source skill is cataloged on AgentStack and links to its original source — we do not rehost the code.

• Author: GPTomics
• Source: GPTomics/bioSkills
• License: MIT

Install and usage instructions live in the source repository linked above.