c++ - 查找最常见的k-mers及其在C++中的出现次数

关于kmer https://en.wikipedia.org/wiki/K-mer

我试图在一个大型的fastq文件中找到最常使用的k-mers。我的计划是使用misra-gries算法查找最频繁的k-mer，然后第二遍搜索文件中每个频繁的k-mer的计数。但是我认为我的算法不够高效。这是我的第一稿。我尝试尽可能提高内存使用效率(程序一定不能耗尽内存)

我也找到了这种DSK算法，但是对于我来说，如果没有简单的实现就很难理解它。 http://minia.genouest.org/dsk/

注意:每个计数器的ID也是整数而不是字符串，我将在稍后的最终草案中对其进行更改。

#include <iostream>
#include <fstream>
#include <sstream>
#include <vector>
using namespace std;

struct node {
    string id;
    int count;
};

void searchCount(vector <node>&memory, string line,int k) {
    int count = 0;
    string kMer;
    for (int i = 0; i < memory.size(); i++) {
        if (memory[i].id != "") {
            for (int j = 0; j < line.length() - k + 1; j++) {
                kMer = line.substr(j, k);
                if (kMer == memory[i].id) {
                    count++;
                }
            }
        }
        memory[i].count = count;
        count = 0;
    }
}

int doWeHaveSpace(vector <node> memory) {
    for (int i = 0; i < memory.size(); i++) {
        if (memory[i].id == "") {
            return i;
        }
    }
    return -1;


}

void MisraGries(string element, vector <node> &memory) {
    bool isHere = false;
    int index;

    for (int i = 0; i < memory.size(); i++) {
        if (memory[i].id == element) {
            isHere = true;
            index = i;
        }
    }
    if (isHere) {
        memory[index].count++;
    }
    else {
        int freeIndex = doWeHaveSpace(memory);
        if (freeIndex+1) {
            memory[freeIndex].count++;
            memory[freeIndex].id = element;
        }
        else {
            for (int i = 0; i < memory.size(); i++) {
                if (memory[i].count != 0) {
                    memory[i].count--;
                    if (memory[i].count == 0) {
                        memory[i].id = "";
                    }
                }
            }
        }
    }
}
void filecheck(ifstream & input, string prompt)  // this function checks and opens input files
{
    string filename;
    cout << "Please enter file directory and name for " << prompt << ": ";
    do
    {
        getline(cin, filename);
        input.open(filename.c_str());
        if (input.fail())
            cout << " wrong file directory. Please enter real directory. ";
    } while (input.fail());
}

int main() {
    int line = 1;
    string filename;
    ifstream input;
    ofstream output;
    string search;
    vector <node> frequent(1000);
    for (int i = 0; i < frequent.size(); i++) {
        frequent[i].id = "";
        frequent[i].count = 0;
    }
    int k = 30;
    string kMer;
    filecheck(input, "input file");

    while (!input.eof())
    {
        getline(input, search); // it gets infos line by line to count lines
        line++;
        if (line == 3) {
            for (int i = 0; i < search.length() - k + 1; i++) {
                kMer = search.substr(i, k);
                MisraGries(kMer, frequent);
            }
            line = -1;
        }

    }

    return 0;
}

最佳答案

您可以通过将最频繁的k-mers存储在散列表而不是数组中来加快代码的速度。这样，如果已经在缓存中，则可以在O(1)时间内处理一个k-mer(假设长度是恒定的)(如果不是，则仍然需要线性传递，但是可能会很大)平均改善)。

如果将很多信息遗漏在某种辅助数据结构(例如优先级队列)中，那么即使有很多遗漏，您也可以使其更快，以便您可以使用count = 0查找该元素并将其删除，而无需检查所有其他元素。

考虑到您的示例中k很小，您可以增加内存中高速缓存的大小(典型的计算机应该很容易在内存中保留几百万个这样的字符串)，以减少遗漏。

您可以在第一遍中通过散列k-mers来存储更多数据(这样，您只需要在内存中存储整数而不是字符串即可)。

总结起来，我建议增大缓存(只要它适合内存)，并使用更合适的数据结构来支持快速查找，例如哈希表(C++中的std::unordered_map)。