Sam:LabNotes/Microbiome-new/2011-2-15: Difference between revisions

From ZhangLabWiki
Jump to navigation Jump to search
>Sam Chiang
No edit summary
>Sam Chiang
 
(24 intermediate revisions by the same user not shown)
Line 1: Line 1:
='''Bowtie mapping(whole reads) and contamination analysis(failed to aligned reads)- HC-MDA-Ecoli-SAGs(NX1-NX8)-Jan17'''=
='''Bowtie mapping on E.coli SAG Nextera libraries HC-MDA-Ecoli-SAGs(NX1-NX8)-Jan09'''=


==Background==
==Background==
*Samples were prepared from -(short and long amplified E.coli SAGs -3rd test 11-07-'10)
*Ecoli SAGs were prepared from flow-sorted single E.coli cells.
*Eight E.coli SAGs samples were made to Illumina sequencing libraries using (N2)Primer-based barcodes ID1 to ID8.
*Eight E.coli SAGs samples were made to Nextera-libraries.
**http://genome-tech.ucsd.edu/LabNotes/index.php/Sam:LabNotes/Microbiome-new/2010-12-29
**http://genome-tech.ucsd.edu/LabNotes/index.php/Sam:LabNotes/Microbiome-new/2011-01-06_exp2
 
*These libraries were equally pulled (20ng per library) and sequenced in pair-end(PE) mode.
*These libraries were equally pulled (20ng per library) and sequenced in pair-end(PE) mode.
*The sequenced data processed by Dr. Zhang for base-calling and saved into FASTQ files.
*The sequenced data processed by Dr. Zhang for base-calling and saved into FASTQ files.
*Since I don't have genome-tech server permission, I asked Dinh's help to transfer data into Genemapster server.
 
Wiki(12-29-'10)
Sample no.  Sample name          Total DNA(ng)
-----------------------------------------------
NX01        WGA102710-(1-4)- 3hr    288.50ng
NX02        WGA102710-(1-4)-10hr    500.00ng
NX03        WGA110310-(5)-  3hr    114.50ng
NX04        WGA110310-(5)-  10hr    500.00ng
 
Wiki(01-06-'11)
Sample no.  Sample name          Total DNA(ng)
-----------------------------------------------
NX05        WGA122010-(1-2)-3hr      181.50ng
NX06        WGA122010-(1-2)-10hr      200.00ng
NX07        WGA122110-(2-3)-3hr        50.00ng
NX08        WGA122110-(2-3)-10hr      200.00ng


==Analysis approaches==
==Analysis approaches==
===Single-end data analysis===
===Data pooling===
*Whole reads alignment using Bowtie => genome coverage calculate
*Using "gunzip" command to unzip SE files.
**Script: Read2coverage-kz-3
*Combine two SE data(read1 + read2) into one PE data
**Data resource: HC-MDA-Ecoli-SAGs-Ind1to8
e.g. samchiang@genemapster:~/Raw-Read-backup/110120_HL087$ cat s_4_1_ID8.txt s_4_2_ID8.txt > s_4_pe_ID8.txt &
**Template:E.coli-k12-mg1655 template
**Bowtie setting: -k 1 -l 28 -n 2 -m 1 --best --best --phred64-quals
**Downstream: Data was used for genome coverage calculation.
[[File:sam120610-wholereadbotie-2.jpg|800px]]


*Using the same alignment processing, I extract the failed-to-aligned and non-specific reads into different folder.
===Genome coverage analysis (whole reads and 1000000 down-sampled reads)===
**Script: Read2coverage-kz4
*Script: Read2coverage-kz-3
**Data resource: HC-MDA-Ecoli-SAGs-Ind1to8
*Data resource: pooled pe reads data from HC-MDA-Ecoli-SAGs(NX1-NX8)-Jan09
**Template:E.coli-k12-mg1655 template
*Command
**Bowtie setting: -k 1 -l 28 -n 2 -m 1 --best --best --phred64-quals --un [failed-to-aligned file output] --max [nonspecific file output]
perl Read2coverage-kz-3 /media/disk-2/samchiang/Raw-Read-backup/110120_HL087/pooled_pe/s_4_pe_ID8_3074085.txt /media/disk-2/samchiang/Tools/bowtie-0.12.5/indexes/e_coli_k12_mg1655/e_coli_k12_mg1655 &
*Bowtie setting
samchiang@genemapster:~/Scripts$ /media/disk-2/samchiang/Tools/bowtie-0.12.5/bowtie -k 1 -l 28 -n 2 -m 1 --best --best --phred64-quals


*ABySS de novo assembly of failed-to-align reads with (using single-end assembly with k-mer 27)
==Results==
**Using abyss-fac function to calculate statistics
*Original excel table [[Media:021511_whole-read-bowtie-NX1-8.xls]]
  '''ABySS-fac report'''
[[File:sam120610-abyss assembly of fail-to-aligne reads.jpg|700px]]
*Leng filtering and keep only >300bp contigs.
**Script: seqLenFilger.pl


*Perform BLAST on filtered contigs
**BLAST setting: -K 1 -m 8 -a 2 -v 1 –b 1(using -m 8 output format to extract accession number)
**BLAST setting: -K 1 -m 8 -a 2 -v 1 –b 1(using -m 9 output format to get BLAST statistics)
***'''NOTE''': Use -K 1 -v 1 -b 1 for reporting only one hit 
**Taxonomic analysis by NCBI Entrez (against nucletide database)


===Results===
{| {{table}}
[[File:sam120810-failed to aligned ID2-4.jpg|px300]]
| align="center" style="background:#f0f0f0;"|'''Sequence Read Name'''
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID1_NX'''
[[File:sam120810-failed to aligned ID5-8.jpg|px300]]
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID2_NX'''
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID3_NX'''
'''Details of blast hits (PDF)'''
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID4_NX'''
[[Media:ID2-failed-to-align-L300-accession.pdf]] [[Media:ID3-failed-to-align-L300-accession.pdf]]
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID5_NX'''
[[Media:ID4-failed-to-align-L300-accession.pdf]] [[Media:ID5-failed-to-align-L300-accession.pdf]]
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID6_NX'''
[[Media:ID6-failed-to-align-L300-accession.pdf]] [[Media:ID7-failed-to-align-L300-accession.pdf]]
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID7_NX'''
[[Media:ID8-failed-to-align-L300-accession.pdf]]
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID8_NX'''
|-
| Library Nick Name||NX01||NX02||NX03||NX04||NX05||NX06||NX07||NX08
|-
| Library Sample Name||WGA102710_(1_4)_3hr||WGA102710_(1_4)_10hr||WGA102710_(5)_3hr||WGA102710_(5)_10hr||WGA122010_(1-2)_3hr||WGA122010_(1-2)_10hr||WGA122010_(2-3)_3hr||WGA122010_(2-3)_10hr
|-
| DNA starting amount (ng)||288.5 ||500.0 ||114.5 ||500.0 ||181.5 ||200.0 ||50.0 ||200.0
|-
| Proccessed reads||13395058||8843072||6451756||15366518||5071132||7353034||6647566||4888674
|-
| Mapped read, no. (unique)||2707193 (20.21%)||5734751 (64.85%)||1988199 (30.82%)||8571871 (55.78%)||742103 (14.63%)||2963238 (40.30%)||1048466 (15.77%)||1590418 (32.53%)
|-
| Failed to align, no.||10515640 (78.50%)||2512119 (28.41%)||4285220 (66.42%)||4736204 (30.82%)||4311962 (85.03%)||4338695 (59.01%)||5549873 (83.49%)||3176843 (64.98%)
|-
| Mapped read, no (non-unique)||172225 (1.29%)||596202 (6.74%)||178337 (2.76%)||2058443 (13.40%)||17067 (0.34%)||51101 (0.69%)||49227 (0.74%)||121413 (2.48%)
|-
| Total bps coverred, bp||2146219||2001793||797791||835527||914841||767267||1282002||738538
|-
| Genome coverage||46.26%||43.15%||17.20%||18.01%||19.72%||16.54%||27.63%||15.92%
|-
|
|}




'''E.coli taxonomic ranks'''
-------------------------
Domain:  Eubacteria;
Phylum:    Proteobacteria;
Class:        Gammaproteobacteria;
Order:          Enterobacteriales;
Family:            Enterobacteriaceae;
Genus:                Escherichia;
Species:                E. coli


'''Yeast taxonomic ranks''' - (to detect if Yeast DNA contamination was introduced through Yeast tRNA carrier)
{| {{table}}
------------------------
| align="center" style="background:#f0f0f0;"|'''Sequence Read Name'''
Super kindom: Eukaryota;
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID1_NX'''
Domain:         Fungi;
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID2_NX'''
Phylum:           Ascomycota;
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID3_NX'''
Class:             Saccharomycetes;
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID4_NX'''
Order:               Saccharomycetales;
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID5_NX'''
Family:                 Saccharomycetaceae;
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID6_NX'''
Genus:                   Saccharomyces;
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID7_NX'''
Species:                    Saccharomyces cerevisiae
| align="center" style="background:#f0f0f0;"|'''s_4_pe_ID8_NX'''
|-
| ||||||||||||||||
|-
| 1000000 uniqe reads||4948046 ||1542020 ||3244646 ||1792757 ||6835270 ||2481390 ||6341154 ||3074085
|-
|  700000 uniqe reads||||||||||4784689 ||1736973 ||||
|-
| Proccessed reads||4948046||1542020||3244646||1792757||4784689 ||1736973 ||6341154||3074085
|-
| Mapped read, no. (unique)||1028262 (20.78%)||1047304 (67.92%)||1032729 (31.83%)||1060752 (59.17%)||701069 (14.65%)||721899 (41.56%)||1001553 (15.79%)||1029838 (33.50%)
|-
| Failed to align, no.||3853931 (77.89%)||381914 (24.77%)||2118408 (65.29%)||467033 (26.05%)||4067507 (85.01%)||1002826 (57.73%)||5292536 (83.46%)||1965747 (63.95%)
|-
| Mapped read, no (non-unique)||65853 (1.33%)||112802 (7.32%)||93509 (2.88%)||264972 (14.78%)||16113 (0.34%)||12248 (0.71%)||47065 (0.74%)||78500 (2.55%)
|-
| Total bps coverred, bp||1890835||1445489||652120||310201||897507||456971||1267625||676415
|-
| Genome coverage||40.76%||31.16%||14.06%||6.69%||19.35%||9.85%||27.32%||14.58%
|-
|
|}


===Discussions===
==Discussions==
*Most of fail-to-align reads can only be assembled to short contigs (under k-mer 27)
===Mappalbe rates===
*The contigs assembled from fail-to-aligned reads(against E.coli template) '''contain a significant amount of human DNA contigs'''.
*The mappalbe rates are very low in these Nextera libraries (NX_ID_1 ~ NX_ID8), even worse then previous standard Illumina libraries (E.coli SAG) http://genome-tech.ucsd.edu/LabNotes/index.php/Sam:LabNotes/Microbiome-new/2010-12-8.
           ID1  ID2  ID3  ID4  ID5    ID6  ID7    ID8
*The low mappable rates are due to high failed-to-align rates, which could be caused by:
-----------------------------------------------------------
**1. Exogeneous DNA contamination -> Validate failed-to-align reads (or de novo assembled contigs) by BLASTn or MEGAN.
Homo hits  N/A    22    0    25  155    86    35    83
**2. Low quality of reads or containing sequencing tail -> Manually check the read structure
All  hits  N/A    75    15  106  200    145    56    120
'''If it is true that exogeneous DNA contamination is the major issue of libraries quality, then using Nexteral-libraries construction were not albe to reduce DNA contamination effectively.'''
(%)              29.3%      23.6% 77.5%  59.3%  62.5% 69.1%
'''It is possible that the primary exogeneous DNA contamination came from prior MDA amplificaiont, so changing library construction won't help much.'''
All failed-to-align contig hit
*'''A certain level fail-to-align contigs also match E.coli genome (e.g. ID-3 and ID-8)'''. It is unknown why some E.coli reads are failed to map to E.coli genome in bowtie.
***Sequencing error => mis-match >3bp => fail to meet mapping criteria in bowtie?
***Chimeric junction in these "single-end" reads?


*ID-1 showed no BLAST hit using length>300 contigs. I tried using length>200 contigs for BLAST but it still came out no hit.
===Genome coverage===
*I am thinking that maybe it's a better idea to '''use MEGAN metagenomics tool for analyzing these fail-to-align reads''' directly, which it can draw a taxonomic distribution plot easily.
*Genome coverage were compared under the same number of unique reads (down-sampled to 1,000,000 unique reads).
**Libraries '''NX-ID5 and NX-ID6''' were down-sampled to 700,000 unique reads due to the insufficient available mappalbe reads.
*In the short-vs-long paired comparisons (NX-ID1 vs NX-ID2, NX-ID3 vs NX-ID4, NX-ID5 vs NX-ID6, NX-ID7 vs NX-ID8), short-amplified libraries always showed higher (1.5~2 fold) genome coverage than long-amplified libraries.

Latest revision as of 06:19, 15 April 2011

Bowtie mapping on E.coli SAG Nextera libraries HC-MDA-Ecoli-SAGs(NX1-NX8)-Jan09[edit]

Background[edit]

  • These libraries were equally pulled (20ng per library) and sequenced in pair-end(PE) mode.
  • The sequenced data processed by Dr. Zhang for base-calling and saved into FASTQ files.
Wiki(12-29-'10)

Sample no.   Sample name          Total DNA(ng)
-----------------------------------------------
NX01         WGA102710-(1-4)- 3hr     288.50ng
NX02         WGA102710-(1-4)-10hr     500.00ng
NX03         WGA110310-(5)-   3hr     114.50ng
NX04         WGA110310-(5)-  10hr     500.00ng
Wiki(01-06-'11)

Sample no.   Sample name          Total DNA(ng)
-----------------------------------------------
NX05         WGA122010-(1-2)-3hr       181.50ng
NX06         WGA122010-(1-2)-10hr      200.00ng
NX07         WGA122110-(2-3)-3hr        50.00ng
NX08         WGA122110-(2-3)-10hr      200.00ng

Analysis approaches[edit]

Data pooling[edit]

  • Using "gunzip" command to unzip SE files.
  • Combine two SE data(read1 + read2) into one PE data
e.g. samchiang@genemapster:~/Raw-Read-backup/110120_HL087$ cat s_4_1_ID8.txt s_4_2_ID8.txt > s_4_pe_ID8.txt &

Genome coverage analysis (whole reads and 1000000 down-sampled reads)[edit]

  • Script: Read2coverage-kz-3
  • Data resource: pooled pe reads data from HC-MDA-Ecoli-SAGs(NX1-NX8)-Jan09
  • Command
perl Read2coverage-kz-3 /media/disk-2/samchiang/Raw-Read-backup/110120_HL087/pooled_pe/s_4_pe_ID8_3074085.txt /media/disk-2/samchiang/Tools/bowtie-0.12.5/indexes/e_coli_k12_mg1655/e_coli_k12_mg1655 &
  • Bowtie setting
samchiang@genemapster:~/Scripts$ /media/disk-2/samchiang/Tools/bowtie-0.12.5/bowtie -k 1 -l 28 -n 2 -m 1 --best --best --phred64-quals

Results[edit]


Sequence Read Name s_4_pe_ID1_NX s_4_pe_ID2_NX s_4_pe_ID3_NX s_4_pe_ID4_NX s_4_pe_ID5_NX s_4_pe_ID6_NX s_4_pe_ID7_NX s_4_pe_ID8_NX
Library Nick Name NX01 NX02 NX03 NX04 NX05 NX06 NX07 NX08
Library Sample Name WGA102710_(1_4)_3hr WGA102710_(1_4)_10hr WGA102710_(5)_3hr WGA102710_(5)_10hr WGA122010_(1-2)_3hr WGA122010_(1-2)_10hr WGA122010_(2-3)_3hr WGA122010_(2-3)_10hr
DNA starting amount (ng) 288.5 500.0 114.5 500.0 181.5 200.0 50.0 200.0
Proccessed reads 13395058 8843072 6451756 15366518 5071132 7353034 6647566 4888674
Mapped read, no. (unique) 2707193 (20.21%) 5734751 (64.85%) 1988199 (30.82%) 8571871 (55.78%) 742103 (14.63%) 2963238 (40.30%) 1048466 (15.77%) 1590418 (32.53%)
Failed to align, no. 10515640 (78.50%) 2512119 (28.41%) 4285220 (66.42%) 4736204 (30.82%) 4311962 (85.03%) 4338695 (59.01%) 5549873 (83.49%) 3176843 (64.98%)
Mapped read, no (non-unique) 172225 (1.29%) 596202 (6.74%) 178337 (2.76%) 2058443 (13.40%) 17067 (0.34%) 51101 (0.69%) 49227 (0.74%) 121413 (2.48%)
Total bps coverred, bp 2146219 2001793 797791 835527 914841 767267 1282002 738538
Genome coverage 46.26% 43.15% 17.20% 18.01% 19.72% 16.54% 27.63% 15.92%


Sequence Read Name s_4_pe_ID1_NX s_4_pe_ID2_NX s_4_pe_ID3_NX s_4_pe_ID4_NX s_4_pe_ID5_NX s_4_pe_ID6_NX s_4_pe_ID7_NX s_4_pe_ID8_NX
1000000 uniqe reads 4948046 1542020 3244646 1792757 6835270 2481390 6341154 3074085
700000 uniqe reads 4784689 1736973
Proccessed reads 4948046 1542020 3244646 1792757 4784689 1736973 6341154 3074085
Mapped read, no. (unique) 1028262 (20.78%) 1047304 (67.92%) 1032729 (31.83%) 1060752 (59.17%) 701069 (14.65%) 721899 (41.56%) 1001553 (15.79%) 1029838 (33.50%)
Failed to align, no. 3853931 (77.89%) 381914 (24.77%) 2118408 (65.29%) 467033 (26.05%) 4067507 (85.01%) 1002826 (57.73%) 5292536 (83.46%) 1965747 (63.95%)
Mapped read, no (non-unique) 65853 (1.33%) 112802 (7.32%) 93509 (2.88%) 264972 (14.78%) 16113 (0.34%) 12248 (0.71%) 47065 (0.74%) 78500 (2.55%)
Total bps coverred, bp 1890835 1445489 652120 310201 897507 456971 1267625 676415
Genome coverage 40.76% 31.16% 14.06% 6.69% 19.35% 9.85% 27.32% 14.58%

Discussions[edit]

Mappalbe rates[edit]

  • The mappalbe rates are very low in these Nextera libraries (NX_ID_1 ~ NX_ID8), even worse then previous standard Illumina libraries (E.coli SAG) http://genome-tech.ucsd.edu/LabNotes/index.php/Sam:LabNotes/Microbiome-new/2010-12-8.
  • The low mappable rates are due to high failed-to-align rates, which could be caused by:
    • 1. Exogeneous DNA contamination -> Validate failed-to-align reads (or de novo assembled contigs) by BLASTn or MEGAN.
    • 2. Low quality of reads or containing sequencing tail -> Manually check the read structure

If it is true that exogeneous DNA contamination is the major issue of libraries quality, then using Nexteral-libraries construction were not albe to reduce DNA contamination effectively. It is possible that the primary exogeneous DNA contamination came from prior MDA amplificaiont, so changing library construction won't help much.

Genome coverage[edit]

  • Genome coverage were compared under the same number of unique reads (down-sampled to 1,000,000 unique reads).
    • Libraries NX-ID5 and NX-ID6 were down-sampled to 700,000 unique reads due to the insufficient available mappalbe reads.
  • In the short-vs-long paired comparisons (NX-ID1 vs NX-ID2, NX-ID3 vs NX-ID4, NX-ID5 vs NX-ID6, NX-ID7 vs NX-ID8), short-amplified libraries always showed higher (1.5~2 fold) genome coverage than long-amplified libraries.