Diagnostic for fine-mapping with summary
statistics
Yuxin Zou
2023-02-17
This vignette demonstrates the use of the diagnostic plot for
assessing consistency of the summary statistics and the reference LD
matrix.
The susie_rss assumes the LD matrix accurately estimate
the correlations among SNPs from the original GWAS genotype data.
Typically, the LD matrix comes from some public database of genotypes in
a suitable reference population. The inaccurate LD information leads to
unreliable fine-mapping result.
The diagnostic for consistency between summary statistics and
refenrence LD matrix is based on the RSS model under the null with
regularized LD matrix. \[
\hat{z} | R, \lambda \sim N(0, (1-\lambda)R + \lambda I),
0<\lambda<1
\] The parameter \(\lambda\) is
estimated by maximum likelihood. A larger \(\lambda\) means a greater inconsistency
between summary statistics and the LD matrix. The expected z score is
computed for each SNP, \(E(\hat{z}_j |
\hat{z}_{-j})\), and plotted against the observed z scores.
library(susieR)
library(curl)
LD information from the reference panel
We use another simulated data where the LD matrix is estimated from a
reference panel. In this example data set, there is one association
signal in the simulated data (red point), and there is one SNP with
mismatched reference and alternative allele between summary statistics
and the reference panel (yellow point).
Note: In some versions of PLINK, these mismatches
can occur when PLINK
automatically flips the alleles to make the minor allele be the effect
allele, leading to different allele encodings in the z scores and LD
matrix. Adding the flag --keep-allele-order will disable
this behaviour in PLINK.
data_file <- tempfile(fileext = ".RData")
data_url <- paste0("https://raw.githubusercontent.com/stephenslab/susieR/",
"master/inst/datafiles/SummaryConsistency1k.RData")
curl_download(data_url,data_file)
load(data_file)
zflip = SummaryConsistency$z
ld = SummaryConsistency$ldref
n=10000
b = numeric(length(zflip))
b[SummaryConsistency$signal_id] = zflip[SummaryConsistency$signal_id]
plot(zflip, pch = 16, col = "#767676", main = "Marginal Associations",
xlab="SNP", ylab = "z scores")
points(SummaryConsistency$signal_id, zflip[SummaryConsistency$signal_id], col=2, pch=16)
points(SummaryConsistency$flip_id, zflip[SummaryConsistency$flip_id], col=7, pch=16)
Using the data with misaligned allele, SuSiE-RSS identifies a true
positive CS containing the true effect SNP; and a false positive CS that
incorrectly contains the mismatched SNP.
fit = susie_rss(zflip, ld, n=n)
susie_plot(fit, y='PIP', b=b)
points(SummaryConsistency$flip_id, fit$pip[SummaryConsistency$flip_id], col=7, pch=16)
The estimated \(\lambda\) is
lambda = estimate_s_rss(zflip, ld, n=n)
lambda
# [1] 0.005005011
In the diagnostic plot, the mismatched SNP shows the largest
difference between observed and expected z-scores, and therefore appears
furthest away from the diagonal.
condz = kriging_rss(zflip, ld, n=n)
condz$plot
After fixing the allele encoding, SuSiE-RSS identifies a single true
positive CS containing the true-effect SNP, and the formerly mismatched
SNP is (correctly) not included in a CS.
z = zflip
z[SummaryConsistency$flip_id] = -z[SummaryConsistency$flip_id]
fit = susie_rss(z, ld, n=n)
susie_plot(fit, y='PIP', b=b)
