Dinh/Dinh 2014/NOTES/2014-11-17

From ZhangLabWiki
Revision as of 18:59, 24 November 2014 by >Dinh
Jump to navigation Jump to search

141112 HiSeqRapidRun

  • WGBS libraries were captured with biotin probes by Noi

Mapping results

Library ID Total sequenced bases Total trimmed bases Total mapped bases %trimmed %mapped
SeqCap-pCancerSetA-Nov4_12 1,557,796,200 1,490,698,729 1,183,722,483 4.31% 79.41%
SeqCap-pCancerSetA-Nov4_13 1,956,384,000 1,875,643,861 1,441,463,344 4.13% 76.85%
SeqCap-pCancerSetA-Nov4_14 1,604,724,200 1,533,952,230 1,236,147,785 4.41% 80.59%
SeqCap-pCancerSetA-Nov4_15 1,656,779,800 1,582,578,941 1,278,421,106 4.48% 80.78%
SeqCap-pCancerSetA-Nov4_16 1,551,344,600 1,482,228,220 1,172,067,971 4.46% 79.07%
SeqCap-pCancerSetA-Nov4_18 1,544,614,600 1,475,609,077 1,188,780,824 4.47% 80.56%
SeqCap-pCancerSetA-Nov4_19 1,937,986,600 1,868,857,625 1,466,874,474 3.57% 78.49%
SeqCap-pCancerSetA-Nov4_2 1,926,683,000 1,840,697,493 1,484,110,129 4.46% 80.63%
SeqCap-pCancerSetA-Nov4_4 2,527,810,600 2,414,697,940 1,897,794,294 4.47% 78.59%
SeqCap-pCancerSetA-Nov4_5 1,685,759,800 1,617,249,202 1,293,881,116 4.06% 80.01%
SeqCap-pCancerSetA-Nov4_6 1,974,733,600 1,891,689,068 1,507,964,542 4.21% 79.72%
SeqCap-pCancerSetA-Nov4_7 1,896,482,800 1,815,259,906 1,391,697,180 4.28% 76.67%
SeqCap-pNormalSetA-Nov4_12 2,074,019,600 1,975,408,680 1,571,662,041 4.75% 79.56%
SeqCap-pNormalSetA-Nov4_13 1,740,938,000 1,658,485,066 1,318,975,132 4.74% 79.53%
SeqCap-pNormalSetA-Nov4_14 1,867,843,400 1,781,077,564 1,442,362,287 4.65% 80.98%
SeqCap-pNormalSetA-Nov4_15 2,428,342,400 2,310,196,699 1,966,339,089 4.87% 85.12%
SeqCap-pNormalSetA-Nov4_16 1,922,011,000 1,827,491,168 1,448,437,222 4.92% 79.26%
SeqCap-pNormalSetA-Nov4_18 1,750,909,200 1,663,264,442 1,302,849,113 5.01% 78.33%
SeqCap-pNormalSetA-Nov4_19 1,714,868,600 1,634,401,979 1,295,380,595 4.69% 79.26%
SeqCap-pNormalSetA-Nov4_2 1,723,048,400 1,647,260,004 1,359,456,644 4.40% 82.53%
SeqCap-pNormalSetA-Nov4_4 2,111,874,400 2,013,677,234 1,563,440,813 4.65% 77.64%
SeqCap-pNormalSetA-Nov4_5 1,505,598,000 1,436,744,421 1,126,814,306 4.57% 78.43%
SeqCap-pNormalSetA-Nov4_6 1,731,991,400 1,653,646,386 1,314,038,490 4.52% 79.46%
SeqCap-pNormalSetA-Nov4_7 1,953,762,600 1,863,167,423 1,495,812,111 4.64% 80.28%

Scripts

  • Go.trimReads.sh
reads_dir="/oasis/tscc/scratch/ddiep/Working/141117_RRBS/Reads"
cur_dir=`pwd`
cope="/home/ddiep/softwares/cope-src-v1.1.3/src/cope -o connect.fq -2 left1.fq -3 left2.fq -m 0"
#THIS ONE for PE 
trim="/home/ddiep/softwares/trim_galore_latest/trim_galore --phred33 --paired --dont_gzip -a AGATCGGAAGAGCACACGTCTGAACTCCAGTCAC -a2 AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTAGATCTCGGTGGTCGCCGTATCATT"

#cd $reads_dir
#FILES=`ls *R1_001.fastq`
cd $cur_dir

for f in Sample_NP-SeqCap-pCancerSetA-Nov4_2 Sample_NP-SeqCap-pNormalSetA-Nov4_12  Sample_NP-SeqCap-pNormalSetA-Nov4_13 
do
	n=`echo $f | sed 's/Sample_NP-//g'`
	echo "#!/bin/csh" > $n.job
        echo "#PBS -l nodes=1:ppn=2" >> $n.job
        echo "#PBS -l walltime=3:00:00" >> $n.job
        echo "#PBS -o $n.log" >> $n.job
        echo "#PBS -e $n.err" >> $n.job
        echo "#PBS -V" >> $n.job
        echo "#PBS -M diep.hue.dinh@gmail.com" >> $n.job
        echo "#PBS -m abe" >> $n.job
        echo "#PBS -A k4zhang-group" >> $n.job
        echo "cd /state/partition1/\$USER/\$PBS_JOBID" >> $n.job
	echo "zcat $reads_dir/$f/*_R1_*gz > $n.R1.fastq" >> $n.job
	echo "zcat $reads_dir/$f/*_R2_*gz > $n.R2.fastq" >> $n.job
	echo "$trim $n.R1.fastq $n.R2.fastq" >> $n.job
	echo "cp *fq $cur_dir/" >> $n.job
	echo "cp *txt $cur_dir/" >> $n.job
	qsub -q hotel $n.job
done

  • Go.Map.sh
#===Change the following paths===#
cur_dir=`pwd`
scripts_dir="/oasis/tscc/scratch/ddiep/BisReadMapper/src"
reads_dir="$cur_dir/Trimmed_reads"

# reference files
ref_dir="/oasis/tscc/scratch/ddiep"
ref_fai="$ref_dir/bisHg19_plusLambda/hg19_lambda.fa.fai"
template_fwd="$ref_dir/bisHg19_plusLambda/hg19_lambda.bis.CT"
template_rev="$ref_dir/bisHg19_plusLambda/hg19_lambda.bis.GA"

# softwares
bwa="/home/ddiep/softwares/bwa-0.7.5a/bwa"

#=== List FASTQ to process ===#
cd $reads_dir
FILES=`ls *val_1.fq`
#FILES="s_1_1_ILMN_Indx01.cope.trimmed.fq"
cd $cur_dir

#===Begin===#
for f in ${FILES}
do
	n=`echo $f | sed 's/.R1_val_1.fq//g'`
	mkdir $cur_dir/$n
	cd $cur_dir/$n
	#1) Run mapper:
        echo "#!/bin/csh" > $n.job
        echo "#PBS -l nodes=1:ppn=4" >> $n.job
        echo "#PBS -l walltime=14:00:00" >> $n.job
        echo "#PBS -o $n.log" >> $n.job
        echo "#PBS -e $n.err" >> $n.job
        echo "#PBS -V" >> $n.job
        echo "#PBS -M diep.hue.dinh@gmail.com" >> $n.job
        echo "#PBS -m abe" >> $n.job
        echo "#PBS -A k4zhang-group" >> $n.job
        echo "cd /state/partition1/\$USER/\$PBS_JOBID" >> $n.job
	echo "cat $reads_dir/${n}.R1_val_1.fq $reads_dir/${n}.R2_val_2.fq > reads.fq" >> $n.job
        echo "$scripts_dir/BisReadMapper.pl -r reads.fq -W $template_fwd -C $template_rev -g $ref_fai -a $bwa -b 64 -p 4 -n $n > $n.status" >> $n.job
	echo "cp *sam $cur_dir/$n" >> $n.job
	echo "cp *status $cur_dir/$n" >> $n.job
	echo "ls *sorted.sam > $cur_dir/list_sams" >> $n.job
        echo "$scripts_dir/BamExtractor.pl -i $cur_dir/list_sams -s $cur_dir/list_paths_Hg19 -o $n -r none -v no -b no -p no -d 5" >> $n.job
	echo "cp *bam $cur_dir/$n" >> $n.job
        echo "cp *methylFreq $cur_dir/$n" >> $n.job
	qsub -q hotel $n.job
done
#===End===#

Assessment of capturing performance

Sample total_bases_coverred total_bases on_target_bases_covered on_target_bases Specificity target_coverred Enrichment_factor
SeqCap-pCancerSetA-Nov4_12 138,931,528 1,168,877,279 5,883,919 450,558,243 38.5% 20.0% 63
SeqCap-pCancerSetA-Nov4_13 369,134,874 1,424,314,151 21,209,194 614,835,055 43.2% 72.2% 77
SeqCap-pCancerSetA-Nov4_14 302,700,391 1,203,488,680 13,952,014 489,036,528 40.6% 47.5% 69
SeqCap-pCancerSetA-Nov4_15 277,562,818 1,250,329,814 12,577,032 511,943,961 40.9% 42.8% 70
SeqCap-pCancerSetA-Nov4_16 274,256,170 1,164,394,023 12,418,863 474,550,198 40.8% 42.3% 70
SeqCap-pCancerSetA-Nov4_18 403,128,602 1,170,296,929 23,073,184 461,126,140 39.4% 78.5% 66
SeqCap-pCancerSetA-Nov4_19 385,167,143 1,417,434,873 22,614,122 614,480,251 43.4% 77.0% 77
SeqCap-pCancerSetA-Nov4_2 354,784,810 1,514,708,612 19,490,111 751,025,646 49.6% 66.3% 99
SeqCap-pCancerSetA-Nov4_4 352,108,468 2,000,885,721 23,609,719 1,141,652,594 57.1% 80.3% 134
SeqCap-pCancerSetA-Nov4_5 196,108,966 1,300,234,129 8,246,274 563,004,318 43.3% 28.1% 77
SeqCap-pCancerSetA-Nov4_6 450,502,450 1,502,032,451 28,135,196 645,300,745 43.0% 95.7% 76
SeqCap-pCancerSetA-Nov4_7 226,755,475 1,432,933,860 10,016,683 629,536,922 43.9% 34.1% 79
SeqCap-pNormalSetA-Nov4_12 368,272,742 1,450,833,379 32,791,421 853,947,465 58.9% 111.6% 145
SeqCap-pNormalSetA-Nov4_13 330,611,335 1,291,478,405 29,764,564 756,289,960 58.6% 101.3% 143
SeqCap-pNormalSetA-Nov4_14 360,000,998 1,364,998,815 32,082,346 801,469,734 58.7% 109.2% 144
SeqCap-pNormalSetA-Nov4_15 375,193,438 1,970,841,697 30,241,737 1,295,408,500 65.7% 102.9% 194
SeqCap-pNormalSetA-Nov4_16 365,138,307 1,417,821,628 31,787,247 825,520,233 58.2% 108.2% 141
SeqCap-pNormalSetA-Nov4_18 326,952,891 1,305,603,890 28,937,598 760,701,025 58.3% 98.5% 141
SeqCap-pNormalSetA-Nov4_19 312,924,926 1,271,134,775 27,505,559 749,645,411 59.0% 93.6% 145
SeqCap-pNormalSetA-Nov4_2 353,552,883 1,275,196,941 31,331,422 731,603,678 57.4% 106.6% 136
SeqCap-pNormalSetA-Nov4_4 354,621,865 1,456,339,179 31,829,571 867,560,025 59.6% 108.3% 149
SeqCap-pNormalSetA-Nov4_5 291,891,632 1,066,253,299 27,061,797 618,610,388 58.0% 92.1% 140
SeqCap-pNormalSetA-Nov4_6 333,830,242 1,234,362,498 31,357,436 723,278,371 58.6% 106.7% 143
SeqCap-pNormalSetA-Nov4_7 377,475,888 1,451,628,039 32,325,421 843,343,133 58.1% 110.0% 140

Assessing the complexity of each library

  • Count the number of starting positions in each BAM file in non-overlapping 50 bp bins
  • Normalize each bin count by per million mapped reads
  • Plot the density for each file. Y-axis is the frequency of each count, X-axis is the log2(normalized counts).
  • Normal samples have greater complexity because most bins have fewer counts whereas cancer have more bins with higher counts (more reads covering same bins).
  • Cancer samples

File:141112 HiSeqRapidRun cancer samples 50bp coverage.png

  • Normal samples

File:141112 HiSeqRapidRun normal samples 50bp coverage.png