Kun:LabNotes/MONOD/2014-7-6
Jump to navigation
Jump to search
Analysis of 140702_RRBS data[edit]
Introduction[edit]
- Noi performed RRBS experiment on Pancreatic tumor tissue and matched plasma from Moores Cancer center, plus two Colon cancer samples from Biochain and normal plasma controls from Kang Zhang's group.
- These samples were sequenced in one HiSeq Rapid run (140702_RRBS) and one MiSeq run (20140613_MiSeq).
- Dinh did the read trimming and mapping: [5]; Mapping summary: Media:140702_RRBS_v2.xlsx
- I started with the bam files and BED files Dinh generated: /media/Ext12T/DD_Ext12T/140702_RRBS
- Some quick checks of the data.
- Numbers of CpG sites called (>=10x)
CTT-FFPE-100ng_1 | 1,405,378 |
CTT-FFPE-100ng_2 | 1,387,013 |
CTT-FFPE-5ng_1 | 867,737 |
CTT-FFPE-5ng_2 | 1,226,656 |
CTT-frozen-100ng_1 | 1,356,236 |
CTT-frozen-100ng_2 | 1,350,437 |
CTT-frozen-5ng_1 | 941,422 |
CTT-frozen-5ng_2 | 1,341,451 |
PC-T-1_1 | 1,222,248 |
PC-T-1_2 | 1,101,232 |
PC-T-2_1 | 1,520,701 |
PC-T-2_2 | 1,181,493 |
PC-T-4_1 | 1,203,834 |
PC-T-4_2 | 1,153,934 |
PC-T-6_1 | 1,427,700 |
PC-T-6_2 | 1,334,831 |
PC-T-7_1 | 1,186,899 |
PC-T-7_2 | 1,278,534 |
PC-P-1 | 75,392 |
PC-P-2 | 227,194 |
PC-P-3 | 24,801 |
PC-P-4 | 76,362 |
PC-P-5 | 62,956 |
PC-P-6 | 88,236 |
PC-P-7 | 78,418 |
PC-P-8 | 337,582 |
PC-P-9 | 172,890 |
PC-P-10 | 102,832 |
NC-P-1 | 155,577 |
NC-P-2 | 290,084 |
NC-P-3 | 527,864 |
NC-P-5 | 142,381 |
NC-P-6 | 205,255 |
NC-P-7 | 165,611 |
NC-P-8 | 45,906 |
NC-P-9 | 35,724 |
Hierarchical clustering of all samples File:140702 RRBS all samples clustering.png CTT: fresh frozen versus FFPE (100ng) File:140702 RRBS CTT FreshFrozen vs FFPE 100ng.pngFile:140702 RRBS CTT FreshFrozen vs FFPE 100ng 2.png 5ng vs 100ng File:140702 RRBS CTT FreshFrozen 5ng vs 100ng.pngFile:140702 RRBS CTT FFPE 5ng vs 100ng.png
Step 0: target identification.[edit]
- The idea here is to extensively search the entire genome and identify regions in which there is no or little methylation in whole blood.
- For this I used a combination of WGBS, RRBS, BSPP and ILMN450k data on whole blood.
- I first identified individual CpG sites (LMS) that have low methylation level (initially use mC<0.05; then relaxed the cutoff to 0.1).
- Then I group these LMSs into clusters.
/home/kunzhang/CpgMIP/MONOD/Public_data/find_LMS.pl > whole_blood_LMS.txt ./extract_clusters.pl whole_blood_LMS.txt > whole_blood_LMS_clusters.txt
mC<0.05 mC<0.1 #LMS 252,923 288,169 #LMS clusters 16,208 18,272
- Finally I manually converted whole_blood_LMS_clusters.txt into the BED format using Excel.
whole_blood_LMS_2_clusters_BED.txt (mC < 0.1) whole_blood_LMS_clusters_BED.txt (mC < 0.05)
Step 1: extract methylation haplotypes.[edit]
- I wrote a script to extract methylation haplotypes within a list of target regions from a merged bam file: mergedBam2hapInfo.pl
- I created a batch shell script to process all data sets at once: step1_batch_command_3.sh
Step 2: identify high-methylation haplotypes (HMHs)[edit]
- Because whole blood has very little methylation in the target regions, I considered a haplotype that contain at least three Cs as a high methylation haplotype (HMH) for the sake of simplicity. Such HMHs should come from the cells of different tissues, potentially from cancer cells.
- I wrote a script methHapClassfier.pl to separate all haplotypes into LMHs and HMHs. Note that in the near future we might need to use a more sophisticated approach to classify HMHs.
- Again I created a batch shell script to process all samples at once: step2_batch_get_methHapCounts.sh
Step 3: trim target list based on HMHs[edit]
- I found all regions that contain at least one HMH in the tumor tissues (colon and pancreatic).
- Then I reported and compared the HMH counts in these regions for all samples.
sh step3_batch_get_HMH_regions.sh
Step 4: compare HMH loads across all samples[edit]
sh step4_batch_report_HMH_regions_counts.sh ./step5_report_summary.pl
Fraction of HMH in all haplotypes. File:140702 RRBS HMH abundance.png
Enrichment of HMH compared with normal plasma (based on Z-scores) File:140702 RRBS HMH enrichment.png
- PC-P-2(UCSD-004-05), PC-P-4(UCSD-004-07), PC-P-10(UCSD-004-13) showed significant enrichment of methylation haplotypes.
Inspect individual regions[edit]
- Prepare two sampleInfo files for UCSD-004-05 and UCSD-004-07.
- Report all regions in which HMHs were present in both tumor and matched plasma, but absence in normal plasma.
./report_tumor_HMH_regions.pl UCSD-004-05_sampleInfo.txt | sort -k 3nr > UCSD-004-05_HMH_regions.txt ./report_tumor_HMH_regions.pl UCSD-004-07_sampleInfo.txt | sort -k 3nr > UCSD-004-07_HMH_regions.txt
- Then I picked individual regions and plotted methylation haplotypes in such regions.
1: Combined multiple hapInfo files cat s_1_1_ILMN_Indx04.hapInfo.txt s_1_1_ILMN_Indx05.hapInfo.txt s_1_2_ILMN_Indx04.hapInfo.txt s_1_2_ILMN_Indx05.hapInfo.txt > PC-T-2.pooled.hapInfo.txt cat s_2_1_ILMN_Indx02.hapInfo.txt s_2_2_ILMN_Indx02.hapInfo.txt > PC-P-2.pooled.hapInfo.txt cat s_2_1_ILMN_Indx14.hapInfo.txt s_2_2_ILMN_Indx14.hapInfo.txt s_2_1_ILMN_Indx15.hapInfo.txt s_2_2_ILMN_Indx15.hapInfo.txt s_2_1_ILMN_Indx16.hapInfo.txt s_2_2_ILMN_Indx16.hapInfo.txt s_2_1_ILMN_Indx27.hapInfo.txt s_2_2_ILMN_Indx27.hapInfo.txt NP-RRBS-NC-P-1ng-p2-Jun20-A013.hapInfo.txt > NC-plasma.pooled.hapInfo.txt 2: Prepare data files for plotting. ../prepare_plotting_files.pl NC-plasma.pooled.hapInfo.txt PC-T-2.pooled.hapInfo.txt PC-P-2.pooled.hapInfo.txt chr6:152128536-152129155 ../prepare_plotting_files.pl NC-plasma.pooled.hapInfo.txt PC-T-2.pooled.hapInfo.txt PC-P-2.pooled.hapInfo.txt chr5:176543913-176544076 3: Upload each "hapMatrix" file to Methylation-Plotter to make a haplotype figure.
chr6:152128536-152129155 (ESR1 promoter) Normal Plasma Primary Tumor Patient Plasma File:Chr6.152128536-152129155 hapMatrix A.txt main plot Mon Jul 07 19 50 59 2014.pngFile:Chr6.152128536-152129155 hapMatrix B.txt main plot Mon Jul 07 19 49 45 2014.pngFile:Chr6.152128536-152129155 hapMatrix Mix.txt main plot Mon Jul 07 19 50 29 2014.png
2nd analysis updates[edit]
- Dinh mapped Read 1 and Read 2 separately, and report the mapping data in two separate bam files. In my first analysis, I extracted haplotypes from the two files and treat them separately. In fact each pair of reads came from one DNA molecule and should be combined for one single haplotype. Treating the two reads separately leads to double counting and also some short haplotypes that might be filtered out.
- First I merged two bam files from each library, sorted the merged files and re-created the bam indexes: batch_merge_N_sort.sh
- Then I modified the mergedBam2hapInfo.pl script in order to combine each pair of reads into one haplotype. Note that the fastq Read IDs are different between MiSeq and HiSeq data. It took me some debugging in order to process both types of data correctly.
- Next I re-ran all the batch processing scripts, with two more changes: (i) I required that a cancer haplotype should have at least four methylated CpG sites (increase the stringency); (ii) Target regions that have at least 1 HMH are included (previously I required 2 or more).
- With these changes, I managed to detect significant (Z score > 3) cancer load in five patient serum samples.
Fraction of HMH in all haplotypes. File:140702 RRBS HMH abundance PE-haps.png
Enrichment of HMH compared with normal plasma (based on Z-scores) File:140702 RRBS HMH enrichment PE-haps.png