CARMA: Novel Bayesian model for fine-mapping in meta-analysis studies
Zikun Yang1 , Chen Wang1,2 , Linxi Liu3 , Atlas Khan2 , Annie Lee4 , Badri Vardarajan4 , Richard Mayeux4 ,Krzysztof Kiryluk2 , Iuliana Ionita-Laza1,#

1

Department of Biostatistics, Columbia University, New YorkDivision of Nephrology, Department of Medicine, College of Physicians and Surgeons, Columbia University, New York3
Department of Statistics, University of Pittsburgh, Pittsburgh, PA
4
Department of Neurology, College of Physicians and Surgeons, Columbia University, New York#
Correspondance: ii2135@columbia.edu
2

Abstract
We propose a novel Bayesian model for fine-mapping in order to identify putative causal variants at GWAS
loci. Relative to existing fine-mapping methods, the proposed model has several appealing features, such as
flexible specification of the prior distribution of effect sizes, joint modeling of summary statistics and large
number of functional annotations, and accounting for discrepancies between summary statistics and external
linkage disequilibrium values in meta-analysis settings. Using simulations, we compare performance with
commonly used fine-mapping methods, including SuSiE and fastPAINTOR, and show that the proposedmodel has higher power and lower FDR when including functional annotations, and higher power, lower
FDR and higher coverage for credible sets in meta-analysis settings. We further illustrate our approachby applying it to a meta-analysis of Alzheimer’s Disease GWAS data where we prioritize putatively causal
variants and genes, including TREM2, ABI3, CASS4, SORL1, APH1B, BIN1, KAT8, ABCA7, INPP5D,SLC24A4, NDUFS2, MPO, PTK2B, EPHA1, CD33, RABEP1, ADAM10, CD2AP, PICALM, PILRB, and
CLNK. Finally, our model can be used more generally in conjunction with other fine-mapping methods suchas SuSiE to minimize false-positive findings in meta-analyses.

Introduction
Meta-analyses of GWAS studies have identified a large number of significant loci. Fine-mapping is the natural
next step in order to identify putative causal genetic variants at these loci. Meta-analyses however poseseveral challenges that can invalidate the results from existing fine-mapping methods as also pointed out in.1
For example, using linkage disequilibrium (LD) from external panels can create inconsistencies with GWAS
summary statistics which can lead fine-mapping methods to prioritize non-causal variants.2 Similarly, unevensample size coverage at different variants can lead to biased posterior inclusion probability (PIP) values. To
illustrate this point we show the example of a GWAS locus SCIMP (SLP adaptor and CSK interacting membraneprotein) for Alzheimer’s disease (AD). We use summary statistics from a large meta-analysis GWAS of clinicallydiagnosed AD and AD-by-proxy with 71,880 cases and 383,378 controls of European ancestry from three
consortia.3 The LD matrix is estimated using individuals of European descent in UK Biobank (UKBB).4 Two
fine-mapping models, SuSiE5 and fastPAINTOR,6 prioritize variants with low Z-scores due to discrepancies
between summary statistics and LD values, whereas the results of the proposed model appear as expected(Figure 1(a-c)). In this paper we introduce a new fine-mapping method that improves the power and reduces
false positives in such situations. Moreover, the proposed method can be used to improve results from existing
fine-mapping methods, such as SuSiE (Figure 1(d)).
There are many statistical fine-mapping methods in the literature.5–13 Most of the existing methods canwork with summary statistics and LD information from relevant reference panels, and make certain assumptions
on the number of causal variants, i.e. some methods restrict the number of possible causal variants,7, 8, 11, 13 whileothers relax this assumption by introducing a prior distribution on model space, and implementing stochastic
1

algorithms such as Markov Chain Monte Carlo (MCMC) to reduce the computational cost.6, 9, 10 Multiplicity
control is another important aspect of any Bayesian fine-mapping method in order to control the false discovery
rate in the context of multiple testing,14 and several existing methods address this issue formally by introducing
prior probabilities on model space.9, 10Here we propose a new Bayesian model, CARMA (CAusal Robust mapping method in Meta-Analysis
studies), that attempts to improve upon existing methods especially in complex meta-analysis settings as
described above. Our proposed model has several technical innovations: (1) it allows in addition to the usual
Normal-Gamma prior family for the effect size distribution the option of a heavy-tail Cauchy distribution,
which leads to a prior distribution of effect sizes that is more adaptive to the possibly different signal-to-noise
ratios across loci; (2) it jointly models summary statistics and high-dimensional functional annotations - thisis different from other recent models such as PolyFun12 that first estimate prior causal probabilities based onfunctional annotations and summary statistics, and then apply fine-mapping methods such as SuSiE with theseprior probabilities; (3) it introduces a novel Bayesian hypothesis testing approach to account for discrepancies
between summary statistics and LD from external reference panels in order to avoid an increase in false positives.
We illustrate the proposed method using simulations and applications to an AD GWAS meta-analysis.

Results
Overview of the proposed model
For a given locus, we assume a standard linear model y = Xβ + ,  ∼ MVN(0, σy2 In ), where y is a n × 1
standardized vector of quantitative phenotype values, X is a standardized n × p genotype matrix, β is a p− 1dimensional vector of effect sizes of variants, and  is a Gaussian noise vector. Let Z = nσy2 2 X 0 y denote− 1
the vector of single-SNP Z-scores, then E [Z|β, X] = nσy2 2 X 0 Xβ and Var [Z|β, X] = X 0 X/n. We propose
a method that can work with summary statistics, i.e. single-SNP Z-scores from a standard GWAS analysis,
and LD from an external reference panel as previously done in GCTA,15 fastPAINTOR,6 and FINEMAP.10√
Let Σ = X 0 X/n denote the LD matrix, and λ = nβ. Assuming σy2 = 1, then the sampling distribution of Z
can be written as:
Z|λ, Σ ∼ MVN(Σλ, Σ).
Let γ 0 = {0, 1}p denote an indicator vector, such that γi = 1 iff λi 6= 0 (βi 6= 0). Given any γ, we assume aspike-and-slab prior:1λγ |τ, γ ∼ MVN(0, Iγ ),τwhere λγ represents the vector of non-zero entries of λ, i.e. γi = 1. The value of the prior precision parameterτ may be tuned using empirical evidence (i.e. with 95% probability, the causal variants at a locus explain lessthan 1% of the trait variation). We want to identify the true model that generated the summary statistics
through posterior inference within a Bayesian paradigm.
We note that the precision parameter τ can be assigned a prior distribution to yield a mixture of normal
distributions on λ. This option provides more flexibility and robustness in situations when tuning τ to beconsistent with the underlying model at a given locus may be difficult (e.g. for loci very different effects).
Specifically, integrating out τ against a prior distribution Gamma(0.5, 0.5) yields a heavy-tailed Cauchy priordistribution on λ, which provides an adaptive Bayesian inference under different signal-to-noise ratios. Although
we focus our description in the main text on the slab-and-spike prior, we derive the corresponding results for the
Cauchy prior distribution in the Supplemental Material, and show simulation results (both the spike-and-slab
and Cauchy priors are provided as options in the software package implementing CARMA).
Pp
We
Ppadditionally assume a truncated Poisson prior on the size of the model, i.e. i=1 γi ∼ Truncated Poisson (η)
and i=1 γi ∈ {0, 1, . . . , p}, to control the total number of causal variants assumed by a given model, and hence
provide multiplicity control (or control of the false discovery rate). We implement a Shotgun stochastic search2algorithm16 for sampling from the posterior distribution over the model space, which has advantages over the
more commonly used MCMC-based methods in that it leads to a more complete exploration of the areas withhigh posterior probability in the model space. This semi-exhaustive search feature alleviates the problem ofunequal PIPs for perfectly correlated variants that can happen with other MCMC algorithms.
High-dimensional functional annotations. CARMA seeks to distinguish among highly correlated variants
by leveraging functional annotations. Let W be a p × (q + 1) matrix for p SNPs and q functional annotations
and a p-dimensional 1 vector. We maximize the joint likelihood of the summary statistics and functionalannotations
L(θ; Z, γ, W ) = m(Z|γ)Pr (γ|W , θ) ,where θ is a (q + 1)-dimensional vector of the corresponding coefficients, through an expectation–maximization(EM) algorithm with the prior probability Pr (γ|W , θ) being linked to functional annotations via a penalized
logistic model.17
Bayesian hypothesis testing for meta-analysis settings. A principal advantage of the proposed modelover existing fine-mapping models is in the context of meta-analyses when discrepancies between GWAS summary statistics and LD from external reference panels can lead to false positive results when using existing
methods as also noted in.1 A first type of discrepancy is due to the external LD panel itself; a second type isdue to sample size heterogeneity across SNPs in the meta-analysis. To account for these possible discrepancies,
we propose a Bayesian hypothesis testing procedure for detecting and removing inconsistencies, leading to a
reduction in false discoveries relative to existing fine-mapping models. We also show that this procedure canbe combined with genotype imputation to reduce the effect of the second type of discrepancy due to sample
size inconsistencies across SNPs. More details are available in the Methods section.
Credible sets and credible models. In5 the authors define a credible set as the smallest subset of correlated
variants (with correlation within the set greater than some threshold r) that has probability ρ or greater ofcontaining at least one causal variant, e.g. the ith and jth SNPs can form a credible set if Pr (γi = 1|Z) +Pr (γj = 1|Z) ≥ ρ and cor (Zi , Zj ) ≥ r. The credible set provides a small set of SNPs for follow-up studies.
However, as also mentioned in,5 a SNP should not be considered as non-causal if it is not included inany credible set. For example, there are circumstances where a group of highly correlated SNPs containing
one causal SNP might not result in a sufficiently large PIP to form a credible set. Hence, we introduce theconcept of credible model at a locus, as a complement to credible sets, based on the top candidate models interms of posterior probability. This is similar to FINEMAP10 which provides configuration-specific posterior
probabilities. Let γ (1) denote the leading model that receives the largest posterior probability among the visited
candidate models; then we define as credible model at a locus the set of candidate models such that the posterior
odds between the leading model and the models in the set is smaller than a pre-determined threshold, such as
3.2 or 10 as recommended in.18 We view the credible model and the variants identified by a credible model ascomplementary to credible sets. Compared to credible sets, credible models typically involve less variants while
identifying a higher proportion of causal variants as shown in simulations below.

Simulations
We perform simulations to investigate the performance of CARMA and competitor methods, fastPAINTOR
and SuSiE, first in a homogeneous setting with in-sample LD and then in more realistic heterogeneous settings
for meta-analyses with LD estimated from external reference panels.
Generating simulated datasets (in-sample LD)
Genotype simulation. We use the R package ‘sim1000G’19 to simulate genotypes based on the 1000
Genomes Project data (phase 3, European population). To select regions representative of GWAS loci in3

terms of size and LD structure, we focus on 94 loci identified as risk regions in a recent GWAS on breast
cancer.20 Within each region, we filter out rare variants (MAF< 0.01); the number of variants in each region
ranges between ∼ 1, 500 − 4, 000. We simulate genotype data for n = 10, 000 individuals.
Causal SNP selection. To generate prior causal probabilities we use data on 200 DeepSEA chromatinfeatures.21 Let θ denote a 201-dimensional coefficient vector including the intercept term; we randomly select100 chromatin features as being related to the causal status of variants, with the coefficients being sampledfrom N (0, 0.12 ). Hence the impact of each individual annotation is weak. We used a simulation setting similar
to that in PolyFun, i.e. the prior causal probability for each SNP is proportional to a linear combination ofannotations, such as Pr (γi = 1) ∝ wi θ, where wi is the vector of annotations. The causal SNPs are selectedamong those with highest prior probabilities in such a way that the absolute value for the pairwise Pearson
correlation between any two causal variants is at most 0.1. For a given locus, let T denote the index set of thetrue causal SNPs, i.e., if i ∈ T , then γi = 1. In simulations we assume |T | ∈ {1, 2, 3} at each locus (results for|T | = 2 are shown in the main manuscript; results for |T | ∈ {1, 3} are shown in Supplemental Figures).
Phenotype generation. For each i ∈ T , the effect sizes of the causal SNPs are drawn independently fromβi ∼ N (0, 0.52 ), and for all i ∈/ T we set βi = 0. The phenotypic variance σy2 is computed such that φ = 0.0075
at each locus, where φ = σ2Var(Xβ). Then we sample y such that y = Xβ + ;  ∼ N(0, σy2 In×n ).
+Var(Xβ)yRemarks on some implementation details:
1. We assume two scenarios (1) no functional annotation and (2) with functional annotations (200 DeepSEA
chromatin features21 ). SuSiE can only include one annotation in the form of prior probabilities. We
estimated PolyFun prior probabilities using the PolyFun package in R, then we provided the estimatedprior probabilities to SuSiE. fastPAINTOR cannot handle high-dimensional annotations therefore weselect a subset of the ten most informative annotations as follows: we first compute correlations between
the true linear predictor wi θ and each functional annotation, then we iteratively select top correlated
annotations such that no annotation has absolute correlation > 0.3 with previously selected annotations.
Then, fastPAINTOR is run based on the selected annotations.
2. For CARMA, the value of the hyper-parameter η is set at the default setting 1 throughout the simulations.
The prior precision parameter τ (in the case of the spike-and-slab prior) is set at 0.0512 corresponding
to an assumed explained variance at the locus φ = 0.0075 and sample size n = 10, 000. More simulation
studies for different values of η and τ can be found in the Supplemental Material (Figures SI2, SI3, andSI5).
3. We run fastPAINTOR and CARMA models at a chromosome level, while SuSiE is run one locus at atime (since SuSiE does not have the option to aggregate multiple loci).
4. We run each model with the default number of maximum causal variants per locus assumed by eachmodel, i.e. CARMA and SuSiE assume 10 causal variants and fastPAINTOR assumes 2 causal variants
per locus.
Credible sets. For SuSiE, we use the credible sets reported by SuSiE. For CARMA and fastPAINTOR wecompute credible sets as in5 (Methods). We assess different measures of performance with respect to credible
sets from each model, including:Power: The overall proportion of simulated causal variants included in any credible set.
Coverage: The proportion of credible sets that contain a causal variant.
Size: The number of variants included in a credible set.
Purity: The mean squared correlation of variants in a credible set.
4

Results for ρ = 0.99 are shown in Figure 2(a) (see Figure S1 for ρ = 0.95). Without annotations, CARMAand SuSiE have similar performance. In contrast, fastPAINTOR has lower power, coverage and purity. The
substantially lower coverage of the credible sets from fastPAINTOR suggest the need of including penalizationof the size of the selected model to control the false positives. When functional annotations are included,
the power of CARMA increases relative to the scenario with no functional annotation, and relative to thecompetitor models, SuSiE and fastPAINTOR. These results suggest that the inclusion of functional annotations has a greater impact for the models that jointly integrate summary statistics and functional annotations(such as CARMA and fastPAINTOR) compared to SuSiE+PolyFun which is based on two separate modelingstages. Therefore algorithms that maximize the marginal likelihood conditional on both summary statistics
and annotations as proposed here may be preferable.
Credible models. Credible models provide complementary information to credible sets. Therefore, we compare results based on the variants included in the credible model in CARMA and the variants included in
all the credible sets for CARMA, SuSiE, and fastPAINTOR at a locus. Compared to credible sets, credible
models generally achieve higher power with smaller sets of selected variants, which may facilitate follow-up
investigations (Figure 2(b)).
Power and FDR based on positive predictions. We also compare methods in terms of power and false
discovery rate (FDR) using positive predictions, i.e. PIP is greater than a given threshold. With no functionalannotation, the performance of CARMA and SuSiE are comparable whereas fastPAINTOR has higher FDR
and lower power (Figure 2(c)). With functional annotations, CARMA performs best with higher power andlower FDR. As also noted above, there is greater improvement in performance (higher power and lower FDR)
for CARMA and fastPAINTOR when introducing functional annotations relative to SuSiE.
PIPs for perfectly correlated SNPs. Methods based on stochastic algorithms such as MCMC often fail toevenly explore the area of the posterior model space which creates problems when there are perfectly correlated
SNPs. Due to the semi-exhaustive search feature of the Shotgun algorithm, CARMA examines the neighborhoodof the currently selected model, including candidate models that exchange perfectly correlated SNPs. This way,
the resulting PIPs tend to be more consistent between highly correlated SNPs. To illustrate this point, we
examine the standard deviation of the PIPs within groups of perfectly correlated SNPs based on the resultsof the three models when no functional annotation is included. Across the 94 loci, there are 30,131 groups ofperfectly correlated SNPs with an average of 317 groups at each locus, and 4.23 SNPs per group. Since thevalues of the PIPs returned by each model can vary across different models, we standardized the PIPs withineach group by dividing the PIPs by the maximum PIP in each group. CARMA has a standard deviation of
0.06 for PIPs of perfectly correlated SNPs. This is in contrast to fastPAINTOR which has a standard deviationof 0.493. Note that SuSiE has identical PIPs for perfectly correlated SNPs by definition due to its algorithm
that runs univariate regression individually.
Effect of LD-summary statistics inconsistencies in meta-analyses
In complex meta-analyses studies, with LD estimated from external reference panels, it is not uncommon to havediscrepancies between Z-scores and LD values, which leads to biased PIPs for existing fine-mapping models, andfalse prioritization of non-causal variants as also illustrated in Figure 1. Furthermore, more complications canarise due to sample size heterogeneity across different SNPs. Allele flipping (i.e. alleles are encoded differently
in the study and reference panel) can also create problems for fine-mapping.22 Here we show the robustness of
CARMA in such scenarios.
Specifically, we start with the same simulated datasets as above, i.e. 94 GWAS loci with 10,000 simulated
individuals, but switch from in-sample LD to UKBB LD. We randomly separate the individuals into three
non-overlapping groups with 1,000, 1,500, and 7,500 individuals. We compute summary statistics for eachgroup, i.e. Z 1 , Z 2 , and Z 3 . We meta-analyze results across groups using the inverse-variance based method in
METAL.23 We first compute Z-scores using a consistent meta-analysis by meta-analyzing all three groups for5each SNP. We then consider an inconsistent meta-analysis, where the final dataset is constructed in such a waythat the Z-scores of 75% of the SNPs are generated by meta-analyzing all three groups, the Z-scores for 15%of the SNPs are based on the two groups with 1,000 and 1,500 individuals, and the remaining 10% are only
based on the group with the 1,000 individuals. Note that in the consistent meta-analysis setting, all Z-scores
are computed based on the complete set of 10,000 individuals, and therefore the use of external UKBB LD in
the only source of discrepancies. In the inconsistent meta-analysis setting, both the use of external LD andheterogeneity in sample size across SNPs leads to discrepancies between Z-scores and LD values.
Power, Coverage and FDR. CARMA’s performance is significantly less affected by discrepancies between
Z-scores and LD values compared with SuSiE and fastPAINTOR in terms of both power and coverage (Figure 3(a)). In the inconsistent setting, the coverage for CARMA is also reduced but still much higher than for
SuSiE and fastPAINTOR (61% vs. 10% and 22%) due to sample size heterogeneity across SNPs. Using imputation before fine-mapping largely alleviates this issue, with the performance of CARMA being similar to thatin the consistent setting. Without imputation, the outlier detection method will drop a significant proportionof SNPs especially due to the sample size heterogeneity across SNPs, whereas imputation greatly lowers the
dropping rate (Table SI1). Similar results hold when looking at credible models (Figure 3(b)). In terms of
FDR, CARMA performs substantially better than SuSiE and fastPAINTOR (Figure 3(c)). In particular, for
SuSiE, there are many non-causal SNPs that are mistakenly assigned PIPs close to 1 (see below).
Credible sets with only one SNP. A particular situation arises when a credible set contains only one
SNP (i.e. its PIP exceeds ρ). Such scenarios require that a SNP has a large summary statistic and be
relatively independent of all the other SNPs. Discordant LD/Z-score values when using external LD may
create the appearance that a SNP with large summary statistic is weakly correlated to surrounding SNPs.
Figure 4(a) shows the total number of credible sets with one SNP for the three models across the 94 loci, andthe corresponding coverage of these credible sets. As shown, with external UKBB LD, SuSiE has increased
number of credible sets with one SNP, and corresponding poor coverage. Similarly, fastPAINTOR has poor
coverage in such settings. In contrast, CARMA’s performance, especially with the imputation method, appears
robust.
Using CARMA and SLALOM to identify outliers. We note that recently another method, SLALOM,1
has been described for identifying discrepancies between Z-scores and LD values in meta-analyses as discussed
here and flag suspicious loci. SLALOM assumes a single causal variant per locus. For a given region with p
SNPs, SLALOM utilizes the leading SNP (the highest PIP based on approximate Bayes factor (ABF) finemapping) as the presumed causal variant and the corresponding LD correlations to all other SNPs in theregion, and conducts p − 1 χ2 tests for detecting suspicious loci and outlier variants at such loci. However, theassumption of only one causal variant residing in the testing region may be too strong, and may lead to low
recall if there are multiple independent signals in the region.
We ran SLALOM and CARMA on the simulated datasets for the inconsistent meta-analysis scenario, andrecorded the outliers identified by each method. Then, we ran SuSiE on the datasets with outliers removed by
SLALOM, and CARMA respectively. Figure 4(b) show that the outlier removal by SLALOM does not make asignificant impact on the performance of SuSiE. On the other hand, the performance of SuSiE has been greatly
improved by removing outliers identified by CARMA. These results suggest that CARMA can successfully
identify the problematic outliers and other fine-mapping methods can benefit as well.
In our analyses SLALOM only identified a small proportion of outliers, and those outliers are highly correlated to the SNP with the largest absolute Z-score, as the minimum of the absolute value of pair-wise correlations
between the SNP with the largest Z-score and the outliers identified by SLALOM is ∼ 0.76. Therefore, given
the fact that the pairwise LD correlation between causal variants in the simulation is at most 0.1, SLALOM
has missed those outliers/inconsistencies associated with secondary or other independent signals in the region
besides the signal containing the leading SNP.

6

Fine-mapping Alzheimer’s disease GWAS lociWe present fine-mapping results (primarily from CARMA) at 28 GWAS loci identified in a large meta-analysisof clinically diagnosed AD and AD-by-proxy with 71,880 cases and 383,378 controls of European ancestry.3 The
clinically diagnosed AD case-control data are from three consortia: Alzheimer’s disease working group of thePsychiatric Genomics Consortium (PGC-ALZ), the International Genomics of Alzheimer’s Project (IGAP), and
the Alzheimer’s Disease Sequencing Project (ADSP). The AD-by-proxy data are based on 376,113 individualsof European ancestry from UKBB. Note that this is a complex meta-analysis dataset, and the sample sizes
available at different genetic variants can vary from 9,703 to 444,006 depending on which data are available
for each SNP. We use the leading SNP at each locus from the meta-analysis (phase 3) and for the purposesof fine-mapping define the locus as ± 500kb centered around the leading SNP (Figures S7-S34). We do notinclude the HLA, APOE and AC074212.3 loci due to long-range LD in these regions and extreme values ofthe summary statistics, which cause numerous false findings when fine-mapping these regions. We also exclude
the SUZ12P1 locus for the lack of globally significant variant (−log (p-value) ≥ 7.3) in the phase 3 data. In
total we perform fine-mapping at 28 loci. We use the LD matrix from the UKBB provided by PolyFun. We
also perform fine-mapping by including functional annotations provided by PolyFun (including 187 individual
annotations and the PolyFun prior causal probability made available by PolyFun based on a meta-analysisof several UKBB traits12 ). Therefore, in total for CARMA we consider three scenarios: (1) no functionalannotations, (2) PolyFun prior causal probability, and (3) 187 annotations provided by PolyFun plus PolyFun
prior causal probability.12
SuSiE fine-mapping with outliers removed by CARMA. For 12 loci (CR1, BIN1, CLNK, HS3ST1,
CD2AP, ZCWPW1, MS4A6A, PICALM, ADAM10, SCIMP, ABI3, and ABCA7), SuSiE reports multiple
SNPs with PIPs equal or close to 1 due to discrepancies between Z-scores and LD values (Figures S35-S46),
e.g. multiple credible sets with single SNPs with weak Z-scores. We have removed the outliers identified by
CARMA and provided the resulting filtered data to SuSiE. By filtering out the outliers identified by CARMA,SuSiE can be properly executed (Table 1 and Figures S35-S46).
Next we focus on CARMA fine-mapping results.
Enrichment of fine-mapped SNPs in functional categories. To assess the enrichment of fine-mapped
SNPs in several primary functional categories in the genome, we estimated the functional enrichment of finemapped SNPs in credible models for six selected binary categories (including non-synonymous, conserved, 3’
UTR flanking and 5’ UTR flanking, H3K4me3, repressed) where the enrichment is defined as the proportion of
SNPs included in the credible models of CARMA lying in a functional category divided by the proportion of
genome-wide SNPs lying in the same functional category. Even when not using annotations we see evidence forenrichment for variants in credible models for the five active functional categories (non-synonymous, conserved,
3’ UTR flanking and 5’ UTR flanking, H3K4me3,), and depletion in the repressed category (Figure 5). By
integrating PolyFun annotations, the functional enrichments for SNPs in credible models become stronger
especially for non-synonymous, conserved and 3’ UTR/5’ UTR flanking categories.
Fine-mapping results at 28 GWAS loci. We show fine-mapping results for individual loci using the three
CARMA scenarios in Figures S47-S74. For a quick glance at the results we highlight in Table 2 those SNPs
with the largest PIP in the analysis with the PolyFun annotations (even though these SNPs may not necessarily be the causal SNPs). Furthermore, we have used several functional annotation tools including CADD,24
Eigen/EigenPC,25 data on eQTL and sQTL from GTEx26 and mQTL from The Religious Orders Study andMemory and Aging Project (ROSMAP)27 to functionally annotate these SNPs. We also report for each SNP
the gene with the highest V2G score in the Open Targets Genetics Portal. There is a rather high concordance
between the gene with highest score in Open Targets and the gene identified by mQTLs in Brain Dorsolateral
Prefrontal Cortex in ROSMAP (19/28; including BIN1, SORL1, APH1B, ABCA7, NDUFS2, PTK2B, KAT8,
INPP5D, CD2AP, PILRB, EPHA1, CD33, PICALM, RABEP1, CASS4, SLTM, CLNK, RABEP1 and MPO).
Many of these top PIP SNPs are also eQTLs in brain or whole blood tissues in GTEx. Additionally, we have7

identified two rare missense variants (rs616338 in ABI3 and rs143332484 in TREM2) both with PIP> 0.99 in
the analysis including PolyFun annotations. Both of these genes are highly expressed in microglia cells andthe same two variants were identified in a previous AD association study.28 These fine-mapped SNPs and theputative effector genes could be interesting for follow-up investigations.

Discussion
We have proposed here a novel Bayesian fine-mapping method, CARMA, which is designed to prioritize potentially causal variants within GWAS risk loci by leveraging the LD structure and functional annotations available
for variants at the locus under investigation. CARMA allows flexible specification of the prior distribution ofeffect sizes, including a spike-and-slab prior and a heavy-tail Cauchy distribution, and jointly maximizes thelikelihood of summary statistics and functional annotations in a unified EM algorithm with multiplicity control.
Importantly, we introduce a novel Bayesian hypothesis testing approach to account for mismatches between
summary statistics and LD values from external reference panels in order to avoid prioritization of non-causal
variants. Through extensive simulations, we demonstrate that CARMA has higher precision and lower FDR
relative to SuSiE and fastPAINTOR, especially in complex meta-analyses settings with LD from external reference panels. Furthermore, in such settings the outlier observations identified by CARMA can be used to
improve the performance of other fine-mapping methods, including SuSiE.
Our procedure for outlier detection and removal may lead to true causal variants being removed from theanalysis. In particular, when summary statistics at SNPs in a region are based on different sample sizes, apossible solution would be to perform imputation before running CARMA. Although imputation is also affectedby the reference LD being used, and can therefore introduce outliers on its own, we have found that it can beof great help in retaining more SNPs including potentially causal variants in the analyses.
We further illustrate the use of CARMA in a large fine-mapping analysis of 28 GWAS loci for Alzheimer’s
disease identified in.3 The results show that CARMA has the ability to better handle possible discrepancies
between the LD matrix from the UKBB and the summary statistics in the meta-analysis relative to competitor
methods such as SuSiE. In particular, for 21/28 loci we highlight a putative effector gene that has multiplelines of evidence including from CARMA, external eQTL/mQTL databases and effector gene mapping tools
such as V2G in OpenTargets.
Very recently, a new version of SuSiE has been released that is applicable to summary statistics and LD
matrix extracted from reference panels.22 In this context, the authors developed a likelihood ratio test for
identifying a particular type of inconsistency, namely the “allele flip” scenario. Note that our proposed outlier
detection method is more general and deals with broad types of inconsistencies, including those generated bydifferent sample size coverage at different SNPs in the data. We also note that the diagnostic procedure in
SuSiE is a separate step, i.e. not integrated into the main algorithm, and requires the inversion of the entire
LD matrix for identifying one outlier, which is computationally intensive. In our applications to the motivating
example in Figure 1 and other loci with biased PIPs in our analyses, this new implementation has failed tosolve the problem (more details are in the Supplemental Material).
CARMA has been implemented in a computationally efficient R package.

Methods
Basic notations and assumptions
For a given locus, we assume a standard linear model:
y = Xβ + ,  ∼ MVN(0, σy2 In ),where y is a standardized n × 1 vector of quantitative phenotype values, X is a standardized n × p genotype
matrix, β is a p-dimensional vector of effect sizes of SNPs, and  is a Gaussian noise vector. GWAS are usually
performed in a univariate fashion, so that for each variant i with i = 1, . . . , p, we compute the marginal Z-score
8

x0 y
σy n

0

Zi = √ i 2 , as in other fine-mapping methods10 .6 Let the SNP correlation (LD) matrix be Σ = XnX ; notethat Σ can be approximated using the appropriate, population-matched reference panel. Assuming σy2 = 1, a√
common practice in fine-mapping, and λ = nβ (i.e., λi 6= 0 if and only if βi 6= 0), the sampling distributionof Z can be represented as
Z|λ, Σ ∼ MVN(Σλ, Σ).
Details are in the Supplemental Material.
Indicator vector γ and true model. Let γ 0 = {0, 1}p denote an indicator vector, such that λi 6= 0 iff
γi = 1. Also, let S = {i; λi 6= 0} denote the index set of causal SNPs assumed by a particular model, such that
if i ∈ S, thenPp γi = 1 and λi 6= 0. Then, each indicator vector γ S uniquely defines a model with dimension equal
to |S| = i=1 γi . For a particular model γ S , we denote the p-dimensional coefficient vector by λS and denotethe |S|-dimensional subvector of the non-zero entries of λS by λγ S , such that λi 6= 0 if i ∈ S.

Review of existing fine-mapping models
We briefly review several representative fine-mapping models, including JAM, fastPAINTOR and SuSiE.
JAM.9 Building upon prior Bayesian fine-mapping methods such as CAVIARBF13 and FINEMAP,10 JAM
assumes the popular g-prior29 to model the effect sizes for a given model γ S :λγ S |σy2 ∼ MVN 0, gσy2 Σ−1γS ,σy2 ∼ Inv-Gamma (0.01, 0.01) ,where g is a pre-determined constant. JAM assigns a beta-binomial distribution as the prior distribution of
the model size, i.e., the number of the assumed causal SNPs in a model:
Pr ({γ; 1 · γ = s}) ∼ Beta-binomial(s + a, p − s + b),
Pr ({γ; 1 · γ = s})Pr (γ) =,ps

where {γ; 1 · γ = s} is the set of models having the same dimension as s, and all models of the same size are
equally likely. a = 1 and b = 9 are used in applications, corresponding to a mean proportion of truly casual
SNPs of 10%.
Relative to previous fine-mapping models, JAM does not restrict the number of causal SNPs in the region
and provides a computationally efficient method based on Reversible Jump MCMC to explore a wide rangeof candidate models. Also, unlike previous fine-mapping methods which assume that σy2 is a plug-in estimate,
JAM marginalizes out both λ and σy2 to improve the power. One potential issue is that the correlation matrix
Σ is not necessarily positive definite, and JAM cannot work properly in such cases.
fastPAINTOR6 .

fastPAINTOR uses as input Z-scores for individual variants, Zi =β̂i, and the LD matrixse(β̂i )Σ, usually estimated from a reference panel. fastPAINTOR assumes the following sampling distribution for Z:
Z|λS , Σ ∼ MVN(ΣλS , Σ),for a given model γ. The assumed prior distribution of λS can be written as
λS |γ, σλ2 ∼ MVN(0, Σγ ),Σγ = σλ2 Diag(γ) + Diag(σy2 ).

9The distribution of Z after integrating out λS is
Z|γ, Σ ∼ MVN(0, Σ + ΣΣγ Σ)P (γ).
The prior distribution of γi is a Bernoulli distribution, where the probability of being causal can be related tofunctional annotations as follows:
exp {w0i θ}),
γi ∼ Bern(
1 + exp {w0i θ}
where wi is the vector of functional annotations for SNP i. Posterior probabilities for each SNP are calculatedusing an importance sampling algorithm. Note that fastPAINTOR cannot handle large number of functionalannotations, and a pre-selection of a small number of annotations may be necessary.12
SuSiE.5 The SuSiE model works particularly well in situations where variables are highly correlated andthe effects are sparse. It is conceptually different from the other fine-mapping models. The inference in the
SuSiE model is based on multiple basic models, the so-called Single-Effect Regression (SER) models. Under
the assumption of no prior information, the SER model is defined as:
y = Xβ + ,  ∼ MVN(0, σy2 In )β = bγb ∼ N (0, σ02 )γ ∼ Mult(1, π),
where Mult(m, π) denotes the multinomial distribution on class counts that is obtained when m samples are
drawn with class probability π. For the SER model, there is only one non-zero effect, and hence the indicator
vector γ has only one non-zero element. SuSiE assumes by default that π = (1/p, . . . , 1/p), but one could
define π based on prior probabilities derived
 from functional annotations. Calculating the posterior inclusion
probabilitles (PIP) Pr γi = 1|y, X, σy2 , σ02 involves fitting p univariate regression of y on the columns of xi of
X.
The SER model assumes only one causal SNP. For extensions to multiple causal variants within a locus,
the final SuSiE model is built based on L SER models. The idea is to introduce multiple single-effect vectors
β 1 , . . . , β L and construct the overall effect vector β as the sum of these single effects. The model is as follows:y = Xβ + ,  ∼ MVN(0, σy2 In )
β=

L
X

βl

l=1

β l = γ l bl
bl ∼ N (0, σ02 )γ l ∼ Mult(1, π).
SuSiE performs an iterative Bayesian stepwise selection (IBSS) algorithm to fit this model; at each iteration it
0
uses the SER model to estimate β l given
P current estimate of β l0 for l 6= l. Specifically, it fits the SER modelfor β l using the residual r̄ ← y − X l0 6=l β l0 . The result of the SuSiE model consists of L fitted β̂ l , and L 0corresponding PIP vectors αl = Pr β l,1 6= 0|X, y , . . . , Pr β l,p 6= 0|X, y . Then the final PIP is definedasYPr β̂j 6= 0|X, y ≈ 1 −(1 − αl,j ),l∈{1,...,L}

assuming that the βl,j are independent across l = 1 . . . L. SuSiE naturally produces credible sets; a level ρcredible set is defined as a subset of variables that has probability ρ or more to contain at least one causal
variable.

10

SuSiE+PolyFun.12 To fully utilize the potential of the fine-mapping methods that can only take as input
univariate prior causal probability, Weissbrod et al.12 proposed to estimate prior probabilities through the
S-LDSC model,30 which essentially estimates global, genome-wide functional genomic enrichments using only
summary data. Specifically, PolyFun estimates the prior probability for each SNP in proportion to per-SNP
0
heritability estimates: Pr (βi 6= 0|wi ) ∝ varˆ [βi |wi ] = θ̂ wi whereθ̂ :=argmin
θ∈Rq

X
i

χ2i − nqX2θj l(i, j) − nb − 1 + α||θ||2 ;

j=1

2
P
(x0 y)
q is the number of functional annotations, χ2i = in is the χ2 statistic of SNP i, l(i, j) = pk=1 cor (xi , xk ) wj,kis the LD-score for SNP i weighted by functional annotation j, and b measures the contribution of confounding
biases.

In applications, SuSiE+PolyFun and FINEMAP+PolyFun were shown to lead to more causal variant discoveries relative to SuSiE without functional annotations. The advantage of PolyFun is that the contributionof a functional annotation to the heritability is estimated using genome-wide SNPs.
We summarize the main features and assumptions of commonly used fine-mapping methods in the literature
in Table 3.

Proposed Bayesian fine-mapping model: CARMA
We introduce here the details of our model; we refer to it as CARMA (CAusal Robust Mapping method with
Annotations). Given any model γ S , we assume that the summary statistics Z follow a multivariate normal
distribution:
Z|λS , Σ ∼ MVN(ΣλS , Σ).
We want to identify the true model that generated the summary statistics through posterior inference withina Bayesian paradigm.
Spike-and-slab prior distribution on the scaled effect sizes. We assume a spike-and-slab prior for theprior distribution of λi , i.e.,
λi ∼ Normal(0, τ −1 ) if γi = 1,λi = 0 if γi = 0.
Specifically, given an index set S, the prior distribution of the assumed non-zero effect sizes λγ S corresponding
to γ S is
1
λγ S |τ, γ S ∼ MVN(0, Iγ S ),
τwhere Iγ S is a |S| × |S| identity matrix. We may tune the value of the precision parameter τ according to empirical evidence (i.e. such that with 95% probability, a causal SNP explains less than 1% of the trait variation).
Alternatively, we can assign a prior distribution f (τ ) on the precision parameter τ , resulting in more flexibility and robustness to situations where different loci may have different effect sizes and therefore explain different
proportions of trait variance. By assigning a prior distribution on the mixing parameter τ , we intrinsically assign a mixture of normal distribution on the scaled effect size λγ S . The choice of f (τ ) directly impacts thecomputation of posterior probabilities for candidate models. Here, for several reasons including computational
efficiency, we explore one specific form of f (τ ), the Zellner-Siow’s Cauchy prior:31 τ ∼ Gamma( 21 , 21 ). We focus
our derivations below on the spike-and-slab prior, and include the corresponding results for the Cauchy priorin the Supplemental Material.

11

Prior distribution on model space. A fine-mapping locus typically contains a large number of variants
which can be a challenge for model selection. Without an appropriate prior distribution, all candidate models
including the finite-dimensional true model would receive very small posterior probabilities, as also noted in.7
A penalization of the size of the selected model, through the prior distribution, is required for achieving model
selection consistency in such high-dimensional regime as shown in.32, 33
We introduce a prior distribution on model space to control the total number of causal SNPsP
that any
candidate model assumes, which is analogous to controlling the false discovery rate. Let |S| =γi ∈γ S γidenote the total number of causal SNPs for a given γ S . We first place a discrete prior distribution on the
random variable |S| and let Pr (|S||η) denote the p.m.f. with η as a hyperparameter. Given |S|, we assume
that all those models have the same prior probability. Hence, the prior probability of γ S is
Pr (γ S |η) =Pr (|S||η)
.
p
|S|

We propose to use the truncated Poisson prior (|S| ∈ {0, . . . , p}), which has been introduced before andshown to enjoy model selection consistency in.33 Let F (·|η) be the cumulative distribution function of a Poisson
distribution with mean η, the truncated Poisson prior is defined as follows:
Pr (|S||η) =η |S| exp {−η}
η |S| exp {−η}∝.
|S|!F (p|η)|S|!

Since the total number of SNPs p is usually a large number and η is chosen to be small in order to reflect thesparse scenario for the true causal SNPs, F (p|η) → 1. Then given any specific model γ S , the prior probabilityof this model under Poisson distribution is
Pr (γ S |η) ∝η |S| exp {−η} (p − |S|)!
.
p!

The Poisson distribution on model space imposes necessary sparsity for multiplicity control. The hyperparameter η plays a critical role in this sparsity encouraging mechanism. Note that the computation of
the PIPs is based on the unnormalized posterior probabilities for the candidate models, i.e., Pr (γ S |Z) ∝f (Z|γ S )Pr (γ S |η). Given that η can be thought of as the prior expectation of the number of causal SNPs at agiven locus, we recommend η = 1 as a default choice and all our analyses are performed under this assumption.
Marginal likelihood and posterior probability. Given the Z-scores Z, LD correlation matrix Σ, and a
non-empty index set S, the marginal likelihood conditional on γ S isZm(Z|γ S , Σ) =f (Z|λS , Σ)f (λS |γ S , τ )dλS ,λS

where f (Z|λS , Σ) is the density function of MVN(ΣλS , Σ), and f (λS |γ S , τ ) is the prior density function of
the scaled effect size, which is the product of p − |S| point mass distribution at 0 and the density function of
MVN(0, τ1 IS ). The marginal likelihood after integrating out λγ S is
(
)
0 −1
0
−1 Z
p |S|
Z
Σ
Z
−
Z
(Σ
+
τ
I
)
1
1
γ
γ
S
S
S
S
m(Z|γ S , Σ) = |Σγ S + τ Iγ S |− 2 |Σ|− 2 (2π)− 2 τ 2 exp −,(1)
2
and the ratio between the marginal likelihood of γ S and the null model γ 0 is
 0

|S|
Z S (ΣγS + τ Iγ S )−1 Z Sm(Z|γ S )− 12 2= |Σγ S + τ Iγ S | τ exp.
m(Z|γ 0 )2

12The details can be found in the Supplemental Material. Let M denote the model set that contains all candidate
models. Then the posterior probability of any non-null model γ S and the posterior probability of γi being equal
to 1 (the PIP) can be computed as
P Oγ S :γ 0,
γ A ∈M P Oγ A :γ 0X
Pr (γi = 1|Z) =Pr (γ S |Z) ,Pr (γ S |Z) = P

γ S :i∈Sf (Z|γ S )f (Z|γ 0 )

Z 0S (ΣγS + τ Iγ S )−1 Z S
2where the posterior odds (P Oγ S :γ 0 ) is defined as the product of the Bayes factorPr(γ S |η)Pr(γ |η) :and the prior odds

0

P Oγ S :γ 0 =m(Z|γ S ) Pr (γ S |η)·m(Z|γ 0 ) Pr (γ 0 |η)1 |S|
η |S| (p − |S|)!
|Σγ S + τ Iγ S |− 2 τ 2 exp
=
p!



.

Shotgun stochastic search algorithm. Most fine mapping methods adopt MCMC-based algorithms to
explore the posterior model space. However such algorithms can be ineffective in high dimension situations
because of slow convergence and inefficient proposal distribution.16 Furthermore, due to high correlations
among SNPs at a locus, MCMC algorithms may not visit perfectly correlated SNPs evenly, which leads tounequal PIPs for such SNPs. Therefore, we adopt here the Shotgun stochastic search (Shotgun) algorithm16
for exploring the posterior distribution over model space, which has been used before in the context of finemapping.10 The main benefit of using the Shotgun algorithm is that it not only records the visited candidate
models, but also the neighborhood of the visited candidate models, hence the Shotgun algorithm has a more
complete exploration of the areas with high marginal likelihoods in the posterior model space. This semiexhaustive search feature alleviates the problem of unequal PIPs for perfectly correlated SNPs. This feature
will also play an important role when integrating functional annotations into the computation of the proposedmodel. We note that running the Shotgun algorithm in this setting is feasible due to the strong dimensional
penalization introduced by the prior Poisson distribution that constrains the search in the posterior model
space.
Note that for a specific model denoted by an index set S, the unnormalized posterior distribution of theindicator vector γ S is proportional to a product of the marginal likelihood and the prior distribution:
Pr (γ S |Z) =f (Z|γ S )f (γ S )∝ f (Z|γ S )f (γ S ).
f (Z)The Shotgun algorithm is an iterative procedure that exhaustively examines the following three neighborhood
sets of the current model:
Γ− (S) := {A : A ⊂ S, |S| − |A| = 1} (one less SNP than S),
Γ+ (S) := {A : A ⊃ S, |A| − |S| = 1} (one more SNP than S),
Γ⇔ (S) := {A : |S| − |A ∩ S| = 1, |A| = |S|} (models that replace one SNP in S).
Then, all the unnormalized posterior probabilities of the neighborhood models, i.e., {Γ− (S) ∪ Γ+ (S) ∪ Γ⇔ (S)},will be computed. To update the current model, the algorithm first randomly selects one candidate model fromeach neighborhood set according to the unnormalized posterior probabilities, then randomly selects the next
current model from the three selected models according to the corresponding posterior probabilities. By doing so, the algorithm stochastically moves towards the high posterior area in the model space. The Shotgun
algorithm is stopped when the sum of the absolute difference of the PIPs between iterations is smaller thana pre-determined threshold, and we keep the top B visited candidate models with the largest unnormalized
posterior probabilities.
13Remarks
The inversion of Σ. Due to high correlations among SNPs at a locus, the LD matrix is not full rank mostof the time. We use the Moore–Penrose inverse of the LD matrix instead.
Posterior probability of candidate models. Most of the fine-mapping methods that have been proposed
focus on the goal of variable selection, accounting for LD. One advantage of our proposed model is that inaddition to PIPs, we provide the posterior probability of candidate models which can be informative about theputative causal variants at the test locus. Note that FINEMAP also provides similar posterior probabilities ofcandidate models.
Specifically, CARMA returns a vast library of visited candidate models that can be ranked according to thecorresponding marginal likelihoods. Let γ (b) , b = 1, . . . , B, denote the ranked candidate models, such as γ (1)receives the largest marginal likelihood. We use γ (1) as the reference model to select all other candidate models
that are not significantly different from γ (1) through a Bayesian hypothesis testing procedure as follows. The
posterior odds of the hypothesis test for comparing γ (1) to all other candidate models is



 

Pr γ (1) |Z, η
Pr Z|γ (1) Pr γ (1) |η

 =

 
,
POγ (1) :γ (b) =Pr γ (b) |Z, η
Pr Z|γ (b) Pr γ (b) |ηfor b = 1, . . . , B. The posterior odds quantifies the strength of evidence in favor of the leading causal configuration, represented by γ (1) relative to other candidate models. If all candidate models have the same prior,then the posterior odds is equal to the Bayes factor. We use commonly accepted thresholds for BFγ (1) :γ (b) as
the threshold for P Oγ (1) :γ (b) to determine if a candidate model γ(b) is significantly different from the reference
model γ(1) .18 We only keep those models with the corresponding Bayes factor smaller than the threshold (e.g.
3.2 or 10), and refer to the corresponding SNPs identified by these selected models as the credible model.
Incorporating functional annotations. Leveraging functional annotation data on genetic variants at a
locus of interest may improve the results of fine-mapping studies.6, 12 Most fine-mapping methods allow inclusionof only one functional annotation in the form of the prior probability for each SNP. However, it is difficult tochoose one best annotation. fastPAINTOR cannot deal with high-dimensional functional annotation data and apre-selection of a small number of functional annotations needs to be performed before running fastPAINTOR.
SuSiE+PolyFun is essentially a two-step procedure: first, it estimates prior probabilities using the S-LDSC
regression framework, and then those priors are provided as input to SuSiE.
Here we propose to incorporate large number of functional annotations using a logistic distribution to model
the relationship between the annotations and the prior probability, such as logit (Pr (γi = 1|wi , θ)) = w0i θ, where
wi is the q + 1 dimensional vector (for q annotations) of the ith SNP. Since the indicator vector γ is unknown
to us, an EM algorithm is typically used to replace the unknown indicator by its expectation conditional on
summary statistics, functional annotations, and the estimated coefficients during the computation. Here we
propose an EM algorithm with regularization on the coefficients to integrate information on large number ofannotations. By adding the functional annotations, the log-likelihood can be written as:`(θ; Z, γ, W ) = log (m(Z|γ)) + log (Pr (γ|W , θ))= log (m(Z|γ)) +

p
X


γi w0i θ − (1 − γi ) log 1 + exp(w0i θ) ,i=1

where m(Z|γ) is given by Equation (1). The EM algorithm is as follows:

14

(2)

Input: Summary statistics Z, functional annotations W , hyperparameter η of the Poisson priordistribution (by default η = 1).
Initialization: Run Shotgun algorithm without annotations to generate the initial posterior probability
Pr (γ|Z, η).
for s = 1, . . . until convergence doE-step- Run Shotgun algorithm with initial prior probability Pr γi |wi , θ (s) = logit−1 w0i θ (s)estimated at step (s − 1).
- Compute the posterior probability Pr γi = 1|Z, W , θ (s) :m(Z|γ S )Pr γ S |W , θ (s),Pr γ S |Z, W , θ (s) = P
(s)
m(Z|γ
)Pr
γ
|W
,
θ
A
A
γ A ∈M
X
Pr γi = 1|Z, W , θ (s) =Pr γ S |Z, W , θ (s) .
γ S :i∈S- Compute the Q function, i.e. the expectation of the log-likelihood with respect to the posterior
probability of γ as follows:Q(θ|γ, θ (s) ) = Eγ|Z,W ,θ(s) [`(θ; Z, γ, W )]= Eγ|Z,W ,θ(s) [log (m(Z|γ))] +

p h
X



i
pi (θ (s) )w0i θ − 1 − pi (θ (s) ) log (1 + exp(w0i θ)) , (3)i=1where we denote bypi (θ (s) ) = Pr γi = 1|Z, W , θ (s) .
M-step We wish to estimate θ by maximizing the Q function with regularization on θ (except for theintercept) to prevent overfitting with large number of annotations. We use cyclical coordinate
descent algorithm34 for finding the maximizer. We maximize the penalized Q function as follows:θ (s+1) :=argmax

p h


i
X
pi (θ (s) )w0i θ − 1 − pi (θ (s) ) log (1 + exp {w0i θ}) − κPa (θ −0 ),θ∈Rq+1 i=1

where(1 − a)||θ −0 ||2`2 + a||θ −0 ||`1 ,2and θ −0 is the coefficient vector without intercept, i.e. the coefficient vector of the q annotations.
The tuning parameter κ is chosen by a standard cross-validation procedure based on the deviance.
Notice that in this step we can drop the first term in Equation (3) because there is no θ involved in
the maximization. We assume a = 0 in our implementation (i.e. ridge penalty).
Pa (θ −0 ) =Then, compute the estimated prior probability for the (s + 1) step:
Pr γi = 1|wi , θ (s+1) =exp(w0i θ (s+1) )1 + exp(w0i θ (s+1) )for i = 1, . . . , p.
end
Algorithm 1: EM algorithm with functional annotations.

15,Remarks
Feature selection on functional annotations. As in PolyFun, we use the ridge penalty (a = 0) to regularize the strength of the coefficients of annotations. If the dimension of annotations is very high, an alternative
penalty is to use the elastic net penalty (a = 0.5).
Identifying SNPs with Z-score/LD discrepancies. A challenge when doing fine-mapping with summary
statistics and LD from external reference panels is the possibility of discrepancies between LD and Z-scorevalues. This can happen for several reasons. For example, when using a reference panel to estimate the LD, orwhen summary statistics at SNPs in a region are based on studies with different sample sizes. Such scenarios
can result in unrealistically large values for the marginal likelihood and high PIP values, especially given thelarge sample sizes for meta-analysis studies. Here, we propose a Bayesian hypothesis testing procedure foroutlier detection embedded in the Shotgun algorithm.
Bayesian hypothesis testing. Suppose that MS is the currently selected model in the Shotgun algorithmwith a non-null index set S corresponding to γ S . Let Σγ S denote the corresponding LD matrix. When a−1discrepancy is identified in the index set S, the quadratic term Z 0S Σγ S + τ Iγ S
Z S will be substantiallyinflated, which biases the corresponding marginal likelihood m(Z|γ S , Σ) and leads to biased posterior analysis.
The proposed outlier detection procedure is based on the work in,35 where the authors proposed a ridge
regularization of the correlation matrix that leads to a stable likelihood-based test. Specifically, for a given LDb γ such as
matrix Σγ S , we construct another LD matrix ΣSb γ = ξΣγ + (1 − ξ)I, for ξ ∈ [0, 1],ΣSSb γ ). We
where I is the identity matrix of same dimension as Σγ S and for which we can compute m(Z S |γ S , ΣSconduct hypothesis testing as follows:
H0 : ξ = 1(there is no discrepancy),
H1 : ξ 6= 1(there is discrepancy).
Given an equal prior probability on the null and alternative hypothesis, we derive the Bayes factor of the test
as a function of ξ following the robust Bayesian approach:36
BF0:1 = inf B (ξ) =ξ∈[0,1]m(Z S |γ S , Σγ S ),bγ )sup m(Z S |γ S , Σξ∈[0,1]S

whereb γ ).
ξˆ :=argmax m(Z S |γ S , ΣSξ∈[0,1]

We reject the null hypothesis, i.e. the marginal likelihood is indeed inflated by the presence of discrepancies, if
BF0:1 < δ. Possible values for the threshold δ are suggested in,18 and we assess sensitivity to different choices.
We find the maximizer ξˆ using the R function optimize. (More discussion is in the Supplemental Material.)
Imputation before fine-mapping. One way to mitigate the effect of differential sample sizes across SNPs
in a meta-analysis study is to impute the Z-scores at SNPs with small sample sizes using data from SNPs withlarger sample sizes. Here, we use a simple imputation method similar to the one used in.1 Let A and B denotethe index set of two cohorts with smaller and larger sample sizes, respectively; then we can impute the Z-scoresof the SNPs in set A, denoted as Z A , conditional on the LD matrix Σ and Z B . In,1 the authors suggest to
approximate E [Zi ], for i ∈ A, by Σi,c Zc , where Zc is the leading Z-score at the locus and Σi,c is the Pearson
correlation between Zi and Zc . Similarly, for a given SNP i ∈ A, we propose to impute Zi by
bi = Σi,m Zm ,Z
where m =argmax {|Σi,j Zj |}. Note that much more advanced imputation methods can be used instead, andj∈B

our simple implementation is used here only as a proof-of-principle.
16

Implementation in the Shotgun algorithm.
Suppose that MS is the currently selected model in the
Shotgun algorithm with a non-null index set S (note that MS has already passed the outlier hypothesis testin the previous step of the Shotgun algorithm, therefore we assume that there are no discrepancies included in
MS for the hypothesis test at the current step). There are three sets of candidate models in the neighborhood
space of MS , i.e. Γ− (S), Γ+ (S), and Γ⇔ (S) as defined in the section on the Shotgun algorithm. The Shotgun
algorithm selects one candidate model from each set according to the posterior probabilities, which could be
inflated by the outliers. Therefore, for each selected model from the model sets Γ+ (S) and Γ⇔ (S), we conduct
the Bayesian hypothesis test described above to examine whether the selection is due to inflated marginal
likelihood caused by possible outliers included in the selected models. Note that the candidate models includedin Γ− (S) have already been tested when the algorithm reached the currently selected model γ S , therefore we
do not need to test the models in set Γ− (S). We drop the included/exchanged SNP from the algorithm if thecorresponding Bayes factor is smaller than a threshold δ. Algorithm 2 shows the exact procedure of the outlier
detection.
Note that, for a current model S, the Shotgun algorithm is a semi-exhaustive searching algorithm thatconsiders all configurations associated with the SNPs in S. Furthermore, the Shotgun algorithm stochasticallymoves towards the area of large posterior probabilities in the model space, therefore the outlier detection
procedure is more efficient at detecting outliers that could bias the fine-mapping analysis.
At any step of the Shotgun algorithm, suppose that the current model is γ S and S 6= ∅.
Input: The index set S for the current model γ S , the threshold on the Bayes factor δ.
Construct two neighborhood model sets Γ+ (S) and Γ⇔ (S):Γ+ (S) := {A : A ⊃ S, |A| − |S| = 1} (one more SNP than S),
Γ⇔ (S) := {A : |S| − |A ∩ S| = 1, |A| = |S|} (models that replaces one SNP in S).
for Γ+ (S) and Γ⇔ (S) do
repeat- Randomly select one candidate model A from the neighborhood model set according to theunnormalized posterior probabilities.
T- i = A \ (A S).
- Define the hypothesis test for the marginal likelihood of A, such that
H0 : ξ = 1;
H1 : ξ 6= 1;

ith SNP is not an outlier
ith SNP is an outlier.

- Compute the corresponding Bayes factor
B̂A = inf B (ξ) =ξ∈[0,1]m(Z|γ A , Σγ A ),bγ )sup m(Z|γ A , Σξ∈[0,1]Ab γ = ξΣγ + (1 − ξ)I
where ΣAA
if B̂A ≤ δ then- Drop ith SNP from consideration for the fine-mapping algorithm.
until B̂A > δ;end
Algorithm 2: The outlier detection procedure implemented in Shotgun algorithm.

17

Multi-locus extension
As shown before,6, 12, 37 a fine-mapping method may benefit from inference based on multiple loci. The multiple
locus model assumes that functional annotations are similarly associated with the probability to be causal. In
our applications we fit the proposed model at the chromosome level.

Credible sets
In5 the authors define a credible set as follows.
Definition 1. A level ρ credible set is defined to be a subset of correlated variables (with correlation within theset greater than some threshold r) that has probability ρ or greater of containing at least one effect variable (i.e.
causal SNP).
Constructing credible sets. Given a correlation threshold r, we can define multiple candidate credible sets
as follows
S := {i ∈ {1, . . . , p} : min {cor(i, j) ≥ r} , for all i, j ∈ S} .
Given an index set S, let s = s1 , . . . , s|S| denote the indices of the selected SNPs in the set S ranked in orderof decreasing PIPs, such as Pr (γs1 |Z) > Pr (γs2 |Z) > . . . > Pr γs|S| |Z , and let Pl denote the cumulative sumof the l largest PIPs:
Pl =

l
X


Pr γsj |Z .

j=1

Then, the credible set is defined as Sl0 = {s1 , . . . , sl0 }, where l0 = min {l : Pl ≥ ρ}. This procedure makes surethat the selected credible sets have the smallest size. Then, given a specific level ρ, we run this procedure on apredetermined sequence of candidate correlation values, such as r ∈ {0.60, 0.65, 0.70, 0.75, 0.80, 0.85, 0.90, 0.95},
and compute the corresponding number of credible sets for each candidate value. The largest correlationthreshold r with level ρ credible sets (those sets with posterior probability greater than ρ) is selected alongwith the corresponding credible sets.
Data Availability. LD matrices for UK Biobank British-ancestry samples are precomputed by the PolyFun
study and accessed from https://data.broadinstitute.org/alkesgroup/UKBB LD. Baseline-LF v2.2.UKB annotations for 19 million imputed UK Biobank SNPs are available fromhttps://data.broadinstitute.org/alkesgroup/LDSCORE/baselineLF v2.2.UKB.tar.gz. Summary statistics of AD
GWAS are available at https://ctg.cncr.nl/software/summary statistics. The 1000 Genomes Project phase 3
genotype data are available through ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/.
Code Availability. CARMA has been implemented in a computationally efficient R package available athttps://github.com/Iuliana-Ionita-Laza/CARMA.

References
1 Masahiro Kanai, Roy Elzur, Wei Zhou, Mark J Daly, Hilary K Finucane, Global Biobank Meta analysis

Initiative, et al. Meta-analysis fine-mapping is often miscalibrated at single-variant resolution. medRxiv,
2022.
2 Christian Benner, Aki Havulinna, Marjo-Riitta Järvelin, Veikko Salomaa, Samuli Ripatti, and Matti Pirinen.

Prospects of fine-mapping trait-associated genomic regions by using summary statistics from genome-wide
association studies. Am J Hum Genet., (101):539–551, 2017.
18

3 Iris E Jansen, Jeanne E Savage, Kyoko Watanabe, Julien Bryois, Dylan M Williams, Stacy Steinberg, Julia

Sealock, Ida K Karlsson, Sara Hägg, Lavinia Athanasiu, et al. Genome-wide meta-analysis identifies new lociand functional pathways influencing alzheimer’s disease risk. Nature genetics, 51(3):404–413, 2019.
4 Clare Bycroft, Colin Freeman, Desislava Petkova, Gavin Band, Lloyd T Elliott, Kevin Sharp, Allan Motyer,

Damjan Vukcevic, Olivier Delaneau, Jared O’Connell, et al. The uk biobank resource with deep phenotyping
and genomic data. Nature, 562(7726):203–209, 2018.
5 Gao Wang, Abhishek K Sarkar, Peter Carbonetto, and Matthew Stephens. A simple new approach to variableselection in regression, with application to genetic fine-mapping. Journal of the Royal Statistical Society:
Series B (Statistical Methodology), 82(5):1273–1300, 2020.
6 Gleb Kichaev, Megan Roytman, Ruth Johnson, Eleazar Eskin, Sara Lindstroem, Peter Kraft, and Bogdan

Pasaniuc. Improved methods for multi-trait fine mapping of pleiotropic risk loci. Bioinformatics, 33(2):248–
255, 2017.
7 Yongtao Guan and Matthew Stephens.

Bayesian variable selection regression for genome-wide association
studies and other large-scale problems. The Annals of Applied Statistics, pages 1780–1815, 2011.

8 Julian B Maller,

Gilean McVean, Jake Byrnes, Damjan Vukcevic, Kimmo Palin, Zhan Su, Joanna MM
Howson, Adam Auton, Simon Myers, Andrew Morris, et al. Bayesian refinement of association signals for 14
loci in 3 common diseases. Nature genetics, 44(12):1294, 2012.

9 Paul J Newcombe, David V Conti, and Sylvia Richardson.

Jam: a scalable bayesian framework for joint
analysis of marginal snp effects. Genetic epidemiology, 40(3):188–201, 2016.

10 Christian Benner, Chris CA Spencer, Aki S Havulinna, Veikko Salomaa, Samuli Ripatti, and Matti Piri-nen. Finemap: efficient variable selection using summary data from genome-wide association studies.
Bioinformatics, 32(10):1493–1501, 2016.
11 Gleb Kichaev, Wen-Yun Yang, Sara Lindstrom, Farhad Hormozdiari, Eleazar Eskin, Alkes L Price, Peter

Kraft, and Bogdan Pasaniuc. Integrating functional data to prioritize causal variants in statistical finemapping studies. PLoS Genet, 10(10):e1004722, 2014.
12 Omer Weissbrod, Farhad Hormozdiari, Christian Benner, Ran Cui, Jacob Ulirsch, Steven Gazal, Armin P

Schoech, Bryce Van De Geijn, Yakir Reshef, Carla Márquez-Luna, et al. Functionally informed fine-mapping
and polygenic localization of complex trait heritability. Nature Genetics, pages 1–9, 2020.
13 Wenan Chen, Beth R Larrabee, Inna G Ovsyannikova, Richard B Kennedy, Iana H Haralambieva, Gregory A

Poland, and Daniel J Schaid. Fine mapping causal variants with an approximate bayesian method using
marginal test statistics. Genetics, 200(3):719–736, 2015.
14 James G Scott and James O Berger.

Bayes and empirical-bayes multiplicity adjustment in the variableselection problem. The Annals of Statistics, pages 2587–2619, 2010.

15 Jian Yang, Teresa Ferreira, Andrew P Morris, Sarah E Medland, Pamela AF Madden, Andrew C Heath,

Nicholas G Martin, Grant W Montgomery, Michael N Weedon, Ruth J Loos, et al. Conditional and joint
multiple-snp analysis of gwas summary statistics identifies additional variants influencing complex traits.
Nature genetics, 44(4):369–375, 2012.
16 Chris Hans, Adrian Dobra, and Mike West.

Shotgun stochastic search for “large p” regression. Journal of
the American Statistical Association, 102(478):507–516, 2007.

17 Hui Zou and Trevor Hastie.

Regularization and variable selection via the elastic net. Journal of the royalstatistical society: series B (statistical methodology), 67(2):301–320, 2005.

18 Robert

E Kass and Adrian E Raftery. Bayes factors. Journal of the american statistical association,
90(430):773–795, 1995.
19

19 Apostolos Dimitromanolakis, Jingxiong Xu, Agnieszka Krol, and Laurent Briollais. sim1000g: a user-friendly

genetic variant simulator in r for unrelated individuals and family-based designs. BMC bioinformatics,
20(1):26, 2019.
20 Laura Fachal, Hugues Aschard, Jonathan Beesley, Daniel R Barnes, Jamie Allen, Siddhartha Kar, Karen A

Pooley, Joe Dennis, Kyriaki Michailidou, Constance Turman, et al. Fine-mapping of 150 breast cancer risk
regions identifies 191 likely target genes. Nature genetics, 52(1):56–73, 2020.
21 Jian Zhou and Olga G Troyanskaya.

Predicting effects of noncoding variants with deep learning–based
sequence model. Nature methods, 12(10):931–934, 2015.

22 Yuxin Zou, Peter Carbonetto, Gao Wang, and Matthew Stephens.

Fine-mapping from summary data withthe “sum of single effects” model. bioRxiv, 2021.
23 Cristen J Willer, Yun Li, and Gonalo R Abecasis.

Metal: fast and efficient meta-analysis of genomewide
association scans. Bioinformatics, 26(17):2190–2191, 2010.

24 Martin Kircher, Daniela M Witten, Preti Jain, Brian J O’roak, Gregory M Cooper, and Jay Shendure.

A
general framework for estimating the relative pathogenicity of human genetic variants. Nature genetics,
46(3):310–315, 2014.

25 Iuliana Ionita-Laza, Kenneth McCallum, Bin Xu, and Joseph D Buxbaum. A spectral approach integrating

functional genomic annotations for coding and noncoding variants. Nature genetics, 48(2):214, 2016.
26 GTEx Consortium et al. The gtex consortium atlas of genetic regulatory effects across human tissues. Science,

369(6509):1318–1330, 2020.
27 Bernard Ng, Charles C White, Hans-Ulrich Klein, Solveig K Sieberts, Cristin McCabe, Ellis Patrick, Jishu

Xu, Lei Yu, Chris Gaiteri, David A Bennett, et al. An xqtl map integrates the genetic architecture of thehuman brain’s transcriptome and epigenome. Nature neuroscience, 20(10):1418–1426, 2017.
28 Rebecca Sims, Sven J Van Der Lee, Adam C Naj, Céline Bellenguez, Nandini Badarinarayan, Johanna

Jakobsdottir, Brian W Kunkle, Anne Boland, Rachel Raybould, Joshua C Bis, et al. Rare coding variants inplcg2, abi3, and trem2 implicate microglial-mediated innate immunity in alzheimer’s disease. Nature genetics,
49(9):1373–1384, 2017.
29 Arnold Zellner. On assessing prior distributions and bayesian regression analysis with g-prior distributions.

Bayesian inference and decision techniques, 1986.
30 Hilary K Finucane, Brendan Bulik-Sullivan, Alexander Gusev, Gosia Trynka, Yakir Reshef, Po-Ru Loh,

Verneri Anttila, Han Xu, Chongzhi Zang, Kyle Farh, et al. Partitioning heritability by functional annotation
using genome-wide association summary statistics. Nature genetics, 47(11):1228, 2015.
31 Arnold Zellner and Aloysius Siow.

Posterior odds ratios for selected regression hypotheses. Trabajos deestadı́stica y de investigación operativa, 31(1):585–603, 1980.

32 Ismaël Castillo, Aad van der Vaart, et al.

Needles and straw in a haystack: Posterior concentration forpossibly sparse sequences. The Annals of Statistics, 40(4):2069–2101, 2012.

33 Andrew J Womack, Claudio Fuentes, and Daniel Taylor-Rodriguez. Model space priors for objective sparsebayesian regression. arXiv preprint arXiv:1511.04745, 2015.
34 Jerome Friedman, Trevor Hastie, and Rob Tibshirani. Regularization paths for generalized linear models via

coordinate descent. Journal of statistical software, 33(1):1, 2010.
35 David I Warton. Penalized normal likelihood and ridge regularization of correlation and covariance matrices.

Journal of the American Statistical Association, 103(481):340–349, 2008.

20

36 MJ Bayarri and James O Berger. Robust Bayesian bounds for outlier detection. De Gruyter, 1992.
37 Wenan Chen, Shannon K McDonnell, Stephen N Thibodeau, Lori S Tillmans, and Daniel J Schaid.

Incorporating functional annotations for fine-mapping causal variants in a bayesian framework using summary
statistics. Genetics, 204(3):933–958, 2016.

21

SCIMP(a)

rs113260531

CARMA

6

1.0

rs4456560

rs4456560

rs75511804

0.6

PIP
0

rs77972827
rs7225151

0.4

2

rs77972827

rs7225151

rs61481506
rs113806568

0.0

−4

rs61481506

rs113260531

0.2

−2

Z−scores40.8

rs75511804

4.94Mb

5.01Mb

5.11Mb

5.18Mb

5.26Mb

5.32Mb

4.94Mb

5.01Mb

Position

5.11Mb

5.18Mb

5.26Mb

5.32Mb

0

0.1

0.2

0.3

0.4

(b)

0.5

0.6

0.7

0.8

0.9

1

r2

Position
rs3887943

SuSiE

6

1.0

rs113215752

rs3887943
rs73976325
rs113215752

rs78534596
rs113061798
rs58124010 rs1046368

rs113215752

rs78534596

5.11Mb

5.18Mb

5.26Mb

0.6
0.4

PIP
5.01Mbrs73976325

0.2

0
−4
4.94Mbrs113061798
rs73976325

rs58124010

0.02

rs113061798−2

Z−scores40.8

rs78534596

rs1046368

5.32Mbrs58124010

4.94Mb

5.01Mb

Position

5.11Mb

5.18Mb

5.26Mb

5.32Mb

0

0.1

0.2

0.3

0.4

(c)
6

1.0

4

0.8
0.6

0.6

0.7

0.8

0.9

1

rs113061798

fastPAINTOR

2

0.5

r2

Position

rs78534596

rs78534596

rs401742

PIP

0.4

0

rs62072535

rs401742

rs62072535

rs78534596

rs62072535

rs12602829

rs7502862

rs12602829

rs12602829

rs7502862

0.0

−4

0.2

−2

Z−scoresrs401742

rs113061798

4.94Mb

5.01Mb

5.11Mb

5.18Mb

5.26Mb

5.32Mb

4.94Mb

5.01Mb

Position

5.11Mb

5.18Mb

5.26Mb

5.32Mb

0

0.1

0.2

0.3

0.4

(d)

SuSiE (outliers removed by CARMA)

0.6

0.7

0.8

0.9

1

rs77972827
rs7225151

rs7225151

1.0

6

0.5

r2

Position

rs75511804

40.8

rs75511804

0.6

rs58124010

0

0.4
0.2
0.0

rs58124010

−2

PIP

2

rs113260531−4Z−scoresrs113260531

rs4456560

rs4456560

rs113260531

rs113806568

rs113806568

4.94Mb

5.01Mb

5.12Mb

Position

5.19Mb

5.27Mb

5.33Mbrs61481506

4.94Mb

5.01Mb

5.12Mb

Position

5.19Mb

5.27Mb

5.33Mb

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

r2

Figure 1: Motivating example. (a-c) GWAS Z-scores along with PIPs from three models (CARMA, SuSiEand fastPAINTOR) for one Alzheimer’s disease risk locus SCIMP (GRCh37/hg19). For each model, theheatmap depicts the LD (r2 based on the UKBB data) between SNPs highlighted in the PIP panel in the
middle. The highlighted top SNPs of CARMA belong to one credible set; the highlighted top SNPs of SuSiE
belong to seven different credible sets; the highlighted top SNPs of fastPAINTOR belong to three different
credible sets. (d) Results of SuSiE applied to a dataset where outliers identified by CARMA are removed. Notethat for better visualization, only a segment of the 1 Mb window analyzed at this locus is shown.

22In−sample LD(a) Credible set
Power

Coverage

1.00

1

1.00

0.9

1

0.99

0.97

0.91
0.81

0.75

0.75

0.74

0.73

0.72

0.75

Value

Value

0.61

0.50

0.25

0.50

0.25

0.00

0.00
CARMA

SuSiE (SuSiE+PolyFun)

fastPAINTOR

CARMA

Credible size

SuSiE (SuSiE+PolyFun)

0.97

15

0.95

0.94

0.91

14

0.9

13

10

8.8
8

0.9

0.86

0.6

Value

Value

fastPAINTOR

Purity

6.2

5

0.3

3.7

0

0.0
CARMA

SuSiE (SuSiE+PolyFun)

fastPAINTOR

Scenario

CARMA

No Annotations

SuSiE (SuSiE+PolyFun)

fastPAINTOR

Including Annotations(b) Credible model
Power

Size
25

1.00
0.93

22

0.9

0.86

21

20

0.81
0.75

0.75

0.73

18

0.74
14

15

Value

Value

0.61

0.50

12
11

10
6.7

0.25

5

0.00

3.8

0
CARMA
(credible model)

CARMA(credible set)

SuSiE (SuSiE+PolyFun)(credible set)

fastPAINTOR(credible set)

Scenario

CARMA(credible model)

No Annotations

CARMA(credible set)

SuSiE (SuSiE+PolyFun)(credible set)

fastPAINTOR(credible set)

Including Annotations(c) FDR vs. Power
Including Annotations
1.00

0.75

0.75

FDR

FDR

No Annotations
1.00

0.50

0.25

0.50

0.25

0.00

0.00
0.00

0.25

0.50

0.75

1.00

0.00

0.25

0.50

Power

Model

0.75

1.00

Power

CARMA

SuSiE (SuSiE+PolyFun)

fastPAINTOR

PIP

PIP>0.5

PIP>0.95

Figure 2: Credible sets, credible models, and FDR vs. Power (in-sample LD, with and without
functional annotations). (a) Performance of credible sets (ρ = 0.99). Power: the proportion of simulatedcausal variants included in any credible set; Coverage: the proportion of credible sets containing a causal variant;Size: the number of variants included in a credible set; Purity: the mean squared correlation of variants ina credible set. (b) Performance of credible models. Power is the percentage of the simulated causal variants
identified by the credible model or the credible sets. Size is the average number of variants included in thecredible model (CARMA) or all the credible sets (CARMA, SuSiE, fastPAINTOR) at a locus. The credible
model is computed based on a threshold of 10 for the posterior odds. (c) FDR vs. power using positiveFP
predictions as the PIP threshold varies from 0.1 to 1. These quantities are calculated as FDR := TP+FP
andTP
power := TP+FN, where FP, TP, FN, TN denote the number of false positives, true positives, false negatives
and true negatives respectively given a certain PIP threshold. Open circles denote the results at PIP threshold
0.5, and solid circles denote the results at PIP threshold 0.95. Two scenarios are shown: no annotation andincluding functional annotations. The true number of causal variants is |T | = 2.

23

(a) Credible set

UKBB LDPower

Coverage
0.88

0.86

0.82

0.81

0.75

0.75
0.61

0.59

Value

Value

0.53

0.50
0.44

0.41

0.52

0.50
0.34

0.31

0.29

0.28

0.25

0.23

0.25

0.22
0.13

0.1

0.00

0.00
CARMA

SuSiE

fastPAINTOR

CARMA

Credible size

SuSiE

fastPAINTOR

Purity

20
1

0.99
0.95

0.95

0.99

0.93

0.91

0.9

15

15

0.9

0.86

10

Value

Value

14

9.6

0.6

7.6
6.3

0.3

5.1

5

4.4
3.7
1.3

0.0

0

CARMA
CARMA

SuSiE

SuSiE
Number of causal variants

fastPAINTOR

Scenario

Consistent

Inconsistent

fastPAINTOR

Inconsistent + Imputation(b) Credible model
Power

Size

1.00
31

0.86
0.82

30

0.82

0.81

27

26

0.75
23

0.64

20

0.53

0.50

0.44

0.41

0.29

0.28

0.23

0.25

20

18

Value

Value

0.59

18

18

18
16
14

10

0.00

8.5

0
CARMA
(credible model)

CARMA(credible set)

SuSiE(credible set)

fastPAINTOR(credible set)

Scenario

Consistent

CARMA(credible model)

Inconsistent

CARMA(credible set)

SuSiE(credible set)

fastPAINTOR(credible set)

Inconsistent + Imputation(c) FDR vs. Power
Inconsistent

Inconsistent + Imputation
1.00

0.75

0.75

0.75

0.50

0.25

FDR

1.00

FDR

FDR

Consistent
1.00

0.50

0.25

0.00

0.25

0.00
0.0

0.2

0.4

0.6

0.50

0.00
0.0

0.2

Power

0.4

0.6

0.0

Power

Model

CARMA

SuSiE

fastPAINTOR

0.2

0.4

0.6

Power

PIP

PIP>0.5

PIP>0.95

Figure 3: Credible sets, credible models, and FDR vs. Power (UKBB LD) (a) Performance of credible
sets (ρ = 0.99). Power: the proportion of simulated causal variants included in any credible set; Coverage:
the proportion of credible sets containing a causal variant; Size: the number of variants included in a credible
set; Purity: the mean squared correlation of variants in a credible set. (b) Performance of credible models.
Power is the percentage of the simulated causal variants identified by the credible model or the credible sets.
Size is the average number of variants included in the credible model (CARMA) or all the credible sets at a
locus. The credible model is computed based on a threshold of 10 for the posterior odds. (c) FDR vs. powerusing positive predictions as the PIP threshold varies from 0.1 to 1. These quantities are calculated as FDR :=FP
TP
TP+FP and power := TP+FN , where FP, TP, FN, TN denote the number of false positives, true positives, false
negatives and true negatives respectively given a certain PIP threshold. Open circles denote the results at PIP
threshold 0.5, and solid circles denote the results at PIP threshold 0.95. Two scenarios are shown: consistentand inconsistent meta-analyses. The true number of causal variants is |T | = 2.

24

(a) Credible sets with only one SNPUKBB LD

Total Number

Coverage

1.00
540
0.86
0.77

440

0.75

400

Value

Value

320

0.50

0.45
0.39

200
0.25

33
13

11

7

23

0.22

0.07

0.041

8

0

0.025

0

0.00
CARMA

SuSiE

fastPAINTOR

Scenario

Consistent

CARMA

Inconsistent

SuSiE

fastPAINTOR

Inconsistent + Imputation(b) Credible set: SLALOM vs. CARMA
Power

Coverage

0.52

0.6

0.53
0.4

Value

Value

0.4

0.28
0.24

0.2

0.2

0.0

0.1

0.094

SuSiE

SuSiE + SLALOM

SuSiE + CARMA

1

1

0.98

SuSiE

SuSiE + SLALOM

SuSiE + CARMA

0.0
SuSiE

SuSiE + SLALOM

SuSiE + CARMA

Credible size

Purity

0.9
10

Value

Value

9.3
0.6

5
0.3

1.3

1.6

SuSiE

SuSiE + SLALOM

0.00
SuSiE + CARMA

Figure 4: (a) Increased single SNP credible sets in the presence of Z-scores/LD inconsistencies.
Results on credible sets (ρ = 0.99) with only one SNP, i.e., the corresponding SNP receives a PIP larger than
0.99. For each model we report the total number of credible sets that contain only one SNP across all 94 lociand three different scenarios, and the coverage of these sets, i.e. the proportion of these credible sets that
contain a causal SNP. (b) SuSiE with CARMA vs. SLALOM. Performance of credible sets (ρ = 0.99).
Power: the proportion of simulated causal variants included in any credible set; Coverage: the proportion of
credible sets containing a causal variant; Size: the number of variants included in a credible set; Purity: themean squared correlation of variants in a credible set. The true number of causal variants is |T | = 2.

25

Figure 5: Functional enrichments of prioritized SNPs at AD GWAS loci. Functional enrichment of SNPs in credible models prioritized by CARMA, CARMA +PolyFun prior causal probability, and
CARMA+Polyfun annotations in AD meta-analysis. The proportion of common SNPs lying in each functional
category is reported above its name.

26

27

0
1

# PIPs > 0.99
# credible sets

1
1

8
9

BIN1
7
7

HS3ST1

AD loci
CD2AP ZCWPW1 MS4A6A
SuSiE on original AD datasets
5
7
5
6
10
5
5
5

PICALM

01

5
6

ADAM10SuSiE on filtered AD datasets based on CARMA outlier detection
0
0
0
0
0
0
1
1
1
1
1
1

5
8

CLNK

0
0

5
7

SCIMP

0
0

4
8

ABI3

1
2

9
10

ABCA7

Table 1: SuSiE on the original AD data and filtered AD data based on CARMA outliers. For each locus, the number of SNPs with PIP
larger than 0.99 and the number of credible sets before and after CARMA outlier filtering are shown.

3
4# PIPs > 0.99
# credible sets

CR1

28

SNP
rs76726049
rs6014724
rs74685827
rs117618017
rs4663105
rs11865499
rs3752241
rs10933431
rs11623883
rs616338
rs4575098
rs143332484
rs180896475
rs77145198
rs2526380
rs184384746
rs3935067
rs28834970
rs10792832
rs75511804
rs13134197
rs11257238
rs3865444
rs2405442
rs9381563
rs116155464
rs6701713
rs6448453

NoAnnot.
0.961
0.292
0.009
0.485
1.000
0.099
0.970
0.327
0.260
0.156
0.716
0.965
0.089
0.008
0.902
0.826
0.421
0.347
0.073
0.097
0.241
0.000
0.430
0.282
0.416
0.217
0.173
0.000

PolyFunProb.
0.961
0.268
0.010
0.442
1.000
0.113
0.970
0.339
0.233
0.161
0.716
0.965
0.062
0.008
0.902
0.826
0.051
0.305
0.074
0.101
0.274
0.000
0.432
0.279
0.405
0.220
0.157
0.140

PolyFunAnnot.
1.000
1.000
1.000
1.000
1.000
1.000
1.000
1.000
0.998
0.995
0.994
0.991
0.990
0.987
0.933
0.905
0.903
0.867
0.764
0.697
0.632
0.534
0.498
0.472
0.461
0.310
0.298
0.129

Z
5.521
-6.2083.919
5.518
14.005-5.499-6.279-6.153-5.160
4.940
6.370
5.540
4.583
-5.059-5.568
5.695
6.562
6.830
-8.666
6.049
5.439
5.724
-5.842-7.906
6.362
5.954
8.837
6.000

N
43322743385843649844400642513443539742885543649843484536485943432638233636390044360871639
36448042696043477043623743574771639
430330436498440597429204364834443992436498

EAF
0.011
0.095
0.018
0.125
0.411
0.299
0.169
0.244
0.466
0.006
0.228
0.007
0.003
0.065
0.445
0.001
0.350
0.366
0.359
0.120
0.291
0.361
0.299
0.323
0.356
0.002
0.192
0.262

CADD
3.554
9.371
6.325
9.705
0.348
4.182
9.364
2.817
0.68
15.9
7.876
11.42
9.641
9.019
0.264
5.481
7.186
5.444
4.29
2.753
4.291
4.933
0.345
8.402
0.27
4.425
3.564
3.22

EigenPC
2.938
8.716
17.254
0.602
0.904
4.315
14.957
4.214
0.287
0.424
14.425
0.744
4.134
8.361
8.529
1.940
6.709
6.171
8.520
22.211
0.823
3.247
7.971
16.60512.254
2.244
8.013
3.222

Eigen
4.090
9.369
4.850
0.663
0.026
5.328
16.883
3.631
0.369
0.240
12.5550.630
4.850
15.064
0.269
2.055
6.781
4.850
4.850
12.297
1.091
4.850
1.550
14.470
2.614
2.765
0.157
5.329

V2G
MALT1
CASS4
SORL1
APH1B
BIN1
KAT8
ABCA7
INPP5D
SLC24A4ABI3NDUFS2
TREM2
CBLIF
SLTM
MPO
IL17RD
EPHA1
PTK2B
PICALM
RABEP1HS3ST1USP6NL
CD33
PILRB
CD2AP
CNTNAP2
CR1
CLNK
X
X

X

X

X
X
X

X
X
X

X

X

X

X
X
X
X

X
X

X

X
X
X
X
X
X
X
X
X
X
X
X

eQTL

meQTL

Table 2: Top SNP (with largest PIP in CARMA with PolyFun Annotations) at each GWAS locus. PIPs from CARMA without andwith functional annotations (PolyFun prior probability or PolyFun annotations). Z: meta-analysis Z-score; N: number of samples; EAF: effect allele
frequency; CADD/EigenPC/Eigen: PHRED-scaled scores for CADD/EigenPC/Eigen; V2G: top scoring gene from Open Targets; meQTL: indicates
if the top gene from V2G is pointed to by meQTL data from ROSMAP Brain Dorsolateral Prefrontal Cortex; eQTL: indicates if the top gene from
V2G is an eGene GTEx in a brain or blood tissue in GTEx.

Locus
ALPK2
CASS4
SORL1
APH1B
BIN1
KAT8
ABCA7
INPPD5
SLC24A4ABI3ADAMTS4
TREM2
MS4A6A
ADAM10BZRAP1-AS1
HESX1
EPHA1
CLU-PTK2B
PICALM
SCIMP
HS3ST1ECHDC3CD33
ZCWPW1CD2AP
CNTNAP2
CR1
CLNKMethod
FINEMAPCAVIARBF
JAM
fastPAINTOR
SuSiE(+PolyFun)
CARMA (proposed model)

Prior on effect size

Prior on model space

Mixture of Normal-GammaNormal-Gamma
g-priorNormal-Gamma

discrete Uniform
or discrete probability
discrete probability
Beta-Binomial
discrete probability
Multinomial(discrete probability)
PoissonNormal-Gamma
Spike-and-slab/Cauchy

Functional annot

Computation

Ref.

Yes

Shotgun

10

Yes
No
Yes

Exhaustive
Exhaustive and MCMC
Exhaustive and MCMC

37

Yes

Variational Bayes

5, 12

Yes

ShotgunTable 3: Summary of commonly used Bayesian fine-mapping methods.

29

96