Latest revision as of 15:07, 14 October 2010

Strain PA-b1

Data

CBCB

Data: 33 bp Solexa reads from 2 libraries
  lib     #reads      #reads(1+N's)   #reads(2+N's)
  s_1    4,105,993    37,933          18,771
  s_7    4,521,907    69,669          59,589
  total  8,627,900 
  ~ 43.55X cvg

Files:

 s_1_sequence.txt  # contains seq & qual; seq names: HWI% ; all reads are 33 bp
 s_1_0001_prb.txt  # 27679 lines of 4*36 values in the -40..40 range  
 s_1_tag.txt       # unique seq count in s_1_sequence.txt
 
 s_7_sequence.txt  # contains seq & qual; seq names HWU%  ; all reads are 33 bp
 s_7_0300_prb.txt
 s_7_tag.txt    
 
 s_1 reads have significantly fewer N's than s_7 reads !!!
 s_1 reads align better to ref than s_2 reads !!!
  
 Conflicts with avg qualities: 28.26 for s_1, 31.05 for s_2 !!!
 
 Base qualities:
 .       #elem           #elem0  #elem<0 min     mean                    max     sum
 s_1     135497769       1863585 3231899 -5      28.2624918053079        40      3829504586
 s_7     149222931       1357384 3113806 -5      31.0508663175903        40      4633501282

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/

NCBI

Same Strains

Ranging from 6.2M to 6.9M

 PA14        CP000438       6,537,648 bp, 66.29 %GC : 5,892 protein coding genes  most similar to PAb1         
 PAO1        AE004091       6,264,404 bp, 66.56 %GC : 5,568                genes  rearrangement vs PA14        
 PACS2       AAQW01000001.1 6,492,423 bp, 66.33% GC : 5,676                genes  rearrangement vs PA14
 PA2192      .              6,905,121 bp, 65.99% GC : 6,157                genes  rearrangement vs PA14
 PAC3719     .              6,222,097 bp, 66.30% GC : 5,650                genes  no rearrangement vs PA14
 PA7         CP000744       6,588,339 bp            : 6,286                genes  no rearrangement vs PA14
 PALESB58    FM209186       6,601,757 bp  66% GC    : 5,925                                           (new)

All vs all nucmer alignments:

 #ref-qry                count   sum     avg%id  comment
 PA14-PAO1               404     6296632 98.66   ~510,582 bp (filter-q) ; ~ 921,919bp (filter-r) in PA14 not present in PAO1
 PA14-PACS2              462     6303773 98.56
 PA14-PA2192             532     6248297 98.47
 PA14-PAC3719            604     6129432 98.44
 PA14-PA7                822     5899920 93.72

Incomplete strains details(Broad):

 Contig length stats:
 desc  #contigs   min     max     mean    stdev   sum
 C3719 124        2079    242903  49572   53770   6146998
 2192  82         2087    398738  83246   88681   682625

5.7Kbp repeat: contains 16S rRNA gene; identical in all the strains

 Coordinates:
 
 PA14:
 732540   738302  +
 4951956  4957714 -  
 5535975  5541728 -
 6312434  6318199 -

 PAO1:
 721550   727325  +
 4788516  4794273 -
 5264042  5269801 -
 6039515  6045289 -

 PA7:
 806558   812299  +
 4982874  4988600 -
 5566182  5571924 -
 6353418  6359148 -

  Pa's vs Pa's (nucmer)

Same Genus, same Group, different Species: Pseudomonas mendocina ymp (Pm)

 Name           Length    %GC    Info
 NC_009439.1    5072807   64.68  4594 genes
 
 PA14 vs Pm (nucmer)
 
                                  #elem   min     max     mean    median  n50     sum      avg%id
 Alignments(1con-1con  : nucmer)  1281    75      19618   1256    1011    1548    1609296  85.26
 Alignments(genes-genes: nucmer)  ...                                 1326 PA14 genes align to 1306 Pm genes (85.36% ID)                                                        
 Alignments(genes-genes: promer)  ...                                 2924 PA14 genes align to 2957 Pm genes (68% ID,75% SYM)

Same Genus, different Species: Pseudomonas syringae pv. tomato str. DC3000 (Ps) : chomosome+2plasmids

 Name           Length    %GC    Info
 NC_004578.1    6,397,126 58.40  chromosome       5478 genes
 NC_004633.1    73,661    55.15  plasmid pDC3000A
 NC_004632.1    67,473    56.17  plasmid pDC3000B
 total          6,538,260        
 
 PA14 vs Ps (nucmer) 
 
                                   #elem   min     max     mean    median  n50     sum      avg%id
 Alignments(1con-1con  : nucmer)   499     75      9435    1296    1161    1619    646614   84.79
 Alignments(genes-genes: nucmer)   438     161     4446    1090    996     1281    477229   83.38      425  PA14 genes align to 418  Ps genes (83.38% ID)
 Alignments(genes-genes: promer)   6193    114     13026   1159    1011    1359    7180527  67.36      2556 PA14 genes align to 2644 Ps genes (67% ID, 75% SYM)

Publications

Dynamics of Pseudomonas aeruginosa genome evolution PNAS Feb 2008

 * 5,021 genes that are conserved across all five genomes analyzed, with at least 70% sequence identity (core genome)
 * Of these, >90% of the genes share at least 98% sequence identity
 
 * The genome for each strain carries a relatively modest number of unique sequences, 10% or less of the number in the core genome, 
   ranging from 79 in C3719, to 507 in PA2192. 
 * Less than 100 genes are unique to pairs of genomes, whereas the number decreases significantly when comparison is carried out among 
   sets of three or four genomes.
 
 * The overall architecture of the pangenome of P. aeruginosa an be represented as a circular chromosome with dispersed polymorphic 
   strain-specific segments, flanked by conserved genes referred to as anchors 
 * These strain-specific segments are regions of genomic plasticity (RGP) and include any region of at least four contiguous ORFs 
   that are missing in at least one of the genomes analyzed.
 * The number of RGPs in individual genomes varies from 27 to 37. A total of 52 RGPs were identified in the five genomes analyzed

Gene-Boosted Assembly of a Novel Bacterial Genome from Very Short Reads Steven L. Salzberg1*, Daniel D. Sommer1, Daniela Puiu1, Vincent T. Lee2, PLoS Sept 2008

Read trimming

No trimming

If no trimming is done and the default AMOScmp pipeline is used (ALIGNWIGGLE=15), the contigs end up containing reference duplications Make-consensus will also take a very long time to run

Duplication causes:

Only the reference nucmer alignment coordinates are stored in the bank and used by casm-layout and make-consensus; the read alignment coordinates are not stored, just the read clear ranges (the ones given in the .afg file); read positions in the assembly layout are approximate
reads can be shifted from their original alignment position by make-consensus by up to 2^3*15=120 bp (constants ALIGN_WIGGLE=15; MAX_ALIGN_ATTEMPTS=3)
The Pseudomonas_aeruginosa reference contains multiple 2 copy 5-10 bp kmers, adjacent within a few dozen bp of each other; reads starting with these kmers align in 2 ways to the reference; if the reads contain errors or SNP's, the 2nd (shorter) alignment to the greedyly built consensus is chosen over the correct one and a duplication is introduced

Quality trimming

Art's script: qual-trim.awk

Min_Qual=15

 .               #elem   #elem0  #elem<0 min     median  max     sum             mean    stdev   n50
 qualTrim        8627900 1266    0       0       23      33      193999902       22.49   7.43    26

!!! Avg clr drops from 33 to 22 bp !!! Duplications still exist though at a lower rate

Alignment based trimming (nucmer)

PA14 is the most similar strain: more PAb1 reads aligned to PA14 than to PAO1 or PA7

 grep -c "^>" PA14-Pa.l-17_c-17.delta PAO1-Pa.l-17_c-17.delta PA7-Pa.l-17_c-17.delta
 PA14-Pa.l-17_c-17.delta:7261629
 PAO1-Pa.l-17_c-17.delta:7036243
 PA7-Pa.l-17_c-17.delta:5393144

Read mapping

PA14 reference

Commands:

blat -noHead -t=dna  -q=dna  -tileSize=10 -stepSize=3 Pa.1con    Pa.seq    Pa.blat
nucmer.pl    -l 20,16,12 -c 20,16,14                  Pa.1con    Pa.seq -p Pa.nucmer
rmap         -m 3 -w 33                            -c Pa.1con    Pa.seq -o Pa.rmap
soap         -v 5                                  -d Pa.1con -a Pa.seq -o Pa.soap
#soap        -v 2 -s 12 -g 3                       -d Pa.1con -a Pa.seq -o Pa.soap

Results:

         runTime                                 hits(all)       hits(>=30bp)
 blat:   1038.28u 3.696s 17:49.00 97.4%          7004925         6353293
 nucmer: .                                       7553734         6832230         # 3 itterations were run for nucmer -l 20,16,12 -c 20,16,14; delta-filter -l 20                        
 rmap:   7318.77u 5.770s 2:02:11.29 99.9%        6958043         6958043
 soap:   685.208u 3.971s 12:45.68 90.0%          7310025         7310025         # fastest, most hits
 all:    .                                       7638420         7201727

Soap seems the "best" read mapping program at current moment.

 #number of hits/read
 cat Pa.soap.5 | awk '{print $4}' | count.pl
 #       total
 1       7185054  # most reads align in only one place
 2       83418
 3       8706
 4       32277
 5       291
 ...
 Total   7310025

 #number of misamtches
 awk '{print $10}' Pa.soap.5 | count.pl
 #       total
 0       3708583   # that many reads align exactly to the reference
 1       2025346
 2       908998
 3       395002
 4       181326
 5       90770

 # mismatch % (Pa's,PA14,Staphylococcus aureus,Strep suis)
 Media:Soap.mismatch.percent.png
 Media:Soap.mismatch_total.percent.png Pa's & Ss profiles are very similar; Sa has the best mappings

Pm reference

 Reads aligned: 1182734

 cat Pm.soap | awk '{print $10}' | count.pl | sort -n
 #number of misamatches
 0       48396
 1       117192
 2       207369
 3       267275
 4       279223
 5       263279
 Total   1182734

Ps reference

 Reads aligned: 632459

 cat Pm.soap | awk '{print $10}' | count.pl | sort -n
 #number of misamatches
 0       16442
 1       36440
 2       81147
 3       136472
 4       176779
 5       185179

Sequence Assemblies

Comparative

PA14 ref assembly

 Location: 2008_0109_AMOSCmp-PA14-relaxed-17-nucmer-redo2
 Command:   AMOScmp Pa -D MINCLUSTER=17 -D MINMATCH=17 -D MINOVL=5 -D MAJORITY=50 -D ALIGNWIGGLE=2

a align all reads to the reference using nucmer. I initially used minmatch=17, mincluster=17 (nucmer -c 17 -l 17)

   8.62M reads, 4.10M HWI 5.32M HWU
   8.30M (...) total reads align to the reference (nucmer -c 14 -l 14)
   7.59M (...) total reads align to the reference (nucmer -c 16 -l 16)
   7.25M (...) total reads align to the reference (nucmer -c 16 -l 16; delta-filter -l 20)
   7.26M (84.16%, 83.50% HWI, 84.76% HWU) total reads align to the reference (nucmer -c 17 -l 17)
   6.54M (75.82%, 74.49% HWI, 77.02% HWU) total reads align to the reference (nucmer -c 20 -l 20)
   5.17M (...) total reads align to the reference (nucmer -c 28 -l 14)
   5.32M (73.27%) reads aligned on their full length (as opposed to 0.94M 33bp quality untrimmed reads)  (nucmer -c 33 -l 17)
   3.69M (42.79%, 38.24% HWI, 46.92% HWU) align exactly (33 bp, 100%id)

   7.00M align using blat -t=dna -q=dna -stepSize=5 -tileSize=10 (20bp alignments, max 3 mismatches)
   6.35M align using blat -t=dna -q=dna -stepSize=5 -tileSize=10 (30bp alignments filtered from above results)
  
   7.11M align using blat -t=dna -q=dna -stepSize=3 -tileSize=10 (20bp alignments, 3 mismatches) 
   500K the unaligned ones were aligned by blat -t=dnax -q=dnax -stepSize=2 -tileSize=4 (20bp alignments, 6 mismatches; only 15% of the aligned reads are unique) 
   => 7.61M aligned by blat/blatx

   6.95M aligned by rmap (-m 3)

   7.30M align by soap (-v 5 -s 12 -g 3) to PA14 chromosome
   6.64M align by soap (-v 2 -s 12 -g 3) to PA14 chromosome
 
   7.24M align by soap (-v 2 -s 12 -g 3 -c 42) to PA14 chromosome (iterative trimming)
   6.72M align by soap (-v 1 -s 12 -g 3 -c 41) to PA14 chromosome (iterative trimming)

   6.36M align by soap (-v 5 -s 12 -g 3) to PA14 genes

   7.43M align by maq

b align all unaligned reads to the reference using nucmer. minmatch=14, mincluster=14 (-c 14 -l 14)
Combine 1a & 1b
trim reads according to their nucmer alignment coordinates; don't trim the ones adjacent to zero cvg regions
assemble them using the AMOScmp pipeline(ALIGN_WIGGLE=2 instead of 15)

Contigs stats:
 desc    #elem   min     max     mean            stdev           sum   singletons
 all     2053    17      170485  3011.84         11917.53        6183320   1127399
 200     428     203     170485  14262.09        22852.74        6104175
 10K     157     10240   170485  35468.89        26531.33        5568616
 
 Data accuracy
 Get all assembled bases with coverage>=20
 count_bases=5926977
 sum_bases=235619001
 sum_errors=2455670

 sum_errors/sum_bases=2455670/235619001=0.01042=1.042% error

SNP's in the 157 >10K contigs: T<->C & A<-G> the most common

  PA14 PAb1    count
   T   C       6153
   C   T       6100
   G   A       6075
   A   G       6021
   G   C       1357
   C   G       1309
   G   T       859
   C   A       846
   T   G       840
   A   C       836
   A   T       338
   T   A       308
   G   .       241
   C   .       196
   .   C       182
   .   G       152
   A   .       101
   T   .       88
   .   T       74
   .   A       70
  Total        32146

Problems:

Some of the short ctgs are not "real"
Some misassemblies due to loose assembly params (MINOVL=5 MAJORITY=50 ): Example: pos 6700

PA14 ref assembly (2) *

run 3 nucmer iterations

Parameters:

 MINCLUSTER=20,16,14   # nucmer -c
 MINMATCH=  20,16,12   # nucmer -l
 MINLEN=20       # delta-filter -l
 MINOVL=10       # casm-layout -o ; previously 5, allowed for missasemblies
 MAJORITY=70     # casm-layout -m ; previously 50 to reduce the number of negative gaps

Read CLR's:

 33      5510573
 32      471289
 31      456959
 30      412738
 29      136423
 28      120148
 27      104213
 26      60338
 25      54828
 24      44854
 23      38574 # smallest value !!!
 22      43302
 21      49923
 20      49572
 Total   7553734

Contig stats:
  desc     #elem   min     max     mean    median  n50     sum     singl   breaks
  all      1843    20      127984  3350    63      30995   6174765 1107200 9
  200bp+   542     200     127984  11245   3263    31106   6094904

Alignment of contigs >=33bp to PA14 reference; qry_hits

 count   totalLen  %id
 940     6130875   99.33       !!! 99.33%id to PA14 vs Pab1

Alignments of PA14 genes to the contigs

 show-coords ctg-genes.filter-1.delta | awk '{print $20}' | count.pl
 #               total
 [CONTAINS]      5318
 [END]           198
 [BEGIN]         191
 [CONTAINED]     190
 [PARTIAL]       88
 [IDENTITY]      3

Location:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0424_AMOSCmp-PA14

PA14 ref assembly (used blat instead of nucmer for read alignment)

run 2 blat iteration:

 blat -t=dna  -q=dna  -stepSize=3 -tileSize=10 Pa.1con Pa.seq 
 blat -t=dnax -q=dnax -stepSize=2 -tileSize=4  Pa.1con Pa.singletons.seq

Contig stats:

         #ctgs   min     max     mean    median  n50     sum     singletons
 all     1372    20      176323  4517    64      47563   6197771 1049302
 200bp+  370     201     176323  16584   3929    47563   6136095
 10Kbp+  146     10097   176323  39176   30457   51554   5719752

Location:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0613_AMOSCmp-PA14-blat

 !!! Best AMOScmp assembly so far : 
 * fewer contig than when using nucmer; 
 * largest contig
 * largest contig sum

PA14 ref assembly (used soap instead of nucmer for read alignment)

soap parameters: -v 5 -g 3 -s 12

 Contig stats:
         #elem   min     max     mean    median  n50     sum
 all     1228    33      132817  5022    106     33779   6167051   !!! smaller assembly  6.16M vs 6.19M with blat
 200bp+  478     200     132817  12789   3824    33988   6113109
 10Kbp+  303     1000    132817  19927   12470   34129   6037852

Location: /fs/szasmg3/dpuiu/Pseudomonas_aeruginosa/Assembly/2008_0625_AMOSCmp-PA14-soap/

soap parameters: -v 2 -g 3 -s 12 -c 42 :
iterative trimming of 2 bp from read ends

 Contig stats:
          #elem   min     max     mean    median  n50     sum
 all      1932    33      79489   3190    153     16122   6162793
 200bp+   874     200     79489   6960    3195    16145   6082694
 10Kbp+   202     10048   79489   20745   17422   21640   4190588

  params                #ctgs  #breaks
  v2.g3.s12.OVL16       2655   34
  v2.g3.s12.OVL20       3446   27
  v2.s12.g3.c42.OVL10   1932   64
  v5.g3.s12.OVL10       1228   83
  v5.g3.s12.OVL3        953    89

v5.g3.s12.OVL10

 all
         #elem   min     max     mean    median  n50     sum       breaks
 PA14    1228    33      132817  5022    106     33779   6167051   83
 PA2192  1784    33      127510  3477    191     16665   6203818
 PAC3719 1335    33      46247   4424    968     14667   5905499
 PACS2   1285    33      87356   4746    165     25837   6098413
 PAO1    969     33      103523  6210    293     25135   6017855
 
 200bp+
         #elem   min     max     mean    median  n50     sum
 PA14    478     200     132817  12789   3824    33988   6113109
 PA2192  880     200     127510  6974    3011    16816   6137484
 PAC3719 838     200     46247   7002    3789    14875   5868046
 PACS2   606     200     87356   9976    3842    25961   6045393
 PAO1    520     200     103523  11510   5601    25135   5985190

PAO1 ref assembly

Same method except that 1.b was not used

 Contigs stats:
 
 desc    #elem   min     max     mean            stdev   sum
 all     2797    17      75626   2161.19         5812.2  6044851
 200     865     200     75626   6893.96         8766.63 5963278
 10K     204     10016   75626   19016.22        10368   3879309
 
 Singletons: 1592525

Location:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0124_AMOSCmp-PAO1-relaxed-17-nucmer-redo2

PAO1 ref assembly (2)

Parameters:

 MINCLUSTER=16   # nucmer -l
 MINMATCH=16     # nucmer -c
 MINLEN=20       # delta-filter -l
 MINOVL=10       # casm-layout -o ; previously 5, allowed for missasemblies
 MAJORITY=70     # casm-layout -m ; previously 50 to reduce the number of negative gaps

Contig stats:

 desc    #elem   min     max     mean    stdev   sum
 all     1961    20      46414   3065    5251    6010800
 200bp+  1166    203     46414   5108    6007    5955960
 1758777 singletons

Alignment of contigs >=33bp to PAO1 reference; qry_hits

 count   totalLen  %id
 1471    5984905   99.06       !!! 99.06%id PAO1 vs Pab1

Location:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0424_AMOSCmp-PAO1

PA14 & PAO1 merge

Use minimus to merge all contigs in the PA14 & PAO1 reference assemblies

Filter contigs that contain at least 2 adjacent PA14 merged by PAO1

 desc       #elem   min     max     mean            stdev           sum
 all        1850    17      236472  3400.31         16863.79        6290586
 200        306     204     236472  20318.3         37147.91        6217401
 10K        113     10520   236472  52647.45        45602.76        5949162
 
 Singletons: 1066226

De novo assembly of singletons using SSAKE Input: 1066226 singletons

 Ssake version 3.0 run with default parameters:
 -m  Minimum number of reads needed to call a base during overhang consensus build up (default -m 16)
 -o  Minimum number of reads needed to call a base during an extension (default -o 2)
 -r  Minimum base ratio used to accept a overhang consensus base (default -r 0.6)

 Contigs stats:
 
 desc    #elem   min     max     mean    stdev   sum
 all     19795   34      2866    67.78   102.73  1341825
 200     879     200     2866    416.72  308.92  366304

 200 bp contigs stats
 
 desc.           #elem   min     mean    max     sum
 contig_len      879     200     416.72  2866    366304
 contig_reads    879     52      414.21  3157    364091
 contig_cvg      879     7.91    31.58   146.52  27764.51 
 
 Singletons:  702007

 Find contigs overlapping 12 bp or more and merge them using EMBOSS merger program:
 
 desc            #elem   min     max     mean    stdev   sum
 original        879     200     2866    416.72  308.92  366304
 new             670     200     4826    539.32  579.26  361350

Pa-b1 AMOScmp reference assembly

 Contig stats:
  
 desc                           #elem   min     max     mean      stdev           sum
 input(reference|final)         770     200     512638  8484.21   37210.58        6532844
 output(AMOScmp contigs)        936     21      260827  6988.98   22288.49        6541689
 
 Singletons: 1257963
 Snp's:      1941

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0215_AMOScmp-PAb1

Use ALL strains as reference

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0228_AMOSCmp-ALL

Steps:

 * Identify related strains (completed & draft assemblies): PA14, PAO1, PA7, PACS2, PA2192, PAC3719
 * Download & merge related strains into a multi-FASTA file
 * Align PAb1 reads to multi-FASTA file using nucmer (Ex: -l 17 -c 17)
 * Count the number of PAb1 reads aligned to each strain

    strain   #reads_aligned
    ======================
    PA14    7261629
    PA2192  7193960
    PACS2   7089455
    PAO1    7036243
    PAC3719 6928744
    PA7     5393144
    ======================

 * Sort the strains in descending order according to the number of reads aligned to each one
 * The strain with most hits is PA14 => closest relative

 * Assuming the order is PA14, PAO1, PA7, PACS2, PA2192, PAC3719 ...
 * Align strains to one other; identify unique regions in each genome (regions which are not present in the previously considered ones)
 
   Unique regions length summary:

   Table 1:
   reference                             #regions  min     max     mean    stdev   sum
   =======================================================================================
   PA14                                  1       6537648 6537648 6537648 0       6537648
 
   PAO1-PA14                             115     1       19263   2038.58 3843.93 234437
   PA7-PA14                              511     1       46020   1876.96 4872.51 959129
   PACS2-PA14                            147     1       90591   3054.19 8766.3  448967
   PA2192-PA14                           231     1       92553   3746.3  11475   865396
   PAC3719-PA14                          272     1       14058   1164.78 2089.85 316822
   =======================================================================================
   ALL                                   1277                                     9362399

   Table 2:
   reference                             #regions  min     max     mean    stdev   sum
   =======================================================================================
   PA14                                    1       6537648 6537648 6537648 0       6537648
 
   PAO1-PA14                               115     1       19263   2038.58 3843.93 234437
   PA7-PAO1-PA14                           487     1       46020   1830.82 4910.86 891614
   PACS2-PA7-PAO1-PA14                     71      2       79031   4112.94 10618.2 292019
   PA2192-PACS2-PA7-PAO1-PA14              181     1       92066   3439.9  9500.07 622622 #some N's maybe?
   PAC3719-PA2192-PACS2-PA7-PAO1-PA14      173     1       11425   998.16  1888.94 172682 #some N's maybe?
   =======================================================================================
   ALL                                     1028    2       6537648 8513.66 203929  8752049

   Table 3:
   reference                             #regions  min     max     mean    stdev   sum
   =======================================================================================
   PA14                                    1       6537648 6537648 6537648 0       6537648
   
   PA14-PAO1                               123     1       43690   4009.28 7903.16 493142
   PA14-PA7                                516     1       31464   1768.16 3534.88 912375
   PA14-PACS2                              146     2       43690   3259.15 6129.12 475836
   PA14-PA2192                             222     1       34854   2189.76 4338.6  486127
   PA14-PAC3719                            277     1       43690   2206.38 5214.55 611170
   =======================================================================================

 * Extract unique regions(Table 2) and add them to the most similar genome

 * Assemble using all these regions as reference => 10,201 contigs; 1,052,178 singletons

   reference                        #ctgs   min     max     mean    stdev           sum
   =======================================================================================
   PA14                             2056    17      57476   2992.27 7078.25         6152122 
   PAO1                             877     17      8159    95.16   482.74          83456
   PA7                              4723    17      6822    41.7    128.87          196996
   PACS2                            914     17      3346    66.96   148.47          61207
   PA2192                           1528    17      49826   134.88  1734.19         206109
   PAC3719                          103     17      693     42.11   85.89           4338
   =======================================================================================
   ALL                              10201   17      57476   657.21  3457.1          6704228
 
 * Compare ALL assembly contigs to PA14 ref assembly :

   reference                        #ctgs   min     max     mean    stdev           sum        #singletons
   ========================================================================================================
   ALL                              10201   17      57476   657.21  3457.1          6704228    1052178 
   ALL(200bp+)                      1020    200     57476   6235.17 9220.76         6359881
   ALL(PA14)                        2056    17      57476   2992.27 7078.25         6152122 
   PA14                             1984    17      82196   3103.21 8274.53         6156788    1375004
   ========================================================================================================
   ALL-PA14                         8217     ......                                 547440                 # 0.5M extra bp
   ALL(PA14)-PA14                   72       ......                                 -4666                  # 4666 bp which were originally mapped to PA14 are now mapped to other strains

  * Compare ALL assembly contigs to ssake 200bp+ contigs :
 
   ssake 200bp+ contigs:
     strain    #ctgs   min     max     mean    stdev   sum
     PAb1      879     200     2866    416.72  308.92  366304

   ALL vs ssake alignments:
    desc    #elem   min     max     mean    stdev   sum
    ALL     265     23      49826   912.9   4187.17 241919
    ssake   306     200     2358    379.85  257.04  116235

  ALL vs ssake CONTAINED alignments:
    desc    #elem   min     max     mean    stdev   sum
    ALL     164     23      49826   1264.71 5268.16 207413
    ssake   260     200     2358    368.09  259.02  95704

  ~300 out of 879 ssake contigs were assembled using this method.
  The ALL assembly contigs ssake contigs align to are longer than the ssake contigs!

  * Compare ALL assembly contigs to final de novo contigs :

       id                          len         #de novo ctgs it CONTAINS
  ======================================================================
    1  2361_PA2192_2249428_2341494 49826       52
    2  2360_PA2192_2249428_2341494 42258       49
    3  2364_PA2192_2346558_2361947 15388       17
    4  10190_PAO1_747565_761293    8159        7
    5  2429_PA2192_2405271_2412668 6103        7
    6  3259_PA2192_6656583_6693732 2594        6
    7  2367_PA2192_2374235_2379148 4914        5
    8  6917_PA7_4563267_4570077    6822        5
    9  2958_PA2192_4496915_4500004 3090        4
   10  2375_PA2192_2383706_2386485 2789        3
   11  2366_PA2192_2366515_2369415 2901        2
   12  6919_PA7_4579629_4582081    1085        2
   13  5364_PA7_3138605_3164438    1931        2
   14  8663_PACS2_3494057_3573088  449         2
   15  9322_PACS2_982458_983976    1519        2
   16  2368_PA2192_2379432_2380661 1230        2
   17  2716_PA2192_3268239_3292117 567         2
   18  9587_PAO1_2439014_2458239   3513        1
   ...
   49  8384_PAC3719_5517906_5518496 591        1
  ======================================================================
      Total                                    201

PA14 Sanger maq on all reads

 Contig stats:
 min_len    #contigs   min     max     mean     stdev           sum
 0          991        33      155551  6199.79  17445.05        6143996  !!!Smaller assembly size than with AMOScmp
 200        368        200     155551  16581.7  25475.92        6102067
 10K        149        10048   155551  38164.35 28513.51        5686489

1197385 singletons; SNP's in the contigs: T<->C & A<-G> the most common

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0206_PA14_maq

PA14 gene ref assembly

Use the 5892 PA14 genes as reference
Use "soap -v 5 -s 12 -g 3" for alignments

=> 5671 scaffolds (that many genes are shared by the 2 strains)

contig stats:

         #elem   min     max     mean    median  n50     sum   
 all     6487    33      13070   854     752     1164    5541032
 200bp+  5666    200     13070   965     840     1174    5470492

Comparative-ABBA

ABBA

Steps

Annotate draft assemby
Find fragmented genes at ends of contigs. Example shown is of putative secretion protein from PAO1
Fill in missing amino acid sequence using reference protein
Merge contigs together with output sequence from ABBA

Pipeline

Input reference amino acid sequence and database of short reads
Search read database (tblastn)‏
Create layout using alignment coordinates
Validate good read coverage across entire sequence
Make consensus from a multiple alignment

PA14:PO1

5892 genes in PA14 vs 5568 in PAO1
5244 genes in PA14 IDENTITY (99.03%) to 5239 genes in PAO1
5379 genes in PA14 IDENTITY|CONTAIN|BEGIN|END (99.003%) to 5367 genes in PAO1
5389 genes in PA14 align to 5372 genes in PAO1 (98.69%id)
501 genes in PA14 & 196 genes in PAO1 don't align to opposite strand genes

!!! The genes that align are very conserved

PA14 gene ref assembly

 desc            #elem   min     max     mean    stdev   sum
 genes(nt)       5892    72      15639   992     752     5846175  # out of a total of 6537648bp genome => 89.42 coding
 genes(aa)       5892    23      5212    330     251     1942833

PAO1 gene ref assembly

 desc            #elem   min     max     mean    stdev   sum
 genes(nt)       5568    72      16884   1006    751     5601468  # out of a total of 6264404bp genome => 89.41 coding
 genes(aa)       5568    23      5627    334     250     1861588

De-novo

SSAKE

version 3.0 with default parameters

 Contigs stats:
 
 desc    #elem   min     max     mean    stdev   sum
 length  185030  34      5490    77.21   141.23  14287079
 reads   185030  2       13352   29.52   127.65  5463405 
 
 Singletons: 3,164,495  (singletons file + ambiguous reads)

Run times:

 18459.472u 31.849s 5:08:59.82 99.7%	0+0k 0+0io 0pf+0w

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0214_ssake/

Latest version is 3.2: run on the same dataset generated fewer contigs and more singletons; also took longer to run

SSAKE (redo)

version 3.2.1 with default parameters

 Contigs stats:
 
 desc    #elem   min     max     mean    median  n50     sum
 length  107046  34      3994    96.39   42      201     10317932
 reads   ...
 
 Singlets: 3,309,974

Velvet *

Version 0.5.05 Ovl len=23

Commands:

 host: sycamore (Opteron)
 velveth . 23 Pa.seq 
 velvetg . -read_trkg yes

Run times:

 velveth : 88.614u  22.865s 2:16.35  81.7%  0+0k 0+0io 0pf+0w
 velvetg : 439.939u 19.258s 8:19.77  91.8%  0+0k 0+0io 0pf+0w
 total   : 528.553u 42.123s 10:36.12

Output files (used)

 contigs.fa : contig FASTA file

Contigs stats:
                #elem   min     max     mean    median  n50     sum
  all           10684   45      16239   640     363     1184    6841458
  200bp+        7382    200     16239   877     586     1252    6474426
  breaks        257  (-b 20) 
  breaks(200+)  251  (-b 20) coming from 170 ctgs
  1,273,164 singletons; 879,840 of them aligned by soap (-v5,-s12,-g3) to PA14

Distance between breakpoints:

 #elem   min     max     mean    median  n50     sum
 251     0       314620  24946   3776    99406   6261404

Align velvet ctgs to PA14 ref (nucmer -c 40; delta-filter -q):

 alignemnts  length  %id
 9927        6164312 99.20  (slightly lower %id than the AMOS PA14 ref assembly: 99.33%)

PA14 0cvg regions

                      #elem   min     max     mean    median  n50     sum
original              2235    1       31503   197     22      7529    439633  # 6098014 of PA14 covered
delta-filter -q       2378    1       31676   220     24      6434    522246  # 6015401 of PA14 covered

rRNA repeat is collapsed; Contigs aligned to the 4 copy 5.7Kbp PA14 rRNA repeat:

     [S1]     [E1]  |     [S2]     [E2]  |  [LEN 1]  [LEN 2]  |  [% IDY]  |  [LEN R]  [LEN Q]  |  [COV R]  [COV Q]  | [TAGS]
 ===============================================================================================================================
       23      409  |        1      388  |      387      388  |    97.94  |     5763      388  |     6.72   100.00  | rRNA        856     [CONTAINS]
      388     2132  |        1     1745  |     1745     1745  |    99.94  |     5763     1745  |    30.28   100.00  | rRNA        642     [CONTAINS]
     2111     2279  |      169        1  |      169      169  |   100.00  |     5763      169  |     2.93   100.00  | rRNA        3433    [CONTAINS]
     2258     5763  |        1     3506  |     3506     3506  |    99.89  |     5763     3506  |    60.84   100.00  | rRNA        260     [CONTAINS]

Question: Are there some velvet contigs contained in each other? Answer: Very few: 302 (out of 10K) ; avg %id is 86 and avg lenght is 106 bp; nucmer -c 20 has been used for alignments

Many contigs seem to overlap by ~ 20bp

Some contigs share end reads labeled as D(uplicates)

 listReadPlacedStatus Pa.bnk | awk '{print $3}' | count.pl
 #       total
 P       7322684
 S       1273164
 D       32052     # 33 of them belong to 3 contigs
 Total   8627900

About ~4607 contig pairs could be merged:

 grep D Pa.status.detailed | awk '{print $5,$6}' | sort -u | wc -l
 4607

Most of the overlaps are shorter than 10bp; if "nucmer -c 20 -l 20" is run on the contig ends; only 627 of them are found

Example:

 Contigs 1915 & 3551 share 15 Duplicate reads
 cat Pa.status | grep 1915 | grep 3551$ | wc -l
 15 

They align to adjacent pos in PA14 :
 3743092 3745646 |  2555  1  | 2555  2555  |    99.49  |  6537648     2555  |     0.04   100.00  | 1      1915    [CONTAINS]
 3745625 3746382 |   758  1  |  758   758  |    99.60  |  6537648      758  |     0.01   100.00  | 1      3551    [CONTAINS]

411870 out of 1273164 singletons (~33%) align to PA14 (nucmer -maxmatch -c 20) 304162 out of 1273164 singletons (~25%) are also singletons in the AMOScmp-PA14 assembly; probably the common ones are low quality sequences that just don't assemble

PA14 reference alignments: nucmer -maxmatch -c 40):

Aligned:

 #elem   min     max     mean    stdev   sum
 9859    45      15342   642     785     6331357

Completely aligned:

 #elem   min     max     mean    stdev   sum
 9681    45      11120   617     712     5969126

 9208 are aligned 1 time
 408  are aligned 2 times  #repeats
 39   are aligned 3 times  #repeats
 25   are aligned 4 times  #repeats
 1    is aligned  6 times  #repeat

Partially aligned:

 #elem   min     max     mean    stdev   sum
 178     78      15342   2035    2139    362231 
 
 120 are aligned once:           Pa.1con-velvet.uniq.not_contained.filter-q.1.coords 
 58  are aligned more than once: Pa.1con-velvet.uniq.not_contained.filter-q.2.coords

 All these contigs have AMOScmp assembly singletons aligned to them!!! 
 These velvet contigs seem correct.

 Example: There is a large insertion in PA14 at position 6700
 
 PA-b1 reads vs PA14 ref  
  find-query-breaks -B Pa.filter-q.delta
  => 3 reads end at base 6700
   
 velvet ctgs vs PA14 ref:
 ...
   6177     6700  |    1   524  |  524    524  | 99.43  |  6537648  875  | 0.01  59.89  | 1  900     # 431 out of 1054 seqs in this contigs are singletons in AMOScmp assembly
   6694     8045  |  743   2094 |  1352   1352 | 99.63  |  6537648  2097 | 0.02  64.47  | 1  1970    # 895 out of 2354 seqs in this contigs are singletons in AMOScmp assembly; unaligned regions aligns at low %id to a Pseudomonas fluorescens hypothetical protein

 Nucmer/Soap alignments of the reads adjacent to position 6700
 Media:PA14.coord.6700 : ~20 reads end at pos 6700, ~12 reads start at pos 6694 (6bp ovl only)
 Media:PA14.soap.v5.s12.g3.6700
 Media:PA14.soap.v2.s12.g3.6700 10 bp ovl between reads that start at pos 6671(2 mismatches) & 6694(0 mismatches)
 Media:PA14.soap.v2.s12.g3.c42.6700 11 bp ovl between reads that start at pos 6676(2 mismatches) & 6694(0 mismatches)

 Other breaks: 
 ...
  66764    67073  |    1   310  |  310    310  | 100.00 |  6537648  1282 | 0.00  24.18  | 1 320
  67181    67322  | 1136  1277  |  142    142  |  98.59 |  6537648  1282 | 0.00  11.08  | 1 320
 ...
 116429   117249  |    1   821  |  821    821  |  99.03 |  6537648  1457 | 0.01  56.35  | 1 1955
 120490   120573  |   84     1  |   84     84  | 100.00 |  6537648  1693 | 0.00  4.96   | 1 1827
...
 234641   235980  | 3294  1955  | 1340   1340  |  99.40 |  6537648  3294 | 0.02  40.68  | 1 139
 235981   237442  | 1466     5  | 1462   1462  |  99.32 |  6537648  3294 | 0.02  44.38  | 1 139
...
 286003   287915  | 1216  3128  | 1913   1913  |  99.22 |  6537648  3128 | 0.03  61.16  | 1 3068
...

Not aligned:

 #elem   min     max     mean    stdev   sum
 825     45      16239   618.3   1209.37 510101        # the longest velvet ctg does not align

Not CONTAINED:

 #elem   min     max     mean    stdev   sum
 1002    45      16239   868.82  1518.56 870566

Velvet vs all Pseudomonas strains:

 strain          aligned         CONTAINED       partially-aligned    
 PA14            9859            9681            178                # smallest number
 PA2192          9940            9607            333
 PACS2           9688            9469            219
 PAO1            9587            9402            185
 PAC3719         9427            9057            370
 PA7             8514            7651            863
 all             10378           10243           135

 not aligned stats: 
 #elem   min     max     mean    median  n50     sum
 306     45      9681    549     119     2310    168094

 strain(s)                               CONTAINED
 PA14                                    9681
 PA14|PA2192                             10081
 PA14|PA2192|PACS2                       10186
 PA14|PA2192|PACS2|PAO1                  10198
 PA14|PA2192|PACS2|PAO1|PAC3719          10210
 PA14|PA2192|PACS2|PAO1|PAC3719|PA7      10243

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24

Velvet version 0.5.07 : does not "seem" to be an improvement over Version 0.5.05 Ovl len=23

 desc              #elem   min     max     mean    stdev   sum
 contig            11230   45      16239   611.64  814.09  6868784
 contig(200bp+)    7353    200     16239   878.61  897.06  6460467

Some of the 306 contigs should probably assemble together:

Assembled with minimus2 REFCOUNT=0 OVL=20 => 54 ctgs that contain 132 reads & 174 singlets
=> reduced the number from 306 to 54+174= 228

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/no_hits

!!! If AMOScmp is run only on the velvet singlets (mostly low quality), some of the resulting contigs will have "artificial" gaps

         #elem   min     max     mean    median  n50     sum
 all     43926   33      963     122     93      157     5367766
 200bp+  6757    200     963     293     263     289     1980273

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/AMOScmp-PA14.singletons/1con-contigs.filter-q.na.coords

Edena

Version 2.1.1 , default parameters

 Contigs stats:  
               #elem   min     max     mean    median  n50     sum
 all           11180   100     11300   552     356     837     6175460
 200bp+        8316    200     11300   693     488     902     5759209
 breaks(PA14)        245 (-b 20)   # 78 of them also present in velvet
 breaks(velvet)      210 (-b 20)
 
 3955865 singletons

Contig don't seem to overlap much

Run times:

 ovl step:      1655.561u 7.560s  27:54.62 99.3%   0+0k 0+0io 0pf+0w
 assembly step: 56.141u   3.791s   1:03.82 93.9%   0+0k 0+0io 0pf+0w
 total:         1711.702u 11.351s 28:58.44

A significant number of edena contigs CONTAIN one(several) velvet contigs:

 $ show-coords -I 94 velvet-edena.delta | egrep 'CONTAINED|BEGIN|END' | awk '{print $19}' | count.pl -m 1 | wc -l # 2269
 $ show-coords -I 94 velvet-edena.delta | egrep 'CONTAINED|BEGIN|END' | awk '{print $19}' | count.pl -m 2 | wc -l # 1455

Align edena ctgs to PA14 ref (nucmer -c 40; delta-filter -q):

       #elem   min     max     mean    median  n50     sum
 0cvg  8809    1       34854   119     38      438     1050135

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0312_edena/

SSAKE/Edena/Velvet summary

host: sycamore.umiacs.umd.edu (Linux 64)

 program  version   run-time
 SSAKE    3.0       8459.472u 31.849s 5:08:59.82  
 Velvet   0.5.05    528.553u  42.123s   10:36.12
 Edena    2.11      1711.702u 11.351s   28:58.44

Input: 8627900 33 bp Solexa reads(2 lanes)

Output: contigs+singletons

Contig statistics:

 program #seqs   min  max    mean  med  n50   sum
 SSAKE   185030  34   5490   77    44   87    14287079
 Velvet  10684   45   16239  640   363  1184  6841458
 Edena   11180   100  11300  552   356  837   6175460
 Euler   7090    22   22671  953   312  2622  6759979

200bp+ contig statistics:

 program #seqs   min  max    mean  med  n50   sum
 SSAKE   12532   200  5490   486   385  549   6090567
 Edena   8316    200  11300  692   488  902   5759209
 Velvet  7382    200  16239  877   586  1252  6474426

Singleton statistics:

 program #seqs
 SSAKE   3164495
 Velvet  1273164 
 Edena   3955865

Euler-SR

 cmd:     Assemble.pl Pa.seq vertexSize
 version: 12/4/08

Ctg len stats:

 vertexSize       #elem   min     max     mean    median  n50     sum
 21               7090    22      22671   953     312     2622    6759979
 25               11130   26      9619    564     359     912     6276548

Gene assemblies

Dan Sommer reduced number of contigs >=200bp from 306 to 120 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Annotation/PAb1.fasta

 Contig stats:
 desc    #elem   min     max     mean    stdev   sum
 contigs 120     212     512638  51438   81999   6172675

 Average Gap: 105 nt bases
 Median Gap: 14 nt bases
 Largest Gap: 1095 nt bases

927 singletons assembled

Contig Assemblies

Comparative

AMOScmp on the velvet contigs (PA14 ref) *

Use only velvet ctgs
This assembly will have some of the repeats collapsed

All:

       #elem   min     max     mean    median  n50     sum
 ctg   3046    45      27957   1928    1101    3744    5873032
 singl 989     45      16239   836     292     2221    826462
 all   4035    45      27957   1660    850     3573    6699494
 breaks 270

200 bp+:

       #elem   min     max     mean    median  n50     sum
 ctg   2632    200     27957   2213    1390    3798    5824679
 singl 586     200     16239   1344    721     2352    787686
 all   3218    200     27957   2055    1243    3631    6612365
 breaks 265

Reduced the number of contigs from 10684 to 4035 !!!

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/AMOScmp-PA14

AMOScmp on the velvet contigs & velvet singletons (PA14 ref) *

Use singletones from velvet ctgs as well
This assembly will have some of the repeats collapsed

                   #elem   min     max     mean    median  n50     sum
  all              2383    33      177630  2592    72      45071   6177155
  200bp+           513     200     177630  11778   738     45351   6042047
  all-velvet       334     59      177630  17920   7295    45351   5985407  # ctgs that contain at least 1 velvet ctg
  singl-velvet     984     45      16239   837     292     2221    823862   # velvet ctgs left as singletons
  singl-read       393416  33      45      33      33      33      12983724 # original Solexa reads not assembled

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/AMOScmp-PA14.singletons

AMOScmp on the velvet contigs (PAO1 ref)

All:

       #elem   min     max     mean    median  n50     sum
 ctg   2982    45      27957   1928    1124    3652    5749909
 singl 1255    45      16239   764     247     2275    958885
 all   4237    45      27957   1583    773     3440    6708794
 breaks(PAO1) 271

200 bp+:

       #elem   min     max     mean    median  n50     sum
 ctg   2610    200     27957   2186    1380    3675    5705769
 singl 707     200     16239   1283    674     2398    906818
 all   3317    200     27957   1994    1205    3494    6612587
 breaks(PAO1) 266

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/AMOScmp-PAO1

AMOScmp on the velvet contigs (PA14 & PA2192 ref)

Input: above contigs & singletons

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     2841    45      34374   2022    2570    5744186
 singl   751     45      17475   1251    2264    939443   # 593 are velvet original ctgs
 all     3592    45      34374   1861    2529    6683629

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     2489    200     34374   2291    2638    5701588
 singl   456     200     17475   2001    2648    912503
 all     2945    200     34374   2246    2641    6614091

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/AMOScmp-PA14-PA2192

AMOScmp on the velvet contigs (PA14, PA2192 & PACS2 ref)

Input: above contigs & singletons

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     2759    45      28263   2029    2554    5597067
 singl   725     45      34374   1495    2744    1084142  # 493 are velvet original ctgs
 all     3484    45      34374   1918    2603    6681209

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     2419    200     28263   2297    2618    5556096
 singl   463     200     34374   2290    3170    1060442
 all     2882    200     34374   2296    2714    6616538

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/AMOScmp-PA14-PA2192-PACS2

AMOScmp on the velvet contigs (PA14, PA2192, PACS2 & PAO1 ref)

Input: above contigs & singletons

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     2696    45      24046   2028    2510    5468683
 singl   753     45      34374   1609    3050    1211809  # 483 are velvet original ctgs
 all     3449    45      34374   1937    2643    6680492

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     2369    200     24046   2292    2569    5429268
 singl   487     200     34374   2439    3527    1187586
 all     2856    200     34374   2317    2756    6616854

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/AMOScmp-PA14-PA2192-PACS2-PAO1

AMOScmp/minimus on the AMOScmp (PA14 ref) & velvet "unique" contigs

Input:

 AMOScmp PA14 ref assembly (2) contigs
 velvet "unique: contigs
 !!! No singletons

AMOScmp2 steps:

Align AMOScmp & velvet contigs (nucmer -c 40) to each other
Identify velvet contigs CONTAINED in AMOScmp contigs(max trimming of 20)
Remove CONTAINED contigs from the velvet set
Create alignment file based to the AMOS.scaff file (if nucmer is repeat contigs are collapsed/no placed exactly)
Align velvet contigs to reference using nucmer ; filter only the CONTAINED ones
Merge AMOS & velvet alignments
Run "casm-layout -S -r ..." ; if -S is not used, many AMOS contigs are not placed
Run "make-consensus ..."

All AMOS input contigs should be part of the new assembly
If the velvet CONTAINED alignments are not filtered, some end conflicts might cause some of the AMOS contigs end up as singletons

Input:

 All:
 desc            #elem   min     max     mean    stdev   sum
 AMOS            1843    20      127984  3350.38 10776   6174765
 velvet-uniq     1189    45      16239   997.85  1485.34 1186452  # 462 aligned to PA14 ref (nucmer -c 40); 
                                                                  # 284 CONTAIEND in the ref
                                                                  # 178 (462-284) partially aligned : 
                                                                    # 58  have 2+ alignments 
                                                                    # 120 have 1  alignment
                                                                  # 727 (1189-462) not aligned to the ref
 all             3032    20      127984  2427.84 8530.33 7361217

 200bp+:
 desc            #elem   min     max     mean    stdev   sum
 AMOS            542     200     127984  11245   17520   6094904
 velvet-uniq     848     200     16239   1354.36 1627.86 1148499  # 422 aligned to PA14 ref (nucmer -c 40); 249 CONTAIEND in the ref; 426 not aligned to the ref
 all             1390    200     127984  5211.08 12019.6 7243403

Output:

 All:
 desc            #elem   min     max     mean    stdev   sum
 ctg             1327    20      221554  4644.9  19079.  6163795    # two of them(530,531) contain only velvet seqs
                                                                    # 1210 contain only AMOS seqs ; 11 contain 2 AMOS seqs (negative gap got merged)
                                                                    # 115 contain both AMOS & velvet seqs        
 singl           817     45      16274   1053.94 1719.12 861076     # all but one(807) are velvet : see .conflict file
                                                                    # 170 out of 178 partially aligned (see above) is a singleton
 all             2144    20      221554  3276.52 15146.  7024871

 200bp+:
 desc            #elem   min     max     mean    stdev   sum
 ctg             282     203     221554  21651   36731   6105733
 singl           571     200     16274   1461.55 1917.62 834546
 all             853     200     221554  8136.31 23189   6940279

minimus steps:

Run minimus2 pipeline on the above contigs(1327) and singletons(817) ; Only uses contig:singletion overlaps (not all:all)

 $ ~dpuiu/bin/minimus2 -D BREAKID=1327 -D OVERLAP=40 Pa

min OVL=40 : use contig:singleton alignments

Output:

 All:
 desc            #elem   min     max     mean    stdev   sum
 ctg             180     135     253489  22817.2 45339   4107096   # 218 of the PA14 unaligned velvet contigs assembled
 singl           1377    20      153036  1996.79 10534   2749587
 all             1557    20      253489  4403.77 19467   6856683

 200bp+:
 desc            #elem   min     max     mean    stdev   sum
 ctg             174     210     253489  23598.3 45919   4106116
 singl           488     200     153036  5519.51 17154   2693524
 all             662     200     253489  10271.3 28846   6799640

min OVL=20 : use contig:singleton alignments (not significantlty better than ,on OVL=40)

Output:

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     279     63      253489  15557   37925   4340305 # 321 of the PA14 unaligned velvet contigs assembled
 singl   972     20      153036  2576    12252   2504000
 all     1251    20      253489  5471    21580   6844305

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     263     214     253489  16495   38868   4338061
 singl   393     200     153036  6261    18680   2460755
 all     656     200     253489  10364   28954   6798816

min OVL=20 : use contig:singleton & singleton:singleton alignments

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     277     103     253489  15983   37991   4427292  # 556 of the PA14 unaligned velvet contigs assembled
 singl   808     20      153036  2986    13439   2412439  # 260 of them are velvet original contigs
 all     1085    20      253489  6304    23111   6839731

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     263     222     253489  16826   38812   4425275
 singl   298     201     153036  7971    21242   2375456  # 165 of them are velvet original contigs
 all     561     201     253489  12123   31042   6800731

min OVL=20 : use contig:singleton & singleton:singleton alignments & adjacent contig ovelaps (best)

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     317     47      253489  16450   40806   5214579
 singl   703     20      112908  2324    11332   1633528  # 335 are original velvet contigs
                                                          # 29 are velvet contigs partially aligned
 all     1020    20      253489  6714    25448   6848107

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     291     200     253489  17909   42289   5211547
 singl   247     201     112908  6480    18431   1600539  # 186 are original velvet contigs
 all     538     200     253489  12662   33969   6812086

Locations:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-velvet
 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-velvet/AMOScmp2
 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-velvet/AMOScmp2.not_CONTAINED

AMOScmp/minimus on the AMOScmp (PA14 ref) & velvet contigs *

Use all velvet contigs
Very similar with the previous case (the velvet CONTAINED contigs were discarded at the beginning)

AMOScmp2

Input:

 All:
 desc    #elem   min     max     mean    stdev   sum
 AMOS    1843    20      127984  3350    10777   6174765
 velvet  10684   45      16239   640     825     6841458
 all     12527   20      127984  1039    4311    13016223
 
 200bp+: 
 desc    #elem   min     max     mean    stdev   sum
 AMOS    542     200     127984  11245   17521   6094904
 velvet  7382    200     16239   877     896     6474426
 all     7924    200     127984  1586    5344    12569330

Output:

 All: 
 desc    #elem   min     max     mean    stdev   sum
 ctg     1361    20      221563  4536    19428   6173665
 singl   838     45      22884   1064    1868    891463
 all     2199    20      221563  3213    15418   7065128

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     299**   203     221563  20448   37378   6113919
 singl   576     200     22884   1499    2114    863431
 all     875     200     221563  7974    23668   6977350

minimus2 on AMOScmp2 output

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     284     47      274913  18279   46193   5191268
 singl   1001    20      133259  1666    9754    1667778
 all     1285    20      274913  5338    24330   6859046

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     246     210     274913  21088   49046   5187663
 singl   285     201     133259  5696    17670   1623459
 all     531     201     274913  12827   36583   6811122

Locations:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-velvet/AMOScmp2
 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-velvet/AMOScmp2/minimus2/

AMOScmp/minimus on the AMOScmp PA14 & PAO1 contigs

AMOScmp

Input:

 All:
 desc    #elem   min     max     mean    stdev   sum
 PA14    1843    20      127984  3350    10777   6174765
 PAO1    1961    20      46414   3065    5251    6010800   # 1655 aligned   
                                                           # 1312 CONTAINED in PA14 contigs
                                                           # 343 aligned but not CONTAINED in PA14 ctgs; 36 have 2+ alignments
                                                           # 306 not aligned
 all     3804    20      127984  3203    8396    12185565

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 PA14    542     200     127984  11245   17521   6094904
 PAO1    1166    203     46414   5108    6007    5955960
 all     1708    200     127984  7056    11405   12050864

Output:

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     1585    20      159598  3896    15496   6174745
 singl   266     20      26709   1900    4486    505395      # all singletons are PAO1
 all     1851    20      159598  3609    14456   6680140

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     356     200     159598  17137   29063   6100717
 singl   75      215     26709   6631    6361    497316
 all     431     200     159598  15309   26837   6598033

minimus2 on the above contigs and singletons:

Output:

 All:
 desc    #elem   min     max     mean    stdev   sum
 ctg     246     30      266798  18632   44344   4583528
 singl   1138    19      126780  1569    9194    1785442  # 161 are original PAO1 contigs
 all     1384    19      266798  4602    21458   6368970

 200bp+:
 desc    #elem   min     max     mean    stdev   sum
 ctg     139     202     266798  32900   54951   4573157
 singl   152     201     126780  11405   22894   1733533 # 23 are original PAO1 contigs
 all     291     201     266798  21672   42726   6306690

AMOScmp2 on the AMOScmp PA14 contigs & PA14 genes *

INPUT:

                 #elem   min     max     mean    median  n50     sum
 AMOScmp         1843    20      127984  3350    63      30995   6174765
 AMOScmp(200bp+) 542     200     127984  11245   3263    31106   6094904
 genes           5892    72      15639   992     861     1206    5846175

OUTPUT:

         #elem   min     max     mean    median  n50     sum
 all     650     20      187557  9954    846     60335   6470391
 200bp+  462     201     187557  13977   1578    60335   6457375
 0 singletons !!!

Contigs that contain at least 1 AMOScmp contig & 1 gene: 356

                 #elem   min     max     mean    median  n50     sum       "singletons"
 ctg-all         356**   111     187557  17861   2692    63102   6358519
 ctg-200bp+      353     217     187557  18012   2723    63102   6358095
 ctg-AMOS        356     1       85      5       2       9       1653       190  
 ctg-genes       356     1       160     16      2       58      5781       111  **  that many genes from PA14 seem to miss from PAb1

Contigs that contain at least 1 AMOScmp contig(>100bp) & 1 gene: 218

                 #elem   min     max     mean    median  n50     sum
 ctg-all         218**   158     187557  28371   13008   66709   6184790
 ctg-200bp+      217     305     187557  28501   13008   66709   6184632
 ctg-AMOS        218     1       19      3       2       4       685
 ctg-genes       218     1       160     26      12      59      5632

Contigs that contain at least 1 AMOScmp contig(>200bp) & 1 gene: 206

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-gene/AMOScmp2

AMOScmp2 on the AMOScmp PA14 contigs , PA14 genes * & velvet contigs

INPUT:

                 #elem   min     max     mean    median  n50     sum
 AMOScmp         1843    20      127984  3350    63      30995   6174765
 AMOScmp(200bp+) 542     200     127984  11245   3263    31106   6094904
 genes           5892    72      15639   992     861     1206    5846175
 velvet          10684   45      16239   640     363     1184    6841458

OUTPUT:

                 #elem   min     max     mean    median  n50     sum
 all             552     20      337673  11718   663     133802  6468490
 200bp+          377     201     337673  17126   1155    133802  6456628
 singl           1006    45      16239   869     306     2308    874312       # all but 2 are velvet ctgs; 
                                                                              # 232 velvet singletons have annotated alignments to the AMOScmp ctgs

About half of the contigs end up in genes

Contigs that contain at least 1 gene & (1 AMOS ctg >=100bp or 1 velvet ctg): 135

                 #elem   min     max     mean    median  n50     sum
 ctg-all         135**   163     337673  45819   14232   141346  6185574
 ctg-200bp+      134     305     337673  46160   16317   141346  6185411
 ctg-AMOS100-vel 135     1       508     77      28      215     10356
 ctg-genes       135     1       287     42      12      126     5631

Contigs that contain at least 1 gene & 1 velvet ctg: 104 Contigs that contain at least 1 gene , 1 velvet ctg & 1 AMOScmp ctg : 104 (same number)

                 #elem   min     max     mean    median  n50     sum
 ctg-all         104**   163     337673  58949   32961   141346  6130650
 ctg-200bp+      103     747     337673  59519   32961   141346  6130487
 ctg-AMOS                                                        1085
 ctg-velvet                                                      9667
 ctg-genes                                                       5572

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-velvet-gene/AMOScmp2

AMOScmp on velvet & edena contigs (PA14 ref)

Input:

               #elem   min     max     mean    median  n50     sum
 velvet        10684   45      16239   640     363     1184    6841458
 edena         11180   100     11300   552     356     837     6175460
 all           21864   45      16239   595     359     992     13016918

PA14 0cvg regions (nucmer -c 40; delta-filter -q):

              #elem   min     max     mean    median  n50     sum
 velvet        2388    1       31832   219     24      6434    522048
 edena         8809    1       34854   119     38      438     1050135
 all           2129    1       31832   234     22      7529    498077

Output:

 All:
       #elem   min     max     mean    median  n50     sum
 ctg   2782    45      27957   2141    1207    4246    5955924
 singl 1874    45      16239   839     343     1972    1572938
 all   4656    45      27957   1617    727     3749    7528862
 breaks 492
 
 200 bp+:
       #elem   min     max     mean    median  n50     sum
 ctg   2424    200     27957   2439    1505    4288    5912995
 singl 1233    200     16239   1217    677     2121    1500822
 all   3657    200     27957   2027    1122    3816    7413817
 breaks 478

AMOScmp on velvet contigs (PA14 gene ref) *

Contig stats:

                             #elem   min     max     mean    median  n50     sum
 nucmer (-l 20,-c 65)        3137    65      3615    490     370     705     1535656
 promer (-l 6, -c 10)        3622    61      6125    623     448     922     2255378

Only the contigs CONTAINED in the genes get assembled; BEGIN/END are discarded

Locations:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/genes/AMOScmp-nucmer
 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/genes/AMOScmp-promer

De-novo

minimus on velvet contigs (OVL=40)

All:

       #elem   min     max     mean    median  n50     sum
 ctg   179     45      1495    75      45      45      13336
 singl 10324   45      16239   661     384     1188    6819550
 all   10503   45      16239   651     374     1185    6832886
 breaks 257

200 bp+:

       #elem   min     max     mean    median  n50     sum
 ctg   6       227     1495    836     1129    1495    5015
 singl 7375    200     16239   877     586     1252    6469371
 all   7381    200     16239   877     586     1252    6474386
 breaks 251

minimus on velvet contigs (OVL=20)

All:

       #elem   min     max     mean    median  n50     sum
 ctg   1812    47      23446   1621    1164    2366    2936650
 singl 6135    45      15342   626     357     1136    3839779
 all   7947    45      23446   853     471     1574    6776429
 breaks(PA14) 299

200 bp+:

       #elem   min     max     mean    median  n50     sum
 ctg   1746    201     23446   1677    1211    2374    2927255
 singl 4246    200     15342   851     568     1216    3615358
 all   5992    200     23446   1092    702     1651    6542613
 breaks(PA14) 293

minimus1 on velvet contigs (OVL=20)

All:

       #elem   min     max     mean    median  n50     sum
 ctg   2046    73      23446   1628    1157    2396    3330276
 singl 5658    45      15342   609     344     1144    3445004
 all   7704    45      23446   879     489     1660    6775280
 breaks(PA14) 568

200 bp+:

       #elem   min     max     mean    median  n50     sum
 ctg   1978    201     23446   1678    1205    2396    3320032
 singl 3786    200     15342   857     575     1222    3243138
 all   5764    200     23446   1139    737     1729    6563170
 breaks(PA14) 559

If OVL is dropped to 16, fewer contigs get merged => fewer breakpoints

minimus2 on velvet & edena ctgs (OVL=40)

Input: velvet & edene ctgs

          #elem   min     max     mean    median  n50     sum       breaks
 velvet   10684   45      16239   640     363     1184    6841458   240 
 edena    11180   100     11300   552     356     837     6175460   226
 v+e      21864   45      16239   595     359     992     13016918  
 v+e(200+)15698   200     16239   779     528     1069    12233635

Output (ctg+singletons):

          #elem   min     max     mean    median  n50     sum        breaks
 all      8835    45      16296   774     424     1448    6841222    253      
 200bp+   6446    200     16296   1018    663     1514    6563221

Comment: for OVL=20 we get 475 breaks (too many?)

minimus2 on the AMOScmp PA14 & PAO1 contigs

Use PA14-PA14 (5 bp, adj ctg ends) , PA14-PAO1 (20 bp, contig ends) , PAO1-PAO1 (5 bp, adj ctg ends) overlaps

All:

       #elem   min     max     mean    median  n50     sum
 ctg   441     43      223777  13972   294     86090   6161798
 singl 1227    20      18974   237     42      6391    290505
 all   1668    20      223777  3868    52      78110   6452303
 breaks(PA14) 130

200 bp+:

       #elem   min     max     mean    median  n50     sum
 ctg   245     202     223777  25080   5100    86090   6144514
 singl 107     210     18974   2180    547     8193    233281
 all   352     202     223777  18119   1678    78110   6377795
 breaks(PA14) 121

Location:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/AMOScmp-PA14-PAO1/minimus2

minimus2 on AMOScmp PA14 contigs and PA14 genes (OVL=20)

INPUT:

         #elem   min     max     mean    median  n50     sum
 AMOScmp 1843    20      127984  3350    63      30995   6174765
 genes   5892    72      15639   992     861     1206    5846175

OUTPUT:

                 #elem   min     max     mean    median  n50     sum       singl
 ctg-all         278**   158     170499  22136   10332   52297   6153935
 ctg-200bp+      277     315     170499  22216   10332   52297   6153777
 ctg-AMOS        278     1       23      2       1       3       604       1239
 ctg-genes       278     1       157     20      10      49      5582      310

minimus2 on velvet contigs and PA14 genes (OVL=20)

min %id = 94

INPUT

                 #elem   min     max     mean    median  n50     sum
 gene            5892    72      15639   992     861     1206    5846175
 velvet          10684   45      16239   640     363     1184    6841458

OUTPUT (nucmer,OVL=20)

                 #elem   min     max     mean    median  n50     sum  
 ctg             923**   181     66497   6728    4737    10446   6210244
 singl-gene      352     93      7347    979     777     1272    344478
 singl-velvet    1144    45      16239   580     214     1603    663940
 all             2419    45      66497   2984    901     9049    7218662

Best alignments:

 show-coords Pa.filter-1.delta | ~/bin/AMOS/mummerAnnotate.pl -all | awk '{print $20}' | count.pl
 [CONTAINED]     5306
 [END]           3482
 [BEGIN]         3453
 [CONTAINS]      1977  # ??? why only 1977 velvet ctgs contain full genes; genes are not repetitive; even worse for edena 
 [PARTIAL]       168
 [IDENTITY]      1
 Total           14387

0cvg regions in the genes

 #genes      #elem   min     max     mean    median  n50     sum
 1391        1940    1       1762    36      18      67      69336

=>4166 out of 5557 assembled genes are fully covered by velvet ctgs

Location: /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/genes-PA14/minimus2-nucmer.OVL20

minimus2 on velvet contigs and Ps genes (OVL20)

min %id: = 74 (nucmer default)
avg %id: = 84

 min%id   overlaps
 74       1162      # all considered
 85       415
 90       110

INPUT

         #elem   min     max     mean    median  n50     sum
 gene    5478    75      18825   996     852     1227    5456880
 velvet  10684   45      16239   640     363     1184    6841458

OUTPUT

               #elem   min     max     mean    median  n50     sum        "singl"
 ctg           652**   323     18825   2150    1768    2451    1401627                 !!! 1.4M of genome gets covered
 ctg-gene      757     96      18825   1271    1122    1455    962199      4721
 ctg-velvet    998     45      11120   938     558     1767    935636      9686

minimus2 on velvet contigs and Ps genes (OVL36: promer -l 6 -c 12)

 min%id   overlaps  genes  velvet 
 22       8496      2985   4827        #default
 50       7209      2813   4495
 60       5948      2571   4087
 70       4361      2122   3269
 74       3700      1880   2827
 80       2615      1407   2036
 85       1714      956    1311
 90       976       516    659

min%id=74 =>

 * 3700 overlaps 
 * 2169 Ps genes
 * 3206 velvet contigs

make-consensus error

                #elem   min     max     mean    median  n50     sum      "singl"
ctg             929**   206     16291   2012    1669    2335    1869539            !!! 1.8M of genome gets covered
ctg-gene        991     111     9522    1118    990     1272    1107648   4487     # 2169-991 genes had overlaps but were not assembled
ctg-velvet      1032    45      16239   1199    910     1742    1237453   9652     # 3206-1032  velvet contigs had overlaps but were not assembled

!!! Slightly more genes/velvet contigs get assembled with promer than nucmer

minimusN on PA14,PA2192,PAC3719,PACS2,PAO1 genes

Output:

 strains        #elem   min     max     mean    median  n50     sum
 2+             5702    72      15450   1015    870     1227    5785616
 3+             5469    153     15450   1007    876     1212    5505737
 4+             5292    153     13067   1008    883     1209    5334219
 5              4962    153     13067   1020    891     1212    5061120  (5021 in the article)

Total and strain specific genes:

 strain          total   singletons
 PA14            5892    519
 PA2192          6157    387   (507 in the article)
 PAC3719         5650    132   (79  in the article)
 PACS2           5676    127
 PAO1            5568    68
 total           28943   1233

all vs all genes: from the aligned ones ~ 1556 are CONTAINED|BEGIN|END (426 are PA14)

Final assembly (v3)

Merged some of the 120 AMOScmp contigs from v2 using velvet contigs. Grouped overlapping contigs and used minimus to merge them.

 Contig length stats:
 
  desc            #ctgs   min     max     mean              stdev           sum
  final3          512     200     512638  13099.41          53526.8         6706902
  AMOScmp3        76      212     512638  82763.22          117199.28       6290005
  novo3           436     200     10493   956.18            1336.44         416897
 
  final3-10KB+    49      10493   512638  127086.12         125757.98       6227220

Comment:

 ~ 53-74 (out of 436) the velvet denovo that contigs align to PA14.
 
 Length statistics:
                    #ctgs   min     max     mean    stdev   sum
   CONTAINED        53      204     7228    514.81  978.25  27285
   ALL              74      204     7228    622.85  926.13  46091

Locations:

 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/final3/
 /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/2008_0311_velvet-24/nucmer/merge/

Simulations

First ~ 100KB from PA14, PAO1

All PA14, PAO1

Uniform read error : 1%

Repetitive regions in PA14 (100 bp found my repeat-match)

         #rep.   min     max     mean    median  n50     sum
 PA14    282     100     5571    340     171     520     95919

AMOScmp

 PA14=sample
 PAO1=reference

Conclusions: (results are below)

more 33bp reads align with soap
more 50bp reads align with nucmer
slightly more reads align to the reference if the errors have exponential than uniform distribution ; no much impact on the assembly

AMOScmp 50bp read assemblies have fewer ctgs than the 33bp read assemblies

33bp reads: more errors(snps & breaks) in the AMOScmp soap than nucmer assembly
50bp reads: more errors(snps & breaks) in the AMOScmp nucmer than soap assembly

more snps & breaks in the non-unique regions => many introduced by repeats
more indels (than substitutions) in the unique regions
more deletions from the reference
the indels are actually longer than 1bp; usually 5-20bp; occur where there are larger differences in the 2 genomes; they seem to be assembly errors

there are fewer longer insertions in the qry than in the ref (for all PA14 33bp, PAO1 33bp, PA14 50bp)
decreasing the make-consensus alignment wiggle -w to 1 makes no much difference than -w 2

if we incease the min% id for the alignments, fewer reads will align , more 1X cvg regions => more snps; these snps are mostly substitutions, not indels

Read alignments to PAO1

Uniform read error : 1%

Alignment programs: nucmer & soap

33bp: min alignment=28bp both nucmer & soap
50bp: min alignment=45bp both nucmer & soap

       error    total    aligned(nucmer)   aligned(soap)       aligned(nucmer&soap)
 33bp  uniform  7528200  6631195(88.08%)   6851948(91.01%)     6602210(87.69%)
 33bp  exp.     7528200  6713822           6879256
 50bp  uniform  4968612  4524844(91.06%)   4497468(90.51%)

0cvg regions in PAO1:

                 #elem   min     max     mean    median  n50     sum
 33bp-nucmer     969     1       13696   245     22      1470    237643
 33bp-soap       543     1       18331   434     26      3262    235621
 50bp-nucmer     326     1       18663   716     96      3178    233556
 50bp-soap       453     1       19256   546     29      4840    247205  !!! largest values

PA14 33 & 50bp simulated reads

Uniform read error : 1%

Legend:

uniq:contigs that align uniquely to the reference
snps-I: substitutions only

33bp

                                #ctg    min     max     mean    median  n50     sum     snps    snps-I   snps/1KB  snps-I/1KB    breaks
 AMOScmp-PAO1.nucmer            1298    33      60570   4641    1342    13211   6023652 1655    915      0.274     0.151         57   *
 AMOScmp-PAO1.nucmer(uniq)      851     40      46537   3482    1167    9474    2963269 425     57       0.143     0.019         2
  
 AMOScmp-PAO1.soap              774     33      214591  7791    527     30892   6030156 2037    1209     0.337     0.200         85
 AMOScmp-PAO1.soap(uniq)        395     40      42671   3496    194     14876   1380824 189     29       0.136     0.021         2

50bp

                                #ctg    min     max     mean    median  n50     sum     snps    snps-I   snps/1KB  snps-I/1KB    breaks
 AMOScmp-PAO1.nucmer            366     50      214838  16466   1998    59036   6026608 1966    1190     0.326     0.197         101   *
 AMOScmp-PAO1.nucmer(uniq)      158     50      91686   4954    233     23389   782667  77      32       0.098     0.040         0
 
 AMOScmp-PAO1.soap              634     50      144624  9492    2302    28330   6018087 1927    1112     0.320     0.184         60 
 AMOScmp-PAO1.soap(uniq)        349     50      62962   4705    1149    15406   1641990 133     13       0.080     0.007         0      !!! fewest snps

                snps    snps-I  breaks
 *common        838     441     29

 PA14: 6537648

                                                                   all               common
                                                        Coverage   snps  snps-I      snps   snps-I  
 --------------------------------------------------------------------------------------------------
 sim.PA14/AMOScmp-PAO1/PA14-ctgs.filter-q               6003782    1655  915         1281   761
 sim.PA14.50bp//AMOScmp-PAO1.redo//PA14-ctgs.filter-q   6017813    1966  1190        1198   734
 common(last 2)                                         5996786

The larger number of snps in the 50bp assemblies is an artifact of gap closing in 33bp assembly
The number of snps in the 33bp covered regions goes slightly down if read length is increased

PA14 33bp simulated reads

Uniform read error : 1%

Legend:

uniq : reads that seem to be unique have been filtered out

                         #elem   min     max     mean    median  n50     sum     snps    snps-I   breaks
 AMOScmp-PAO1            1298    33      60570   4641    1342    13211   6023652 1655    915      57* 
 AMOScmp-PAO1.uniq       1383    33      60570   4297    1026    13192   5943384 1207    57       48
  
 AMOScmp-PAO1.soap       774     33      214591  7791    527     30892   6030156 2037    1209     85*
 
 velvet                  4145    45      25746   1577    529     4317    6536518 151     148      0
 
 minimus2                582     33      275737  11306   261     91069   6579905 1610    1056     69   # AMOS & velvet ctg phred score=30
 minimus2.qual           582     33      275730  11305   261     91069   6579772 1240    928      67   # AMOS ctg phred score=20; velvet ctg phred score=30
 minimus2.pre            575     23      275737  11371   238     84812   6538289 1569    1131     34   # some of the AMOS ctg are broken by velvet ctg

 AMOScmp-PACS2           1713    33      55026   3521    423     13041   6031518 1344    831     56    # more ctgs than in the AMOScmp-PAO1 assembly

          snps    snps-I   breaks
 *common  1039    598      34

PA14 33bp simulated reads

Exponential read error : for a 33 bp read avg=~1%

Base 0 error: 1/(1.2**48-1) Multiplied by 1.2 for next base ... Base 48+ error: 1

 $ ~/bin/mutateSeq.pl -exp 48 -test
 0     0.000158259181534141
 1     0.000189911017840969
 2     0.000227893221409163
 ...
 30    0.0375669811375429
 31    0.0450803773650515
 32    0.0540964528380618
 ----------------------------------
 total 0.323

 47    0.833465215984612
 48    1.00015825918153

AMOScmp alignment

                                 #ctgs   min     max     mean    median  n50     sum     snps    snps-I  breaks  #reads-aligned
 AMOScmp-PAO1.exp                1396    33      55370   4315    1188    12277   6023140 1773    877     56      6713822(89.18%)
 AMOScmp-PAO1.exp.uniq           908     40      46548   3394    1054    9533    3081827 481     22      3       .

Annotation

( 5769 genes /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Annotation/merge/gene_list.ptt ) OLD ?
5602 genes /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Annotation/final/gene_list.ptt
512 contigs /fs/szasmg2/Bacteria/Pseudomonas_aeruginosa/Assembly/final3/PAb1.final3.fasta

Submission

PAb1 genome project
Taxonomy
Genbank accession: ABKZ01000000
SRA001120
PLoS aricle
- Our final assembly contains one large scaffold with 76 contigs whose total length is 6,290,005 bp, with the longest contig at 512,638 bp. The 436 unplaced contigs, which should fit into the remaining gaps, represent sequence that is unique to PAb1 ... The final assembly thus consists of 512 contigs covering 6,706,902 bp ...
- Our annotation of the PAb1 genome identified 5,602 protein-coding genes, ...

Pseudomonas aeruginosa: Difference between revisions

Latest revision as of 15:07, 14 October 2010

Data

CBCB

NCBI

Publications

Read trimming

No trimming

Quality trimming

Alignment based trimming (nucmer)

Read mapping

PA14 reference

Pm reference

Ps reference

Sequence Assemblies

Comparative

PA14 ref assembly

PA14 ref assembly (2) *

PA14 ref assembly (used blat instead of nucmer for read alignment)

PA14 ref assembly (used soap instead of nucmer for read alignment)

PAO1 ref assembly

PAO1 ref assembly (2)

PA14 & PAO1 merge

Pa-b1 AMOScmp reference assembly

Use ALL strains as reference

PA14 Sanger maq on all reads

PA14 gene ref assembly

Comparative-ABBA

ABBA

PA14 gene ref assembly

PAO1 gene ref assembly

De-novo

SSAKE

SSAKE (redo)

Velvet *

Edena

SSAKE/Edena/Velvet summary

Euler-SR

Gene assemblies

Contig Assemblies

Comparative

AMOScmp on the velvet contigs (PA14 ref) *

AMOScmp on the velvet contigs & velvet singletons (PA14 ref) *

AMOScmp on the velvet contigs (PAO1 ref)

AMOScmp on the velvet contigs (PA14 & PA2192 ref)

AMOScmp on the velvet contigs (PA14, PA2192 & PACS2 ref)

AMOScmp on the velvet contigs (PA14, PA2192, PACS2 & PAO1 ref)

AMOScmp/minimus on the AMOScmp (PA14 ref) & velvet "unique" contigs

AMOScmp/minimus on the AMOScmp (PA14 ref) & velvet contigs *

AMOScmp/minimus on the AMOScmp PA14 & PAO1 contigs

AMOScmp2 on the AMOScmp PA14 contigs & PA14 genes *

AMOScmp2 on the AMOScmp PA14 contigs , PA14 genes * & velvet contigs

AMOScmp on velvet & edena contigs (PA14 ref)

AMOScmp on velvet contigs (PA14 gene ref) *

De-novo

minimus on velvet contigs (OVL=40)

minimus on velvet contigs (OVL=20)

minimus1 on velvet contigs (OVL=20)

minimus2 on velvet & edena ctgs (OVL=40)

minimus2 on the AMOScmp PA14 & PAO1 contigs

minimus2 on AMOScmp PA14 contigs and PA14 genes (OVL=20)

minimus2 on velvet contigs and PA14 genes (OVL=20)

minimus2 on velvet contigs and Ps genes (OVL20)

minimus2 on velvet contigs and Ps genes (OVL36: promer -l 6 -c 12)

minimusN on PA14,PA2192,PAC3719,PACS2,PAO1 genes

Final assembly (v3)

Simulations

First ~ 100KB from PA14, PAO1

All PA14, PAO1

Read alignments to PAO1

PA14 33 & 50bp simulated reads

PA14 33bp simulated reads

PA14 33bp simulated reads

Annotation

Submission

Navigation menu

Search