Noi/NOTES/2012-4-26
Jump to navigation
Jump to search
- Link to calendar: [[1]]
- Previous analysis of mQTL on UCLA data set: http://genome-tech.ucsd.edu/LabNotes/index.php/Noi/NOTES/2012-4-11
- Continued from: http://genome-tech.ucsd.edu/LabNotes/index.php/Noi/NOTES/2012-4-21
mQTL analysis on UCLA-SZ data set (no CpG-SNP)
- Details will be filled later
- Filter out SNP from methylMatrix (Note: the sample IDs were modified to be compatible with tfam file and plink)
./remove_snps.pl ./UCLA-4batches_methyl_min76_minSTD0.1-edit /t4T/GenomeDB/HsGenome/snp134_snv.txt > UCLA-4batches_methyl-noSNP_min76_minSTD0.1
61664 UCLA-4batches_methyl-noSNP_min76_minSTD0.1
- Split the methylMatrix in to 10,000 CpG sites each
split -10000 UCLA-4batches_methyl-noSNP_min76_minSTD0.1 UCLA_noSNP_ 10000 UCLA_noSNP_aa 10001 UCLA_noSNP_ab 10001 UCLA_noSNP_ac 10001 UCLA_noSNP_ad 10001 UCLA_noSNP_ae 10001 UCLA_noSNP_af 1666 UCLA_noSNP_ag
- Run mQTL (no CpG-SNP)
nohup ../mQTL_2MB_plink_allP.pl ../UCLA.hg19SNP ../UCLA_noSNP_aa > 2MB_noNSP_aa nohup ../mQTL_2MB_plink_allP.pl ../UCLA.hg19SNP ../UCLA_noSNP_ab > 2MB_noNSP_ab nohup ../mQTL_2MB_plink_allP.pl ../UCLA.hg19SNP ../UCLA_noSNP_ac > 2MB_noNSP_ac nohup ../mQTL_2MB_plink_allP.pl ../UCLA.hg19SNP ../UCLA_noSNP_ad > 2MB_noNSP_ad nohup ../mQTL_2MB_plink_allP.pl ../UCLA.hg19SNP ../UCLA_noSNP_ae > 2MB_noNSP_ae nohup ../mQTL_2MB_plink_allP.pl ../UCLA.hg19SNP ../UCLA_noSNP_af > 2MB_noNSP_af nohup ../mQTL_2MB_plink_allP.pl ../UCLA.hg19SNP ../UCLA_noSNP_ag > 2MB_noNSP_ag
- Concatenate the out put files to the same file
cat mQTL_noSNP_a*/2MB_noNSP_a* > 2MB_mQTL-noSNP_rawoutput awk '{if ($1 ~ /chr/) print $0}' 2MB_mQTL-noSNP_rawoutput > 2MB_mQTL-noSNP_output
- Number of test
433,947 2MB_mQTL-noSNP_output
- Sort out adj. p-val and count the number of significant associated CpG-SNP roughly at different p-value
sort column7 (adj p_val) sort -g -k7 2MB_mQTL-noSNP_output | awk '{if ($7<=0.05) print $0}' | wc -l --> 6,418 sort -g -k7 2MB_mQTL-noSNP_output | awk '{if ($7<=0.01) print $0}' | wc -l --> 2,292 sort -g -k7 2MB_mQTL-noSNP_output | awk '{if ($7<=0.001) print $0}' | wc -l --> 1,026 sort -g -k7 2MB_mQTL-noSNP_output | awk '{if ($7<=0.0001) print $0}' | wc -l --> 698 sort -g -k7 2MB_mQTL-noSNP_output | awk '{if ($7<=0.00001) print $0}' | wc -l --> 526
- Print out the site at min p-value 0.05 before multiple test correction
sort -g -k7 2MB_mQTL-noSNP_output | awk '{if ($7<=0.05) print $0}' > sorted_2MB_mQTL-noSNP_minPVal0.05.txt 6418 sorted_2MB_mQTL-noSNP_minPVal0.05.txt
- Concatenate QVal to the mQTL file
paste -d"\t" sorted_2MB_mQTL-noSNP_minPVal0.05.txt sorted_2MB_mQTL-noSNP_QVal.txt.txt > sorted_2MB_mQTL-noSNP_minPVal0.05_QVal.txt.txt awk '{if ($8 <=0.1) print $0}' sorted_2MB_mQTL-noSNP_minPVal0.05_QVal.txt.txt > 10%FDR_2MB_mQTL-noSNP 741 10%FDR_2MB_mQTL-noSNP p-value cut off: 0.0001636 or 1.636E-4 *Look at # of associated CpG
awk '{print $1}' 10%FDR_2MB_mQTL-noSNP | sort | uniq -u | wc -l 295
awk '{print $1}' 10%FDR_2MB_mQTL-noSNP | sort | uniq -d | wc -l 121 294+121X = 741 --> X ~ (741-295)/121 ~3.69 this show >1 CpGs associated with the same SNP Total single CpG: 294+121 = 415 *Look at # of associated SNP awk '{print $2}' 10%FDR_2MB_mQTL-noSNP | sort | uniq -u | wc -l --> 195 (subtract header line) awk '{print $2}' 10%FDR_2MB_mQTL-noSNP | sort | uniq -d | wc -l --> 109 195+109Y = 741 --> Y ~ (741-195)/109 ~ 5 this also suggest >1 SNP associated with the same CpG Total single SNP 195+109 = 304
- Separate the CpG-SNP association based on the distance between CpG and SNP
- Note: Dr. Zhang suggested to break down the CpG-SNP distance with equal length
sed s'/:/\t/g' 10%FDR_2MB_mQTL-noSNP > 10%FDR_2MB_mQTL-noSNP_split_pos Run sh print_distance.sh
- Note again: P-value of 10%FDR is p <= 1.636E−4
Associated CpG-SNP distance | Number of assiciated CpG-SNP | % |
10%FDR_10kb | 378 | 50.94% |
10%FDR_10-20kb | 27 | 3.64% |
10%FDR_20-30kb | 21 | 2.83% |
10%FDR_30-40kb | 36 | 4.85% |
10%FDR_40-50kb | 2 | 0.27% |
10%FDR_50kb-2Mb | 278 | 37.47% |
Total in 2Mb | 742 |
Associated CpG-SNP distance | Number of assiciated CpG-SNP | % |
minpval0.05_10kb | 862 | 13.43% |
minpval0.05_10-20kb | 124 | 1.93% |
minpval0.05_20-30kb | 133 | 2.07% |
minpval0.05_30-40kb | 108 | 1.68% |
minpval0.05_40-50kb | 75 | 1.17% |
minpval0.05_50kb-2Mb | 5,116 | 79.71% |
Total in 2Mb | 6,418 |
- Get sequences for manhattan plots
Run get_manht_seq.sh
Comparing the no CpG-SNP mQTL to the whole CpG and SNP mQTL result
- Analyzing the association of CpG and SNP with 10%FDR derived from each analysis
' | no CpG-SNP | SNP-free CpG and CpG-SNP | % of no CpG-SNP to all association |
# of association (10%FDR) | 742 | 3,202 | 23.17 |
significant p-value | 1.636E-04 | 5.448E-04 |
- From the table above, I could see the different of significant p-value from each test. The reason is that the number of tests were different between the test of the whole methylation matrix and the whose SNPs have been filtered out. Therefore the BH correction which deal with the number of test in calculation provide the q-value with at the different significant p-value. I turned out that the significant p-value of no SNP analysis was lower than the whole CpG. This implied that there might be more no CpG-SNP between 1.64E-04 to 5.45E-04 p-value of the 2MB_UCLA_mQTL_3202_10%FDR.txt file (Total significant CpG-SNP association 3,202 sites, 10%FDR)