Ns126:Xliu2014GB: Difference between revisions

From ZhangLabWiki
Jump to navigation Jump to search
>Shicheng
>Shicheng
Line 330: Line 330:
  cd /oasis/tscc/scratch/shg047/Xliu2014/bam
  cd /oasis/tscc/scratch/shg047/Xliu2014/bam
  bismark_methylation_extractor --paired-end --bedGraph --multicore 3 --ignore 3 --ignore_3prime 3 --ignore_r2 5 --ignore_3prime_r2 5 --gzip --buffer_size 4G --zero_based --comprehensive --output ../methyfreq SRR1232311_1_val_1.fq.gz_bismark_bt2_pe.bamrm CHG* CpG* CHH* *bedGraph.gz *bismark.cov.gz
  bismark_methylation_extractor --paired-end --bedGraph --multicore 3 --ignore 3 --ignore_3prime 3 --ignore_r2 5 --ignore_3prime_r2 5 --gzip --buffer_size 4G --zero_based --comprehensive --output ../methyfreq SRR1232311_1_val_1.fq.gz_bismark_bt2_pe.bamrm CHG* CpG* CHH* *bedGraph.gz *bismark.cov.gz
 
Detail Parameters as Following:
  Parameters used to extract methylation information:
  Parameters used to extract methylation information:
  Bismark Extractor Version: v0.14.4
  Bismark Extractor Version: v0.14.4

Revision as of 19:15, 24 February 2016

Background

  • Construction protocol: Genomic DNA from the tumor and corresponding adjacent tissues was prepared using the QIAamp DNA Blood Mini Kit(Qiagen) following the manufacturer’s instructions. Prior to the library construction, 3µg of genomic DNA from each sample was fragmented using a Covarias sonication system to mean sizes of approximately 200-300bp. After fragmentation, libraries were constructed according to the Illumina Paired-End protocol. Briefly, the purified, randomly fragmented DNA was treated with a mix of T4 DNA polymerase, Klenow fragments, T4 polynucleotide kinase and a nucleotide triphosphate mix to repair the ends by blunting and phosphorylation. The blunted DNA fragments were subsequently 3’-adenylated using the Klenow fragment (3’-5’exo) and ligated by T4 DNA ligase to adapters synthesized with 5’ –methyl-cytosine instead of cytosine. The adaptor-ligated library was purified using AMPure XP beads (Beckman Coulter Genomics). 500 ng of each library was hybridized to Agilent SureSelect Methyl-Seq biotinylated RNA baits (84 Mb) for 24 h at 65 C. Biotinylated target hybrids were captured on Dynal MyOne Streptavidin T1 (Invitrogen), and purified through MinElute PCR column (Qiagen). Bisulfite conversion of the purified captured library was performed using the EZ DNA Methylation Gold Kit (Zymo Research) as per manufacturer’s instructions. The bisulfite converted captured library was amplified by PCR with 14 PCR cycles, then purified by AMPure XP beads and quantified by Agilent 2100 system and qPCR. Methyl-Seq libraries were then sequenced on the HiSeq 2000 platform according to the manufacturer’s instructions, and 2 ×90bp paired-end reads were generated. Methyl-Seq libraries were prepared for sequencing using standard Illumina protocols.
  • We propose a statistical algorithm MethylPurify that uses regions with bisulfite reads showing discordant methylation levels to infer tumor purity from tumor samples alone. With purity estimate, MethylPurify can identify differentially methylated regions (DMRs) from individual tumor samples without genomic variation information or prior knowledge from other datasets. In simulations with mixed bisulfite reads from cancer and normal cell lines, MethylPurify correctly inferred tumor purity and identified over 96% of the DMRs. On real patient data where tumor to normal comparison were used as golden standard, MethylPurify called DMR from tumor samples alone at over 57% sensitivity and 91% specificity.
  • Lung adenocarcinoma cancer and normal tissues from 5 patients were captured by Agilent SureSelect Methyl-Seq system, followed by bisulfite sequencing.

Method

Sample

title geo_accession channel_count tissue patient id smoking Gender Alignment processing sequencing library relation relation relation
75A-2 GSM1367123 1 adjacent nomal lung tissue 75 smoker male bsmap MethylPurify HiSeq 2000 Bisulfite-Seq SAMN02725476 SRX515122 SRR1232302
75B GSM1367124 1 lung adenocarcinoma tissue 75 smoker male bsmap MethylPurify HiSeq 2001 Bisulfite-Seq SAMN02725471 SRX515123 SRR1232303
109A-2 GSM1367125 1 adjacent nomal lung tissue 109 non-smoker male bsmap MethylPurify HiSeq 2002 Bisulfite-Seq SAMN02725468 SRX515124 SRR1232304
109B GSM1367126 1 lung adenocarcinoma tissue 109 non-smoker male bsmap MethylPurify HiSeq 2003 Bisulfite-Seq SAMN02725475 SRX515125 SRR1232305
109C GSM1367127 1 corresponding lymph node metastasis tissue 109 non-smoker male bsmap MethylPurify HiSeq 2004 Bisulfite-Seq SAMN02725470 SRX515126 SRR1232306
137A GSM1367128 1 adjacent nomal lung tissue 137 non-smoker female bsmap MethylPurify HiSeq 2005 Bisulfite-Seq SAMN02725467 SRX515127 SRR1232307
137B GSM1367129 1 lung adenocarcinoma tissue 137 non-smoker female bsmap MethylPurify HiSeq 2006 Bisulfite-Seq SAMN02725477 SRX515128 SRR1232308
137C GSM1367130 1 corresponding lymph node metastasis tissue 137 non-smoker female bsmap MethylPurify HiSeq 2007 Bisulfite-Seq SAMN02725469 SRX515129 SRR1232309
156A-2 GSM1367131 1 adjacent nomal lung tissue 156 non-smoker female bsmap MethylPurify HiSeq 2008 Bisulfite-Seq SAMN02725472 SRX515130 SRR1232310
156B-2 GSM1367132 1 lung adenocarcinoma tissue 156 non-smoker female bsmap MethylPurify HiSeq 2009 Bisulfite-Seq SAMN02725478 SRX515131 SRR1232311
201A GSM1367133 1 adjacent nomal lung tissue 201 non-smoker female bsmap MethylPurify HiSeq 2010 Bisulfite-Seq SAMN02725474 SRX515132 SRR1232312
201B GSM1367134 1 lung adenocarcinoma tissue 201 non-smoker female bsmap MethylPurify HiSeq 2011 Bisulfite-Seq SAMN02725473 SRX515133 SRR1232313

Sample Matching Table

SRR_Normal SRR_Cancer patient id smoking Gender
SRR1232302 SRR1232303 75 smoker male
SRR1232304 SRR1232305 109 non-smoker male
SRR1232307 SRR1232308 137 non-smoker female
SRR1232310 SRR1232311 156 non-smoker female
SRR1232312 SRR1232313 201 non-smoker female

SRA Download

  • sra saved in: /home/shg047/ncbi/public/sra
prefetch -v SRR1232302 &
prefetch -v SRR1232303 &
prefetch -v SRR1232304 &
prefetch -v SRR1232305 &
prefetch -v SRR1232306 &
prefetch -v SRR1232307 &
prefetch -v SRR1232308 &
prefetch -v SRR1232309 &
prefetch -v SRR1232310 &
prefetch -v SRR1232311 &
prefetch -v SRR1232312 &
prefetch -v SRR1232313 &

SRA to Fastq

fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232302.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232303.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232304.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232305.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232306.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232307.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232308.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232309.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232310.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232311.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232312.sra &
fastq-dump --outdir /home/shg047/oasis/Xliu2014 --split-files /home/shg047/ncbi/public/sra/SRR1232313.sra &

Fastq to Bam

trim_glore

#!/bin/csh
#PBS -q glean
#PBS -l nodes=1:ppn=1
#PBS -l walltime=6:00:00
#PBS -o SRR1232309_2.log
#PBS -e SRR1232309_2.err
#PBS -V
#PBS -M shihcheng.guo@gmail.com
#PBS -m abe
#PBS -A k4zhang-group
cd /oasis/tscc/scratch/shg047/Xliu2014/fastq
gzip SRR1232309_2.fastq
trim_galore --phred33 --fastqc --illumina --clip_R1 10 --clip_R2 10 --three_prime_clip_R1 5 --three_prime_clip_R1 5 SRR1232309_2.fastq --output_dir ../fastq_trim
Sample N(reads) P(adaptor) P(pass1) Base(filting) PB(pass2)
SRR1232302_1.fastq.gz_trimming_report.txt 10901814 34.00% 100.00% 2.10% 97.20%
SRR1232302_2.fastq.gz_trimming_report.txt 10896296 55.10% 100.00% 1.80% 97.40%
SRR1232303_1.fastq.gz_trimming_report.txt 20216497 34.20% 100.00% 2.00% 97.40%
SRR1232303_2.fastq.gz_trimming_report.txt 20227949 54.60% 100.00% 1.70% 97.50%
SRR1232304_1.fastq.gz_trimming_report.txt 11083565 32.40% 100.00% 2.40% 97.00%
SRR1232304_2.fastq.gz_trimming_report.txt 11080866 54.30% 100.00% 2.10% 97.20%
SRR1232305_1.fastq.gz_trimming_report.txt 19402456 36.20% 100.00% 1.80% 97.50%
SRR1232305_2.fastq.gz_trimming_report.txt 19481098 54.80% 100.00% 1.60% 97.60%
SRR1232306_1.fastq.gz_trimming_report.txt 51510022 32.10% 100.00% 2.00% 97.30%
SRR1232306_2.fastq.gz_trimming_report.txt 51401414 52.70% 100.00% 1.80% 97.50%
SRR1232307_1.fastq.gz_trimming_report.txt 15046190 37.80% 100.00% 2.60% 96.40%
SRR1232307_2.fastq.gz_trimming_report.txt 15015457 60.50% 100.00% 2.30% 96.50%
SRR1232308_1.fastq.gz_trimming_report.txt 35592228 40.50% 100.00% 2.40% 96.50%
SRR1232308_2.fastq.gz_trimming_report.txt 35440880 61.50% 100.00% 2.20% 96.60%
SRR1232309_1.fastq.gz_trimming_report.txt 34522276 36.20% 100.00% 2.60% 96.40%
SRR1232309_2.fastq.gz_trimming_report.txt 34294567 59.10% 100.00% 2.40% 96.50%
SRR1232310_1.fastq.gz_trimming_report.txt 24123686 28.40% 100.00% 17.70% 81.70%
SRR1232310_2.fastq.gz_trimming_report.txt 23880345 54.90% 100.00% 2.40% 96.90%
SRR1232311_1.fastq.gz_trimming_report.txt 12239714 28.80% 100.00% 17.20% 82.30%
SRR1232311_2.fastq.gz_trimming_report.txt 12157130 54.70% 100.00% 2.20% 97.00%
SRR1232312_1.fastq.gz_trimming_report.txt 9216970 33.50% 100.00% 2.60% 96.70%
SRR1232312_2.fastq.gz_trimming_report.txt 9221739 56.30% 100.00% 2.60% 96.60%
SRR1232313_1.fastq.gz_trimming_report.txt 32811085 34.60% 100.00% 2.60% 96.70%
SRR1232313_2.fastq.gz_trimming_report.txt 32789834 56.70% 100.00% 2.60% 96.60%

Alignment

#!/bin/csh
#PBS -n Xie.LUN.SRR1232310_1.fastq.gz
#PBS -q glean
#PBS -l nodes=1:ppn=1
#PBS -l walltime=72:00:00
#PBS -o SRR1232310_1.fastq.gz.trim.log
#PBS -e SRR1232310_1.fastq.gz.trim.err
#PBS -V
#PBS -M shihcheng.guo@gmail.com
#PBS -m abe
#PBS -A k4zhang-group
cd /oasis/tscc/scratch/shg047/Xliu2014/fastq
bismark --bowtie2 --non_directional --phred33-quals --fastq -L 30 -N 1 -s 10 --multicore 6 /home/shg047/db/hg19/meth/bismark -1 ../fastq2_trim/SRR1232310_1_val_1.fq.gz -2 ../fastq2_trim/SRR1232310_2.fastq.gz_val_2.fq.gz  -o ../bam2
  • check the fastqc and find the first 10 base were not suitable to be aligned in the mapping.

Alignment

Pair End Mode

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)
SRR1232302_1_val_1.fq.gz_bismark_bt2_PE_report.txt 10896020 266017 2.40% 10986559 540111 45475 311108 558390 2521979 7009496 49.20% 1.80% 4.20%
SRR1232303_1_val_1.fq.gz_bismark_bt2_PE_report.txt 20216162 532701 2.60% 21017209 877814 62873 428697 973603 4695560 13978662 47.40% 1.30% 3.00%
SRR1232304_1_val_1.fq.gz_bismark_bt2_PE_report.txt 9233853 52146 0.60% 2154088 104373 8962 61461 113360 495917 1370015 47.90% 1.80% 4.30%
SRR1232305_1_val_1.fq.gz_bismark_bt2_PE_report.txt 19402159 359121 1.90% 14853592 750976 48615 303789 702872 3436588 9610752 51.70% 1.40% 3.10%
SRR1232306_1_val_1.fq.gz_bismark_bt2_PE_report.txt 51400600 1780106 3.50% 71928234 3884870 188607 772020 3360341 16663027 47059369 53.60% 1.10% 1.60%
SRR1232307_1_val_1.fq.gz_bismark_bt2_PE_report.txt 15015151 255663 1.70% 10806322 480441 44420 298374 509555 2351403 7122129 48.50% 1.90% 4.00%
SRR1232308_1_val_1.fq.gz_bismark_bt2_PE_report.txt 32508902 278849 0.90% 11147653 519877 35048 231897 507661 2493803 7359367 50.60% 1.40% 3.10%
SRR1232309_1_val_1.fq.gz_bismark_bt2_PE_report.txt 28578330 270922 0.90% 10742909 496990 27651 125158 386967 2329220 7376923 56.20% 1.20% 1.70%
SRR1232310_1_val_1.fq.gz_bismark_bt2_PE_report.txt 23879809 240025 1.00% 9580472 431160 36022 285577 493439 2136766 6197508 46.60% 1.70% 4.40%
SRR1232311_1_val_1.fq.gz_bismark_bt2_PE_report.txt 6078443 57482 0.90% 2159612 91805 6617 49301 91535 469450 1450904 50.10% 1.40% 3.30%
SRR1232312_1_val_1.fq.gz_bismark_bt2_PE_report.txt 9216625 194894 2.10% 8091546 392407 33197 226776 428295 1871966 5138905 47.80% 1.70% 4.20%

Bsmap Alignment

  • Bsmap can use 16 CPU, however, with only 8G memory. it is different with bismark. in bismark memory usage is proportional with CUP number with 3G.
#!/bin/csh
#PBS -n Xie.LUN.SRR1232302_1.fastq.gz
#PBS -q glean
#PBS -l nodes=1:ppn=16
#PBS -l walltime=72:00:00
#PBS -o SRR1232302_1.fastq.gz.trim.log
#PBS -e SRR1232302_1.fastq.gz.trim.err
#PBS -V
#PBS -M shihcheng.guo@gmail.com
#PBS -m abe
#PBS -A k4zhang-group
cd /oasis/tscc/scratch/shg047/Xliu2014/fastq
bsmap -u -s 12 -v 0.04 -p 16 -a ../fastq2_trim/SRR1232302_1_val_1.fq.gz -b ../fastq2_trim/SRR1232302_2_val_2.fq.gz -d /home/shg047/db/hg19/meth/bismark/hg19.fa  -o ../bsmap/SRR1232302_1.bam
  • two job failed since glean. and this time pdafm were used to re-alignment. Resubmit (Feb/18/2016)
cd /home/shg047/oasis/Xliu2014/fastq
qsub SRR1232310_1.fastq.gz.job
qsub SRR1232312_1.fastq.gz.job
  • two job completed (Feb/16/2016) and all the alignment were succeeded.the updated mapping statisitic
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)
SRR1232302_1_val_1.fq.gz_bismark_bt2_PE_report.txt 10879469 8166399 75.10% 283705221 14059544 396832 1083414 14827795 65987979 187349657 48.70% 0.60% 0.60%
SRR1232303_1_val_1.fq.gz_bismark_bt2_PE_report.txt 20185783 15643793 77.50% 538753835 24467542 731361 1999839 29198598 125132925 357223570 45.60% 0.60% 0.60%
SRR1232304_1_val_1.fq.gz_bismark_bt2_PE_report.txt 11064669 8668123 78.30% 301062107 14716791 435437 1166303 16357108 70648275 197738193 47.40% 0.60% 0.60%
SRR1232305_1_val_1.fq.gz_bismark_bt2_PE_report.txt 19379802 11358953 58.60% 405848758 21963190 606473 1550397 22103090 96959575 262666033 49.80% 0.60% 0.60%
SRR1232306_1_val_1.fq.gz_bismark_bt2_PE_report.txt 51331076 39232137 76.40% 1375770206 70954836 1959550 5129368 73767418 324505783 899453251 49.00% 0.60% 0.60%
SRR1232307_1_val_1.fq.gz_bismark_bt2_PE_report.txt 14987691 8977542 59.90% 312378988 14249586 467917 1334605 15164065 70650657 210512158 48.40% 0.70% 0.60%
SRR1232308_1_val_1.fq.gz_bismark_bt2_PE_report.txt 35378719 20706958 58.50% 723041259 35238860 1069487 2948419 35619931 165087550 483077012 49.70% 0.60% 0.60%
SRR1232309_1_val_1.fq.gz_bismark_bt2_PE_report.txt 34228833 20426178 59.70% 704014483 32162156 1007755 2905755 30199793 156286525 481452499 51.60% 0.60% 0.60%
SRR1232310_1_val_1.fq.gz_bismark_bt2_PE_report.txt 23850149 17995294 75.50% 563331726 27402528 831180 2213192 31796570 132428948 368659308 46.30% 0.60% 0.60%
SRR1232311_1_val_1.fq.gz_bismark_bt2_PE_report.txt 12142320 8808806 72.50% 272573717 13022456 405024 1118769 14270047 63239639 180517782 47.70% 0.60% 0.60%
SRR1232312_1_val_1.fq.gz_bismark_bt2_PE_report.txt 6133332 4094748 66.80% 142984922 6981841 202170 531307 8206242 33563875 93499487 46.00% 0.60% 0.60%
SRR1232313_1_val_1.fq.gz_bismark_bt2_PE_report.txt 32727127 23917117 73.10% 833607030 42703068 1162364 3069766 42946919 194253467 549471446 49.90% 0.60% 0.60%

unSortBam to Methylfreq (Pair-end)

When you save pair-end reads into one bam file, usually the read1 and read2 would be one by one in the order. However, when you use samtools sort, this order will be changed and therefore, some further analysis would be report error when they require pair-wise order reads, such as bismark_methylation_extractor

#!/bin/csh
#PBS -N SRR1232311_1_val_1.fq.gz_bismark_bt2_pe.bam
#PBS -q glean
#PBS -l nodes=1:ppn=6
#PBS -l walltime=6:00:00
#PBS -o SRR1232311_1_val_1.fq.gz_bismark_bt2_pe.bam.log
#PBS -e SRR1232311_1_val_1.fq.gz_bismark_bt2_pe.bam.err
#PBS -V
#PBS -M shihcheng.guo@gmail.com
#PBS -m abe
#PBS -A k4zhang-group
cd /oasis/tscc/scratch/shg047/Xliu2014/bam
bismark_methylation_extractor --paired-end --bedGraph --multicore 3 --ignore 3 --ignore_3prime 3 --ignore_r2 5 --ignore_3prime_r2 5 --gzip --buffer_size 4G --zero_based --comprehensive --output ../methyfreq SRR1232311_1_val_1.fq.gz_bismark_bt2_pe.bamrm CHG* CpG* CHH* *bedGraph.gz *bismark.cov.gz

Detail Parameters as Following:

Parameters used to extract methylation information:
Bismark Extractor Version: v0.14.4
Bismark result file: paired-end (SAM format)
Ignoring first 3 bp of Read 1
Ignoring first 5 bp of Read 2
Ignoring last 3 bp of Read 1
Ignoring last 5 bp of Read 2
Output specified: comprehensive
No overlapping methylation calls specified

Bam to SortBam (Feb 18 2016 Finished)

/home/shg047/oasis/Xliu2014/sortbam/SRR1232310.sort.bam
  • two bam file was not intact and therefore alignment were re-done to SRR1232310_1 and SRR1232312_1

SortBam to Hapinfo

MethylKit

MethylPurify