Dinh/Dinh 2013/NOTES/2013-6-24: Difference between revisions

From ZhangLabWiki
Jump to navigation Jump to search
>Dinh
mNo edit summary
>Dinh
Line 281: Line 281:
  <nowiki>
  <nowiki>


23892 mandatory_13:chr15:100249085-100249289
23892 mandatory_13:chr15:100249085-100249289 +
12876 mandatory_10:chr7:3025508-3025737
12876 mandatory_10:chr7:3025508-3025737 -
5540 mandatory_3:chr4:7526550-7526760
5540 mandatory_3:chr4:7526550-7526760 -
4183 mandatory_15:chr17:75369137-75369338
4183 mandatory_15:chr17:75369137-75369338 -
3944 mandatory_1:chr4:154710370-154710640
3944 mandatory_1:chr4:154710370-154710640 +
2965 mandatory_8:chr7:140218006-140218261
2965 mandatory_8:chr7:140218006-140218261 -
2579 mandatory_5:chr2:9518172-9518402
2579 mandatory_5:chr2:9518172-9518402 +
1839 mandatory_6:chr17:80709116-80709327
1839 mandatory_6:chr17:80709116-80709327 -
1358 mandatory_9:chr7:26206450-26206690
1358 mandatory_9:chr7:26206450-26206690 -
805 mandatory_14:chr4:147557728-147557999
805 mandatory_14:chr4:147557728-147557999 -
633 mandatory_9:chr7:26206446-26206652
633 mandatory_9:chr7:26206446-26206652 +
561 mandatory_12:chr2:42275654-42275925
561 mandatory_12:chr2:42275654-42275925 +
537 mandatory_7:chr3:142837896-142838117
537 mandatory_7:chr3:142837896-142838117 -
484 mandatory_7:chr3:142837897-142838117
484 mandatory_7:chr3:142837897-142838117 +
331 mandatory_12:chr2:42275614-42275844
331 mandatory_12:chr2:42275614-42275844 -
231 mandatory_13:chr15:100249082-100249302
231 mandatory_13:chr15:100249082-100249302 -
230 mandatory_3:chr4:7526557-7526773
230 mandatory_3:chr4:7526557-7526773 +
183 mandatory_6:chr17:80709042-80709272
183 mandatory_6:chr17:80709042-80709272 +
67 mandatory_2:chr1:110052333-110052595
67 mandatory_2:chr1:110052333-110052595 -
32 mandatory_4:chr2:164593111-164593326
32 mandatory_4:chr2:164593111-164593326 +
17 mandatory_4:chr2:164593173-164593378
17 mandatory_4:chr2:164593173-164593378 -
16 mandatory_14:chr4:147557730-147558008
16 mandatory_14:chr4:147557730-147558008 +
9 mandatory_8:chr7:140218000-140218270
9 mandatory_8:chr7:140218000-140218270 +
5 mandatory_2:chr1:110052282-110052542
5 mandatory_2:chr1:110052282-110052542 +
5 mandatory_1:chr4:154710373-154710633
5 mandatory_1:chr4:154710373-154710633 -
5 mandatory_10:chr7:3025518-3025736
5 mandatory_10:chr7:3025518-3025736 +
0 mandatory_15:chr17:75369134-75369334
0 mandatory_15:chr17:75369134-75369334 +
0 mandatory_5:chr2:9518220-9518443
0 mandatory_5:chr2:9518220-9518443 -
0 mandatory_11:chr7:138229771-138230020
0 mandatory_11:chr7:138229771-138230020 +
0 mandatory_11:chr7:138229843-138230113
0 mandatory_11:chr7:138229843-138230113 -
0 mandatory_16:chr7:93520143-93520413
0 mandatory_16:chr7:93520143-93520413 -
</nowiki>
===Recommended===
<nowiki>
65845 recommended_33:chr8:3316771-3316991 +
19194 recommended_34:chr1:25257453-25257723 -
12188 recommended_14:chr17:43339401-43339611 -
12048 recommended_11:chr21:47783987-47784263 +
11668 recommended_24:chr1:115124392-115124598 -
11231 recommended_31:chrX:135333528-135333788 -
11132 recommended_5:chr7:27154944-27155156 -
10747 recommended_19:chr16:1017670-1017944 -
9879 recommended_6:chr10:7451206-7451436 -
7055 recommended_32:chr17:1633600-1633820 -
6674 recommended_10:chr1:161442627-161442877 +
5537 recommended_28:chr17:43044892-43045132 +
4945 recommended_15:chr16:85676275-85676545 +
4350 recommended_23:chr6:32042899-32043159 -
4311 recommended_3:chr19:38085548-38085818 -
2204 recommended_9:chr13:103052745-103053015 -
2117 recommended_13:chr20:36013306-36013546 +
2010 recommended_32:chr17:1633600-1633820 +
1950 recommended_22:chr7:71682175-71682435 -
1822 recommended_1:chr5:54281139-54281408 -
1516 recommended_19:chr16:1017669-1017943 +
1285 recommended_7:chr17:7165161-7165397 -
1185 recommended_13:chr20:36013306-36013536 -
1163 recommended_18:chr11:69197082-69197351 -
981 recommended_11:chr21:47783999-47784254 -
934 recommended_33:chr8:3316763-3316986 -
900 recommended_26:chr1:3567532-3567752 -
773 recommended_28:chr17:43044890-43045150 -
720 recommended_17:chr17:76921742-76921982 +
714 recommended_26:chr1:3567513-3567773 +
676 recommended_10:chr1:161442625-161442877 -
627 recommended_25:chr12:130589110-130589320 -
504 recommended_24:chr1:115124345-115124605 +
448 recommended_18:chr11:69197079-69197355 +
410 recommended_7:chr17:7165139-7165351 +
354 recommended_3:chr19:38085549-38085819 +
297 recommended_2:chr20:20349143-20349401 -
241 recommended_24:chr1:115124345-115124605 -
206 recommended_17:chr17:76921738-76921961 -
123 recommended_15:chr16:85676279-85676544 -
119 recommended_12:chr2:8597284-8597533 +
94 recommended_6:chr10:7451190-7451465 +
91 recommended_8:chr2:9614480-9614708 -
88 recommended_1:chr5:54281139-54281409 +
58 recommended_34:chr1:25257451-25257727 +
53 recommended_31:chrX:135333568-135333778 +
43 recommended_16:chr4:118975389-118975609 +
39 recommended_12:chr2:8597159-8597438 -
16 recommended_27:chr16:1243456-1243666 +
8 recommended_9:chr13:103052712-103052935 +
7 recommended_20:chr8:74878400-74878650 -
4 recommended_27:chr16:1243426-1243663 -
2 recommended_23:chr6:32042893-32043163 +
0 recommended_21:chr10:125853024-125853277 -
0 recommended_25:chr12:130588995-130589226 +
0 recommended_32:chr17:1633599-1633825 +
0 recommended_14:chr17:43339405-43339665 +
0 recommended_2:chr20:20349125-20349399 +
0 recommended_8:chr2:9614470-9614722 +
0 recommended_16:chr4:118975389-118975619 -
0 recommended_5:chr7:27154941-27155143 +
0 recommended_22:chr7:71682176-71682426 +
0 recommended_20:chr8:74878348-74878592 +
</nowiki>
</nowiki>

Revision as of 19:51, 26 June 2013

Analysis of Blueprint MiSeq Test Run

Pre-processing

Standard trimming + UMI

  • Obtain UMI from the first 10bp of read 1, label both reads with UMI
  • Trim 27 bp from 5 prime end of read 1 and read 2

Adapter removal with fastq-mcf

  • Remove adapter sequences using fastq-mcf ?
  • Since the insert = (target_size= 200-280) + (2 arms ~ 64 bp) + (UMI = 10 bp) = 274-354bp, we would not have sequenced the adapters with just 250bp from each end.
  • Try removing adapters because mapping rate was ~40%.
  • Adapters:
>Linker
TTGGAGGCTCATCGTTCCTATTCAGGCAGATGTTATCGAGGTCCGAC
>Linker_rev
GTCGGACCTCGATAACATCTGCCTGAATAGGAACGATGAGCCTCCAA

fastq-mcf results

  • Note that only R1 was end-trimmed using "Linker" and only R2 was end-trimmed using "Linker_rev"
index Total reads Clipped 'end' reads % clipped Too short after clip % too short
1, R1 546,091 265,765 48.67% 16,468 3.02%
1, R2 546,091 330,822 60.58% 16,112 2.95%
2, R1 1,201,482 599,205 49.87% 46,756 3.89%
2, R2 1,201,482 733,958 61.09% 45,880 3.82%
3, R1 2,107,061 1,159,863 55.05% 30,221 1.43%
3, R2 2,107,061 1,348,459 64.00% 29,076 1.38%
4, R1 1,576,107 821,955 52.15% 16,478 1.05%
4, R2 1,576,107 1,000,442 63.48% 16,703 1.06%
5, R1 2,416,376 1,172,370 48.52% 30,785 1.27%
5, R2 2,416,376 1,470,535 60.86% 29,880 1.24%
6, R1 2,025,366 1,044,010 51.55% 17,056 0.84%
6, R2 2,025,366 1,266,104 62.51% 16,306 0.81%
7, R1 1,632,245 737,076 45.16% 12,924 0.79%
7, R2 1,632,245 398,506 24.41% 12,788 0.78%
8, R1 2,699,273 1,252,029 46.38% 24,410 0.90%
8, R2 2,699,273 668,800 24.78% 22,635 0.84%

Merge R1 and R2 with COPE

  • Create kmer_table for COPE
kmerfreq -k 17 -t 4 -c -1 -p kmer_table read.lst >kmerfreq.cout 2>kmerfreq.cerr
  • Use COPE to combine overlapping read 1 and read 2
for f in 1Index1_S1 1Index2_S2 1Index3_S3 2Index4_S6 1Index5_S4 2Index6_S7 2Index7_S8 1Index8_S5
do

	./getUMI.pl 130620_MiSeq/TES1-${f}_L001_R1_001.fastq 130620_MiSeq/TES1-${f}_L001_R2_001.fastq $f
	~/softwares/cope-src-v1.1.3/src/cope -a $f.R1.fq -b $f.R2.fq -o $f.fq -2 $f.leftR1.fq -3 $f.leftR2.fq -m 1  -t kmer_table.freq.cz -f kmer_table.freq.cz.len >cope.$f.log 2>cope.$f.error
	rm $f.R1.fq $f.R2.fq
done

COPE results

  • Without adapter removal:
index total_pairs connected_pairs connect_ratio(%) low_quality_pairs low_quality_ratio(%)
1 546,091 37,624 6.88969 239413 43.8412
2 1,201,482 76,332 6.35315 603335 50.2159
3 2,107,061 120,889 5.73733 1155859 54.8565
5 2,416,376 157,448 6.51587 1124046 46.5178
8 2,699,273 191,332 7.08828 1138412 42.1748
4 1,576,107 100,314 6.36467 805229 51.0897
6 2,025,366 131,396 6.48752 990421 48.9008
7 546,091 35,800 6.55568 239459 43.8497
  • With adapter (linker sequence) removal:
index total_pairs connected_pairs connect_ratio(%) low_quality_pairs low_quality_ratio(%)
1 529623 20770 3.92166 195892 36.9871
2 1154726 42045 3.64112 477352 41.339
3 2076840 65820 3.16924 970053 46.7081
5 2385591 96387 4.04038 964834 40.4442
8 2674863 126170 4.71688 1012342 37.8465
4 1559629 55184 3.53828 684665 43.8992
6 2008310 76903 3.82924 820560 40.8582
7 1619321 80497 4.97103 584325 36.0846

Alignment with Bowtie2

  • The long reads were not compatible with our pipeline using Bowtie2, because the aligner suppresses read name lines which are longer than 256 characters. This causes truncated original reads stored in the read name line.
  • Added code to split long reads into multiple short reads to use Bowtie2
  • File:SmartBisReadMapper.txt
  • Mapping pipeline:
cur_dir="/media/3TB_Dinh/Blueprint"
reads_dir="/media/3TB_Dinh/Blueprint"
email="diep.hue.dinh@gmail.com"

bisReadMapper="/home/ddiep/scripts/MethylationPipeline/scripts/smartBisReadMapper.pl"
template_fwd="/media/2TB_storeA/BisRef/bisHg19/hg19.fa.bis.fwd.bowtie2"
template_rev="/media/2TB_storeA/BisRef/bisHg19/hg19.fa.bis.rev.bowtie2"
template_fa="/media/2TB_storeA/BisRef/bisHg19/hg19.fa"
soap="/home/ddiep/softwares/soap2.21release/soap"
bowtie="bowtie2"
cpg="/media/2TB_storeA/BisRef/bisHg19/C_Pos/hg19.fa.cpg.positions.txt"

INDX="1Index1_S1 1Index2_S2 1Index3_S3 2Index4_S6 1Index5_S4 2Index6_S7 2Index7_S8 1Index8_S5"

# make sure the qual_base variable is set correctly
cd $cur_dir

for s in ${INDX}
do

f="$s.fq"
g="$s.leftR1.fq"
h="$s.leftR2.fq"
n="$s-MS"
mkdir $n
echo "cd $cur_dir/$n" > $n.job
echo "$bisReadMapper -r $reads_dir/$f -W $template_fwd -C $template_rev -g $template_fa -a $bowtie -p 16 -b 33 -n $n.f1 -q 20 -l $cpg > $n.f1.statusMbias 2>$n.f1.err" >> $n.job
echo "$bisReadMapper -r $reads_dir/$g -W $template_fwd -C $template_rev -g $template_fa -a $bowtie -p 16 -b 33 -n $n.f2 -q 20 -l $cpg > $n.f2.statusMbias 2>$n.f2.err" >> $n.job
echo "$bisReadMapper -r $reads_dir/$h -W $template_fwd -C $template_rev -g $template_fa -a $bowtie -p 16 -b 33 -n $n.f3 -q 20 -l $cpg > $n.f3.statusMbias 2>$n.f3.err" >> $n.job
echo "rm *encoded" >> $n.job

done

Mapping statistics

  • f1 = COPE combined PE reads, f2 = leftover R1, f3 = leftover R2.
  • With clipping, we got 2.846 Gbps mapped, without clipping, we got 2.842 Gbps mapped. Clipping did not improve the amount of usable bp by a lot.
  • Libraries Index 7 and Index 8 have the best mapping rates.
index read file total clipped bases total clipped mapped bases %mapped(clipped) total bases total mapped bases %mapped
1 f1 8255745 5487013 66.46% 16078514 9619049 59.83%
1 f2 79964964 53864094 67.36% 104318024 50945748 48.84%
1 f3 89749717 50800000 56.60% 99998101 48797421 48.80%
2 f1 16566696 11031637 66.59% 32606092 19771959 60.64%
2 f2 171105395 112357034 65.67% 230586489 106229965 46.07%
2 f3 189062773 102110571 54.01% 215612917 97902718 45.41%
3 f1 25104735 15365281 61.20% 51734457 30631757 59.21%
3 f2 285378749 171910645 60.24% 405686018 161655992 39.85%
3 f3 330732139 151852920 45.91% 372674273 145065779 38.93%
4 f1 21403635 13378820 62.51% 42833085 24714027 57.70%
4 f2 222330667 138200084 62.16% 300568658 131048023 43.60%
4 f3 250430069 125847956 50.25% 281078201 120620387 42.91%
5 f1 38452894 26019949 67.67% 67355222 42886852 63.67%
5 f2 360645974 247314571 68.58% 463290839 235545215 50.84%
5 f3 401191540 226901252 56.56% 438788886 218980079 49.91%
6 f1 30067978 18916755 62.91% 56145035 34232803 60.97%
6 f2 289462431 191584810 66.19% 389224257 183501240 47.15%
6 f3 328112598 176045539 53.65% 364481484 169262473 46.44%
7 f1 32538699 22178956 68.16% 49348430 34280682 69.47%
7 f2 253380410 192127161 75.83% 313644714 186792389 59.56%
7 f3 279659004 180487623 64.54% 300856044 175441163 58.31%
8 f1 50524286 33292157 65.89% 81821672 53601105 65.51%
8 f2 415435035 301628773 72.61% 516798645 291900532 56.48%
8 f3 463261828 277764207 59.96% 496294120 269456511 54.29%

Remove clonal reads

  • Modified Prof. Zhang's code to remove clonal reads using UMI.

Final Blueprint capture experiment results

  • The capture setup for Index 6 provided the best results.
  • Mapping was performed without adapter end-trimming.
  • capture specificity = # on target bp / total usable bp
  • capture sensitivity = # target bp at 1x / target size
  • capture enrichment = # on target bp / target size
index total PE reads total bps total bps mapped % bps mapped total bps after clonal removal % bps clonal genome bp average genome depth of coverage total CpG depth of coverage number of CpGs (1x) average CpGs depth of coverage capture specificity capture sensitivity capture enrichment f/r correlation f/r correlation (w/o clonal removal) #mandatory CpG (max 16) #recommended CpG (max 32) #optional CpG (max 1024)
1 546,091 273,045,500 109,362,218 40.05% 13,695,568 87.48% 439,410 31 736,350 16,720 44 82.77% 55.12% 35.07 0.8204 0.7376 8 22 675
2 1,201,482 600,741,000 223,904,642 37.27% 31,069,492 86.12% 700,376 44 1,638,154 23,633 69 83.07% 64.59% 79.85 0.7707 0.7222 10 23 786
3 2,107,061 1,053,530,500 337,353,528 32.02% 59,858,799 82.26% 1,448,705 41 3,082,267 42,205 73 80.03% 73.22% 148.22 0.8222 0.7619 10 27 869
4 1,576,107 788,053,500 276,382,437 35.07% 86,119,440 68.84% 2,394,363 36 4,322,851 67,792 64 82.13% 78.36% 218.84 0.8762 0.8343 12 29 920
5 2,416,376 1,208,188,000 497,412,146 41.17% 88,991,063 82.11% 1,618,255 55 4,617,178 47,342 98 84.13% 75.74% 231.64 0.8639 0.8170 11 27 898
6 2,025,366 1,012,683,000 386,996,516 38.21% 127,581,738 67.03% 3,036,498 42 6,328,115 81,327 78 84.30% 81.11% 332.75 0.9089 0.8677 12 30 946
7 1,632,245 816,122,500 396,514,234 48.59% 116,341,160 70.66% 1,499,425 78 6,022,204 44,654 135 90.63% 78.05% 326.23 0.9029 0.8467 12 28 919
8 2,699,273 1,349,636,500 614,958,148 45.56% 79,887,075 87.01% 1,314,474 61 4,271,697 38,539 111 87.39% 72.11% 216.01 0.8372 0.7962 11 27 852

Probe Bias

Mandatory sites


23892	mandatory_13:chr15:100249085-100249289	+
12876	mandatory_10:chr7:3025508-3025737	-
5540	mandatory_3:chr4:7526550-7526760	-
4183	mandatory_15:chr17:75369137-75369338	-
3944	mandatory_1:chr4:154710370-154710640	+
2965	mandatory_8:chr7:140218006-140218261	-
2579	mandatory_5:chr2:9518172-9518402	+
1839	mandatory_6:chr17:80709116-80709327	-
1358	mandatory_9:chr7:26206450-26206690	-
805	mandatory_14:chr4:147557728-147557999	-
633	mandatory_9:chr7:26206446-26206652	+
561	mandatory_12:chr2:42275654-42275925	+
537	mandatory_7:chr3:142837896-142838117	-
484	mandatory_7:chr3:142837897-142838117	+
331	mandatory_12:chr2:42275614-42275844	-
231	mandatory_13:chr15:100249082-100249302	-
230	mandatory_3:chr4:7526557-7526773	+
183	mandatory_6:chr17:80709042-80709272	+
67	mandatory_2:chr1:110052333-110052595	-
32	mandatory_4:chr2:164593111-164593326	+
17	mandatory_4:chr2:164593173-164593378	-
16	mandatory_14:chr4:147557730-147558008	+
9	mandatory_8:chr7:140218000-140218270	+
5	mandatory_2:chr1:110052282-110052542	+
5	mandatory_1:chr4:154710373-154710633	-
5	mandatory_10:chr7:3025518-3025736	+
0	mandatory_15:chr17:75369134-75369334	+
0	mandatory_5:chr2:9518220-9518443	-
0	mandatory_11:chr7:138229771-138230020	+
0	mandatory_11:chr7:138229843-138230113	-
0	mandatory_16:chr7:93520143-93520413	-

Recommended

65845	recommended_33:chr8:3316771-3316991	+
19194	recommended_34:chr1:25257453-25257723	-
12188	recommended_14:chr17:43339401-43339611	-
12048	recommended_11:chr21:47783987-47784263	+
11668	recommended_24:chr1:115124392-115124598	-
11231	recommended_31:chrX:135333528-135333788	-
11132	recommended_5:chr7:27154944-27155156	-
10747	recommended_19:chr16:1017670-1017944	-
9879	recommended_6:chr10:7451206-7451436	-
7055	recommended_32:chr17:1633600-1633820	-
6674	recommended_10:chr1:161442627-161442877	+
5537	recommended_28:chr17:43044892-43045132	+
4945	recommended_15:chr16:85676275-85676545	+
4350	recommended_23:chr6:32042899-32043159	-
4311	recommended_3:chr19:38085548-38085818	-
2204	recommended_9:chr13:103052745-103053015	-
2117	recommended_13:chr20:36013306-36013546	+
2010	recommended_32:chr17:1633600-1633820	+
1950	recommended_22:chr7:71682175-71682435	-
1822	recommended_1:chr5:54281139-54281408	-
1516	recommended_19:chr16:1017669-1017943	+
1285	recommended_7:chr17:7165161-7165397	-
1185	recommended_13:chr20:36013306-36013536	-
1163	recommended_18:chr11:69197082-69197351	-
981	recommended_11:chr21:47783999-47784254	-
934	recommended_33:chr8:3316763-3316986	-
900	recommended_26:chr1:3567532-3567752	-
773	recommended_28:chr17:43044890-43045150	-
720	recommended_17:chr17:76921742-76921982	+
714	recommended_26:chr1:3567513-3567773	+
676	recommended_10:chr1:161442625-161442877	-
627	recommended_25:chr12:130589110-130589320	-
504	recommended_24:chr1:115124345-115124605	+
448	recommended_18:chr11:69197079-69197355	+
410	recommended_7:chr17:7165139-7165351	+
354	recommended_3:chr19:38085549-38085819	+
297	recommended_2:chr20:20349143-20349401	-
241	recommended_24:chr1:115124345-115124605	-
206	recommended_17:chr17:76921738-76921961	-
123	recommended_15:chr16:85676279-85676544	-
119	recommended_12:chr2:8597284-8597533	+
94	recommended_6:chr10:7451190-7451465	+
91	recommended_8:chr2:9614480-9614708	-
88	recommended_1:chr5:54281139-54281409	+
58	recommended_34:chr1:25257451-25257727	+
53	recommended_31:chrX:135333568-135333778	+
43	recommended_16:chr4:118975389-118975609	+
39	recommended_12:chr2:8597159-8597438	-
16	recommended_27:chr16:1243456-1243666	+
8	recommended_9:chr13:103052712-103052935	+
7	recommended_20:chr8:74878400-74878650	-
4	recommended_27:chr16:1243426-1243663	-
2	recommended_23:chr6:32042893-32043163	+
0	recommended_21:chr10:125853024-125853277	-
0	recommended_25:chr12:130588995-130589226	+
0	recommended_32:chr17:1633599-1633825	+
0	recommended_14:chr17:43339405-43339665	+
0	recommended_2:chr20:20349125-20349399	+
0	recommended_8:chr2:9614470-9614722	+
0	recommended_16:chr4:118975389-118975619	-
0	recommended_5:chr7:27154941-27155143	+
0	recommended_22:chr7:71682176-71682426	+
0	recommended_20:chr8:74878348-74878592	+