計(jì)算設(shè)計(jì)改造Thermobifida fusca 5-羧基-2-戊烯酰-輔酶A還原酶促進(jìn)己二酸生產(chǎn)

楊菊,毛銀,黃曉強(qiáng),周勝虎,鄧禹*

1(江南大學(xué) 生物工程學(xué)院,江蘇 無(wú)錫,214122)2(密歇根大學(xué) 計(jì)算醫(yī)學(xué)與生物信息系,密歇根 安娜堡,48109)

己二酸是一種重要的二元羧酸,廣泛應(yīng)用于醫(yī)療衛(wèi)生、食品、化工等行業(yè)[1],據(jù)統(tǒng)計(jì),超過(guò)60%的己二酸被用于生產(chǎn)尼龍類纖維(如尼龍6-6)[2-3]。目前,己二酸的大規(guī)模生產(chǎn)仍依賴化學(xué)合成,主要通過(guò)硝酸對(duì)環(huán)己醇-環(huán)己酮的混合物(醇酮油)進(jìn)行氧化反應(yīng)制取[4-5]。化學(xué)合成法存在著工藝流程長(zhǎng)、副產(chǎn)物(尤其是氮氧化物)排放多、產(chǎn)品收率低等問(wèn)題。因此,人們迫切希望尋找可替代化學(xué)合成己二酸的新方法。目前,全生物法合成己二酸已成為可能[6]。

本課題組于2015年發(fā)現(xiàn)天然菌株ThermobifidafuscaB6中存在一條全生物合成己二酸的代謝途徑(圖1),并且將該途徑可移植到大腸桿菌中[7]。該途徑涉及5個(gè)酶：β-酮硫解酶(β-ketothiolase)、3-羥酰基-輔酶A脫氫酶(3-hydroxyacyl-CoA dehydrogenase)、3-羥基己二酰-輔酶A脫氫酶(3-hydroxyadipyl-CoA dehydrogenase)、5-羰基-2-戊烯酰基-輔酶A還原酶(5-carboxy-2-pentenoyl-CoA reductase)和己二酰輔酶A合成酶(adipyl-CoA synthetase)[7]。總體來(lái)講,該代謝過(guò)程能通過(guò)消耗D-葡萄糖來(lái)合成己二酸,但產(chǎn)量較低(0.2 g/L)[8]。其中,5-羰基-2-戊烯酰基-輔酶A還原酶是該途徑的限速酶,其催化5-羰基-2-戊烯酰基-輔酶A(5-carboxy-2-pentenoyl-CoA)反應(yīng)生成己二酰輔酶A(adipyl-CoA)的活性很低[7]。目前,關(guān)于改造5-羰基-2-戊烯酰基-輔酶A還原酶以促進(jìn)己二酸合成的研究尚未見(jiàn)報(bào)道。

近年來(lái),基于結(jié)構(gòu)的計(jì)算酶設(shè)計(jì)方法被廣泛應(yīng)用于改造天然酶以適應(yīng)工業(yè)需要[9-14],甚至是從頭設(shè)計(jì)具有新功能的酶[15-17]。T.fusca5-羧基-2-戊烯酰-輔酶A還原酶的結(jié)構(gòu)尚未通過(guò)實(shí)驗(yàn)解析出來(lái)。本研究首先通過(guò)同源建模構(gòu)建5-羧基-2-戊烯酰-輔酶A還原酶的結(jié)構(gòu)模型并對(duì)其進(jìn)行分析,進(jìn)而基于酶-配體復(fù)合物的結(jié)構(gòu),通過(guò)計(jì)算酶設(shè)計(jì)方法來(lái)改造5-羧基-2-戊烯酰-輔酶A還原酶的結(jié)合殘基,嘗試引入氫鍵網(wǎng)絡(luò)來(lái)增強(qiáng)突變酶與目的底物5-羧基-2-戊烯酰-輔酶A的結(jié)合能力以提高酶催化效率,并通過(guò)分子動(dòng)力學(xué)模擬來(lái)檢驗(yàn)設(shè)計(jì)的合理性。

圖1 己二酸的生物合成途徑Fig.1 Biosynthetic pathway for producing adipic acid

1 材料和方法

1.1 蛋白質(zhì)序列

5-羧基-2-戊烯酰-輔酶A還原酶包含385個(gè)氨基酸,其蛋白質(zhì)序列為：

MSDFDLYRPTEEHEALREAIRSVAEDKIAPHAADVDEQ-SRFPQEAYEALRASDFHAPHVAEEYGGVGADALATCIVIEEIA-RVCASSSLIPAVNKLGSMPLILSGSDEVKQRYLPELASGEAMFS-YGLSEREAGSDTASMRTRAVRDGDDWILNGQKSWITNAGISK-YYTVMAVTDPDGPRGRNISAFVVHIDDPGFSFGEPERKLGIKG-SPTRELIFDNVRIPGDRLVGKVGEGLRTALRTLDHTRVTIGAQ-AVGIAQGALDYALGYVKERKQFGKAIADFQGIQFMLADMAM-KLEAARQMVYVAAAKSERDDADLSFYGAAAKCFASDVAMEI-TTDAVQLLGGYGYTRDYPVERMMRDAKITQIYEGTNQIQRVV-MARQLLKK

1.2 同源建模

蛋白質(zhì)數(shù)據(jù)庫(kù)(protein data bank,PDB)[18]中沒(méi)有5-羧基-2-戊烯酰-輔酶A還原酶的實(shí)驗(yàn)結(jié)構(gòu),因此本文通過(guò)同源建模來(lái)構(gòu)建其結(jié)構(gòu)模型。對(duì)上述序列進(jìn)行BLASTp[19]搜索(采用默認(rèn)參數(shù)),發(fā)現(xiàn)PDB中存在多個(gè)同源結(jié)構(gòu)。因?yàn)橛?jì)算酶設(shè)計(jì)過(guò)程中需要同時(shí)對(duì)配體分子(FAD和5-羧基-2-戊烯酰-輔酶A)進(jìn)行建模,所以本文選擇包含了FAD和類似于5-羧基-2-戊烯酰-輔酶A配體的實(shí)驗(yàn)結(jié)構(gòu)作為模板。將上述蛋白質(zhì)序列提交到SWISS-MODEL服務(wù)器[20],搜索模板,結(jié)果顯示有超過(guò)50個(gè)結(jié)構(gòu)模板可供使用。本文選擇分辨率1.8?的晶體結(jié)構(gòu)4L1F_A作為模板建模[21],因?yàn)樵摻Y(jié)構(gòu)中包含了配體分子FAD和過(guò)硫化輔酶A(CoA persulfide,PDB文件中對(duì)應(yīng)名字為COS,與目標(biāo)底物5-羧基-2-戊烯酰-輔酶A結(jié)構(gòu)相似),并且與5-羧基-2-戊烯酰-輔酶A還原酶具有46.28%的序列相同性。

1.3 分子動(dòng)力學(xué)模擬

本文使用GROMACS 5.4.0軟件進(jìn)行分子動(dòng)力學(xué)模擬研究[22]。首先使用其中的pdb2 gmx模塊構(gòu)建蛋白模型的拓?fù)湮募?蛋白的力場(chǎng)選擇GROMOS96 43A1聯(lián)合原子力場(chǎng)。然后將蛋白-配體復(fù)合物模型置于一個(gè)充滿TIP3P水模型的十二面體盒子中心,盒子邊緣距離蛋白分子邊緣為1 nm。向盒子中添加適量的Na+以平衡體系電荷。然后利用grompp模塊完成體系的能量最小化。先將體系至于NVT(T=300 K)條件下平衡100 ps,然后置于NPT條件下平衡100 ps。在NVT和NPT過(guò)程中,對(duì)蛋白和配體分子上的重原子加入位置約束(1 000 kJ/mol)[23]。最后,參考文獻(xiàn)[24-25],執(zhí)行無(wú)約束分子動(dòng)力學(xué)模擬,時(shí)長(zhǎng)16 ns。

2 結(jié)果與分析

2.1 構(gòu)建蛋白-配體復(fù)合物結(jié)構(gòu)模型

由圖2可知,本研究構(gòu)建的5-羧基-2-戊烯酰-輔酶A還原酶模型(藍(lán)綠色)具有典型的α/β折疊結(jié)構(gòu)。因?yàn)樵摻Y(jié)構(gòu)模型與模板4L1F_A處于相同的坐標(biāo)系,將4L1F_A中的FAD和COS配體分子復(fù)制到5-羧基-2-戊烯酰-輔酶A還原酶結(jié)構(gòu)模型中。需要指出的是,COS并不能參與反應(yīng),它不是一個(gè)可以參與反應(yīng)的底物。如圖3所示,COS與目的底物5-羧基-2-戊烯酰-輔酶A具有相同的母核結(jié)構(gòu),但其側(cè)鏈不同。COS的側(cè)鏈為巰基,而5-羧基-2-戊烯酰-輔酶A的側(cè)鏈為5-羧基-2-戊烯酰基。本文使用蛋白質(zhì)設(shè)計(jì)工具EvoEF2[26]的配體生成模塊將配體COS修改為底物5-carboxy-2-pentenoyl-CoA(圖2)。考慮到5-羧基-2-戊烯酰基側(cè)鏈具有較大的柔性,配體側(cè)鏈生成過(guò)程中對(duì)可自由旋轉(zhuǎn)的二面角進(jìn)行離散(-180°～180°,間隔30°離散),一共生成20 726個(gè)底物分子構(gòu)象。放置過(guò)程中，同時(shí)計(jì)算出每個(gè)底物構(gòu)象的內(nèi)能及其與蛋白質(zhì)骨架的范德華排斥能(兩者均為正值,越小越好)。接下來(lái)對(duì)所有構(gòu)象分別按照2種能量由低到高排序,共有1 728個(gè)構(gòu)象，在2種能量排序中都處于前25%。最后,按照適當(dāng)?shù)木礁嚯x偏差(root mean squared error,RMSD,單位為?)對(duì)1 728個(gè)構(gòu)象進(jìn)行篩選去重。本文采用的RMSD閾值為0.2 ?,篩選后得到48個(gè)底物分子構(gòu)象；其中與酶骨架排斥能最小的分子如圖2所示。

圖2 5-羧基-2-戊烯酰-輔酶A還原酶的模型Fig.2 Structure model of 5-carboxy-2-pentenoyl-CoA reductase in complex with FAD and 5-carboxy-2-pentenoyl-CoA注：藍(lán)綠色為5-羧基-2-戊烯酰-輔酶A還原酶,粉紅色球棍模型為FAD配體,淺藍(lán)色球棍模型為5-羧基-2-戊烯酰-輔酶A,綠色棍狀模型為結(jié)合殘基

圖3 三種輔酶A配體的比較Fig.3 Comparison of three CoA ligands注:過(guò)硫化輔酶A(CoA persulfide,COS)為PDB結(jié)構(gòu)4L1F中的配體；乙酰輔酶A(acetyl-CoA)為推測(cè)的天然底物；5-羧基-2-戊烯酰-輔酶A(5-carboxy-2-pentenoyl-CoA)為目的底物

同源建模過(guò)程中搜索到的結(jié)構(gòu)模板絕大部分被注釋為乙酰輔酶A脫氫酶(acyl-CoA dehydrogenase)。其中與5-羧基-2-戊烯酰-輔酶A還原酶序列最相近為Mycobacteriumthermoresistibile乙酰輔酶A脫氫酶(PDB ID:3PFD)[27]；兩者的序列相同性達(dá)到73.35%。據(jù)此推斷,T.fusca5-羧基-2-戊烯酰-輔酶A還原酶的天然底物也可能為乙酰輔酶A或其類似物。從圖3可以看出,目的底物5-羧基-2-戊烯酰-輔酶A與乙酰輔酶A的側(cè)鏈有較大差異,5-羧基-2-戊烯酰基較長(zhǎng),有更強(qiáng)的柔性,且包含強(qiáng)極性的羧基。從圖2可知,側(cè)鏈附近的結(jié)合殘基大多數(shù)為疏水殘基。因此,可推測(cè)這些疏水基團(tuán)不能與5-羧基-2-戊烯酰基側(cè)鏈很好地結(jié)合,從而導(dǎo)致酶催化5-羧基-2-戊烯酰-輔酶A反應(yīng)的活性較低。

2.2 計(jì)算設(shè)計(jì)改造5-羧基-2-戊烯酰-輔酶A還原酶

基于上述分析,本研究選擇與5-羧基-2-戊烯酰-輔酶A側(cè)鏈直接接觸(5?范圍內(nèi))的氨基酸殘基作為設(shè)計(jì)位點(diǎn)。這樣的殘基共有15個(gè),分別為Thr 161、Ile 90、Leu 89、Ser 88、Pro 91、Ala 92、Val 93、Lys 95、Leu 96、Ile 250、Gln 253、Thr 249、Tyr 367、Thr 246和Glu 368(圖4)。其中Glu 368為催化殘基,設(shè)計(jì)過(guò)程中不改變其殘基類型。此外,本文選擇距離5-羧基-2-戊烯酰基7?范圍內(nèi)的其他氨基酸殘基作為側(cè)鏈重新安裝的位點(diǎn)(未顯示)；這些殘基的側(cè)鏈構(gòu)象可改變,但氨基酸類型不變。從圖4可以看出,14個(gè)設(shè)計(jì)位點(diǎn)中只有Ser 88、Lys 95、Thr 161、Thr 246、Thr 249 和 Gln 253為親水性殘基,其余9個(gè)位點(diǎn)為疏水性殘基,并且所有殘基均未能與5-羧基-2-戊烯酰基側(cè)鏈的羧基官能團(tuán)形成氫鍵或鹽橋。

圖4 5-羧基-2-戊烯酰-輔酶A側(cè)鏈5?范圍內(nèi)的氨基酸殘基Fig.4 Amino acids within 5? to the side chain of 5-carboxy-2-pentenoyl-CoA

鑒于5-羧基-2-戊烯酰-輔酶A與乙酰輔酶A的差異(圖3),本研究嘗試通過(guò)計(jì)算設(shè)計(jì)引入氫鍵網(wǎng)絡(luò)來(lái)穩(wěn)定極性的5-羧基-2-戊烯酰基。設(shè)計(jì)過(guò)程通過(guò)計(jì)算蛋白質(zhì)設(shè)計(jì)工具EvoEF2[26, 28]完成。每個(gè)蛋白質(zhì)設(shè)計(jì)位點(diǎn)均可在20種天然氨基酸類型中自由選擇,每個(gè)氨基酸的側(cè)鏈構(gòu)象取自SHAPOVALOV等[29]開發(fā)的蛋白質(zhì)骨架依賴旋轉(zhuǎn)異構(gòu)體庫(kù)(backbone-dependent rotamer library)；而5-羧基-2-戊烯酰-輔酶A配體的空間結(jié)構(gòu)從前面生成的48個(gè)構(gòu)象中自由選擇。計(jì)算酶設(shè)計(jì)的過(guò)程相當(dāng)于求解組合優(yōu)化問(wèn)題,即從每個(gè)位點(diǎn)(包括配體分子)選擇一種殘基類型(或構(gòu)象),使得酶設(shè)計(jì)體系的整體能量達(dá)到最低。由于需要搜索的構(gòu)象空間巨大,EvoEF2采用模擬退火優(yōu)化算法來(lái)搜索低能量的氨基酸序列,但并不一定保證能取得全局能量最低序列。由于設(shè)計(jì)過(guò)程存在一定隨機(jī)性,每次設(shè)計(jì)的低能量序列可能并不完全相同,因此可以通過(guò)多次獨(dú)立運(yùn)行EvoEF2計(jì)算程序來(lái)產(chǎn)生許多備選序列進(jìn)行分析。EvoEF2自動(dòng)輸出每個(gè)設(shè)計(jì)的結(jié)構(gòu)和體系能量(單位為EvoEF2 energy unit,簡(jiǎn)記為EEU)。

如表1所示,本研究生成10個(gè)不同的設(shè)計(jì)(記為Des0～Des9)作為對(duì)照,對(duì)天然序列(WT)的設(shè)計(jì)位點(diǎn)的側(cè)鏈構(gòu)象進(jìn)行重新安裝,并計(jì)算其體系能量。不同設(shè)計(jì)的能量值存在一定差異,說(shuō)明EvoEF2每次優(yōu)化求解的過(guò)程可能收斂于不同的解,但能量差異并不大,說(shuō)明EvoEF2優(yōu)化算法收斂較好。10個(gè)設(shè)計(jì)的體系能量(-1462.58～-1452.97 EEU)明顯低于WT的能量(-1247.53 EEU)。這說(shuō)明天然酶與底物5-羧基-2-戊烯酰-輔酶A的結(jié)合作用可能較弱,而通過(guò)設(shè)計(jì)產(chǎn)生突變則可能引入一些有利的結(jié)合作用力。從表1來(lái)看,除了92、95和246位點(diǎn)十分保守(始終選擇對(duì)應(yīng)的WT氨基酸類型)之外,其余位點(diǎn)都是可變的(可能選擇非WT氨基酸類型)。

表1 基于計(jì)算改造5-羧基-2-戊烯酰-輔酶A還原酶產(chǎn)生的10個(gè)設(shè)計(jì)Table 1 Ten designs for computational engineering the 5-carboxy-2-pentenoyl-CoA reductase

仔細(xì)分析設(shè)計(jì)的結(jié)構(gòu)發(fā)現(xiàn)Des0、Des3、Des4和Des9這4個(gè)設(shè)計(jì)中底物分子的羧基可與氨基酸殘基形成氫鍵(或鹽橋)網(wǎng)絡(luò)。如圖5所示,Des0中形成了如下氫鍵：R253的胍基和Q250側(cè)鏈的胺基可與配體羧基的其中一個(gè)氧原子分別形成距離為3.0?和2.8?的鹽橋或氫鍵,S89的羥基可與R253的胍基形成距離為2.6?的氫鍵,并且R253胍基可與T364(非設(shè)計(jì)位點(diǎn))的羥基形成一個(gè)較弱的氫鍵(3.7?)。Des9中引入的氫鍵網(wǎng)絡(luò)與Des0類似；形成氫鍵網(wǎng)絡(luò)的氨基酸殘基完全相同,配體5-carboxy-2-pentenoyl-CoA的羧基側(cè)鏈構(gòu)象稍有不同。Des3和Des4中,位點(diǎn)89選擇為蘇氨酸(T)；與Des0和Des9相比,Des3和Des4中的T89取代S89與R253的胍基形成氫鍵。

圖5 Des0設(shè)計(jì)中的氫鍵網(wǎng)絡(luò)Fig.5 Hydrogen-bonding network in Des0

2.3 基于分子動(dòng)力學(xué)模擬檢驗(yàn)設(shè)計(jì)的氫鍵網(wǎng)絡(luò)

計(jì)算酶設(shè)計(jì)過(guò)程中假定蛋白質(zhì)骨架固定不變,但在實(shí)際酶反應(yīng)過(guò)程中,酶、配體及其相互作用都可動(dòng)態(tài)變化。本研究通過(guò)分子動(dòng)力學(xué)模擬來(lái)檢驗(yàn)設(shè)計(jì)的氫鍵在動(dòng)態(tài)過(guò)程中能否穩(wěn)定存在,以此來(lái)評(píng)估設(shè)計(jì)的優(yōu)劣。統(tǒng)計(jì)分析設(shè)計(jì)的氫鍵在整個(gè)動(dòng)力學(xué)模擬過(guò)程中的變化,氫鍵具有方向性,統(tǒng)計(jì)氫鍵一般考慮3個(gè)參數(shù)：XHY鍵角α,HY距離d1和XY距離d2(圖6)。在一定范圍內(nèi),α越大、d1和d2越小氫鍵越強(qiáng)。通常認(rèn)為α>130°,d1<2.2?,d2<3.2?可算作穩(wěn)定的氫鍵[24-25]。

圖6 氫鍵及相關(guān)幾何參數(shù)圖示Fig.6 Schematic of hydrogen bond and related geometric parameters

對(duì)分子動(dòng)力學(xué)模擬過(guò)程進(jìn)行統(tǒng)計(jì)分析表明,Des0設(shè)計(jì)中構(gòu)建的氫鍵網(wǎng)絡(luò)較穩(wěn)定。Des0中引入的4個(gè)氫鍵統(tǒng)計(jì)結(jié)果如圖7所示。

a～d-不同氫鍵的分析結(jié)果：第一列描述所要分析的氫鍵,第二列描述氫鍵距離隨模擬時(shí)間的變化情況,第三列描述氫鍵距離與角度的關(guān)系圖7 Des0設(shè)計(jì)中引入的氫鍵在分子動(dòng)力學(xué)模擬過(guò)程中的統(tǒng)計(jì)結(jié)果Fig.7 Statistical result of hydrogen bonds introduced in Des0 during 16 ns MD simulation

在16 ns的模擬過(guò)程中,5-carboxy-2-pentenoyl-CoA羧基的其中一個(gè)氧原子可以與R253的胍基以及Q250的胺基形成穩(wěn)定的氫鍵(圖7-a和b,d1平均值約為2?,α平均值約為160°)。相較之下,S89的羥基與R253的胍基之間的氫鍵較弱,d1平均值約為3?,α平均值約為110°(圖7-c)。這可能是因?yàn)門364的羥基可競(jìng)爭(zhēng)性地與R253的胍基形成較強(qiáng)的氫鍵(圖7-d,d1平均值約為2?,α平均值約為170°)。從圖5來(lái)看,Des0設(shè)計(jì)中R253可與S89形成的氫鍵較強(qiáng),而與T364的氫鍵較弱。在動(dòng)力學(xué)模擬過(guò)程中,這2個(gè)氫鍵的強(qiáng)弱雖然發(fā)生了轉(zhuǎn)換,但R253至少可與S89和T364其中之一形成較強(qiáng)的氫鍵。因此,分子動(dòng)力學(xué)模擬結(jié)果表明部分設(shè)計(jì)(如Des0)中的氫鍵網(wǎng)絡(luò)可以起到穩(wěn)定5-羧基-2-戊烯酰-輔酶A柔性側(cè)鏈的作用。這些設(shè)計(jì)可用于進(jìn)一步實(shí)驗(yàn)驗(yàn)證。

3 結(jié)論

菌株ThermobifidafuscaB6中的野生型5-羰基-2-戊烯酰基-輔酶A還原酶(5-carboxy-2-pentenoyl-CoA reductase)為己二酸生物合成途徑中的決速酶。其主要原因是目標(biāo)底物5-羰基-2-戊烯酰基-輔酶A與其天然底物乙酰輔酶A相比柔性更大、極性更強(qiáng)。本文嘗試采用基于結(jié)構(gòu)的計(jì)算酶設(shè)計(jì)方法對(duì)野生型5-羰基-2-戊烯酰基-輔酶A還原酶的結(jié)合位點(diǎn)進(jìn)行改造,試圖通過(guò)設(shè)計(jì)引入氫鍵網(wǎng)絡(luò)來(lái)更好地結(jié)合5-羰基-2-戊烯酰基-輔酶A的極性側(cè)鏈,目的是降低酶促反應(yīng)的活化能以提高催化效率。為檢驗(yàn)設(shè)計(jì)的合理性,本研究通過(guò)分子動(dòng)力學(xué)模擬來(lái)觀察10個(gè)設(shè)計(jì)中氫鍵網(wǎng)絡(luò)的變化情況,結(jié)果發(fā)現(xiàn)Des0設(shè)計(jì)中的氫鍵網(wǎng)絡(luò)很穩(wěn)定,與野生型酶相比,可增強(qiáng)與5-羰基-2-戊烯酰基-輔酶A的結(jié)合作用。因此,Des0突變可用于進(jìn)一步實(shí)驗(yàn)驗(yàn)證。