Dinh/Dinh 2013/NOTES/2013-10-7

From ZhangLabWiki
Revision as of 00:50, 9 October 2013 by >Dinh
Jump to navigation Jump to search

Blueprint Capture, low mapping issue

  • Low mapping can be caused by any of the following problems
  1. Contamination - contamination of library after capture would not be amplified by the last PCR step or sequenced due to lack of adapter sequences. Contamination of reagents during capture may lead to low mapping rates. Contamination may explain why we have such variable mapping rates, since each DNA sample had varying quality, it's possible that the lower quality samples would have more capture from contaminating sequences and thus lower mapping rate.
    1. Do the amount of input converted gDNA correlate with mapping rate?
    2. Do the captured sequences map to other genomes in NT database?
  2. Chimeric capture sequences - When imperfect hybridization of capture arms to genomic sequences which occurs frequently in the genome, leading to off-target & chimeric sequences. Or when the capture arms capture each other.
    1. How many of the sequences map to the adapters, linkers, and primers?
    2. Do the sequences map to human repeats? - No, when mapping with Bowtie, most of the reads either map uniquely or not at all.
    3. Do the sequences map to non-human repeats? - Check by performing Blast to RepBase18.08
  3. Error from sequencing - We did not use a PhiX control for the MiSeq runs, but had a control for the HiSeq RapidRun (control was ~ 89% of flowcell, BP libraries was ~ 11%). The mapping rate of the HiSeq RapidRun was still low for many samples. Also, the majority of mapped reads have low error rates at the 3' end.
  4. Error from oligo synthesis - error from oligo synthesis would still allow us to recover the alignments from the off-target regions. Even so, we should get lower the specificity of capture, but this is not the case.

Reads sampling and trimming

  • Randomly sample 1000 reads from the following libraries for analysis:
 1. Blueprint samples #1 and #2 from the Aug 28 MiSeq run. (compare different BP samples) - BP1, BP2
 2. Technical replicates for Blueprint samples #1 and #2 from the Aug 28 MiSeq run. (compare different captures) - BP1rep, BP2rep
 3. Technical replicates for Blueprint samples #1 and #2 from the Sep 13 HiSeq Rapid run. (compare different runs) - BP1rep-hiseq, BP2rep-hiseq
 4. NA12878 captured with the original probes set and sequenced in HL157 (+ control) - HMrep1
 5. NA12878 captured with the normalized probes set and sequenced in HL157 (+ control) - HMrep2
 6. NA12878 captured with the normalized probes set and sequenced in HiSeq Rapid run - HMrep3
 7. DMR330K capture from the Mar22 MiSeq run as positive control. (+ control) - H1-dmr330k
  • All reads were sampled from the first read from each run.
  • All reads were trimmed 27 bp from the 5' end and quality trimmed using mean quality score of 30.

BLAST sequences to adapters, linkers, and primers

  • Since these were not bisulfite converted, I use the original reads to blast.
  • Library free adapters/linkers/primers
>AmpF6.3Sol
AATGATACGGCGACCACCGACACTCTCAGATGTTATCGAGGTCCGAC
>AmpR6.3Ind1
CAAGCAGAAGACGGCATACGAGATCGTGATGCTAGGAACGATGAGCCTCCAAC
>SolSeqV6.3.2
ACCGACACTCTCAGATGTTATCGAGGTCCGAC
>SolSeqV6.3.2r
GCTAGGAACGATGAGCCTCCAAC
>AmpR6.3IndSeq
GTTGGAGGCTCATCGTTCCTAGC
>AmpF6.4Sol
AATGATACGGCGACCACCGAGATCTACACCACTCTCAGATGTTATCGAGGTCCGAC
>SolSeqV6.3.3
TACACCACTCTCAGATGTTATCGAGGTCCGAC
  • BLAST args:
 $blastn -p blastn -d $db_file -K 100 -e 0.0001 -g F -F "m D" -a 4 -W 7 -m 8 -i $f -o $f.primers.blast
  • Number of reads mapped to library_free_primers:
 (BP - Blueprint sample, HM - HapMap sample NA12878, H1 - Human ESC H1)
 590 BP1
 667 BP2
 458 BP1rep
 812 BP2rep
 222 BP1rep-hiseq
 558 BP2rep-hiseq
 12  HMrep1 (but this library has ~18.70% mapping rate??)
 5   HMrep2 (but this library has ~15.07% mapping rate??)
 505 HMrep3-hiseq
 23  H1-dmr330K
  • Number of hits to each primer sequences
AmpF6.3Sol AmpR6.3Ind1 SolSeqV6.3.2 SolSeqV6.3.2r AmpR6.3IndSeq AmpF6.4Sol SolSeqV6.3.3
BP1 510 555 590 555 555 510 510
BP2 581 639 667 639 639 581 581
BP1rep 363 444 458 444 444 363 363
BP2rep 653 786 811 785 785 653 653
BP1rep-hiseq 8 220 222 220 220 8 8
BP2rep-hiseq 12 556 558 556 556 12 12
H1-dmr330K 12 23 23 23 23 12 12
HMrep1 0 5 5 5 5 0 0
HMrep2 1 12 12 12 12 1 1
HMrep3-hiseq 10 502 505 502 502 10 10

BLAST sequences to NT database

  • If capture reagents were contaminated, we would expect that reads from capture of the converted sample would not map to the unconverted NT database. Reads that map to the unconverted NT database would means that we have contaminated capture reagents.
  • Waiting for mapping results....

BLAST sequences to RepBase

  • If there were off-target capture to repeats, we would get hits to the converted RepBase sequences, but they would be mostly human.
  • Reads were encoded by converting G to A (all read 1 would be on the complement strand)
  • BLAST args:
 $blastn -p blastn -d $db_dir/all_repeats.bis.fwd -K 100 -e 0.0001 -g F -F "m D" -a 4 -W 7 -i $f.encoded_2 -o $f.bis.fwd.all_repeats.blast
' Unconverted repeats Unconverted repeats Unconverted repeats Bis-forward strand Bis-forward strand Bis-forward strand
Unmapped Human repeat hits Non-human repeat hits Unmapped reads Human repeat hits Non-human repeat hits
BP1 977 13 81 396 141 7589
BP2 970 35 249 438 251 7699
BP1rep 930 69 519 455 493 8842
BP2rep 967 19 205 441 319 6028
H1-dmr330K* 980 1 99 549 72 3748
Hmrep2 990 22 45 702 143 1870
Hmrep1* 991 2 40 750 45 2867
BP1rep-hiseq 972 81 218 640 631 6668
BP2rep-hiseq 976 21 86 660 270 3021
Hmrep3 976 47 96 652 385 4286