1 引言

　　代謝通路 (Metabolic Pathways) 是指生物體內(nèi)將源代謝物轉(zhuǎn)化為目標(biāo)代謝物的一系列連續(xù)的酶催化反應(yīng)，如糖酵解通路、三羧酸循環(huán)、梭狀芽胞桿菌固定二氧化碳的Wood-Ljungdahl 通 路、藍藻吸收氨的谷氨酰胺合成酶循環(huán)(GSGOGAT) 等。在合成生物學(xué)領(lǐng)域，常常需要對代謝通路進行設(shè)計與分析。代謝通路設(shè)計的目標(biāo)是給定一個或者若干個起始和目標(biāo)化合物，從代謝數(shù)據(jù)中找到并返回生物相關(guān)的、由酶催化反應(yīng)構(gòu)成的、有特定功能和作用的、從起始化合物到目標(biāo)化合物的可行路徑。在過去的十多年里，隨著實驗技術(shù)的進步， 高通量組學(xué)數(shù)據(jù)日益增長，越來越多的代謝通路得以解析。以KEGG 和 MetaCyc 為代表的代謝數(shù)據(jù)庫發(fā)展迅速，為識別和發(fā)現(xiàn)新的合成替代路徑帶來了新的機遇和挑戰(zhàn)。研究人員根據(jù)不同物種的代謝數(shù)據(jù)特點和生物合成的實際應(yīng)用需求，開發(fā)了多個代謝通路設(shè)計算法和工具，如 MAPPS、路徑搜索系統(tǒng)、novoPathFinder、MRSD、FogLight 和 Metabolic Tinker 等。這些軟件從算法選擇上可以歸為基于圖論、基于化學(xué)計量學(xué)和基于逆合成等不同類別。利用這些工具，研究人員可方便地從代謝數(shù)據(jù)中識別、設(shè)計出具有潛在價值的代謝通路。

　　搜索算法是代謝通路設(shè)計工具的核心步驟。在代謝網(wǎng)絡(luò)較大時，現(xiàn)有的代謝設(shè)計工具存在搜索時間相對較長的問題。為此，現(xiàn)有的代謝設(shè)計工具嘗試從數(shù)據(jù)結(jié)構(gòu)上對算法進行優(yōu)化， 或在路徑選取策略上進行優(yōu)化。此外，在多路徑搜索算法方面，傳統(tǒng)的前 K 條最短路徑算法 (K-Shortest Path，KSP) 是以 Yen’s 算法為代表。該算法是基于 Dijkstra 算法進行演化的版本，在搜索功能上實現(xiàn)了較為有效的前 K 條最短路徑搜索，但其算法過程從第 3 條路徑開始則需要提升。另外，綜合來看 KSP 算法的核心思想都是邊刪除選擇法，而過去眾多算法圍繞數(shù)據(jù)結(jié)構(gòu)及路徑選取策略上進行過優(yōu)化，也取得了一定的性能提升效果。近年來，由于計算機硬件性能的提升，傳統(tǒng)KSP 算法優(yōu)勢變得不再明顯。因此，本文基于 Yen’s 算法進行優(yōu)化改進的 KSP 算法，在對路徑選取策略進行優(yōu)化的同時， 結(jié)合多核心 CPU 計算機硬件特點，使用并行計算方式來進一步提升代謝通路搜索算法的運算性能。

　　2 數(shù)據(jù)建模

　　本文主要是針對合成生物學(xué)的代謝網(wǎng)絡(luò)路徑搜索進行優(yōu)化，代謝網(wǎng)絡(luò)數(shù)據(jù)來源于 KEGG 通路數(shù)據(jù)庫 (https://www. kegg.jp/kegg/pathway.html)。通過 KEGG API 接口批量下載不同物種及參考的代謝反應(yīng)圖的 KGML 文件，并進行解析來提取化合物和反應(yīng)方向。將代謝通路表示為圖，圖中的節(jié)點表示為代謝通路中的化合物 ( 代謝物 )，邊表示化合物之間的轉(zhuǎn)化 ( 化學(xué)反應(yīng)或代謝反應(yīng) )。

　　本文設(shè)計的代謝網(wǎng)絡(luò)路徑搜索模型針對 KEGG 數(shù)據(jù)格式的特點進行了專門優(yōu)化。在不考慮能耗和確保路徑正確的前提下，針對可能存在的代謝通路，盡可能搜索多條路徑，把代謝網(wǎng)絡(luò)數(shù)據(jù)簡化以提高搜索的效率。具體處理方法為：把具有通路的化合物之間的 Cost 值設(shè)為 1，不連通的化合物之間的Cost 值設(shè)為 max( 實驗時，考慮到 KEGG 代謝網(wǎng)絡(luò)參考圖中化合物數(shù)量約為 3 000 個，設(shè) max = 9999)。所生成的代謝網(wǎng)絡(luò)圖的鄰接表如下：

　　3 算法設(shè)計

　　3.1 傳統(tǒng) KSP 算法

　　傳統(tǒng) KSP 算法主要以 Yen’s 算法為代表。其中， Yen’s 算法是 Yen 在 1971 年提出并以其名字命名的算法。Yen’s 算法采用了遞推法中的偏離路徑算法思想，適用于非負(fù)權(quán)邊的有向無環(huán)圖結(jié)構(gòu)。其主算法步驟如下：

　　圖 1 加速選擇法

　　圖 2 相同的候選路徑

　　(1) 在非負(fù)權(quán)邊的有向無環(huán)圖 G 中， 使用 Dijkstra 算法找出第一條從起點 S 到終點 D 的最短路徑 P1，加入集合A( 已選路徑集合);

　　(2) 從集合 A 中選擇最后加入的一條路徑 Pn，在 G 中依次刪除路徑 Pn 中的邊 ( 把邊的 Cost 值設(shè)為最大 )，并使用Dijkstra 算法計算出對應(yīng)的最短路徑 PNi 加入集合 B( 候選路徑集合 ); 

       (3) 從集合 B 中選出一條 Cost 最小的路徑 Pm，加入集合 A 并從集合 B 中刪除 Pm;

       (4) 重復(fù)步驟 (2) ～ (3) 直到集合 A 的數(shù)量達到 K，則集合 A 為所求的前 K 條最短路徑。

　　本文將步驟 (2) ～ (3) 定義為次短路徑搜索 (Second Shortest Path Search，SSPS) 過程。在 KSP 算法中， 最耗時的運算是 Dijkstra 運算，而該算法中的 SSPS 過程存在較多重復(fù)的 Dijkstra 計算。因此，在節(jié)點數(shù)量較多的網(wǎng)絡(luò)中進行KSP 運算將會耗費較長時間。本文主要研究如何對 SSPS 過程進行優(yōu)化以減少不必要的 Dijkstra 運算，從而提高 KSP 算法性能。

　　3.2 加速選擇法

　　在 Dijkstra 路徑搜索算法中，一次只能找出一條最短路徑，但 Cost 相同的最短路徑可能并不唯一。在 SSPS 過程前， 集合 B 中可能已經(jīng)存在所需要的最短路徑。在合成生物代謝網(wǎng)絡(luò)中，所有連通的邊的 Cost 值設(shè)為 1，此時 Cost 值相同的路徑都認(rèn)可為是最短路徑，但并不設(shè)先后順序。因此，若在集合 B 中存在與集合 A 最后一條路徑 Cost 值相同的路徑， 則可認(rèn)為是在集合 A 之后的最短路徑。此時，可直接從集合 B 中選擇出最短路徑 Pn，則省去了 SSPS 的計算過程。

　　如圖 1 所示，當(dāng)前需要搜索第 4 條最短路徑時集合 A 中已有 3 條路徑，集合 B 中 Cost 值最小的路徑與集合 A 的最后一條路徑 Cost 值相等，此時直接從集合 B 中選取 Cost 值為 5 的路徑加入集合 A 中作為第 4 條最短路徑，從而節(jié)省了一輪 SSPS 運算 ( 此例子中可節(jié)省 4 次 Dijkstra 運算 )。

　　圖 3 記錄已選邊

　　圖 4 記錄關(guān)鍵邊

　　3.3 刪除重復(fù)候選路徑

　　在一個 SSPS 過程中，得出的候選路徑可能與上一輪SSPS 得到的候選路徑存在相同的，此時需要判斷，當(dāng)路徑已存在集合 B 中則不加入到集合 B。如圖 2 所示，刪除邊 BC 或CD 所產(chǎn)生的候選路徑是相同的，此時只選擇一條加入集合 B 即可。

　　3.4 記錄已選邊

　　在已找到的路徑集合 A 中，除第一條路徑外，其他路徑都有對應(yīng)刪除的邊，這些邊的集合記為 EA。在進行 SSPS 過程前，先在圖 G 中刪除 EA 里的邊，可減少重復(fù)的 Dijkstra 運算。

　　如圖 3 所示，集合 A 中第二條最短路徑所對應(yīng)的邊為BD，在進行下一輪 SSPS 運算前把邊 BD 刪除以減少 SSPS 的運算時間 ( 此例子中減少了 1 次 Dijkstra 運算 )。

　　3.5 記錄關(guān)鍵邊

　　關(guān)鍵邊 (Critical Edge) 是從起點到終點的必經(jīng)之路，在一個 SSPS 過程中跳過關(guān)鍵邊的 Dijkstra 運算可節(jié)省運算時間。具體做法為：在 SSPS 過程中，若 Dijkstra 運算的結(jié)果Cost 是最大值，則把對應(yīng)的邊加入到關(guān)鍵邊的集合 C 中;隨后在下一輪 SSPS 過程中跳過關(guān)鍵邊的 Dijkstra 運算，可節(jié)省運算時間。

　　如圖 4 所示，當(dāng)經(jīng)過一輪 SSPS 運算后，可以得到邊 AB 和 DH 為關(guān)鍵邊，并在下一輪 SSPS 運算時把關(guān)鍵邊排除，以減少運算量 ( 此例子中節(jié)省 2 次 Dijkstra 運算 )。

　　4 平臺實現(xiàn)

　　圖 5 程序主流程

　　平臺采用 Flask + Vue 框架，實現(xiàn)前后端分離。核心算法代碼采用 C++ 實現(xiàn)以提升運算速率。代謝網(wǎng)絡(luò)數(shù)據(jù)主要來源于 KEGG 數(shù)據(jù)庫，在路徑搜索前把 KEGG 數(shù)據(jù)進行初始化(KEGG 數(shù)據(jù)庫轉(zhuǎn)化為鄰接表 ) 以便進行 KSP 算法運算。系統(tǒng)主要程序流程如圖 5 所示。KEGG 數(shù)據(jù)的初始化在系統(tǒng)啟動前期執(zhí)行，此部分程序運行時間約為 3 ～ 5 s，這對用戶來說并不產(chǎn)生影響。對系統(tǒng)用戶產(chǎn)生影響的過程為 KSP 運算部分， 這也是本文主要的研究工作。

　　圖 6 所示為改進的 KSP 算法流程，其中 SSPS 過程采用了并行方式同時進行多個 Dijkstra 運算，算法時間復(fù)雜度從 O ( n×m ) 減少到 O(n)，算法結(jié)束后會產(chǎn)生 K 條最短路徑。為提高用戶體驗，在實際應(yīng)用系統(tǒng)中每產(chǎn)生一條最短路徑則返回到用戶界面上顯示，這樣可以減少用戶主觀的等待時間。圖 7 所示用戶系統(tǒng)界面，搜索得到的路徑為按短到長依次顯示。

　　5 算法性能測試

　　在實驗數(shù)據(jù)中，隨機選取起點和終點，分別設(shè)定最短路徑數(shù)量為 1 ～ 7 條，觀察算法所調(diào)用 Dijkstra 的次數(shù)來分析算法的性能。具體數(shù)據(jù)如表 1 所示。

　　從表 1 可知，當(dāng)路徑數(shù)量 K 大于 2 時，改進的 KSP 算法在性能上有 2 ～ 3 倍的提升，并行的 KSP 算法在性能上有5 ～ 9 倍的提升，并隨著路徑數(shù)量的增加而不斷提升。本次改進的 KSP 算法性能提升較為顯著。

　　4 條目標(biāo)路徑進行前 7 條最短路徑的搜索，每條路徑的平均搜索時間為 0.27 ～ 0.39 s。其中，所選取的 KEGG 參考代謝網(wǎng)絡(luò)圖節(jié)點總數(shù)為 2828，邊總數(shù)為 4051。即使在 KSP 算法運算時間上增加 0.5 s 的 http 請求及獲取化合物及反應(yīng)信息所需要的時間，在用戶端平均每條路徑的搜索總時間仍然小于 1 s，提升了系統(tǒng)的用戶體驗。

　　圖 6 改進的 KSP 算法流程

　　圖 7 系統(tǒng)界面

　　表 1 改進的 KSP 算法性能對比

　　表 2 改進的 KSP 算法性能實測數(shù)據(jù)

　　6 討論與分析

　　在合成生物學(xué)設(shè)計尤其是細胞工廠 / 底盤細胞設(shè)計中，代謝通路設(shè)計工具可有效輔助研究人員快速找到出發(fā)底物及目標(biāo)產(chǎn)物之間的連接路徑，提高設(shè)計效率。例如，Yim 等對 1,4- 丁二醇 (Butane-1,4-diol，BDO) 的合成代謝通路進行 1 萬多次的計算調(diào)查后，設(shè)計獲得最佳途徑，并在大腸桿菌中獲得高達 18 g/L 的 BDO。隨著下游酶的改善，BDO 滴度提高到110 g/L。該案例成功展示了代謝通路設(shè)計算法在生物合成中的巨大應(yīng)用潛力。

　　自 Dijkstra 算法提出以來，許多學(xué)者都對它進行過不同程度的優(yōu)化以提高其性能。在后來的多路徑搜索 (KSP) 算法上，大部分研究圍繞刪除路徑核心思想進行了許多的改進。這些改進的算法在過去計算機性能有限的條件下取得了較好的效果，算法時間復(fù)雜度從 O(n×n) → O(n×m) → O(n×logm) 不斷提升。但隨著計算機性能的提升，特別是多核心 CPU 以及多CPU 架構(gòu)的計算機系統(tǒng)的廣泛應(yīng)用，傳統(tǒng)以單線程進行算法迭代運算的方式并不能發(fā)揮很好的效果。本文正是利用了多核心 CPU 的硬件特點，在算法優(yōu)化的同時采用并行編程方式，把算法的時間復(fù)雜度提升到 O(n) 水平，大幅提升了運算性能。

　　7 結(jié)論

　　本文針對合成生物學(xué)代謝網(wǎng)絡(luò)中代謝通路非唯一以及傳統(tǒng) KSP 算法效率偏低的問題，對 KSP 算法進行改進優(yōu)化以提升運算性能。在對 KSP 算法策略上優(yōu)化的同時，采用并行計算方式進一步提升算法的性能。使用 Python 實現(xiàn)代謝通路設(shè)計 Web 平臺，并采用 C++ 編寫核心算法的代碼。通過引入KEGG 代謝網(wǎng)絡(luò)圖，驗證改進的 KSP 算法比傳統(tǒng) KSP 算法有較大的性能提升，在合成生物學(xué)代謝通路設(shè)計上具有一定的應(yīng)用價值。然而，由于代謝反應(yīng)在不同物種的代謝網(wǎng)絡(luò)中會有差別，本文基于 KEGG 參考圖上進行了搜索算法的研究，并未對物種的特性進行區(qū)分。因此，后續(xù)工作可在 KEGG 參考圖基礎(chǔ)上針對不同物種的約束條件進行路徑算法研究，以適應(yīng)合成生物學(xué)對不同物種代謝網(wǎng)絡(luò)通路設(shè)計的特定需求。

　　黃祖成 1 沈夢圓 2 侯至丞 1 TOKUYASU Taku Andrew 3 蒙海林 2

　　1 廣州中國科學(xué)院先進技術(shù)研究所 機器人與智能裝備研究中心

　　2 廣州中國科學(xué)院先進技術(shù)研究所 生物工程研究中心

　　3 中國科學(xué)院深圳先進技術(shù)研究院 深圳合成生物學(xué)創(chuàng)新研究院 中國科學(xué)院定量工程生物學(xué)重點實驗室 廣東省合成基因組學(xué)重點實驗室轉(zhuǎn)載自《集成技術(shù)》