Introduction

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.

Which haplotype carries the rare allele?

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.

Local ancestry at a rare variant in three people. Only person C, with one African and one European haplotype, needs phase to tell which ancestry the allele belongs to.

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:

SHAPEIT5 switch-error rate by minor allele count in the UK Biobank (Hofmeister et al. 2023).

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).

What the split gains

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.

Power at a significance level of 0.05 by mean effect size, averaged over 100 genes (Wang et al. 2026). AF, EU: tests of the African and European dosages; Phased, Unphased: 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.

Workflow

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.

Installation

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("lantern")

library(lantern)

Example data

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"

Step 1: split dosages by local ancestry

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

Step 2: write per-ancestry GDS files

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/RtmpFyncRG/lantern-13ad4f42ca4c/AFR.gds"
#> 
#> $EUR
#> [1] "/tmp/RtmpFyncRG/lantern-13ad4f42ca4c/EUR.gds"

Step 3: test each gene in each ancestry

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.

References

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.

Session information

sessionInfo()
#> R version 4.6.2 (2026-10-29)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 26.04.1 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.32.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] 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] digest_0.6.39          rpart_4.1.27           lifecycle_1.0.5       
#> [10] survival_3.8-12        magrittr_2.0.5         compiler_4.6.2        
#> [13] rlang_1.3.0            sass_0.4.10            tools_4.6.2           
#> [16] yaml_2.3.12            data.table_1.18.6.1    knitr_1.52            
#> [19] SeqArray_1.53.6        purrr_1.2.2            BiocGenerics_0.59.12  
#> [22] sys_3.4.3              GWASExactHW_1.2        nnet_7.3-21           
#> [25] grid_4.6.2             stats4_4.6.2           jomo_2.7-6            
#> [28] GMMAT_1.5.0            mice_3.19.0            iterators_1.0.14      
#> [31] MASS_7.3-66            cli_3.6.6              rmarkdown_2.32        
#> [34] crayon_1.5.3           reformulas_0.4.4       generics_0.1.4        
#> [37] otel_0.2.0             minqa_1.2.8            cachem_1.1.0          
#> [40] operator.tools_1.6.3.1 splines_4.6.2          parallel_4.6.2        
#> [43] BiocManager_1.30.27    XVector_0.53.0         vctrs_0.7.3           
#> [46] boot_1.3-32            glmnet_5.2             Matrix_1.7-6          
#> [49] jsonlite_2.0.0         IRanges_2.47.5         S4Vectors_0.51.10     
#> [52] mitml_0.4-5            maketools_1.3.2        foreach_1.5.2         
#> [55] tidyr_1.3.2            jquerylib_0.1.4        glue_1.8.1            
#> [58] nloptr_2.2.1           pan_2.0                codetools_0.2-20      
#> [61] shape_1.4.6.1          GenomicRanges_1.65.4   lme4_2.1-0            
#> [64] tibble_3.3.1           pillar_1.11.1          htmltools_0.5.9       
#> [67] Seqinfo_1.3.2          R6_2.6.1               SeqVarTools_1.51.3    
#> [70] Rdpack_2.6.6           formula.tools_1.7.1    evaluate_1.0.5        
#> [73] lattice_0.23-1         Biobase_2.73.2         rbibutils_2.4.1       
#> [76] backports_1.5.1        broom_1.0.13           bslib_0.12.0          
#> [79] Rcpp_1.1.2             nlme_3.1-172           mgcv_1.9-5            
#> [82] logistf_1.26.1         xfun_0.62              buildtools_1.0.0      
#> [85] pkgconfig_2.0.3

Getting started with lantern