mhtplot.trunc | R Documentation |
Truncated Manhattan plot
mhtplot.trunc(
x,
chr = "CHR",
bp = "BP",
p = NULL,
log10p = NULL,
z = NULL,
snp = "SNP",
col = c("gray10", "gray60"),
chrlabs = NULL,
suggestiveline = -log10(1e-05),
genomewideline = -log10(5e-08),
highlight = NULL,
annotatelog10P = NULL,
annotateTop = FALSE,
cex.mtext = 1.5,
cex.text = 0.7,
mtext.line = 2,
y.ax.space = 5,
y.brk1 = NULL,
y.brk2 = NULL,
trunc.yaxis = TRUE,
cex.axis = 1.2,
delta = 0.05,
...
)
x |
A data.frame. |
chr |
Chromosome. |
bp |
Position. |
p |
p values, e.g., "1.23e-600". |
log10p |
log10(p). |
z |
z statistic, i.e., BETA/SE. |
snp |
SNP. Pending on the setup it could either of variant or gene ID(s). |
col |
Colours. |
chrlabs |
Chromosome labels, 1,2,...22,23,24,25. |
suggestiveline |
Suggestive line. |
genomewideline |
Genomewide line. |
highlight |
A list of SNPs to be highlighted. |
annotatelog10P |
Threshold of -log10(P) to annotate. |
annotateTop |
Annotate top. |
cex.mtext |
axis label extension factor. |
cex.text |
SNP label extension factor. |
mtext.line |
position of the y lab. |
y.ax.space |
interval of ticks of the y axis. |
y.brk1 |
lower -log10(P) break point. |
y.brk2 |
upper -log10(P) break point. |
trunc.yaxis |
do not truncate y-axisx when FALSE. |
cex.axis |
extension factor for x-, y-axis. |
delta |
a value to enable column(s) of red points. |
... |
other options. |
To generate truncated Manhattan plot, e.g., of genomewide significance (P values) or a random variable that is uniformly distributed.
The rationale of this function is to extend mhtplot() to handle extremely small p values as often seen from a protein GWAS. Optionally, the function also draws an ordinary Manhattan plot.
The plot is shown on or saved to the appropriate device.
James Peters, Jing Hua Zhao
mhtplot
.
## Not run:
options(width=120)
require(gap.datasets)
mhtdata <- within(mhtdata, {z=qnorm(p/2, lower.tail=FALSE)})
mhtplot.trunc(mhtdata, chr = "chr", bp = "pos", z = "z", snp = "rsn",
y.brk1=6, y.brk2=10, y.ax.space=1, mtext.line=2.5)
# https://portals.broadinstitute.org/collaboration/
# giant/images/c/c8/Meta-analysis_Locke_et_al%2BUKBiobank_2018_UPDATED.txt.gz
gz <- gzfile("work/Meta-analysis_Locke_et_al+UKBiobank_2018_UPDATED.txt.gz")
BMI <- within(read.delim(gz,as.is=TRUE), {Z <- BETA/SE})
print(subset(BMI[c("CHR","POS","SNP","P")],CHR!=16 & P<=1e-150))
library(Rmpfr)
print(within(subset(BMI, P==0, select=c(CHR,POS,SNP,Z)),
{P <- format(2*pnorm(mpfr(abs(Z),100),lower.tail=FALSE));
Pvalue <- pvalue(Z); log10P <- -log10p(Z)}))
png("BMI.png", res=300, units="in", width=9, height=6)
par(oma=c(0,0,0,0), mar=c(5,6.5,1,1))
mhtplot.trunc(BMI, chr="CHR", bp="POS", z="Z", snp="SNP",
suggestiveline=FALSE, genomewideline=-log10(1e-8),
cex.mtext=1.2, cex.text=1.2,
annotatelog10P=156, annotateTop = FALSE,
highlight=c("rs13021737","rs17817449","rs6567160"),
mtext.line=3, y.brk1=200, y.brk2=280, trunc.yaxis=TRUE,
cex.axis=1.2, cex=0.5,
y.ax.space=20,
col = c("blue4", "skyblue")
)
dev.off()
## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.