People with ancestry from more than one population, such as African Americans, European, are called admixed. In each of their chromosomes, some segments were inherited from African ancestors and others from European ancestors. Software such as RFMix (Maples et al. 2013) can infer, for every segment in each of the two copies of a chromosome, which population it came from. This is so called local ancestry.
Admixed people are often left out of genetic association studies, or analysed as if their whole genome came from one population. Yet the same allele can relate to a trait differently depending on the ancestral background it sits on. Linkage disequilibrium differs between populations, so a variant may tag a nearby causal variant on African haplotypes but not on European ones. Allele frequencies differ, so a variant can be common on one background and almost absent on the other. And an allele’s effect can depend on the other variants carried on the same haplotype. Tractor (Atkinson et al. 2021) uses local ancestry to estimate a separate effect for each ancestry in genome-wide association studies of common variants, and showed that this gives accurate ancestry-specific effect sizes, boosts power, and places association signals closer to the causal variants.
Rare-variant association tests have the same problem. The rare variants in a gene are grouped and tested whether carrying them is associated with a trait. If the variants act on African haplotypes but not on European ones (or the reverse), pooling all carriers mixes carriers with and without the effect, the signal is diluted, and the test loses power.
To bridge the gap, lantern aims to use local ancestry to split each person’s count of alternate alleles (their dosage, 0, 1 or 2) into one dosage per ancestry. A heterozygote whose allele sits on an African haplotype contributes 1 to the African dosage and 0 to the European dosage. lantern then runs the gene-based SMMAT test from the GMMAT package separately for each ancestry, and combines the per-ancestry p-values for each gene into one p-value with the Cauchy combination test.
For most carriers the split is clear. If both copies of the chromosome are African at the variant (person A below), the allele is African; if both are European (person B), it is European. The hard case is a heterozygote with one African and one European copy (person C): the allele is on one of the two, and only phase, the assignment of each allele to a chromosome copy, says which.
Phasing is least reliable for exactly the variants that rare-variant tests use. In the UK Biobank, SHAPEIT5 places a singleton allele on the wrong haplotype about one time in three:
lantern therefore offers two ways to split. mode = "haplotype" uses the
phased genotypes, like Tractor does for common variants. mode = "dosage"
does not use phase: it shares each ambiguous heterozygote’s allele between
the two ancestries, with shares estimated from the variant’s other carriers
(vignette("split-intuition") explains how).
We simulated traits on chromosome 19 of 3,313 participants of the Jackson Heart Study, about 80% of whose ancestry is African, with causal effects carried by alleles on African haplotypes, on European haplotypes, or on both. The ancestry-split tests did not inflate the type I error, down to a significance level of \(2.5 \times 10^{-6}\). When the effect was carried on the European (here the minority-ancestry) haplotypes, testing the European dosage had about twice the power of SMMAT on unsplit genotypes, and the unphased split about 1.5 times. When the effect was on African haplotypes or on both, the African-ancestry test lost almost nothing against SMMAT.
mode = "haplotype" and mode = "dosage".The details are in the preprint: Wang Y, Tuftin B, Raffield LM, Hidalgo B, Kerns SL, DeWan AT, Leal SM, Auer PL (2026). LANTERN: Leveraging local ancestry tracts to enhance rare-variant aggregate association testing. medRxiv. doi: 10.64898/2026.04.24.26351693.
The analysis has three steps as follows:
ancestry_split() reads a phased VCF and an RFMix local-ancestry file
and returns one dosage matrix per ancestry.write_ancestry_gds() saves each matrix as a SeqArray
GDS file, the format GMMAT reads.ancestry_smmat() runs SMMAT on each GDS file and combines the
results.if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("lantern")
library(lantern)
lantern includes a small simulated dataset: 60 admixed individuals, 200
rare variants on chromosome 19 and eight genes. We simulated a continuous
trait y that is raised by alternate alleles in GENE2, but only when the
allele sits on an African haplotype. A good analysis should find GENE2
in the African-ancestry test and not in the European one.
ext <- system.file("extdata", package = "lantern")
vcf <- file.path(ext, "toy.vcf.gz")
msp <- file.path(ext, "toy.msp.tsv")
The VCF holds phased genotypes. A | separates the two copies of the
chromosome, so 0|1 means the alternate allele is on the second copy:
vcf_lines <- readLines(vcf)
data_lines <- vcf_lines[!startsWith(vcf_lines, "##")]
substr(data_lines[1:3], 1, 60)
#> [1] "#CHROM\tPOS\tID\tREF\tALT\tQUAL\tFILTER\tINFO\tFORMAT\tS001\tS002\tS003"
#> [2] "chr19\t1168394\t.\tT\tA\t.\tPASS\t.\tGT\t0|0\t0|0\t0|0\t0|0\t0|0\t0|0\t0|0\t"
#> [3] "chr19\t1240171\t.\tT\tC\t.\tPASS\t.\tGT\t0|0\t0|0\t0|0\t0|1\t0|0\t0|0\t0|0\t"
The MSP file is RFMix’s local-ancestry output. Its first line gives the
population codes (here 0 = AFR and 1 = EUR). Each later row is a segment
of the chromosome, with one ancestry call per sample and chromosome copy
(columns S001.0, S001.1, …):
substr(readLines(msp, n = 4), 1, 70)
#> [1] "#Subpopulation order/codes: AFR=0\tEUR=1"
#> [2] "#chm\tspos\tepos\tsgpos\tegpos\tn snps\tS001.0\tS001.1\tS002.0\tS002.1\tS003.0\tS"
#> [3] "chr19\t1000000\t4999999\t1.00\t5.00\t25\t0\t0\t0\t0\t0\t0\t0\t0\t0\t1\t0\t1\t0\t0\t1\t0\t0\t0"
#> [4] "chr19\t5000000\t8999999\t5.00\t9.00\t25\t0\t0\t1\t0\t0\t0\t0\t1\t0\t0\t0\t0\t0\t1\t0\t0\t0\t0"
ancestry_split() matches samples between the two files, finds the
ancestry segment that contains each variant, and splits the dosages.
split <- ancestry_split(vcf, msp, mode = "dosage", chrom = "chr19",
verbose = FALSE)
The result has one variants x samples dosage matrix per population, named after the populations in the MSP file:
split$AFR[1:4, 1:6]
#> [,1] [,2] [,3] [,4] [,5] [,6]
#> [1,] 0 0 0 0 0 0
#> [2,] 0 0 0 1 0 0
#> [3,] 0 0 0 0 0 0
#> [4,] 0 0 0 0 0 0
split$EUR[1:4, 1:6]
#> [,1] [,2] [,3] [,4] [,5] [,6]
#> [1,] 0 0 0 0 0 0
#> [2,] 0 0 0 0 0 0
#> [3,] 0 0 0 0 0 0
#> [4,] 0 0 0 0 0 0
Most entries are 0 because the variants are rare. mode = "dosage"
ignores phase. For a heterozygote with one African and one European
haplotype, the genotype alone does not say which haplotype carries the
allele, so the allele is shared between the two ancestries, with shares
estimated from the other carriers of the same variant. Here is the first
such heterozygote in the data; its two shares add up to its one allele:
i <- which(split$AFR %% 1 > 0, arr.ind = TRUE)[1, ]
c(AFR = split$AFR[i[1], i[2]], EUR = split$EUR[i[1], i[2]])
#> AFR EUR
#> 0.875 0.125
vignette("split-intuition") explains how the shares are calculated. With
phased data you can instead use mode = "haplotype", which assigns every
allele to the ancestry of its own haplotype and never produces fractions.
The result also lists the variants that were kept. Variants with no alternate allele in these samples are dropped (16 here):
head(split$variant_info)
#> chrom pos ref alt
#> 1 chr19 1168394 T A
#> 2 chr19 1240171 T C
#> 3 chr19 2252809 A G
#> 4 chr19 2330989 C T
#> 5 chr19 2436455 A G
#> 6 chr19 2459485 C G
split$overlap[c("n_common", "n_variants_kept", "n_monomorphic_filtered")]
#> $n_common
#> [1] 60
#>
#> $n_variants_kept
#> [1] 184
#>
#> $n_monomorphic_filtered
#> [1] 16
GMMAT reads genotypes from GDS files, so each dosage matrix is saved to its own file. The function returns the file paths, named by population.
out_dir <- tempfile("lantern-")
gds <- write_ancestry_gds(split, out_dir, verbose = FALSE)
gds
#> $AFR
#> [1] "/tmp/Rtmpl5TObQ/lantern-39e8564132d886/AFR.gds"
#>
#> $EUR
#> [1] "/tmp/Rtmpl5TObQ/lantern-39e8564132d886/EUR.gds"
SMMAT needs three more inputs: a phenotype table, a kinship matrix and a gene-group file. The phenotype table has one row per sample, with the trait and covariates:
pheno <- read.delim(file.path(ext, "toy_pheno.tsv"))
head(pheno, 3)
#> id y age sex
#> 1 S001 1.1948 62 1
#> 2 S002 2.2396 62 1
#> 3 S003 7.6596 64 0
The kinship matrix describes how related the samples are. Our simulated samples are unrelated, so we use the identity matrix. In a real study this would be a genetic relatedness matrix:
kinship <- diag(nrow(pheno))
dimnames(kinship) <- list(pheno$id, pheno$id)
The gene-group file assigns variants to genes. It has no header and six
columns: gene, chromosome, position, reference allele, alternate allele
and a weight. These are the 25 variants of GENE2:
genes <- read.delim(file.path(ext, "toy_genes.tsv"), header = FALSE,
col.names = c("gene", "chr", "pos", "ref", "alt", "weight"))
head(genes[genes$gene == "GENE2", ])
#> gene chr pos ref alt weight
#> 26 GENE2 19 5249792 A T 1
#> 27 GENE2 19 5436416 G C 1
#> 28 GENE2 19 5496995 T A 1
#> 29 GENE2 19 5510001 G C 1
#> 30 GENE2 19 5601633 G T 1
#> 31 GENE2 19 5969726 G T 1
table(genes$gene)
#>
#> GENE1 GENE2 GENE3 GENE4 GENE5 GENE6 GENE7 GENE8
#> 25 25 25 25 25 25 25 25
ancestry_smmat() fits the null model (the trait explained by the
covariates and relatedness, without genotypes) once, then runs SMMAT for
each ancestry. Passing ancestry_counts and variant_info from Step 1
weights the Cauchy combination by how many samples of each ancestry carry
a gene’s variants, so the ancestry with more data counts for more.
res <- ancestry_smmat(gds, pheno, y ~ age + sex, kinship,
gene_group_file = file.path(ext, "toy_genes.tsv"),
ancestry_counts = split$ancestry_counts,
variant_info = split$variant_info,
verbose = FALSE)
res$results
#> gene p_AFR p_EUR w_AFR w_EUR p_cauchy
#> 1 GENE1 8.269778e-01 0.2250896 42 1 8.212171e-01
#> 2 GENE2 8.625395e-09 0.7771531 39 2 9.067722e-09
#> 3 GENE3 5.162240e-01 0.2822809 41 3 4.974372e-01
#> 4 GENE4 3.977574e-01 0.6121023 40 1 4.026767e-01
#> 5 GENE5 3.039969e-01 0.2189521 40 2 2.989101e-01
#> 6 GENE6 6.483093e-01 0.5622577 33 2 6.438543e-01
#> 7 GENE7 1.595176e-01 0.4246312 41 1 1.623335e-01
#> 8 GENE8 5.571825e-01 0.3854922 33 3 5.427531e-01
The table has one row per gene: p_AFR and p_EUR are the SMMAT p-values
from each ancestry, w_AFR and w_EUR the weights, and p_cauchy the
combined p-value. GENE2 stands out:
res$results[res$results$gene == "GENE2", ]
#> gene p_AFR p_EUR w_AFR w_EUR p_cauchy
#> 2 GENE2 8.625395e-09 0.7771531 39 2 9.067722e-09
Its African-ancestry p-value is very small, its European-ancestry p-value is not, and the combined p-value keeps the signal. For comparison, we test the same genes without splitting by ancestry. Adding the two dosage matrices gives back each sample’s total dosage:
pooled_gds <- file.path(out_dir, "pooled.gds")
write_dosage_gds(split$AFR + split$EUR, split$variant_info,
split$sample_ids, pooled_gds)
pooled <- ancestry_smmat(c(POOLED = pooled_gds), pheno, y ~ age + sex,
kinship, file.path(ext, "toy_genes.tsv"),
verbose = FALSE)
pooled$results[pooled$results$gene == "GENE2", c("gene", "p_POOLED")]
#> gene p_POOLED
#> 2 GENE2 2.38154e-08
The pooled test also detects GENE2, with a p-value about three times
larger: carriers whose allele sits on a European haplotype add to the
pooled dosage without adding to the effect.
vignette("split-intuition") works through the split on small
matrices, for two and for more populations.vignette("toy-vcf-msp-example") shows the VCF and MSP formats in more
detail, including a four-population example.vignette("real-data-chr19") analyses chromosome 19 in a real cohort
with binary, continuous and count traits.Atkinson EG, Maihofer AX, Kanai M, Martin AR, Karczewski KJ, Santoro ML, Ulirsch JC, Kamatani Y, Okada Y, Finucane HK, Koenen KC, Nievergelt CM, Daly MJ, Neale BM (2021). Tractor uses local ancestry to enable the inclusion of admixed individuals in GWAS and to boost power. Nature Genetics 53(2):195–204. doi: 10.1038/s41588-020-00766-y.
Maples BK, Gravel S, Kenny EE, Bustamante CD (2013). RFMix: a discriminative modeling approach for rapid and robust local-ancestry inference. American Journal of Human Genetics 93(2):278–288. doi: 10.1016/j.ajhg.2013.06.020.
sessionInfo()
#> R version 4.6.1 (2026-06-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.5 LTS
#>
#> Matrix products: default
#> BLAS: /home/biocbuild/bbs-3.24-bioc/R/lib/libRblas.so
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0 LAPACK version 3.12.0
#>
#> locale:
#> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
#> [3] LC_TIME=en_GB LC_COLLATE=C
#> [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: America/New_York
#> tzcode source: system (glibc)
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] lantern_0.99.1 BiocStyle_2.41.0
#>
#> loaded via a namespace (and not attached):
#> [1] tidyselect_1.2.1 dplyr_1.2.1 Biostrings_2.81.9
#> [4] gdsfmt_1.49.12 fastmap_1.2.0 CompQuadForm_1.4.4
#> [7] TH.data_1.1-5 digest_0.6.39 rpart_4.1.27
#> [10] estimability_2.0.0 lifecycle_1.0.5 survival_3.8-12
#> [13] magrittr_2.0.5 compiler_4.6.1 rlang_1.3.0
#> [16] sass_0.4.10 tools_4.6.1 yaml_2.3.12
#> [19] data.table_1.18.6.1 knitr_1.52 multcomp_1.4-32
#> [22] SeqArray_1.53.6 purrr_1.2.2 GWASExactHW_1.2
#> [25] BiocGenerics_0.59.12 nnet_7.3-21 grid_4.6.1
#> [28] stats4_4.6.1 jomo_2.7-6 GMMAT_1.5.0
#> [31] xtable_1.8-8 mice_3.19.0 emmeans_2.0.4
#> [34] iterators_1.0.14 MASS_7.3-66 cli_3.6.6
#> [37] mvtnorm_1.4-2 rmarkdown_2.32 crayon_1.5.3
#> [40] reformulas_0.4.4 generics_0.1.4 otel_0.2.0
#> [43] minqa_1.2.8 cachem_1.1.0 operator.tools_1.6.3.1
#> [46] splines_4.6.1 parallel_4.6.1 BiocManager_1.30.27
#> [49] XVector_0.53.0 vctrs_0.7.3 boot_1.3-32
#> [52] glmnet_5.2 Matrix_1.7-6 sandwich_3.1-3
#> [55] jsonlite_2.0.0 bookdown_0.48 IRanges_2.47.5
#> [58] S4Vectors_0.51.10 mitml_0.4-5 foreach_1.5.2
#> [61] tidyr_1.3.2 jquerylib_0.1.4 glue_1.8.1
#> [64] nloptr_2.2.1 pan_2.0 codetools_0.2-20
#> [67] shape_1.4.6.1 GenomicRanges_1.65.4 lme4_2.1-0
#> [70] tibble_3.3.1 pillar_1.11.1 htmltools_0.5.9
#> [73] Seqinfo_1.3.2 R6_2.6.1 SeqVarTools_1.51.3
#> [76] Rdpack_2.6.6 formula.tools_1.7.1 evaluate_1.0.5
#> [79] lattice_0.23-1 Biobase_2.73.2 rbibutils_2.4.1
#> [82] backports_1.5.1 broom_1.0.13 bslib_0.12.0
#> [85] Rcpp_1.1.2 coda_0.19-4.1 nlme_3.1-172
#> [88] mgcv_1.9-5 logistf_1.26.1 xfun_0.62
#> [91] zoo_1.9-1 pkgconfig_2.0.3