1 引言

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

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

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

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

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

　　3 算法設(shè)計(jì)

　　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 算法找出第一條從起點(diǎn) S 到終點(diǎn) D 的最短路徑 P1，加入集合A( 已選路徑集合);

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

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

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

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

　　3.2 加速選擇法

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

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

　　圖 3 記錄已選邊

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

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

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

　　3.4 記錄已選邊

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

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

　　3.5 記錄關(guān)鍵邊

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

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

　　4 平臺(tái)實(shí)現(xiàn)

　　圖 5 程序主流程

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

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

　　5 算法性能測(cè)試

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

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

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

　　圖 6 改進(jìn)的 KSP 算法流程

　　圖 7 系統(tǒng)界面

　　表 1 改進(jìn)的 KSP 算法性能對(duì)比

　　表 2 改進(jìn)的 KSP 算法性能實(shí)測(cè)數(shù)據(jù)

　　6 討論與分析

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

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

　　7 結(jié)論

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

　　黃祖成 1 沈夢(mèng)圓 2 侯至丞 1 TOKUYASU Taku Andrew 3 蒙海林 2

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

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

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