求de Bruijn图上与目标串最小编辑距离

求de Bruijn图上与目标串最小编辑距离

题目表述

• 给出一个 kk 阶 de Bruijn 图，以及字符串 a
• 求 de Bruijn 图上的一条路径，使其代表的字符串与 aa 的编辑距离尽可能小，并输出编辑过程
• 其中
  
  • len(a)≤10000
  • k≤30
  • 字符集 |Σ|≤4 （为ATCG）
  • de Bruijn 图中节点数量 ≤≤ 10000

de Bruijn 图

在图论中，一个 k 阶 de Bruijn 图是一个用来表示序列交叠的有向图。

• 图中的每个节点表示一个长度为 k 的字符串
• 若两个节点的字符序列有 k−1的交叠即在这两个节点之间连一条有向边
  
  • 示例 一个 4 阶 de Bruijn*图
    ABCD -> BCDE -> CDEF
• 因此在一个 k 阶 de Bruijn 图上，一个长度为 l 的路径可表示一个长度为 k+l−1 的字符串

思路

大概思路就是在图上跑动态规划

编辑距离

• 定义函数edit(i,j)表示第一个字符串（记为a）前i个字符的子串与第二个字符串（记为b）前j个字符的子串的编辑距离
• 使edit(i, j) == min{ edit(i-1, j) + 1, edit(i, j-1) + 1, edit(i-1, j-1) + f(i, j) }，当第一个字符串的第i个字符不等于第二个字符串的第j个字符时，f(i, j) = 1；否则，f(i, j) = 0。
• 其中 edit(i-1, j) + 1代表最后一步为a串删除最后一个字符
• edit(i, j-1) + 1 代表最后一步为a串增加字符
• edit(i-1, j-1)+ f(i, j)代表替换字符或者不操作

正确性证明

首先串a最小编辑成串b,那么：
1. 插入的字符不能被替换或删除.
2. 替换的字符不能被替换或删除.
否则的话显然不是最小编辑。也就是说最小编辑只能对字符串a已经存在的字符，至多操作一次。要么删除要么替换，要么不操作。现在我们考虑串a转换成串b，考虑b[j]的来源。分为三种情况：
1. 串a的某个字符
2. 串a的某个字符替换得到的
3. 串a增加一个字符得到的
对于情况1，考虑串a的最后字符（记为a[i]）的情况，分为以下几种情况：
1.1 a[i]被删除，那么情况就转化到dp(i-1,j)，
dp(i,j)=dp(i-1,j)+1
1.2 a[i]被替换，不可能。因为不能满足b[j]是串a的字符。
1.3 a[i]==b[j], 情况转化成dp(i-1,j-1)
对于情况2，考虑串a的最后字符（记为a[i]）的情况，分为以下几种情况：
2.1 a[i]被删除，那么情况就转化到dp(i-1,j)
2.2 a[i]被替换b[j], 情况转化成dp(i-1,j-1)
2.3 a[i]被替换但不是b[j], 不可能。因为不能满足b[j]是串a的字符替换这个条件。
2.4 a[i]保留不动，不可能。因为不能满足b[j]是串a的字符替换这个条件。
对于情况3，情况直接转化到dp(i,j-1).
综上，递归求解过程正确。

建图

• 图上的每个节点代表一个k （k<30）长度的字符串, 这边使用long long 类型储存， 也就是说将字符串转成4进制的整数
• 使用的hash函数直接用C++的unordered_map库，将long long类型数据映射到0 ~ n - 1
• 若节点i 有到节点j 的边，代表节点i的后k-1字符与节点j的前k-1字符相同，如节点i（ATCG），节点j（TCGG）

主要过程

以图上的每一个节点为起始节点，使用BFS，获得新的串的编辑距离，其中fdp(i, j)代表从节点start开始到节点i组成的字符串与目标字符串的前j个字符的编辑距离。如start=0, 节点0(AATC), 节点1(ATCC), 目标串(AAABB), 则fdp(1, 3)代表(AATCC)与(AAAB)编辑距离

需要注意

• 从节点start到节点i表述不够准确，因为图中可能有环，所以节点i可能出现多次，但是每次表示的字符串都不一样，所以应该说是上次串长度加1组成的新串。同时，为了防止BFS死循环（因为出现环），我们在BFS时必须记录当前节点是否以前出现过，如果出现过，要考虑这次的编辑距离是否更新过，即fdp(i, j)是否为目前出现过最小的， 否则不加入队列
• 同时，由于数据范围较大，所以动态规划使用循环数组，也就是fdp(i)为循环数组
• 维护队列q,进行BFS, q中储存的是需要更新的新串的编辑距离
• 维护helper队列，储存fdp数组编号信息，如对于节点i,先在helper队列获取编号p, 则fdp[p]代表fdp(i)了，同时满足约束条件helper.front一定是helper.end的前驱。
• 打印转换序列就是再求一次编辑距离，获取状态，依次打出即可

    int findShortID(int start){
        //clearToZero();
        vector<vector<int>> dp(k + 1, vector<int>(dstlen + 1));
        long long a = strMap[start];
        string &b = dstStr;
        for (int i = 0; i <= k; i++){
            dp[i][0] = i;
        }
        for (int j = 0; j <= dstlen; j++){
            dp[0][j] = j;
        }
        for (int i = 0; i < k; i++){
            int c1 = getCharAt(a, i);
            for (int j = 0; j < dstlen; j++){
                char c2 = b[j];
                if (c1 == c2){
                    dp[i + 1][j + 1] = dp[i][j];
                }
                else{
                    int replace = dp[i][j] + 1;
                    int delet = dp[i][j + 1] + 1;
                    int insert = dp[i + 1][j] + 1;
                    int min = replace;
                    if (min > insert){
                        min = insert;
                    }
                    if (min > delet){
                        min = delet;
                    }
                    dp[i + 1][j + 1] = min;
                }
            }
        }
        // 
        tempAns.bestdis = 0.31 * dstlen;
        //tempAns.str = strMap[start];
        queue<int> que;
        Node* p = nodeMap[start];
        int cnum = 0;
        while (p)
        {
            que.push(p->v);
            p = p->next;
            cnum++;
        }
        Trie strTree;
        TrieNode* curChar = strTree.root;
        for (int i = 0; i < k; i++){
            TrieNode* temp = strTree.insert(curChar, getCharAt(a, i));
            curChar = temp;
        }
        FDPHelper helper;
        fdp[helper.push(cnum)].assign(dp[k].begin(), dp[k].end());
        fdpstrings[0] = curChar;
        tempAns.str = strTree.getString(fdpstrings[0]);
        vector<int> checkTimeMap(totalNode + 1);
        checkTimeMap[start] = 1;
        while (!que.empty()){
            int strindex = que.front();
            que.pop();
            long long str = strMap[strindex];
            char c1 = getCharAt(str, k - 1);
            Node* p = nodeMap[strindex];
            int cnum = 0;
            while (p)
            {
                p = p->next;
                cnum++;
            }
            int dpcurindex = helper.push(cnum);
            int dppreindex = helper.addone();
            fdp[dpcurindex][0] = fdp[dppreindex][0] + 1;
            for (int j = 0; j < dstlen; j++){
                char c2 = b[j];
                if (c1 == c2){
                    fdp[dpcurindex][j + 1] = fdp[dppreindex][j];
                }
                else{
                    int replace = fdp[dppreindex][j] + 1;
                    int delet = fdp[dppreindex][j + 1] + 1;
                    int insert = fdp[dpcurindex][j] + 1;
                    int min = replace;
                    if (min > insert){
                        min = insert;
                    }
                    if (min > delet){
                        min = delet;
                    }
                    fdp[dpcurindex][j + 1] = min;
                }
            }
            bool flag = false;
            if (fdp[dpcurindex][dstlen] < tempAns.bestdis){
                tempAns.bestdis = fdp[dpcurindex][dstlen];
                cout << tempAns.bestdis << endl;
                fdpstrings[dpcurindex] = strTree.insert(fdpstrings[dppreindex], c1);
                tempAns.str = strTree.getString(fdpstrings[dpcurindex]);
                return tempAns.bestdis;
                flag = true;
            }
            if (checkTimeMap[strindex] == 0 || tempAns.bestdis == fdp[dpcurindex][dstlen]){
                p = nodeMap[strindex];
                while (p)
                {
                    que.push(p->v);
                    p = p->next;
                }
                if (!flag) fdpstrings[dpcurindex] = strTree.insert(fdpstrings[dppreindex], c1);
            }
            else if (tempAns.bestdis != fdp[dpcurindex][dstlen]){
                helper.reverse();
            }
            checkTimeMap[strindex]++;
        }
        strTree.clearTree(strTree.root);
        return tempAns.bestdis;
    }