Yan Holtz, Zhihong Zhu, Julanne Frater, Perry Bartlett, Jian Yang, John McGrath
This document applies the SMR tool on both the Expression and Methylation data in the same time. It also takes into account the descriptive classification of chromatin.
Methylation and expression dataset are extensively described in the previous tabs.
Run SMR on both expression and methylation datasets in the same time.
# General parameters
# Good folder
cd /shares/compbio/Group-Wray/YanHoltz/VITAMIND_XIA_ET_AL/10_INTEGRATION
# list of option and file:
ref="/gpfs/gpfs01/polaris/Q0286/UKBiobank/v2EURu_HM3/ukbEURu_imp_chr11_v2_HM3_QC"
sum="/shares/compbio/Group-Wray/YanHoltz/VITAMIND_XIA_ET_AL/1_GWAS/GWAS_vitaminD_XiaEtAL.ma"
rug="/shares/compbio/Group-Wray/YanHoltz/DATA/EQTL/eQTL_DATA_EQTLGEN_CONSORTIUM_32K/WITH_GENE_NAME/rug_32k_cis"
genelist="/shares/compbio/Group-Wray/YanHoltz/DATA/GENE_POSITION/glist_hg19_refseq.txt"
window="1000"
psmr="5e-8"
thread="6"
# list of option and file:
METHY="/gpfs/gpfs01/polaris/Q0286/uqywu16/omics/methy/bl_meqtl_cis_std_chr11"
probe="ENSG00000152270"
out="./integration_result_Chr11_15Mb"
# Run the analysis
tmp_command="smr_Linux --bfile $ref --gwas-summary $sum --beqtl-summary $METHY --beqtl-summary $rug --plot --probe $probe --probe-wind $window --gene-list $genelist --psmr $psmr --out $out --thread-num $thread"
qsubshcom "$tmp_command" 1 30G smr_integration 10:00:00 ""
# list of option and file:
METHY="/gpfs/gpfs01/polaris/Q0286/uqywu16/omics/methy/bl_meqtl_cis_std_chr11"
probe="ENSG00000254682"
out="./integration_result_Chr11_71Mb"
# Run the analysis
tmp_command="smr_Linux --bfile $ref --gwas-summary $sum --beqtl-summary $METHY --beqtl-summary $rug --plot --probe $probe --probe-wind $window --gene-list $genelist --psmr $psmr --out $out --thread-num $thread"
qsubshcom "$tmp_command" 1 30G smr_integration 10:00:00 ""
It is not possible to create a figure for the chromosome 4 loci since there is no eSMR significant result.
# list of option and file:
METHY="/gpfs/gpfs01/polaris/Q0286/uqywu16/omics/methy/bl_meqtl_cis_std_chr4"
probe="ENSG00000080493"
out="./integration_result_Chr4"
window="4000"
# Run the analysis
tmp_command="smr_Linux --bfile $ref --gwas-summary $sum --beqtl-summary $METHY --beqtl-summary $rug --plot --probe $probe --probe-wind $window --gene-list $genelist --psmr $psmr --out $out --thread-num $thread"
qsubshcom "$tmp_command" 1 30G smr_integration 10:00:00 ""
How to read the plots (from top to bottom):
IMPORTANT NOTE: I made this plot with a HEIDI threshold of 0 = all significant to avoid a bug in the script.
# good folder
cd /shares/compbio/Group-Wray/YanHoltz/VITAMIND_XIA_ET_AL/10_INTEGRATION
# Create an R script to make the plot for a specific region and put it in a file script.R:
# Read the data file in R:
args <- commandArgs(TRUE)
input <- args[1]
output <- args[2]
source("/shares/compbio/Group-Wray/YanHoltz/SOFT/plot_OmicsSMR_24_April.r")
SMRData = ReadomicSMRData(input)
# Plot
pdf(output, height=15, width=18)
omicSMRLocusPlot(data=SMRData,
# Threshold for SMR tests:
esmr_thresh=6.4*10^-6, msmr_thresh=5.9*10^-07, m2esmr_thresh=5.9*10^-07,
# Threshold for HEIDI tests:
esmr_heidi = 0, msmr_heidi = 0, m2esmr_heidi = 0,
# Other features
max_anno_probe = 6, anno_methyl=TRUE, max_plot_mprobe = 1, window = 700, epi_plot=T )
dev.off()
# Send the script for all the regions I need.
tmp_command="Rscript script.R /shares/compbio/Group-Wray/YanHoltz/VITAMIND_XIA_ET_AL/10_INTEGRATION/plot/integration_result_Chr11_15Mb.ENSG00000152270.txt plot_integ_ENSG00000152270.pdf"
qsubshcom "$tmp_command" 1 10G smr_plot 10:00:00 ""
tmp_command="Rscript script.R /shares/compbio/Group-Wray/YanHoltz/VITAMIND_XIA_ET_AL/10_INTEGRATION/plot/integration_result_Chr11_71Mb.ENSG00000254682.txt plot_integ_ENSG00000254682.pdf"
qsubshcom "$tmp_command" 1 10G smr_plot 10:00:00 ""
# Trasfer result locally
cd /Users/y.holtz/Dropbox/QBI/4_UK_BIOBANK_GWAS_PROJECT/VitaminD-GWAS/IMG/SMR_PLOT
scp y.holtz@delta.imb.uq.edu.au://shares/compbio/Group-Wray/YanHoltz/VITAMIND_XIA_ET_AL/10_INTEGRATION/*pdf .
Observation
rs12785878
, position 71,167,449. Its allele T
is associated with a higher vitamin D concentration.transcript
in every tissue EXCEPT brain where it is Enha = Enhancer
and EnhW = `Weak Enhancer.Litterature
Need to read this publication which goes deep into this loci.
Gene DHCR7: a 7-dehydrocholesterol reductase: it is involved in the cholesterol pathway. Its role is hydroxylation.
Observation
Detail of function arguments:
# Plot caractéristics
data=SMRData
esmr_thresh=0.05
msmr_thresh=0.05
m2esmr_thresh=0.05
esmr_heidi = 0
msmr_heidi = 0
m2esmr_heidi = 0
probeNEARBY=c("ENSG00000152270")
interest_only=T
max_anno_probe = 6
anno_methyl=T
max_plot_mprobe = 1
window = 700
epi_plot=T
esmr_thresh_plot=NULL
msmr_thresh_plot=NULL
highlight=NULL
pointsize=20
anno_self=TRUE
ASO=TRUE
anno_dist=1.05
trait_name=NULL
A work by Yan Holtz
Yan.holtz.data@gmail.com