基于加權基因共表達網絡及分子動力學探討二妙散影響鐵死亡治療類風濕關節(jié)炎的機制*

肖劍偉，蔡旭，郭粉蓮，尹志華，陳新鵬，葉志中

深圳市福田區(qū)風濕病?？漆t(yī)院，廣東 深圳 518000

類風濕關節(jié)炎(rheumatoid arthritis，RA)已成為影響全球0.5%～1.0%人口的常見病之一[1]。作為一種慢性自身免疫性疾病，RA滑膜組織內細胞的死亡刺激可導致炎癥進一步加重，加快疾病進程。脂質過氧化誘發(fā)的鐵死亡是新近發(fā)現的一種非凋亡的細胞死亡方式，近十年來得到了廣泛的研究[2]。大量研究顯示，鐵死亡在多種疾病的發(fā)生發(fā)展過程中發(fā)揮著重要作用，如癌癥、神經系統(tǒng)疾病等[3-4]。淫羊藿苷能夠通過系統(tǒng)Xc-/谷胱甘肽過氧化物酶4(system Xc-/glutathione peroxidase 4，Xc-/GPX4) 軸抑制鐵死亡，從而增強脂多糖誘導的類風濕關節(jié)炎滑膜成纖維細胞活力[5]，提示鐵死亡在RA的發(fā)病中發(fā)揮著重要作用，調控鐵死亡過程成為一個新的治療RA的潛在靶點。

目前，常見的治療RA的藥物包括非甾體類、抗風濕病藥物、生物制劑等，雖然效果明顯，但不良反應也日益凸顯。隨著中藥在疾病中的應用越來越受到重視，從中藥中挖掘潛在的分子用來治療疾病已成為當前研究的熱點之一。二妙散是中醫(yī)藥治療痹癥的經典方劑，由蒼術和黃柏組成，具有治療濕熱下注、筋骨疼痛之作用。研究顯示，二妙散對RA患者具有明顯療效[6]，然而其作用機制尚不明確。本研究以加權基因共表達網絡(weighted gene co-expression network，WGCNA)分析為基礎，結合分子對接、分子動力學技術探索二妙散通過影響鐵死亡治療RA的作用機制。

1 資料與方法

1.1 數據材料以“ rheumatoid arthritis”為關鍵詞，在GEO數據庫(www.ncbi.nlm.nih.gov/geo)中檢索下載數據集。見表1。

表1 基因表達數據譜

1.2 加權基因共表達網絡分析選取方差>50%的基因，在R軟件中使用“WGCNA”包[7]進行分析。選擇合適的軟閾值(β)，將系統(tǒng)聚類樹(MEDissThres)設置為0.25以合并相似的模塊。將矩陣數據轉化為鄰接矩陣并基于拓撲重疊進行聚類，設置模塊至少含50個基因。將GSE77298數據集設置為正常組及RA組，構建臨床信息并計算基因顯著性。選擇相關性最大的模塊中的基因，通過string網站進行基因本體(gene oncology，GO)功能富集分析和京都基因與基因組百科全書(kyoto encyclopedia of genes and genomes，KEGG)信號通路富集分析。

1.3 二妙散活性成分篩選及作用靶點預測在中藥系統(tǒng)藥理學數據庫與分析平臺(traditional Chinese medicine systems pharmacology database and analysis platform，TCMSP)[8]檢索二妙散組方中藥黃柏、蒼術的化學成分，設定口服生物利用度(oral bioavailability，OB)≥30%、類藥性(drug-likeness，DL)≥0.18，篩選得到活性成分及其作用靶點，使用Cytoscape軟件繪制“二妙散-活性成分-活性成分靶點”網絡圖。

1.4 鐵死亡相關靶點及二妙散調控鐵死亡治療RA的潛在靶點獲取通過FerrDb數據庫[9]獲取鐵死亡相關基因。將結果與1.2中與臨床相關性最強的模塊基因及二妙散活性成分靶點進行映射取交集，獲得二妙散調控鐵死亡治療RA的潛在靶點。

1.5 治療靶點多數據集驗證分別在GSE55235和GSE55457數據集中，通過t檢驗分析二妙散調控鐵死亡治療RA的潛在靶點在正常組及RA組滑膜組織中的表達情況。

1.6 分子對接在RCSB PDB(www.rcsb.org/)獲取靶點的晶體結構，采用Pymol[10]去除晶體水后使用AutoDock對蛋白進行加氫及電荷。從Pubchem獲取活性小分子成分結構，使用ChemBio3D進行能量最小化。使用AutoDock Vina進行半柔性分子對接并選取結合能最小的復合物體系作為最佳配體。使用ligplot+計算和繪制復合物的氫鍵和疏水鍵。

1.7 分子動力學分子動力學模擬使用Gromacs 2020.3軟件進行[11]。在ATB網站[12]上傳化合物小分子結構，將其轉換為Gro格式并生成位置限制文件。使用Gromacs內置pdb2gmx命令分別生成靶蛋白坐標文件及拓撲文件。將小分子及蛋白的Gro格式合并為復合物的Gro文件。力場選用Gromos54a7_atb，水分子模型選用SPC，創(chuàng)建大小為 1.0 nm 的十二面體水盒子。根據計算結果加入鈉離子中和電荷，使模擬系統(tǒng)總電荷為0 C。將動力學模擬體系的溫度設定為300 K。采用最速下降法優(yōu)化后進行NVT和NPT 系綜平衡，最后進行總時長為20 ns的動力學模擬，步長設定為2fs。去除體系的周期性邊界并提取20 ns時的最終復合物構象。使用VMD[13]、Pymol進行可視化處理，計算復合物的均方根偏差(root mean squared deviation，RMSD)、溶劑可及表面(solvent accessible surface area，SASA)、蛋白回旋半徑(gyrate，GR)，并通過Gromacs插件 g_mmpbsa計算結合自由能[14]，通過拉氏圖評估動力學模擬后蛋白結構的合理性。

2 結果

2.1 WGCNA結果GSE77298數據集共涉及 10 827 個基因納入網絡的構建。β值一般選擇在平均連接度變化呈平穩(wěn)處，值在16以下，且R方要高于0.8，本研究選擇β=12作為軟閾值，見圖1，圖中縱軸代表對應網絡log(k)與log(p(k))相關系數的平方，相關系數的平方越高，說明該網絡越逼近無網絡尺度的分布；紅線處對應的上方的β值為本次最合適的β值，這意味著特征基因相關性高于0.9的模塊將被合并。因而選擇基于拓撲重疊以及分配的模塊顏色，對具有相似性基因的樹狀圖進行聚類，并對基因相關性高于0.9的模塊進行合并，見圖2，每個彩色行代表一個顏色編碼的模塊，其包含一組高度相關的基因，共確定12個模塊。對各模塊進行臨床信息相關性研究發(fā)現，黃色模塊的相關性為+0.57和-0.57，其中正值表示正相關，負值表示負相關，且黃色模塊的P值為0.004，即黃色模塊相關性最大且P值最小，見圖3及圖4，因而選取黃色模塊進一步進行GO功能富集分析及KEGG信號通路相關分析，見圖5，該黃色模塊共含有744個基因，GO功能富集分析發(fā)現92條分子功能、732條生物學過程及130條細胞成分，其中分子功能(molecular function，MF)涉及腫瘤壞死因子受體超家族結合等，生物學過程(biological process，BP)涉及免疫系統(tǒng)過程等，細胞成分(cellular component，CC)主要定位于囊泡等。KEGG信號通路富集分析發(fā)現75條信號通路，主要涉及腫瘤壞死因子-α(tumor necrosis factor，TNF-α)信號通路、核因子-κB(nuclear factor kappa B，NF-κB)信號通路等，與RA疾病的發(fā)生發(fā)展密切相關。

圖1 確定軟閾值圖

圖2 聚類分析合并后得到的模塊

注：每個方塊中上面的數字表示相關性(紅色方塊為正值表示正相關，綠色方塊為負值表示負相關)，下面括號中數字代表P值圖3 各模塊與RA的相關性

圖4 黃色模塊各基因與臨床性狀相關性

注：5A：GO富集分析氣泡圖；5B：KEGG信號通路富集分析氣泡圖圖5 GO功能富集分析和KEGG信號通路富集分析氣泡圖

2.2 二妙散活性成分及其相關靶點基于TCMSP篩選得到蒼術和黃柏有效活性成分分別為9個和37個，兩者共同作用靶蛋白103個，將蒼術、黃柏有效活性成分及相關靶點導入Cytoscape軟件繪制“二妙散-活性成分-活性成分靶點”網絡圖，見圖6。

注：橘紅色表示藥物；綠色表示有效活性成分；紫色表示靶點圖6 “二妙散-活性成分-活性成分靶點”網絡圖

2.3 二妙散調控鐵死亡治療RA的相關靶點通過FerrDb數據庫獲得鐵死亡相關基因97個，與WGCNA分析中黃色模塊包含的744個基因及103個二妙散活性成分相關靶點取交集，得到二妙散調控鐵死亡治療RA的相關靶點1個，即花生四烯5-脂氧合酶(arachidonate 5-lipoxygenase，ALOX5)。見圖7。

圖7 二妙散調控鐵死亡治療RA的潛在靶點的韋恩圖

2.4 多數據集驗證ALOX5表達結果在GSE55235與GSE55457數據集中分析ALOX5蛋白的表達量，結果顯示，與正常組比較，不同數據集的RA組滑膜組織中ALOX5蛋白表達水平均明顯升高(P<0.05)，見圖8。

圖8 滑膜組織中ALOX5蛋白的差異表達情況

2.5 分子對接結果黃柏、蒼術的有效活性成分中，槲皮素、漢黃芩素共同作用的靶點為ALOX5，因而，選擇槲皮素及漢黃芩素與ALOX5進行分子對接。結果顯示，槲皮素、漢黃芩素均能夠與ALOX5形成多個氫鍵及疏水鍵，其中槲皮素與ALOX5的對接結合能為-9.4 kcal·mol-1，與Leu448、Gln549、Ala453、Ser447、Phe450、Arg457、Tyr470、Val243、Ile365、Thr366形成疏水作用；漢黃芩素與ALOX5的對接結合能為-9.1 kcal·mol-1，與Ile167、Arg401、Glu622、Asp166、Try100、His624、Phe393、Arg101形成疏水作用。見圖9。

圖9 分子對接示意圖

2.6 分子動力學結果槲皮素-ALOX5及漢黃芩素-ALOX5 2個復合體系在長達20 ns的整個模擬過程中，RMSD波動小于0.2 nm并且在5 ns時達到平衡狀態(tài)，見圖10。復合物體系的GR越小，提示復合物在水溶液中的結構越穩(wěn)定，見圖11。SASA結果顯示，隨著時間延長，2個復合物體系的SASA逐漸縮小，提示其隨著時間延長而逐漸趨于穩(wěn)定，見圖12。與模擬前對比，經過20 ns的動力學模擬后，漢黃芩素、槲皮素依然與ALOX5緊密結合。漢黃芩素與ALOX5結合后，小分子在活性口袋內部移動變化幅度較大且其Gyrate變化較明顯，考慮可能為活性小分子使蛋白質的三級結構發(fā)生改變，抑制了靶蛋白的活性，見圖13。拉氏圖顯示，與槲皮素和漢黃芩素結合后，ALOX5蛋白99.7%和99.8%的氨基酸殘基均處在合理范圍內，提示模擬后的蛋白結構是可靠而且合理的，見圖14。結合自由能結果顯示，槲皮素與ALOX5之間主要通過范德華力、靜電勢能結合；漢黃芩素與ALOX5主要通過范德華力結合；極性溶劑化作用則抑制槲皮素、漢黃芩素與ALOX5的結合，見表2。

圖10 小分子蛋白復合物體系的RMSD

圖11 小分子蛋白復合物體系的蛋白回旋半徑

表2 槲皮素、漢黃芩素與ALOX5的結合自由能

表2 槲皮素、漢黃芩素與ALOX5的結合自由能

化合物受體蛋白范德華力靜電勢能極性溶劑化作用表面溶劑化作用結合自由能槲皮素ALXO5-113.91±13.71-139.82±19.95134.90±15.06-16.26±0.86-135.10±16.46漢黃芩素ALXO5-130.81±10.68-17.23±9.5055.59±20.67-16.20±1.01-108.66±21.53

圖12 小分子蛋白復合物體系的溶劑可及表面

圖13 分子動力學模擬示意圖

圖14 分子動力學模擬后ALOX5蛋白拉氏圖

3 討論

RA是最常見的慢性炎癥性疾病，反復的炎癥反應會導致滑膜增厚，并引起關節(jié)腫脹[15]，越早啟動治療，其癥狀緩解的可能性就越大[16]。炎癥通常是由有害刺激引發(fā)的，包括感染、細胞死亡和組織損傷等[17]。細胞凋亡和壞死是細胞死亡的主要形式，而鐵死亡是一種依賴于鐵和脂毒性的非caspase依賴性細胞死亡形式[18]，它是由非酶促(fenton)反應和脂氧合酶(lipoxygenase，LOXs)引發(fā)，并由鐵催化脂質過氧化反應[19]。細胞死亡常伴隨著炎癥反應，而在炎癥發(fā)生的過程中，花生四烯酸作為細胞膜脂質的主要成分，通過3條代謝途徑代謝為生物活性促炎介質前體，其中一條為LOXs途徑，LOXs轉化為白三烯和脂氧素參與炎癥反應[20]。LOXs家族成員之一ALOX5是觸發(fā)脂質過氧化的關鍵酶[21]。研究表明，ALOX5抑制劑Zileuton通過抑制谷氨酸誘導的HT22細胞鐵死亡起到神經保護作用[22]。在脂肪性肝炎小鼠模型中，鐵死亡抑制劑羅格列酮和去鐵螯合劑去鐵酮能夠改善肝細胞死亡并抑制炎癥因子的表達[23]。綜上可見，鐵死亡可能是炎癥的啟動因素，或者至少具有促炎作用，抑制鐵死亡可以抑制炎癥反應。最近的研究發(fā)現，許多鐵死亡抑制劑(包括Fer-1、liproxstatin-1和LOXs抑制劑)可防止在鐵死亡期間fenton反應驅動的多不飽和脂肪酸的自氧化和非酶促反應的破壞[24-25]。本研究發(fā)現，在RA滑膜組織中ALOX5蛋白表達水平明顯升高，由此推測ALOX5可能導致滑膜組織鐵死亡增加，導致炎癥進一步加重，因此抑制滑膜組織的鐵死亡，可能能夠減輕炎癥，改善RA患者癥狀。

在RA疾病活動期，患者以多關節(jié)疼痛為主要臨床表現，中醫(yī)認為邪氣盛則實，其腫、熱屬于“邪實”的表現，導致滑膜循環(huán)血流淤滯，實驗室檢查可見炎癥指標明顯升高，推測局部血管循環(huán)的變化是RA的基本病理因素之一；同時，炎癥過程中產生的各種炎癥介質、自身抗體等可作為毒邪入侵而影響氣血運行，導致脈絡淤滯。上述這類變化與中醫(yī)的“瘀”一致，故中醫(yī)治療RA以攻邪為主，扶正為輔，給予行瘀、清熱、除濕之法，與西醫(yī)的抗炎治療具有相同的含義。二妙散治療RA已取得了明確的療效，可減少IL-1β分泌，對炎癥產生抑制作用[26-27]，提示二妙散能夠通過減輕炎癥反應改善RA患者的癥狀。而本研究亦顯示，二妙散可能通過影響鐵死亡來抑制炎癥反應，以上研究結果與二妙散清熱燥濕化瘀之功效具有異曲同工之處。

本研究結果顯示，漢黃芩素、槲皮素與ALOX5之間的結合自由能較低，經過動力學模擬后仍與ALOX5緊密結合。研究顯示，漢黃芩素通過抑制ALOX5減少乳腺癌細胞的侵襲和炎癥[28]?；⒄忍崛∥镏婚纹に啬軌蛞种颇[瘤細胞中ALOX5的活性[29]。張良等研究顯示，槲皮素能夠通過抑制前列腺組織內ALOX5表達改善自身免疫性慢性前列腺炎[30]。根據以上研究推測，漢黃芩素、槲皮素均能使ALOX5活性降低，從而減少鐵死亡，減輕炎癥反應。二妙散中黃柏苦寒清除濕熱，為治下焦?jié)駸嶂?；蒼術苦溫辛燥，兩者配伍能夠寒溫協(xié)調，陰陽相濟，從而達到清熱燥濕，標本兼顧的治療效果。而槲皮素作為黃柏的主要成分，與ALOX5具有較低的結合自由能，在治療中發(fā)揮著更大的作用，這與二妙散中黃柏作為主藥是相符合的。二妙散為復方中藥湯劑，在治療中具有多組分、多靶點、多途徑的特點，不同的成分也可作用于同一靶點，發(fā)揮協(xié)同作用。本研究從鐵死亡角度闡述了二妙散可能通過影響滑膜組織中ALOX5的活性，減少鐵死亡，進而降低滑膜中的炎癥反應，為進一步驗證和開發(fā)提取中藥有效成分治療RA提供了理論依據。

綜上，基于加權基因共表達網絡分析和分子動力學探討發(fā)現二妙散可通過調控滑膜組織中ALOX5的活性減少鐵死亡進而改善類風濕關節(jié)炎。然而，本研究結果需要進一步通過使用體內外實驗檢測來驗證。