Docstoc

Paul Greenfield

Document Sample
Paul Greenfield Powered By Docstoc
					Analysing Mountains of K-mers

Paul Greenfield
CMIS Transformational Biology
11/06/2009
Claims About K-mer Uniqueness

• K-mer space is big (for reasonable ‘k’)
         • 10^15 25-mers; 10^30 50-mers; even 10^12 20-mers
         • We’ll mostly talk about 25-mers in bacteria
• Collisions are rare in 25-mer space
         • ‘Unrelated’ organisms share few/none 25-mers
         • Shared 25-mers seem to often say something about
           conservation or gene-swapping or …
         • Shared 20-mers still say something but more noise
• Shared 25-mers…
         •   Some highly conserved fundamental genes/RNA
         •   Otherwise unrelated organisms don’t share 25-mers
         •   Closely related organisms share many 25-mers
         •   Can shared 25-mers be used as a relatedness metric?

CSIRO.
Exact K-mer Matching Is Fast

• Exact k-mer matching can be very fast
         • Indexed database tables
            • 100,000 lookups/second
            • Can lookup against thousands of target genomes at once
            • Need to structure DB & lookup code very carefully
         • In-memory hash tables
            • Cache parts of database locally
            • Millions of lookups/second
• Fast enough to do large-scale surveys
         • Compare all bacteria to all bacteria
            • All k-mers in all genomes compared to all k-mers in all genomes
• Fast enough to generate and compare variants
         • Handle errors & SNPs by matching permutated k-mers

CSIRO.
Applications (and Backing Up Claims)

• Cross-comparing genomes of all sequenced bacteria
         • Is k-mer sharing really uncommon?
         • Does k-mer sharing reflect accepted taxonomy?
         • Is k-mer sharing random drift or something else?
• Asking questions about relatedness
         •   How do 2 organisms differ? How can they be distinguished?
         •   Do sets of shared k-mers indicate taxonomic groups?
         •   Do sets of shared k-mers indicate functions?
         •   Applications in metagenomics, homology & taxonomy?
• Healing sequence errors through consensus
         • Fast matching of permutations
• Human gene methylation studies
         • Fast matching of permutations (SNPs, errors, partial methyl)
CSIRO.
Comparing All to All

• Can k-mers from this bacteria be found elsewhere?
         • Generate (overlapping) k-mer ‘tiles’ from a source organism
         • Compare these up against the genomes of all other organisms
         • Count various interesting things when we find a match
• Practicalities
         • About 2.5 billion bases in the genomes
            • Double this because matches could be on other strand
         • 2.5 billion source k-mers & 5 billion target k-mers
            • A mere 1.25E19 comparisons??
• Actualities
         • Indexes allow for 1-many comparisons
            • Retrieve all places the k-mer can be found in 1 lookup
         • 18 hours for database version, 5 hours for hash-table
         • Need pre-computed indexes, sorting, partitioning
CSIRO.
A Few Results




         Part of the Gammaproteobacteria
         (Enterobacteriales in the middle)
CSIRO.
Some Comments

• Colour-coded according to degree of sharing
         • Falls off pretty quickly
            • Oranges & yellows are 0.5% down to 0.01%
• Some taxonomic structures show up
         • Holes indicate possible problems
         • Lots of background low-level matching
• This picture includes ribosomal RNA
         • Includes 16S sequence




CSIRO.
Now Without rRNA




Now considerable ‘white space’
And distinct ‘shadows’ of other taxonomic groups
CSIRO.
Shadows




CSIRO.
rRNA Only??




CSIRO.
Shadows and Uniqueness

• Most of the BxB matrix is uncoloured
         • Less than 0.001% shared 25-mers
• Shadows (higher sharing)
         •   Only significant sharing for close relatives
         •   Sharing highly non-random
         •   Often shared non-coding (RNAs? Transposons?)
         •   Genes with shared k-mers often have same annotations
• Accidental k-mer collisions are rare
• Shared k-mers seems to mean something
         • Even when only a few are shared…
         • Highly-conserved domains?? (TBA)



CSIRO.
                  Exploring the Shadows

                  • Look very regular…
                           • Drift from common ancestors?
                           • Horizontal gene transfer?
                  • Random smear or conservation?
                           • Answered by looking at gene-level counts
                           • Count matches for [source gene, target gene]
Deinococcus_geothermalis_DSM_1
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    7   0.2% Mycobacterium_avium_104          1    4118   118467332 cobaltochelatase   3576
Deinococcus_geothermalis_DSM_1                                                                        Mycobacterium_avium_paratuber
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    9   0.3% culosis                          1    3197    41407903 cobaltochelatase   3576
Deinococcus_geothermalis_DSM_1
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    8   0.2% Mycobacterium_bovis              1    3660    31793245 cobaltochelatase   3585
Deinococcus_geothermalis_DSM_1                                                                        Mycobacterium_bovis_BCG_Paste
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    8   0.2% ur_1173P2                        1    3645   121637947 cobaltochelatase   3585
Deinococcus_geothermalis_DSM_1                                                                        Mycobacterium_smegmatis_MC2
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    8   0.2% _155                             1    6524   118470146 cobaltochelatase   3636
Deinococcus_geothermalis_DSM_1                                                                        Mycobacterium_tuberculosis_CDC
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    8   0.2% 1551                             1    3785    15841550 cobaltochelatase   3588
Deinococcus_geothermalis_DSM_1
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    8   0.2% Mycobacterium_tuberculosis_F11   1    3683   148823277 cobaltochelatase   3585
Deinococcus_geothermalis_DSM_1                                                                        Mycobacterium_tuberculosis_H37
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    8   0.2% Ra                               1    3739   148661875 cobaltochelatase   3585
Deinococcus_geothermalis_DSM_1                                                                        Mycobacterium_tuberculosis_H37
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    8   0.2% Rv                               1    3697    15609199 cobaltochelatase   3585
Deinococcus_geothermalis_DSM_1                                                                        Mycobacterium_vanbaalenii_PYR-
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    2   0.1% 1                                1    5831   120404398 cobaltochelatase   3618
Deinococcus_geothermalis_DSM_1
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335   14   0.4% Nocardia_farcinica_IFM10152      1    5583    54025106 cobaltochelatase   3627
Deinococcus_geothermalis_DSM_1                                                                        Saccharopolyspora_erythraea_NR
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335    2   0.1% RL_2338                          1   10672   134102401 cobaltochelatase   3606
Deinococcus_geothermalis_DSM_1
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335   15   0.4% Streptomyces_avermitilis         1   11778    29832957 cobaltochelatase   3654
Deinococcus_geothermalis_DSM_1
1300                             1   200   94971821 cobaltochelatase. CobN subunit   4335   16   0.4% Streptomyces_coelicolor          1    3324    21220338 cobaltochelatase   3654

                  CSIRO.
            Intersections in more detail
Campylobacter_hominis_ATCC_BAA-381\NC_009714 (1711125) compared to Helicobacter_pylori_Shi470
      geneNo               geneID      strand        length                               product   geneNo             geneID   strand   length                         product   POmers   RCmers   mers
           96 non-coding                    1         4408 rRNA                                      1754 non-coding                1     4464 rRNA                                    0     482    482
           96 non-coding                    1         4408 rRNA                                      2393 non-coding                1     4581 rRNA                                    0     482    482
          121 non-coding                    1         3659 rRNA                                      1869 non-coding                1      529 non-coding                             55        0    55
          121 non-coding                    1         3659 rRNA                                      1967 non-coding                1     2214 rRNA                                    0     315    315
          121 non-coding                    1         3659 rRNA                                      2458 non-coding                1     2062 rRNA                                    0     315    315
          576 non-coding                    1         5098 rRNA                                      1754 non-coding                1     4464 rRNA                                    0     482    482
          576 non-coding                    1         5098 rRNA                                      2393 non-coding                1     4581 rRNA                                    0     482    482
          710 non-coding                    1         5618 rRNA                                      1869 non-coding                1      529 non-coding                             55        0    55
          710 non-coding                    1         5618 rRNA                                      1967 non-coding                1     2214 rRNA                                    0     315    315
          710 non-coding                    1         5618 rRNA                                      2458 non-coding                1     2062 rRNA                                    0     315    315
          712 non-coding                    1         4520 rRNA                                      1754 non-coding                1     4464 rRNA                                    0     482    482
          712 non-coding                    1         4520 rRNA                                      2393 non-coding                1     4581 rRNA                                    0     482    482
                                                           putative cell division protease FtsH-
          756        154149469              1         1932 like protein                               658        188527188          -1    1899 cell division protein (ftsH)            0        2      2

                                                                                                                                               Holliday junction DNA helicase
         1157        154148614              1         1023 Holliday junction DNA helicase B           673        188527198          -1    1011 B                                       0        4      4
         1199 non-coding                    1          382 non-coding                                 400 non-coding                1      239 non-coding                             20        0    20

                                                           enoyl-(acyl carrier protein)                                                        enoyl-(acyl carrier protein)
         1226        154148152              1          822 reductase                                  332        188527000          1      828 reductase                               1        0      1
         1313        154147926             -1         1059 elongation factor Ts                      2658        188528329          -1    1068 elongation factor Ts                    2        0      2
         1439 non-coding                    1          127 non-coding                                1674 non-coding                1      411 non-coding                              1        0      1
                                                           transcription termination factor                                                    transcription termination factor
         2107        154149028              1         1296 Rho                                       1393        188527605          1     1317 Rho                                     5        0      5
         2220 non-coding                    1         1327 non-coding                                 400 non-coding                1      239 non-coding                             20        0    20
         2220 non-coding                    1         1327 non-coding                                1777 non-coding                1      919 non-coding                             13        0    13
         2367 non-coding                    1         3207 rRNA                                      1869 non-coding                1      529 non-coding                              0      55     55
         2367 non-coding                    1         3207 rRNA                                      1967 non-coding                1     2214 rRNA                                  315        0   315
         2367 non-coding                    1         3207 rRNA                                      2458 non-coding                1     2062 rRNA                                  315        0   315
         2449 non-coding                    1          508 non-coding                                2042 non-coding                1      326 non-coding                              9        0      9
         2492 non-coding                    1         1829 non-coding                                  11 non-coding                1     3200 non-coding                              6        0      6
         2593 non-coding                    1         1573 non-coding                                1674 non-coding                1      411 non-coding                              1        0      1
         2593 non-coding                    1         1573 non-coding                                1676 non-coding                1     2349 non-coding                              5        0      5
         2624 non-coding                    1          725 non-coding                                2090 non-coding                1      722 non-coding                             56        0    56


         2758        154148859             -1         1809 GTP-binding protein TypA/BipA              786        188527265          1     1800 hypothetical protein                    0        4      4




            CSIRO.
Asking More Questions…

• Inspiration from Jim Gray’s work on astronomy
         • Can we answer scientific questions using SQL queries?
• Load intersections into database
         • Genome-onto-gene & gene-onto-gene
• Run interesting SQL queries
         •   What genes are shared across families?
         •   How can we distinguish Citrobacter from its relatives?
         •   Do the shadows represent shared genes?
         •   Can we find functions for hypothetical genes?
         •   Find genes unique to an organism
         •   ….



CSIRO.
Sample queries




CSIRO.
            Finding Citrobacter Unique Mers…

                                                                                                                   Salmonella_typhimurium_LT2
            Sal, Esch, Yers onto Citro Gap 1                                                                       Salmonella_typhimurium_LT2
                                                                                                                   Escherichia_coli_536
                                                                                                                   Escherichia_coli_536
                                                                                                                   Yersinia_pestis_Angola
    1
  0.9
  0.8
  0.7
  0.6
  0.5
  0.4
  0.3
  0.2
  0.1
    0
   380000                      390000                    400000                410000                  420000                     430000           440000


                                        Source
                                        Sequen       MersHitsN TargetSeq TargetGen GeneLen
SpeciesID                               ceNo MerHits onrRNA uenceNo eNo            gth     GeneStart GeneEnd gap    %hits Product
Agrobacterium_tumefaciens_C58_Cereon         2      35        35        1      700    1896    412547 414442           1.8%hypothetical protein
Bradyrhizobium_BTAi1                         2       8         8        1      700    1896    412547 414442           0.4%hypothetical protein
Bradyrhizobium_ORS278                        1      11        11        1      700    1896    412547 414442           0.6%hypothetical protein
Chromobacterium_violaceum                    1       2         2        1      700    1896    412547 414442           0.1%hypothetical protein
Hahella_chejuensis_KCTC_2396                 1       3         3        1      713    2544    419013 421556           0.1%cAMP phosphodiesterase
Parvibaculum_lavamentivorans_DS-1            1       2         2        1      703      681   414783 415463           0.3%hypothetical protein
Pseudomonas_putida_GB_1                      1       5         5        1      700    1896    412547 414442           0.3%hypothetical protein
Rhodopseudomonas_palustris_CGA009            1       5         5        1      700    1896    412547 414442           0.3%hypothetical protein
Rhodopseudomonas_palustris_TIE_1             1       8         8        1      700    1896    412547 414442           0.4%hypothetical protein
Thiomicrospira_crunogena_XCL-2               1       2         2        1      713    2544    419013 421556           0.1%cAMP phosphodiesterase
            CSIRO.
    Where are highly shared genes?

                                     Entero onto Citro 60%


                                                                                                       Enterobacter


                                                                                                       E.coli


                                                                                                       Klebsiella


                                                                                                       Salmonella


                                                                                                       Shigella


                                                                                                       Yersinia




0       500000   1000000   1500000   2000000   2500000   3000000   3500000   4000000   4500000   5000000

    CSIRO.
Yet Another Query




CSIRO.
More Applications of Uniqueness

• Map metagenomic samples
         • What families/genera/… are present?
            • Tile the sample , match and see what lights up…
         • What functions/genes are present???
            • Same process but see what ‘genes’ light up
• Allow for small mismatches
         • Do SNPs change the overall picture?
• Healing sequence errors
         • Improve read data without needing reference genomes
• Determining methylation status (CMHT/pHealth)
         • Mapping methylated reads to human reference
         • Partial methylation, sequencing errors, SNPs


CSIRO.
Methylation Study

• Region of interest cut biochemically
         • Reads can only start at known places
         • Only millions of possible targets (50-mers)
• Separate effects of partial methylation & errors/SNPs
         • In-silico fully-methylate reads (reduces uniqueness but…)
         • Match against fully-methylated targets
            • Could use healing to improve read quality first
            • Generate permutations to correct for sequencing errors/SNPs
            • Find reference target sequence
         • Find what bases are methylated
            • Compare reads with unmethylated target sequence
• All done using hash tables (and databases)
         • And relying on very fast exact matches

CSIRO.
Summary

• K-mer spaces are big (for reasonable ‘k’)
         • Almost all 25-mers are unique
         • Sharing 25-mers almost always means something
• Shared mers can be used as a closeness metric
         • Sharing matrix matches taxonomic structure quite well
         • Shared mers seems to mean shared genes
• SQL queries can be used to answer some questions
         • Extensive pre-processing to populate tables first
         • Queries over all bacteria come back in ‘seconds’
• Fast matching of unique k-mers
         • Error correction via ‘healing’
         • Methylation via in-silico full-methylation & matching


CSIRO.

				
DOCUMENT INFO