Multivariate Latent-Factor Genetic Analysis
Summary¶
Multivariate latent-factor genetic analysis compresses many correlated observed traits into a smaller set of statistically uncorrelated latent factors, calculates GWAS summary statistics for those factors, and fine-maps their genetic signals individually and jointly. The approach implemented by flashfmZero can recover pleiotropic associations that are weak in every single observed-trait GWAS and can use cross-factor sharing to shrink credible sets [1].
Why the Method Exists¶
High-dimensional phenotypes such as complete blood counts, molecular measurements and metabolic profiles contain structured correlation. Related measurements may reflect a common cell lineage, shared physiology, measurement process or biological pathway. Treating every trait independently has three costs:
- a pleiotropic variant may remain below the significance threshold in every univariate GWAS even when its combined pattern is strong;
- repeated single-trait tests obscure the lower-dimensional biological processes that drive groups of traits;
- separate fine-mapping analyses cannot borrow evidence from related association patterns and may retain large, redundant credible sets.
The latent-factor workflow addresses these problems by first estimating the shared phenotypic structure and then carrying the factor representation into genetic association analysis and Bayesian fine-mapping.[1]
Conceptual Model¶
Observed traits are modeled as combinations of latent factors plus trait-specific residual variation. The factor-loading matrix describes how strongly each latent factor contributes to each observed trait. In this protocol, maximum-likelihood factor analysis with varimax rotation produces a zero-correlation factor representation, making the latent factors parsimonious and reducing redundant genetic signals.[1]
The method does not treat a latent factor as a directly measured biological quantity. Its meaning is inferred from the pattern of factor contributions across observed traits, so interpretation requires reporting the loading or contribution structure alongside the genetic results.
End-to-End Data Flow¶
Observed quantitative traits
|
+--> phenotypic covariance/correlation matrix
|
Observed-trait GWAS summary statistics
|
v
Trait QC and cross-GWAS allele harmonization
|
v
Maximum-likelihood factor analysis + varimax rotation
|
+--> factor loadings and trait contributions
|
v
Latent-factor GWAS summary statistics
|
+--> genome-wide association discovery
|
v
Select associated latent factors at a locus
|
+--> aligned in-sample or ancestry-matched LD
|
v
Single-factor JAM fine-mapping + joint flashfmZero fine-mapping
|
v
Credible sets, variant MPPs, SNP-group probabilities and interpretation plots
Required Inputs¶
| Input | Role | Critical requirements |
|---|---|---|
| Individual-level observed traits | Estimates missingness and the phenotypic covariance structure | Quantitative traits measured in a single ancestry group for the documented workflow |
| GWAS summary statistics for every retained trait | Supplies variant effects and uncertainty | Common variant identifiers, consistent columns, alleles and genome build |
| Factor-loading matrix | Maps observed traits to latent factors | Derived from the covariance matrix and kept in the same trait order as GWAS inputs |
| Sample size | Used by factor selection and fine-mapping | Protocol uses the median observed-trait sample size for factor analysis and the actual GWAS sample size for fine-mapping |
| Linkage disequilibrium matrix | Separates correlated candidate variants | Preferably calculated in-sample; otherwise ancestry-matched and allele-aligned |
| Regional genotype data, if LD is not precomputed | Generates the locus-specific LD matrix | Variant order and IDs must match the summary statistics |
Trait Quality Control¶
The protocol applies several filters before factor analysis:
- Remove traits with excessive missing measurements. In the INTERVAL example, the threshold is 50%; for samples smaller than 20,000, the authors recommend considering 25% for more stable calculations.
- Calculate pairwise-complete Pearson correlations among retained traits.
- Remove redundant traits until no block of three retained traits has every pairwise absolute correlation above 0.99. Traits are prioritized for removal by the number of high-correlation connections plus their missing-data fraction.
- Restrict each GWAS to a common chromosome or locus and harmonize all GWASs to the same effect allele.
- Retain only variants present in every observed-trait GWAS. The documented defaults are minor-allele frequency at least 0.005 and imputation INFO at least 0.4.[1]
The complete-case requirement is on variants across GWAS files, not on participants across all traits. This lets the latent-GWAS calculation accommodate missing trait measurements while still requiring a common set of variants.[1]
Factor Selection and Interpretation¶
The number of factors is chosen with Horn's parallel analysis, using the trait covariance matrix and median trait sample size. Maximum-likelihood factor analysis is then fit with varimax rotation. The factor_contributions() output provides:
- a contributions matrix, with observed traits in rows and latent factors in columns;
- a reordered factor-loading matrix used by the latent-GWAS transformation;
- groups of traits ordered by their dominant contributing factor.
Contribution heat maps should accompany association results. The protocol displays traits whose contribution from a fine-mapped factor is at least 10%, which reveals block structure without implying that smaller contributions are exactly zero.[1]
For observed trait $i$ and latent factor $j$, the primary paper defines the scaled contribution as:
$$C_{ij} = \frac{L_{ij}^{2}}{\sum_{k=1}^{L} L_{ik}^{2}}$$
where $L_{ij}$ is the factor loading. Contributions therefore sum to one across fitted factors for each observed trait and describe the factor's share of variance explained by the full factor solution.[2]
Biological Structure of the Blood Factors¶
Factor analysis of 99 classical and non-classical blood-cell traits in 18,310 complete-case INTERVAL participants selected 25 orthogonal factors. Traits with the same leading factor usually belonged to the same broad cell class, but the factors separated biological dimensions within a class.[2]
- ML4 primarily represented reticulocyte measurements.
- ML5 was platelet-specific, while ML10 tracked platelet count without parallel changes in volume or other platelet properties.
- ML8 represented basophil traits.
- ML17 tracked reactive lymphocytes.
- ML22 jointly captured volume, granularity and nucleic-acid content in neutrophils and eosinophils, suggesting shared granulocyte biology.
- ML23 represented coordinated lymphocyte and monocyte changes in cell volume, complexity and RNA content.
The structure is many-to-many. A measured trait can receive major contributions from several factors, and traits from the same cell class can have different leading factors. This is biologically useful because it separates cell abundance from morphology, maturity, granularity and nucleic-acid content rather than collapsing an entire lineage into one score.[2]
Latent-Factor GWAS¶
latentGWAS() combines harmonized observed-trait GWAS effects, standard errors, effect-allele frequencies, the observed-trait covariance matrix and the factor loadings. Its output is a GWAS summary-statistics data frame for each latent factor. In the worked example, the transformation produced summary statistics for 25 factors across 891,327 chromosome 1 variants in under five minutes on an Intel Xeon Platinum 8368Q processor.[1]
A latent-factor association is interpreted as evidence that a variant affects the shared trait pattern represented by that factor. It does not by itself establish which observed trait mediates the effect or whether the factor corresponds to a single biological pathway.
Summary-Statistics Transformation¶
Let $\hat{\boldsymbol{\beta}}$ be the vector of observed-trait GWAS effect estimates for one variant and $L$ the factor-loading matrix. The latent-factor effects are estimated as:[2]
$$\hat{\boldsymbol{\beta}}^{*} = (L^{T}L)^{-1}L^{T}\hat{\boldsymbol{\beta}}$$
Their covariance is propagated from the observed-trait effect covariance:
$$\operatorname{Var}(\hat{\boldsymbol{\beta}}^{*}) = (L^{T}L)^{-1}L^{T}\operatorname{Var}(\hat{\boldsymbol{\beta}})L(L^{T}L)^{-1}$$
For observed traits $i$ and $j$, the off-diagonal covariance uses their phenotypic covariance and GWAS standard errors. Association for factor $j$ is tested with:
$$Z_j = \frac{\hat{\beta}^{}_j}{\sqrt{\operatorname{Var}(\hat{\boldsymbol{\beta}}^{})_{jj}}}$$
This transformation is why individual-level complete cases are unnecessary once the trait covariance, factor loadings and harmonized trait-GWAS statistics are available.[2]
Genetic Discovery Gains¶
In 18,310 participants, latent-factor GWAS detected 31 LD clumps with genome-wide significant factor associations but no genome-wide significant association for any connected observed blood trait. Across 3,399 lead variant-factor associations, 211 had only suggestive evidence in connected traits and 152 had no trait association even at the suggestive threshold.[2]
After conditioning on 3,559 lead variants from larger published blood-cell GWASs, six signals remained genome-wide significant for latent factors while remaining sub-threshold in every single blood-trait GWAS from the same participants. The biologically annotated signals included:[2]
| Factor and lineage | Variant / candidate locus | Biological connection |
|---|---|---|
| ML8, basophils | rs9310935 near IL5RA | IL-5 receptor signaling controls eosinophil and basophil differentiation; the factor aggregated moderate evidence across several basophil traits |
| ML18, platelets | two low-frequency variants near C7 | Complement component 7 participates in the terminal complement pathway and has experimental links to platelet aggregation |
| ML8, basophils | rs72878322 near ZNF521 | ZNF521 regulates erythroid and B-lymphoid differentiation programs |
| ML6, mature red cells | rs6064377 near FAM210B | The mitochondrial membrane protein is linked to iron import, heme synthesis and erythropoiesis |
| ML6, mature red cells | rs117617749 in a coding-gene desert | Mechanism unresolved; retained as a genetic association rather than assigned a fabricated target |
These results demonstrate the gain from aggregating a coherent pattern, but the candidate-gene assignments remain hypotheses unless supported by fine-mapping and functional evidence.
Fine-Mapping Strategy¶
Only factors with a genome-wide significant association in the target region are fine-mapped; the protocol uses p < 5 x 10^-8 and stops if no factor passes that screen. The locus workflow then:
- aligns latent-factor summary statistics to the alleles and ordering used for LD;
- hard-calls or imports regional genotypes and filters variants with more than 20% missing genotype data;
- calculates pairwise LD correlations;
- runs
JAMdwithGroups()when one factor is associated; - runs
FLASHFMZEROwithJAMd()when multiple factors are associated, returning both single-factor and joint results; - constructs credible sets at a specified credibility threshold, with 99% as the software default and 95% used in the protocol example.[1]
Variant identifiers must begin with a letter, contain none of :, < or >, and match exactly between GWAS and LD inputs. This formatting constraint is operational rather than biological, but violations are a common cause of failed runs.[1]
Outputs and Metrics¶
For each factor and analysis approach, FMsummary_table() reports:
- credible-set size;
- the variant with the highest marginal posterior probability (MPP);
- the top variant's MPP;
- its high-LD SNP group and group size;
- SNP-group marginal posterior probability (MPPg);
- variants with MPP greater than 0.90;
- SNP groups with MPPg greater than 0.90.
MPP is the probability that a variant is causally associated with the factor after marginalizing over compatible fine-mapping models. It is the flashfmZero terminology corresponding to the causal probability commonly summarized as posterior inclusion probability. MPPg aggregates causal probability over a high-LD group whose members remain statistically interchangeable.[1]
The protocol groups variants with non-negligible MPP when r-squared exceeds 0.8. These groups depend jointly on LD and the fitted fine-mapping results; joint-factor groups can therefore be smaller than groups from single-factor analysis.[1]
SMIM1 Worked Example¶
The example analyzes blood traits from the INTERVAL BioResource and fine-maps a chromosome 1 region spanning 3,441,528 to 3,959,487 around SMIM1. Four latent factors—ML12, ML14, ML20 and ML4—had significant signals.[1]
| Analysis | Factor | 95% credible-set interpretation | Top variant and probability |
|---|---|---|---|
| Single-factor JAM | ML12 | 1 variant | rs1175550 / chr1_3691528, MPP 0.993 |
| Single-factor JAM | ML14 | 37 variants | chr1_3709093, MPP 0.378 |
| Single-factor JAM | ML20 | 4 variants | chr1_3691434, MPP 0.677 |
| Single-factor JAM | ML4 | 5 variants | chr1_3691528, MPP 0.920 |
| Joint flashfmZero | all four factors | 1 variant for every factor | rs1175550 / chr1_3691528, MPP 1.000 |
The joint analysis converted three ambiguous single-factor results into a shared one-variant credible set. This illustrates the central gain: consistent association patterns across latent factors can distinguish a shared candidate variant from neighboring variants that remain plausible when each factor is analyzed alone.[1]
The factors remain biologically distinct despite converging genetically. ML20 primarily contributed to Delta-HGB, ML12 to MCHC, ML14 to red-cell-distribution measures and macrocyte ratio, and ML4 to reticulocyte fluorescence and maturity traits.[1]
The primary 18,310-participant analysis used 99% credible sets. Observed-trait sets containing rs1175550 ranged from 30 to 58 variants; single-factor latent sets ranged from 5 to 27; joint flashfmZero reduced every associated factor to rs1175550 alone with MPP above 0.99. UK Biobank independently prioritized the same variant for nine red-cell traits, and published data identify it as an SMIM1 eQTL and modulator of Vel antigen expression.[2]
Comparative Fine-Mapping Evidence¶
Across 217 blood-cell regions, the paper compared univariate JAMdynamic on observed traits, JAMdynamic on latent factors and joint flashfmZero on associated factors.[2]
| Comparison | Equal or smaller CS99 | Strictly smaller CS99 |
|---|---|---|
| Single latent factor vs observed blood trait | 76% (937/1,238) | 58% (725/1,238) |
| Joint flashfmZero factor vs observed blood trait | 86% (624/725) | 71% (517/725) |
| Joint flashfmZero vs single latent factor | 97% (205/211) | 45% (95/211) |
In 36 regions with high-confidence variants in both INTERVAL and UK Biobank, 25 had a matching prioritized variant. Four additional regions matched UK Biobank only through latent-factor fine-mapping, not through observed-trait fine-mapping. The comparison supports better resolution and external concordance, but it is not a randomized benchmark: regions and matched trait-factor pairs were selected under the study's signal criteria.[2]
The PIEZO1 region illustrates the value of separating shared from distinct signals. Joint analysis prioritized rs551118 for red-cell factors ML12 and ML13 while leaving a basophil factor with a different candidate signal, whereas mvSuSiE needed a much larger sample to recover comparably specific trait-to-variant assignments.[2]
Generalization to Metabolic Traits¶
The method was also applied to 184 Nightingale Health NMR traits measured in up to 40,849 INTERVAL participants. Factor analysis selected 21 factors: some were highly specific, such as acetone or conjugated linoleic acid, while several captured broad lipid and lipoprotein structure.[2]
Across eleven one-megabase regions, latent-factor fine-mapping generally narrowed observed-trait sets and joint flashfmZero further improved resolution. Examples included the PCSK9 p.R46L loss-of-function variant for pro-atherogenic lipids, a CPT1A missense variant for omega-3 fatty acids, and separation of at least three signals near TOMM40L/APOA2.[2]
Relationship to Disease-Subtyping Decomposition Methods¶
A related family of methods decomposes multi-trait genetic architecture into latent components for disease subtyping rather than fine-mapping resolution — see Genomic SEM, DeGAs / dPRS, and Pleiotropic Decomposition Regression, surveyed together in Mechanism-First Polygenic Risk Score Clustering. These methods share the general strategy of replacing many correlated observed quantities with a smaller set of latent components, but target patient/mechanism stratification rather than credible-set narrowing.
Relationship to Other Multitrait Methods¶
| Method | Trait representation | Main modeling stage | Best suited to |
|---|---|---|---|
| flashfmZero latent-factor workflow | Orthogonal latent factors derived from observed-trait covariance | Latent GWAS followed by single and joint Bayesian fine-mapping | Many correlated traits where dimension reduction and interpretable shared structure are desired |
| mvSuSiE | Observed traits retained directly | Joint multivariate fine-mapping with learned effect-sharing patterns | Estimating trait-specific and shared causal effects without first replacing traits by factors |
| Separate univariate fine-mapping | Each observed trait analyzed alone | Independent locus models | Baseline analysis or traits with little shared structure |
The methods answer related but non-identical questions. flashfmZero asks which variants drive latent patterns of phenotypic variation; direct multitrait methods ask which variants affect the original traits and how their effects are shared.
Assumptions and Limitations¶
- The documented protocol is for quantitative traits from one genetic ancestry group.
- LD must come from the analysis sample or a closely ancestry-matched reference.
- GWAS and LD alleles, IDs and variant order must be aligned exactly.
- A variant must be present in every retained observed-trait GWAS to receive latent-factor summary statistics.
- Factor solutions depend on trait selection, covariance estimation, missingness, the number-of-factors rule and rotation.
- Orthogonality is a statistical construction; biological mechanisms can overlap even when fitted factors are uncorrelated.
- Factor labels are not intrinsically biological and require interpretation from contribution patterns.
- Improved credible-set resolution does not prove molecular mechanism or functional causality.
- The protocol's timing estimates are hardware- and dataset-specific.
- The summary-statistics formulation requires every analyzed variant to be available for every retained trait.
- The current method assumes one genetic ancestry; cross-ancestry factors must first be shown to represent equivalent structures before meta-analysis or multi-ancestry fine-mapping.
- Binary traits require linear-model statistics or an appropriate log-odds transformation and correlation model; rare-disease behavior has not been established.
Practical Reporting Checklist¶
A reproducible analysis should report:
- ancestry, sample sizes and trait measurement platform;
- missingness and high-correlation removal rules;
- genome build, allele harmonization and variant QC thresholds;
- factor-number selection method, estimator and rotation;
- the complete loading or contribution matrix;
- latent-factor GWAS sample-size assumptions;
- locus boundaries and genome-wide signal-selection threshold;
- LD source, ancestry, variant ordering and genotype QC;
- single-factor and joint credible sets at the same credibility level;
- variant MPPs, group MPPg values and group LD threshold;
- sensitivity to factor number, trait inclusion and LD reference where feasible.
Citations¶
[1] Astle, Butterworth and Asimit (2025), "Protocol for genetic discovery and fine-mapping of multivariate latent factors from high-dimensional traits" [2] Zhou et al. (2025), "Improved genetic discovery and fine-mapping resolution through multivariate latent factor analysis of high-dimensional traits"