A "kmer" is a DNA sequence of length K. A valid DNA sequence (for my purposes) can only contain the 4 following bases: A,C,T,G. I am looking for a C++ algorithm that simply outputs all possible combinations of these bases in alphabetical order into a string array. For example, if K = 2, The program should generate the following array:
kmers[0] = AA
kmers[1] = AC
kmers[2] = AG
kmers[3] = AT
kmers[4] = CA
kmers[5] = CC
kmers[6] = CG
kmers[7] = CT
kmers[8] = GA
kmers[9] = GC
kmers[10] = GG
kmers[11] = GT
kmers[12] = TA
kmers[13] = TC
kmers[14] = TG
kmers[15] = TT
If I'm thinking about this correctly, the problem really breaks down to converting a decimal integer to base 4 then substituting the appropriate bases. I thought I could use itoa for this, but itoa is not C standard, and my compiler did not support it. I welcome any clever ideas. Here is my sample code:
#include <iostream>
#include <string>
#include <math.h>
#define K 3
using namespace std;
int main() {
int num_kmers = pow(4,K);
string* kmers = NULL;
/* Allocate memory for kmers array */
kmers = new string[num_kmers];
/* Populate kmers array */
for (int i=0; i< pow(4,K); i++) {
// POPULATE THE kmers ARRAY HERE
}
/* Display all possible kmers */
for (int i=0; i< pow(4,K); i++)
cout << kmers[i] << "\n";
delete [] kmers;
}
You would need to use recursion to be flexible (i.e. so that you could change K easily).
void populate(int depth, string base, string* kmers, int* kmers_offset)
{
if(depth == K)
{
kmers[*kmers_offset].assign(base);
(*kmers_offset)++;
}
else
{
static char bases[] = { 'A', 'C', 'G', 'T' };
for(int i = 0; i < 4; ++i)
populate(depth + 1, base + bases[i], kmers, kmers_offset);
}
}
And then call it like this:
int kmers_offset = 0;
populate(0, "", kmers, &kmers_offset);
Cheers.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With