前提条件：后缀数组。
什么是k-mers?
术语k-mer通常是指字符串中包含的所有可能的长度为k的子字符串。对DNA / RNA测序读数中的所有k-mers进行计数是许多生物信息学应用的第一步。
什么是后缀数组?
后缀数组是字符串的所有后缀的排序数组。它是一种数据结构，尤其是在全文索引，数据压缩算法中使用。更多信息可以在这里找到。
问题：我们给了一个字符串str和一个整数k。我们必须找到所有对(substr，i)，使得substr是一个长度–恰好发生i次的str的k个子串。
该方法涉及的步骤：
让我们以“ banana $”一词为例。
步骤1：计算给定文本的后缀数组。
6 $
5 a$
3 ana$
1 anana$
0 banana$
4 na$
2 nana$
步骤2：迭代后缀数组，保持“ curr_count” 。
1.如果当前后缀的长度小于k，则跳过迭代。也就是说，如果k = 2 ，则当当前后缀为$时，将跳过迭代。
2.如果当前后缀以与上一个后缀相同的长度-k子字符串开头，则递增curr_count。例如，在第四次迭代期间，当前后缀“ anana $”以长度k “ an”相同的子串开头，与前一个后缀“ ana $”相同。因此，在这种情况下，我们将增加curr_count。
3.如果不满足条件2，则如果前一个后缀的长度等于k，则它是有效对，并且将其与当前计数一起输出，否则，将跳过该迭代。
curr_count Valid Pair
6 $ 1
5 a$ 1
3 ana$ 1 (a$, 1)
1 anana$ 1
0 banana$ 2 (an, 2)
4 na$ 1 (ba, 1)
2 nana$ 1 (na, 2)
例子：
Input : banana$ // Input text
Output : (a$, 1) // k- mers
(an, 2)
(ba, 1)
(na, 2)
Input : geeksforgeeks
Output : (ee, 2)
(ek, 2)
(fo, 1)
(ge, 2)
(ks, 2)
(or, 1)
(sf, 1)
以下是上述方法的C代码：
// C program to solve K-mer counting problem
#include
#include
#include
// Structure to store data of a rotation
struct rotation {
int index;
char* suffix;
};
// Compares the rotations and
// sorts the rotations alphabetically
int cmpfunc(const void* x, const void* y)
{
struct rotation* rx = (struct rotation*)x;
struct rotation* ry = (struct rotation*)y;
return strcmp(rx->suffix, ry->suffix);
}
// Takes input_text and its length as arguments
// and returns the corresponding suffix array
char** computeSuffixArray(char* input_text,
int len_text)
{
int i;
// Array of structures to store rotations
// and their indexes
struct rotation suff[len_text];
// Structure is needed to maintain old
// indexes of rotations after sorting them
for (i = 0; i < len_text; i++) {
suff[i].index = i;
suff[i].suffix = (input_text + i);
}
// Sorts rotations using comparison function
// defined above
qsort(suff, len_text, sizeof(struct rotation), cmpfunc);
// Stores the suffixes of sorted rotations
char** suffix_arr =
(char**)malloc(len_text * sizeof(char*));
for (i = 0; i < len_text; i++) {
suffix_arr[i] =
(char*)malloc((len_text + 1) * sizeof(char));
strcpy(suffix_arr[i], suff[i].suffix);
}
// Returns the computed suffix array
return suffix_arr;
}
// Takes suffix array, its size and valid length as
// arguments and outputs the valid pairs of k - mers
void findValidPairs(char** suffix_arr, int n, int k)
{
int curr_count = 1, i;
char* prev_suff = (char*)malloc(n * sizeof(char));
// Iterates over the suffix array,
// keeping a current count
for (i = 0; i < n; i++) {
// Skipping the current suffix
// if it has length < valid length
if (strlen(suffix_arr[i]) < k) {
if (i != 0 && strlen(prev_suff) == k) {
printf("(%s, %d)\n", prev_suff, curr_count);
curr_count = 1;}
strcpy(prev_suff, suffix_arr[i]);
continue;
}
// Incrementing the curr_count if first
// k chars of prev_suff and current suffix
// are same
if (!(memcmp(prev_suff, suffix_arr[i], k))) {
curr_count++;
}
else {
// Pair is valid when i!=0 (as there is
// no prev_suff for i = 0) and when strlen
// of prev_suff is k
if (i != 0 && strlen(prev_suff) == k) {
printf("(%s, %d)\n", prev_suff, curr_count);
curr_count = 1;
memcpy(prev_suff, suffix_arr[i], k);
prev_suff[k] = '\0';
}
else {
memcpy(prev_suff, suffix_arr[i], k);
prev_suff[k] = '\0';
continue;
}
}
// Modifying prev_suff[i] to current suffix
memcpy(prev_suff, suffix_arr[i], k);
prev_suff[k] = '\0';
}
// Printing the last valid pair
printf("(%s, %d)\n", prev_suff, curr_count);
}
// Driver program to test functions above
int main()
{
char input_text[] = "geeksforgeeks";
int k = 2;
int len_text = strlen(input_text);
// Computes the suffix array of our text
printf("Input Text: %s\n", input_text);
char** suffix_arr =
computeSuffixArray(input_text, len_text);
// Finds and outputs all valid pairs
printf("k-mers: \n");
findValidPairs(suffix_arr, len_text, k);
return 0;
}
输出：
Input Text: banana$
k-mers:
(a$, 1)
(an, 2)
(ba, 1)
(na, 2)
时间复杂度： O(s * len_text * log(len_text))，假设s是最长后缀的长度。
资料来源：
1.后缀数组Wikipedia
2.后缀阵列CMU