volcano plot

### Not to run ###

# read in EWIS results
res <- read.table('results/ewis_adultast_bmi_nonatopic_wo_neu.txt',
                  header = TRUE, sep = '\t')
res <- res[order(res$p.value), ]
res <- res[!is.na(res$Estimate), ]

# remove the unreliable probes reported by Chen et al 2013
res <- res[!res$ID %in% chen, ]

# highlight the enriched pathways
idx.ppar <- which(res$ID %in% set[['Glucocorticoid/PPAR signaling']])
idx.mapk <- which(res$ID %in% set[['MAPK signaling']])
idx.nfkb <- which(res$ID %in% set[['NF-kB signaling']])
idx.pi3k <- which(res$ID %in% set[['PI3K/AKT Signaling']])
idx.infl <- which(res$ID %in% set[['Global Inflammation']])

# set up colors for the pathways ever enriched
col.pathway <- c('#8fd87a', '#ac004e', '#009c88', '#ffa6ba', '#474d8b', '#e2c467')
names(col.pathway) <- c('nlrp', 'infl', 'ppar', 'mapk', 'nfkb', 'pi3k')

# volcano plot
jpeg(file = 'documents/manuscripts/IJERPH_submission/figure1.jpg',
    width = 1600, height = 1600, res = 200)
plot(res$Estimate, -log10(res$p.value), pch = 16, cex = .7, col = 'darkgrey',
     xlim = c(-4,4), ylab = '-log10(p)',
     xlab = 'Effect modification per 1 SD increase in BMI by 1 SD increase in residuals of beta',
     main = 'EWIS of DNA methylation and BMI on adult-onset asthma')
points(res[idx.infl, 'Estimate'], -log10(res[idx.infl, 'p.value']),
       pch = 16, cex = .7, col = col.pathway['infl'])
points(res[idx.mapk, 'Estimate'], -log10(res[idx.mapk, 'p.value']),
       pch = 16, cex = .7, col = col.pathway['mapk'])
points(res[idx.pi3k, 'Estimate'], -log10(res[idx.pi3k, 'p.value']),
       pch = 16, cex = .7, col = col.pathway['pi3k'])
points(res[idx.nfkb, 'Estimate'], -log10(res[idx.nfkb, 'p.value']),
       pch = 16, cex = .7, col = col.pathway['nfkb'])
points(res[idx.ppar, 'Estimate'], -log10(res[idx.ppar, 'p.value']),
       pch = 16, cex = .7, col = col.pathway['ppar'])
legend('bottomright', bty = 'n', pch = 16,
       col = col.pathway[c('ppar', 'mapk', 'nfkb', 'pi3k', 'infl')],
       legend = c('Glucocorticoid/PPAR signaling', 'MAPK signaling',
                  'NF-kB signaling', 'PI3K/AKT Signaling',
                  'Global inflammation'))
dev.off()