nievergeltlab
7/13/2017 - 9:04 PM

Correlate and plot test statistics (pvalues) on a -log10 scale. Useful for comparing two similar analyses.

Correlate and plot test statistics (pvalues) on a -log10 scale. Useful for comparing two similar analyses.

correlate_stats.r
content_copyfile_download
/home/laramie/CLEANED_DATA/EXTERNAL_DATA/MIRECC6_Duke/duke_gwas_WH/12_15_15/DUKE_WH_model3_w_A2_danerFORMAT

grep -v NA MIRE_eur_analysis1_mf > MIRE_eur_analysis1_mf_nona

adam <- read.table('MIRE_eur_analysis1_mf_nona', nr=10000000,stringsAsFactors=F,header=T)
lara <- read.table('/home/laramie/CLEANED_DATA/EXTERNAL_DATA/MIRECC6_Duke/duke_gwas_WH/12_15_15/DUKE_WH_model3_w_A2_danerFORMAT', nr=8300000,stringsAsFactors=F,header=T)

lara$p <- 2*pnorm(abs(log(lara$OR)/lara$SE),lower.tail=F)
rescor <- merge(adam,lara,by="SNP")


cor.test(-log10(rescor$P),-log10(rescor$p))

png("correlations_me_larame.png")
plot(-log10(rescor$P),-log10(rescor$p))
dev.off()

cor.test(-log10(rescor$OR.x),-log10(rescor$OR.y))

rescor2 <- rescor[order(rescor$P),]




    unadj_filtered <- sort(adam$P)

    UNADJ <- -log(unadj_filtered,10)
    QQ <- -log(ppoints(length(UNADJ)),10)

    png('mirecc_output_qqplot.png')
    par(bty='l')

    chisqs <- qchisq(unadj_filtered,1,lower.tail=F)
    median_chisq <- median(chisqs,na.rm=T)
    GCfactor=median_chisq/.456

    
    plot(QQ, UNADJ, xlab='', ylab='', col='blue', cex=1.3, cex.axis=1.2, cex.lab=1.5,ylim=c(0,8.5),pch=20)
    abline(0,1,col='red', lwd=2)
    title(xlab='Expected -log10(p)', ylab='Observed -log10(p)', cex.lab=1.5)
    legend('topleft', paste('GC Lambda =', GCfactor),  bty='n', cex=1.5, xjust=1)
    dev.off()