K-mer

For broader coverage of this topic, see n-gram.

The sequence ATGG has two 3-mers: ATG and TGG.

In bioinformatics, k-mers are subsequences of length $k$ contained within a biological sequence. Primarily used within the context of computational genomics and sequence analysis, in which k-mers are composed of nucleotides (i.e. A, T, G, and C), k-mers are capitalized upon to assemble DNA sequences, improve heterologous gene expression, identify species in metagenomic samples, and create attenuated vaccines. Usually, the term k-mer refers to all of a sequence's subsequences of length $k$, such that the sequence AGAT would have four monomers (A, G, A, and T), three 2-mers (AG, GA, AT), two 3-mers (AGA and GAT) and one 4-mer (AGAT). More generally, a sequence of length $L$ will have $L-k+1$ k-mers and $n^{k}$ total possible k-mers, where $n$ is number of possible monomers (e.g. four in the case of DNA).

Introduction

k-mers are simply length $k$ subsequences. For example, all the possible k-mers of a DNA sequence are shown below:

An example 8-mer spectrum for E. coli comparing 8-mers' frequency (i.e. multiplicities) with their number of occurrences.

k-mers for GTAGAGCTGT

k	k-mers
1	G, T, A, G, A, G, C, T, G, T
2	GT, TA, AG, GA, AG, GC, CT, TG, GT
3	GTA, TAG, AGA, GAG, AGC, GCT, CTG, TGT
4	GTAG, TAGA, AGAG, GAGC, AGCT, GCTG, CTGT
5	GTAGA, TAGAG, AGAGC, GAGCT, AGCTG, GCTGT
6	GTAGAG, TAGAGC, AGAGCT, GAGCTG, AGCTGT
7	GTAGAGC, TAGAGCT, AGAGCTG, GAGCTGT
8	GTAGAGCT, TAGAGCTG, AGAGCTGT
9	GTAGAGCTG, TAGAGCTGT
10	GTAGAGCTGT

A method of visualizing k-mers, the k-mer spectrum, shows the multiplicity of each k-mer in a sequence versus the number of k-mers with that multiplicity. The number of modes in a k-mer spectrum for a species's genome varies, with most species having a unimodal distribution. However, all mammals have a multimodal distribution. The number of modes within a k-mer spectrum can vary between regions of genomes as well: humans have unimodal k-mer spectra in 5' UTRs and exons but multimodal spectra in 3' UTRs and introns.

Forces Affecting DNA k-mer Frequency

The frequency of k-mer usage is affected by numerous forces, working at multiple levels, which are often in conflict. It is important to note that k-mers for higher values of k are affected by the forces affecting lower values of k as well. For example, if the 1-mer A does not occur in a sequence, none of the 2-mers containing A (AA, AT, AG, and AC) will occur either, thereby linking the effects of the different forces.

k = 1

When k = 1, there are four DNA k-mers, i.e., A, T, G, and C. At the molecular level, there are three hydrogen bonds between G and C, whereas there are only two between A and T. GC bonds, as a result of the extra hydrogen bond (and stronger stacking interactions), are more thermally stable than AT bonds. Mammals and birds have a higher ratio of Gs and Cs to As and Ts (GC-content), which led to the hypothesis that thermal stability was a driving factor of GC-content variation. However, while promising, this hypothesis did not hold up under scrutiny: analysis among a variety of prokaryotes showed no evidence of GC-content correlating with temperature as the thermal adaptation hypothesis would predict. Indeed, if natural selection were to be the driving force behind GC-content variation, that would require that single nucleotide changes, which are often silent, to alter the fitness of an organism.

Rather, current evidence suggests that GC‐biased gene conversion (gBGC) is a driving factor behind variation in GC content. gBGC is a process that occurs during recombination which replaces Gs and Cs with As and Ts. This process, though distinct from natural selection, can nevertheless exert selective pressure on DNA biased towards GC replacements being fixed in the genome. gBGC can therefore be seen as an "impostor" of natural selection. As would be expected, GC content is greater at sites experiencing greater recombination. Furthermore, organisms with higher rates of recombination exhibit higher GC content, in keeping with the gBGC hypothesis's predicted effects. Interestingly, gBGC does not appear to be limited to eukaryotes. Asexual organisms such as bacteria and archaea also experience recombination by means of gene conversion, a process of homologous sequence replacement resulting in multiple identical sequences throughout the genome. That recombination is able to drive up GC content in all domains of life suggests that gBGC is universally conserved. Whether gBGC is a (mostly) neutral byproduct of the molecular machinery of life or is itself under selection remains to be determined. The exact mechanism and evolutionary advantage or disadvantage of gBGC is currently unknown.

k = 2

Despite the comparatively large body of literature discussing GC-content biases, relatively little has been written about dinucleotide biases. What is known is that these dinucleotide biases are relatively constant throughout the genome, unlike GC-content, which, as seen above, can vary considerably. This is an important insight that must not be overlooked. If dinucleotide bias were subject to pressures resulting from translation, then we would expect to see differing patterns of dinucleotide bias in coding and noncoding regions driven by some dinucelotides' reduced translational efficiency. Because we do not, we can therefore infer that the forces modulating dinucleotide bias are independent of translation. Further evidence against translational pressures affecting dinucleotide bias is the fact that the dinucleotide biases of viruses, which rely heavily on translational efficiency, are shaped by their viral family more than by their hosts, whose translational machinery the viruses hijack.

Counter to gBGC's increasing GC-content is CG suppression, which reduces the frequency of CG 2-mers due to deamination of methylated CG dinucleotides, resulting in substitutions of CGs with TGs, thereby reducing the GC-content. This interaction highlights the interrelationship between the forces affecting k-mers for varying values of k.

One interesting fact about dinucleotide bias is that it can serve as a "distance" measurement between phylogenetically similar genomes. The genomes of pairs of organisms that are closely related share more similar dinucleotide biases than between pairs of more distantly related organisms.

k = 3

There are twenty natural amino acids that are used to build the proteins that DNA encodes. However, there are only four nucleotides. Therefore, there cannot be a one-to-one correspondence between nucleotides and amino acids. Similarly, there are 16 2-mers, which is also not enough to unambiguously represent every amino acid. However, there are 64 distinct 3-mers in DNA, which is enough to uniquely represent each amino acid. These non-overlapping 3-mers are called codons. While each codon only maps to one amino acid, each amino acid can be represented by multiple codons. Thus, the same amino acid sequence can have multiple DNA representations. Interestingly, each codon for an amino acid is not used in equal proportions. This is called codon-usage bias (CUB). As a parenthetical, when k = 3 we must distinguish between true 3-mer frequency and CUB. For example, the sequence ATGGCA has four k-mer words within it (ATG, TGG, GGC, and GCA) while only containing two codons (ATG and GCA). However, CUB is a major driving factor of 3-mer usage bias (accounting for up to ⅓ of it, since ⅓ of the k-mers in a coding region are codons) and will be the main focus of this section.

The exact cause of variation between the frequencies of various codons is not fully understood. It is known that codon preference is correlated with tRNA abundances, with codons matching more abundant tRNAs being correspondingly more frequent and that more highly expressed proteins exhibit greater CUB. This suggests that selection for translational efficiency or accuracy is the driving force behind CUB variation.

k = 4

Similar to the effect seen in dinucleotide bias, the tetranucleotide biases of phylogenetically similar organisms are more similar than between less closely related organisms. The exact cause of variation in tetranucleotide bias is not well understood, but it has been hypothesized to be the result of the maintenance of genetic stability at the molecular level.

Sequence assembly

Overview

In sequence assembly, k-mers are typically used during the construction of De Bruijn graphs. In order to create a De Bruijn Graph, the strings stored in each edge with length, $L$, must overlap another string in another edge by $L-1$ in order to create a vertex. Reads generated from next-generation sequencing will typically have different read lengths being generated. For example, reads by Illumina’s sequencing technology capture reads of 100-mers. However, the problem with the sequencing is that only small fractions out of all the possible 100-mers that are present in the genome are actually generated. This is due to read errors, but more importantly, just simple coverage holes that occur during sequencing. The problem is that these small fractions of the possible k-mers violate the key assumption of de Bruijn graphs that all the k-mer reads must overlap its adjoining k-mer in the genome by $k-1$ (which can’t occur when all the possible k-mers aren’t present). The solution to this problem is to break these k-mer sized reads into smaller k-mers, such that the resulting smaller k-mers will represent all the possible k-mers of that smaller size that are present in the genome. Furthermore, splitting the k-mers into smaller sizes also helps alleviate the problem of different initial read lengths. An example of the solution of splitting the reads into smaller k-mers is shown in figure 1. In this example the 5 reads do not account for all the possible 7-mers of the genome, and as such, a de Bruijn graph cannot be created. But when they are split into 4-mers, the resultant subsequences are enough to reconstruct the genome using a de Bruijn graph.

This figure shows the process of splitting reads into smaller k-mers (4-mers in this case) in order to be able to be used in a de Bruijn graph. (A) Shows the initial segment of DNA being sequenced. (B) Shows the reads that were made output from sequencing and also shows how they align. The problem with this alignment though is that they overlap by k-2 not k-1 (which is needed in de Bruijn graphs). (C) Shows the reads being split into smaller 4-mers. (D) Discards the repeated 4-mers and then shows the alignment of them. Note that these k-mers overlap by k-1 and can then be used in a de Bruijn graph.

Choice of k-mer

The choice of the k-mer size has many different effects on the sequence assembly. These effects vary greatly between lower sized and larger sized k-mers. Therefore, an understanding of the different k-mer sizes must be achieved in order to choose a suitable size that balances the effects. The effects of the sizes are outlined below.

Lower k-mer sizes

• A lower k-mer size will decrease the amount of edges stored in the graph, and as such, will help decrease the amount of space required to store DNA sequence.
• Having smaller sizes will increase the chance for all the k-mers to overlap, and as such, have the required subsequences in order to construct the de Bruijn graph.
• However, by having smaller sized k-mers, you also risk having many vertices in the graph leading into a single k-mer. Therefore, this will make the reconstruction of the genome more difficult as there is a higher level of path ambiguities due to the larger amount of vertices that will need to be traversed.
• Information is lost as the k-mers become smaller.
  • E.g. The possibility of AGTCGTAGATGCTG is lower than ACGT, and as such, holds a greater amount of information (refer to entropy (information theory) for more information).
• Smaller k-mers also have the problem of not being able to resolve areas in the DNA where small microsatellites or repeats occur. This is because smaller k-mers will tend to sit entirely within the repeat region and is therefore hard to determine the amount of repetition that has actually taken place.
  • E.g. For the subsequence ATGTGTGTGTGTGTACG, the amount of repetitions of TG will be lost if a k-mer size less than 16 is chosen. This is because most of the k-mers will sit in the repeated region and may just be discarded as repeats of the same k-mer instead of referring the amount of repeats.

Higher k-mer sizes

• Having larger sized k-mers will increase the amount of edges in the graph, which in turn, will increase the amount of memory needed to store the DNA sequence.
• By increasing the size of the k-mers, the number of vertices will also decrease. This will help with the construction of the genome as there will be fewer paths to traverse in the graph.
• Larger k-mers also run a higher risk of not having outward vertices from every k-mer. This is due to larger k-mers increasing the risk that it will not overlap with another k-mer by $k-1$. Therefore, this can lead to disjoints in the reads, and as such, can lead to a higher amount of smaller contigs.
• Larger k-mer sizes help alleviate the problem of small repeat regions. This is due to the fact that the k-mer will contain a balance of the repeat region and the adjoining DNA sequences (given it are a large enough size) that can help to resolve the amount of repetition in that particular area.

Applications of k-mer in bioinformatics analysis

The frequency of a set of k-mers in a species' genome, in a genomic region, or in a class of sequences, can be used as a "signature" of the underlying sequence. Comparing these frequencies are computationally easier than sequence alignment, and is an important method in alignment-free sequence analysis. It can also be used as a first stage analysis before an alignment.

• separating different species in a mixture of genetic material (metagenomics, microbiome); phase/frame information can be added
• DNA barcoding of species
• de novo sequence assembly
• human mitochondrial haplogroup classification
• detect genome mis-assembly
• de novo detection of repeated sequence such as transposable element
• characterize a protein-binding sequence motif. Besides k-mer, gapped k-mers (also named gapped q-grams or spaced seeds) can also be used
• identification of mutation or polymorphism using next generation sequencing data
• characterizing CpG island by the flanking regions
• detect horizontal transfer
• detect bacterial contamination in eukaryotic genome assembly
• detect recombination site
• using the k-mer frequency vs k-mer depth to estimate the genome size
• estimate RNA-seq level

Pseudocode

Determining the possible k-mers of a read can be done by simply cycling over the string length by one and taking out each substring of length, k. The pseudocode to achieve this is as follows:

procedure k-mer(String, k : length of each k-mer)

     n = length(String)
     
     /* cycle over the length of String till k-mers of length, k, can still be made */
     for i = 1 to  n-k+1 inclusive do
          /* output each k-mer of length k, from i to i+k in String*/
          output String[i:i+k]
     end for

end procedure

Examples

Here are some examples showing the possible k-mers (given a specified k value) from DNA sequences:

Read:     AGATCGAGTG
3-mers: AGA GAT ATC TCG CGA GAG AGT GTG

Read:     GTAGAGCTGT
5-mers: GTAGA TAGAG AGAGC GAGCT AGCTG GCTGT

External links

• Sciencedirect topic:k-mer
• bioXriv:k-mer
• arXiv: k-mer