Ns126:Encode Methylation: Difference between revisions
Jump to navigation
Jump to search
>Shicheng No edit summary |
>Shicheng No edit summary |
||
Line 188: | Line 188: | ||
bedtools intersect -wa -u /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.5.bed | bedtools intersect -wa -u /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.5.bed | ||
=RRBS Data Analysis to Encode Project (Bed/MethylFreq)= | |||
# Data Download | |||
Fastq and bed files can be downloaded from Encode Project. 101 RRBS data (bed files) were downloaded. 2,646,999 CpG loci were covered by 101 RRBS data while 866,979 CpG loci (32.8%) were detected in at least 80% samples. | |||
* | # bed11 to bedGraph | ||
for i in `ls *bed.gz` | |||
do | |||
zcat $i | grep -v "^track" |sort -k1,1 -k2,2n | awk '$5>9 {print $1"\t"$2"\t"$3"\t"$11}'> $i.bedGraph | |||
done | |||
# bedGraph to bigwig | |||
# Download file | |||
wget http://hgdownload.cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeHaibMethylRrbs/files.txt | |||
# md5sum file | |||
wget http://hgdownload.cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeHaibMethylRrbs/md5sum.txt | |||
# Download fastq file with perl script | |||
perl fastqDownloadRRBSEncode.pl files.txt | |||
# Fastq to BAM | |||
Revision as of 08:21, 21 February 2016
RRBS Data Analysis to Encode Project (Fastq)
Aim 1: Haib39bioChain
- Aim 1: methylation haplotype block calling with RRBS dataset
Background
- 1, Encode: http://genome.ucsc.edu/ENCODE/downloads.html
- 2, Encode|RRBS1: http://hgdownload.cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeHaibMethylRrbs/
- 3, Encode|RRBS2: http://genome.ucsc.edu/cgi-bin/hgTrackUi?hgsid=437674359_aUhx08DjchWwtBjyCYv61EB7Yy8S&c=chr1&g=wgEncodeHaibMethylRrbs
- big problem for the fastq files from the above link. 1) not fastq but rather RAS file 2) reads is not true reads.
- 4, Encode|Methy450K: http://hgdownload.cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeHaibMethyl450/
- 5, Encode|software: https://www.encodeproject.org/software/
Method
- Download 39 biochain RRBS dataset from Encode Project (Fastq). Alignment with Bismark and merge all the BAM files
- make the haplotype calling with previous perl script. Finally MHB calling were down with Dr. Zhang perl script.
- Methyfreq based MHB calling were conducted with MethBed files download from UCSC.
Result
- Summary Excel Haib Dataset manifest
- Alignment: Maybe walltime time is too short , not all the samples were aligned completely, so I extend the walltime to 72 hours.
#PBS -q glean #PBS -l nodes=1:ppn=8 #PBS -l walltime=72:00:00 bismark --bowtie2 --phred64-quals --fastq -L 30 -N 1 --multicore 2 /home/shg047/db/hg19/meth/bismark ../fastq_trim/HOT197_trimmed.fq.gz -o ../bam2
- With above setting, the alignment time usage is about 16 hours
- The time is propotional to the size of the fastq file
Sample | N(reads) | N(mapped) | P(mapping) | N(C) | N(MCPG) | N(MCHG) | N(MCHH) | N(UCPG) | N(UCHG) | N(UCHH) | P(MCPG) | P(MCHG) | P(MCHH) |
ENCFF000LVN | 42665329 | 31372051 | 73.50% | 227374577 | 9979920 | 566174 | 2568709 | 16516656 | 50531408 | 147211710 | 37.70% | 1.10% | 1.70% |
ENCFF000LWL | 37938401 | 25076517 | 66.10% | 232580318 | 11488185 | 403775 | 1382726 | 44058434 | 51038500 | 124208698 | 20.70% | 0.80% | 1.10% |
ENCFF000LVR | 36779630 | 25226017 | 68.60% | 216621540 | 12232832 | 438455 | 1846136 | 33041219 | 46480917 | 122581981 | 27.00% | 0.90% | 1.50% |
ENCFF000LVW | 36169867 | 22446797 | 62.10% | 224690704 | 12356273 | 372366 | 871608 | 41484525 | 51598506 | 118007426 | 22.90% | 0.70% | 0.70% |
ENCFF000LVB | 35375019 | 22874896 | 64.70% | 211937214 | 13517127 | 359167 | 1121727 | 31366464 | 47860817 | 117711912 | 30.10% | 0.70% | 0.90% |
ENCFF000LUP | 34214415 | 20212548 | 59.10% | 221862905 | 16119238 | 340947 | 580628 | 49637185 | 51949984 | 103234923 | 24.50% | 0.70% | 0.60% |
ENCFF000LWA | 33924878 | 20895342 | 61.60% | 226298397 | 14274674 | 319602 | 583229 | 52752860 | 53011846 | 105356186 | 21.30% | 0.60% | 0.60% |
ENCFF000LUV | 33798173 | 23173663 | 68.60% | 249630261 | 10489904 | 225941 | 428583 | 62372541 | 57063134 | 119050158 | 14.40% | 0.40% | 0.40% |
ENCFF000LWY | 33715116 | 19999502 | 59.30% | 213287381 | 13105506 | 316152 | 583822 | 46433586 | 48748868 | 104099447 | 22.00% | 0.60% | 0.60% |
ENCFF000LWW | 33642793 | 20282649 | 60.30% | 216267059 | 15823449 | 340126 | 638093 | 43292443 | 49869268 | 106303680 | 26.80% | 0.70% | 0.60% |
ENCFF000LWP | 32858930 | 20654452 | 62.90% | 219544115 | 13223974 | 319601 | 599593 | 48455363 | 50129553 | 106816031 | 21.40% | 0.60% | 0.60% |
ENCFF000LUU | 32815340 | 17866807 | 54.40% | 196534662 | 12582324 | 421909 | 799252 | 48097172 | 45287260 | 89346745 | 20.70% | 0.90% | 0.90% |
ENCFF000LVF | 31276946 | 22130843 | 70.80% | 241339808 | 11057578 | 348695 | 664774 | 64178792 | 55276052 | 109813917 | 14.70% | 0.60% | 0.60% |
ENCFF000LXB | 29681695 | 14796658 | 49.90% | 160388755 | 8182288 | 279840 | 480680 | 41788931 | 37024968 | 72632048 | 16.40% | 0.80% | 0.70% |
ENCFF000LWK | 29349446 | 12860133 | 43.80% | 137008343 | 8665157 | 302283 | 567461 | 29580229 | 31441673 | 66451540 | 22.70% | 1.00% | 0.80% |
ENCFF000LWE | 29211774 | 11162094 | 38.20% | 115410335 | 8071486 | 277005 | 511615 | 23863725 | 26859483 | 55827021 | 25.30% | 1.00% | 0.90% |
ENCFF000LVK | 27335110 | 16938233 | 62.00% | 193589572 | 12327753 | 322954 | 495345 | 51995628 | 43692466 | 84755426 | 19.20% | 0.70% | 0.60% |
ENCFF000LVO | 26190075 | 19276422 | 73.60% | 150096025 | 7012551 | 300245 | 1332640 | 19959469 | 31819734 | 89671386 | 26.00% | 0.90% | 1.50% |
ENCFF000LVA | 25660341 | 16009728 | 62.40% | 152217865 | 12748317 | 546338 | 1806399 | 26835875 | 36223859 | 74057077 | 32.20% | 1.50% | 2.40% |
ENCFF000LWD | 25467744 | 15459206 | 60.70% | 149205383 | 11335949 | 535669 | 1829440 | 30226946 | 34383379 | 70894000 | 27.30% | 1.50% | 2.50% |
ENCFF000LVU | 23511285 | 15078397 | 64.10% | 141659239 | 8631088 | 272002 | 734661 | 21915485 | 33328412 | 76777591 | 28.30% | 0.80% | 0.90% |
ENCFF000LWO | 22844877 | 14206467 | 62.20% | 153005363 | 9813445 | 247506 | 461786 | 33745829 | 35223721 | 73513076 | 22.50% | 0.70% | 0.60% |
ENCFF000LVI | 22656883 | 13335953 | 58.90% | 143634615 | 13714133 | 389619 | 650547 | 32241351 | 34451974 | 62186991 | 29.80% | 1.10% | 1.00% |
ENCFF000LUQ | 22247066 | 13650824 | 61.40% | 150112632 | 9777330 | 287617 | 504659 | 38167142 | 34898963 | 66476921 | 20.40% | 0.80% | 0.80% |
ENCFF000LUT | 22240097 | 14849726 | 66.80% | 162779228 | 9543204 | 285378 | 485629 | 41824570 | 38682843 | 71957604 | 18.60% | 0.70% | 0.70% |
ENCFF000LWH | 21943620 | 12872848 | 58.70% | 137662143 | 11585659 | 414164 | 781429 | 29526002 | 32473470 | 62881419 | 28.20% | 1.30% | 1.20% |
ENCFF000LVE | 21473578 | 14519161 | 67.60% | 157450743 | 7444292 | 221900 | 417095 | 40943018 | 36372568 | 72051870 | 15.40% | 0.60% | 0.60% |
ENCFF000LWX | 21208447 | 13218301 | 62.30% | 145215571 | 9668358 | 223882 | 356827 | 34200783 | 33840207 | 66925514 | 22.00% | 0.70% | 0.50% |
ENCFF000LVV | 21080863 | 13197949 | 62.60% | 142640926 | 10137634 | 255425 | 440789 | 32877645 | 33833167 | 65096266 | 23.60% | 0.70% | 0.70% |
ENCFF000MLE | 20338545 | 12207095 | 60.00% | 110465581 | 11201397 | 457519 | 1827464 | 19226614 | 25845153 | 51907434 | 36.80% | 1.70% | 3.40% |
ENCFF000LWS | 20241908 | 11560661 | 57.10% | 126245157 | 11387405 | 245297 | 424016 | 29669508 | 30429870 | 54089061 | 27.70% | 0.80% | 0.80% |
ENCFF000LVZ | 20058311 | 11004008 | 54.90% | 120506134 | 9791631 | 242794 | 448332 | 25470114 | 28547035 | 56006228 | 27.80% | 0.80% | 0.80% |
ENCFF000LWT | 19407909 | 10887950 | 56.10% | 119134111 | 9343415 | 203527 | 331921 | 25331949 | 28315688 | 55607611 | 26.90% | 0.70% | 0.60% |
ENCFF000MLD | 19184685 | 14418868 | 75.20% | 162604971 | 6716952 | 305536 | 587608 | 54429999 | 38467374 | 62097502 | 11.00% | 0.80% | 0.90% |
ENCFF000MLP | 18102015 | 10032177 | 55.40% | 109233809 | 10521402 | 215592 | 370895 | 25497472 | 26468867 | 46159581 | 29.20% | 0.80% | 0.80% |
ENCFF000LVJ | 17856591 | 11324104 | 63.40% | 102317536 | 8591993 | 307833 | 1097735 | 17653908 | 24511494 | 50154573 | 32.70% | 1.20% | 2.10% |
ENCFF000MLJ | 16177211 | 8675428 | 53.60% | 94678187 | 9430914 | 200777 | 318563 | 22135878 | 23265164 | 39326891 | 29.90% | 0.90% | 0.80% |
ENCFF000LUN | 10913594 | 3643074 | 33.40% | 40450477 | 3011262 | 76408 | 126647 | 9141529 | 9549119 | 18545512 | 24.80% | 0.80% | 0.70% |
ENCFF000MLM | 5127779 | 2931330 | 57.20% | 36675584 | 1281379 | 40172 | 75178 | 10725366 | 8687156 | 15866333 | 10.70% | 0.50% | 0.50% |
- MHB regions identification
#!/bin/csh #PBS -N bam2MHB #PBS -q pdafm #PBS -l nodes=1:ppn=16 #PBS -l walltime=72:00:00 #PBS -o bam2MHB.log #PBS -e bam2MHB.err #PBS -V #PBS -M shihcheng.guo@gmail.com #PBS -m abe #PBS -A k4zhang-group cd /home/shg047/oasis/Haib/sortBam # samtools cat -h header.sam -o haib.merge.bam *sort.bam samtools sort -@ 16 haib.encode.merge.bam -o haib.merge.sort.bam samtools index haib.merge.sort.bam bedtools genomecov -bg -split -ibam haib.merge.sort.bam > haib.merge.bam.pool.bed awk '$4>9 { print $1"\t"$2"\t"$3}' haib.merge.bam.pool.bed | bedtools merge -d 10 -i - > haib.RD10.genomecov.bed awk '$3-$2>80 {print $1"\t"$2"\t"$3"\t"$3-$2+1}' haib.RD10.genomecov.bed > haib.RD10_80up.genomecov.bed
- Statistic: haib.RD10_80up.genomecov.bed
cat haib.RD10_80up.genomecov.bed|awk '{sum+=$4} END { print "N = ", NR, "Sum = ", sum, " Average = ",sum/NR}' N = 120994 Sum = 17702479 Average = 146.309
- haploinfo to MHB
cd /home/shg047/oasis/Haib/sortBam /home/shg047/oasis/Haib/mhb/haib.RD10_80up.genomecov.bed /home/shg047/oasis/Haib/mhb/haib.merge.sort.bam /home/shg047/oasis/Haib/hapInfo2mld_blocks.pl ../mergedBam2hapInfo.pl ./haib.RD10_80up.genomecov.bed haib.merge.sort.bam > Haib.merge.RD10_80up.hapinfo.txt # get hapinfo ../hapInfo2mld_block.pl ./Haib.merge.RD10_80up.hapinfo.txt 0.5 > Haib.merge_RD10_80up.mld_blocks_r2-0.5.bed
- MHB identified with different threshold: R-square from 0.1-0.9
/home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.1.bed /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.2.bed /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.3.bed /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.4.bed /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.5.bed /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.6.bed /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.7.bed
R-square threshold | MHB counts |
0.1 | 14933 |
0.2 | 13367 |
0.3 | 11667 |
0.4 | 9754 |
0.5 | 8155 |
0.6 | 7683 |
0.7 | 7445 |
bedtools intersect -wa -u /home/shg047/oasis/Haib/mhb/Haib.merge_RD10_80up.mld_blocks_r2-0.5.bed
RRBS Data Analysis to Encode Project (Bed/MethylFreq)
- Data Download
Fastq and bed files can be downloaded from Encode Project. 101 RRBS data (bed files) were downloaded. 2,646,999 CpG loci were covered by 101 RRBS data while 866,979 CpG loci (32.8%) were detected in at least 80% samples.
- bed11 to bedGraph
for i in `ls *bed.gz` do zcat $i | grep -v "^track" |sort -k1,1 -k2,2n | awk '$5>9 {print $1"\t"$2"\t"$3"\t"$11}'> $i.bedGraph done
- bedGraph to bigwig
- Download file
wget http://hgdownload.cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeHaibMethylRrbs/files.txt
- md5sum file
wget http://hgdownload.cse.ucsc.edu/goldenPath/hg19/encodeDCC/wgEncodeHaibMethylRrbs/md5sum.txt
- Download fastq file with perl script
perl fastqDownloadRRBSEncode.pl files.txt
- Fastq to BAM