雞胸肌條紋/木質表型相關基因的篩選

樊慶燦，鄭路程，王倩倩，劉 犇,2，楊 雪，胡 威,2，鄭文亞

（1.宜春學院生命科學與資源環境學院，江西宜春 336000；2.江西綠科農牧科技有限公司，江西宜春 336000）

過去幾年中，快大型肉雞的胸肌中白條紋肉（White Striping，WS）和木質肉（Woody Breast，WB）異常表型不斷上升。WS 存在數量不等的與肌纖維平行的白色條紋[1]。WB 呈木質化，發病嚴重的肌肉表面可同時出現白色條紋和彌漫分布的脊狀突起[2]，嚴重影響雞肉的屠宰外觀和肉品質。WS 和WB 表型有所不同，但肌肉組織均有所損傷。組織學評估表明，異常樣品均呈現肌肉表層的退化，包括纖維壞死和溶解、異常纖維化和脂質增多[3]。為了闡明WS/WB 表型發展的分子機制，Zambonelli 等[4]比較了的Ross 708 肉雞WS/WB胸肌和正常胸肌基因表達模式，鑒定出的差異表達基因（DEG）與肌肉發育、炎癥、多糖代謝和鈣信號傳導有關。Malila 等[5]則篩選出鈣信號傳導和Ras 相關蛋白1信號傳導等8 個信號通路與WS/WB 的產生有關，并發現細胞內Ca2+等離子、氧化應激和細胞程序性死亡引起WS/WB 表型。本研究利用GEO2R 軟件的相同標準分別分析2 組表達譜數據，統計快大型肉雞正常肌肉表型和WS/WB 肌肉表型的差異表達基因，以進行共同基因分析，并對篩選到的共同基因進行功能注釋和功能富集分析，探討WS/WB 肌肉表型的發生機理，為雞屠宰和肉質性狀的研究提供參考。

1 材料和方法

1.1 數據及樣本 本研究數據來自于美國國家生物信息數據庫中的GEO 數據庫，編號分別為GSE79276 和GSE107362，2 個數據系列的實驗動物分別為羅斯708和羅斯308 的公雞，屠宰日齡均為49 d，樣品均為胸大肌，樣品表達量檢測分別采用GPL21385 和GPL24307平臺。具體信息見表1。

1.2 軟件與方法 打開GEO 數據庫中的分析工具GEO2R（https://www.ncbi.nlm.nih.gov/geo/geo2r/），輸入GSE 79276，按照表1 設定正常表型組（NM）和WS/WB 表型組（WS/WB），選擇Benjamini &Hochberg（False discovery rate）矯正，數據進行log2轉換，進行差異表達分析。下載差異分析結果、數據表達量圖和差異分析火山圖，將Padj＜0.05 且|log2FC|＞1 的轉錄本設為差異表達轉錄本（Differential Expression Transcriptions，DET）。根據GPL21385 平臺信息，查找DET 對應的基因，設為差異表達基因（Differential Expression Gene，DEG）。刪除無基因信息的DETs。在DAVID V6.8（https://david.ncifcrf.gov/tools.jsp）中進行官方基因名轉換。GSE107362 數據分析步驟如上。將GSE79276 和GSE107362 中的DEGs 輸入到Bioinformatics &Evolut ionary Genomics（http://bioinformatics.psb.ugent.be/webtools/Venn/）中繪制韋恩圖，下載共同差異表達基因（Identical Differential Expression Gene，IDEG）。利用David v6.8 進行IDEGs 的信號通路分析。將IDEGs 導入到String（https://string-db.org/cgi/input.pl?sessionId=NPdugdyy2MoP&input_page_show_search=on）中進行互作分析，并下載蛋白互作文件，導入Cytoscape v3.7.2 中繪制互作圖，cytoHubba 插件MHC 算法計算核心基因，BINGO 插件進行生物過程的富集分析。

表1 數據及樣本信息統計表

2 結果

2.1 2 個數據系篩選到的DEGs 2 個數據系的DETs 結果見圖1。Log2FC＜0 表示WS/WB 胸肌中基因表達量下調。表2 中為2 個數據系的差異DEGs。GSE79276數據系共篩選出104 個DEGs（Padj＜0.05），其中17個上調，87 個下調。下調基因白介素3 調節核因子基因（Nuclear Factor,Interleukin 3 Regulated，NFIL3）是P值最小的基因，下調基因半胱氨酸和甘氨酸富集蛋白質3 基因（Cysteine and Glycine-Rich Protein 3，CSRP3）是差異倍數最大的基因。GSE107362 數據系共 篩 選 出953 個DEGs（Padj＜0.05），其中179 個 上調，774 個下調。上調基因細胞因子信號傳導抑制因子1（Suppressor Of Cytokine Signaling 1，SOCS1）是P值最小的基因，上調基因B 細胞表面抗原CD40（CD40 Molecule，CD40）基因是差異倍數最大的基因。

表2 2 個數據系差異基因統計表

圖1 2 個數據系差異表達轉錄本火山圖

2.2 2 個數據系篩選到的IDEGs 2 個數據系IDEG 分析見圖2。2 組DEGs 共篩選出23 個IDEGs，且所有基因均下調，基因相關信息見表3。NFIL3和CSRP3分別為GSE79276 數據系差異分析的P值最小和差異倍數最大的基因。這些基因中，纖連蛋白1（Fibronectin，FN1）、VI 型膠原蛋白α3 鏈編碼基因（Collagen Type VI Alpha 3 Chain，COL6A3）和血小板反應蛋白1 基因（Thrombospondin 1，THBS1）為篩選到的核心基因，3 個基因間存在互作關系（圖3）。

表3 2 個數據系的IDEGs 統計信息

圖3 部分有互作關系的IDEGs

2.3 IDEGs 的聚類和信號通路分析 23 個IDEGs 的聚類分析結果見圖4a。GO 分析共篩選出6 條“P＜0.05”的條目，其中生物黏附（Biological Adhesion）和細胞黏附（Cell Adhesion）是P值最小的條目，參與2個條目的基因均為FN1、COL6A3、THBS1和血小板反應蛋白2（Thrombospondin 2，THBS2）。剩余的4 個條目中，3 個條目為信號傳導相關條目，最后一個為生物調節條目。23 個IDEGs 的信號通路分析結果見圖4b。信號通路分析共篩選出3 個“P＜0.05”的信號通路，其中細胞外基質受體互作通路（ECM-receptor interaction）是P值最小的信號通路，其次為黏附斑信號通路（Focal Adhesion），參與2 個信號通路的基因均為軟骨寡聚基質蛋白（Cartilage Oligomeric Matrix Protein，COMP），I 型膠原蛋白α3 鏈編碼基因（Collagen Type VI Alpha 2 Chain，COL1A2）、FN1、COL6A3、THBS2和THBS1。最后一個通路為吞噬體信號（Phagosome）通路，參與該通路的基因為COMP、THBS2、THBS1基因。

圖4 IDEGs 的GO 分析和信號通路分析結果

3 討論

3.1 差異表達分析方法 本文使用GEO2R 軟件進行分析，與Zambonelli 等[4]和Malila 等[5]分析方法有所差異。生物信息學盡管有多種方法，但并沒有公認的最準確的方法，只是從不同角度盡量挖掘轉錄組數據中的更多信息，以便進行后續的功能富集分析。利用GEO2R 軟件同樣的標準和方法對2 組數據進行差異表達分析，一方面利用GEO2R 挖掘轉錄組數據中的更多信息，另一方面盡量減少分析方法等對結果的影響，使篩選到的共同基因更為準確可靠。

3.2 差異表達基因分析結果 本研究對2 個數據系進行轉錄本差異表達分析，發現WS/WB 表型與正常表型胸肌肉相比，下調的DEGs 的數目均遠大于上調基因的數目，說明在WS/WB 表型雞中，胸肌部分正常基因的表達降低。韋恩圖分析篩選到的基因均為下調，也進一步證明了上述內容，即營養、氧化應激和缺氧等因素導致快大型雞胸肌中正常的基因表達受到抑制[6]，從而導致雞胸肌不能維持正常的生理狀態和功能，出現肌纖維壞死和脂質增多等表型。

3.3 IDEGs 的聚類分析 IDEGs 的GO 聚類分析結果中，篩選出細胞黏附和信號傳導相關的5 個條目。信號通路分析中發現基質受體互作和黏附斑信號通路，二者結果一致。細胞外基質（ECM）由結構和功能大分子的復雜混合物組成，細胞和ECM 之間的特異性相互作用是由跨膜分子整合素介導的[7]。黏附斑信號通路是一種依賴于黏附斑激酶（Focal Adhesion Kinase 1，FAK）的信號通路，該通路利用整合素介導ECM 和細胞之間的黏附作用，并刺激幾種細胞內信號傳導途徑的活性[8]。研究表明，黏附斑信號通路與成肌細胞分化、肌肉纖維形成和肌肉大小的發育狀態相關[9]，并可以通過整合素翻譯跨細胞質膜傳遞的細胞骨架應力和應變信號，引起肌動蛋白細胞骨架的重組[10]。基質受體互作和黏附斑信號通路可能介導了胞外異常信號進入胞內，引入肌細胞的異常代謝。

3.4 IDEGs 功能分析 23 個IDEGs 中，COMP、FN1、COL1A2、COL6A3、THBS2、THBS1基因與細胞黏附和信號傳導有關。研究發現這些基因還可能參與了肌肉發育和調控。COMP 蛋白是一種在肌肉骨骼和心血管系統中都存在的細胞外基質非膠原糖蛋白，研究發現，COMP

可維持線粒體穩態并控制平滑肌細胞的表型特征[11]。FN1纖連蛋白可以結合肌細胞和肌動蛋白表面，參與細胞運動和維持細胞形態[12]。COL1A2和COL6A3均編碼膠原蛋白多肽鏈。膠原蛋白家族成員是在肌肉發育中具有重要調控作用的ECM 蛋白。XIII 型膠原蛋白能夠調節神經肌肉接頭正確組裝[13]。XXV 膠原蛋白對于成肌細胞融合到肌纖維中是必需的[14]。血小板反應蛋白（THBS）是參與多種生物學過程的基質細胞蛋白。THBS 家族每個成員受組織生長、愈合和重塑等誘導。研究發現，肌纖維中THBS 的表達可穩定肌膜[15]。這些信號蛋白的下調可能會影響肌細胞的正常形態和功能。

IDEGs 中，膜聯蛋白A1（Annexin A1，ANXA1）、CSRP3、肌球蛋白2（Myozenin 2，MYOZ2）、ATPase胞質/ 內質網Ca2+轉運2（ATPase Sarcoplasmic/Endop lasmic Reticulum Ca2+，ATP2A2）直接參與肌肉發育與功能的調節。ANXA1 蛋白參與和調節糖皮質激素介導的炎癥反應，促進肌動蛋白細胞骨架的重排[16]，并介導吞噬體和肌動蛋白細胞骨架之間的相互作用[17]。CSRP3 蛋白是肌發生的正調節劑，用作肌源性bHLH（Basic Helix-Loop-Helix）轉錄因子，直接結合肌動蛋白絲，將肌動蛋白絲交聯成束，防止絲切蛋白介導的解聚反應[18]。MYOZ2 是一種細胞內結合蛋白，在肌原纖維形成中起重要作用，參與連接α-肌動蛋白，并將鈣調神經磷酸酶信號轉導至肌小節，參與肌肉活動的調節[19]。ATP2A2基因編碼心肌肌網Ca2+-ATP 酶，它是位于骨骼肌的肌質或內質網中的細胞內泵。該酶催化ATP 的水解，并促進鈣從胞質向肌漿網腔的轉運，參與收縮/松弛周期的調節[20]。這4 個基因的下調，可能直接影響肌動蛋白的形態和正常功能，導致肌肉表型異常。

IDEGs 中，肌 醇5-磷 酸 酶2（Synaptojanin 2，SYNJ2）參與不同的膜運輸和信號轉導途徑，對網格蛋白介導的內吞具有抑制作用[21-22]，可能會導致攜帶膽固醇的低密度脂蛋白（LDL）大量進入細胞。載脂蛋白A1（Apolipoprotein A1，APOA1）是血漿中高密度脂蛋白（HDL）的主要蛋白成分，HDL 具有抗氧化和保護LDL 的作用，并促進胞內膽固醇和甘油三酯等脂類外排[23]。ANXA1介導膽固醇從內質網到多囊泡體（MVB）的轉運，并刺激MVB 分泌具有保護和運輸脂質功能的外泌體[24]。推斷SYNJ2、APOA1、ANXA1表達量的下降可能引起膽固醇等脂質在肌肉細胞的積累和破壞，誘導炎性反應[25]，產生大量吞噬LDL 和氧化型低密度脂蛋白（oxLDL）的巨噬細胞，最終形成泡沫細胞并破裂，在這個過程中將脂蛋白降解為各種脂質[26]，這可能是發生炎性細胞浸潤和脂質沉積的因素之一。微管相關蛋白1B 基因（Microtubule Associated Protein 1B，MAP1B）編碼微管相關蛋白。該家族的蛋白質被認為與微管組裝有關，微管與其他纖維一起構成細胞骨架，微管雙螺旋結構支撐著細胞生理形態[27]。該基因表達量的下調可能引起肌肉細胞形態的變化。

綜上所述，WS/WB 表型的產生可能是由于各種因素導致肌肉組織接受到異常信號，經整合素介導后，通過黏附斑信號通路引起胞內相關基因表達的異常，最終導致胸肌細胞功能異常，肌肉纖維和肌動蛋白的結構和功能受到破壞，引起異常表型。

4 結論

本研究利用GEO 數據的2 個數據系篩選WS/WB表型的差異表達基因，篩選到23 個共同表達基因，功能富集分析發現部分基因富集在信號傳導和細胞黏附通路，參與ECM 受體互作通路和黏附斑信號通路。此外，23 個IDEGs 中，ANXA1、CSRP3、MYOZ2和ATP2A2基因直接參與肌動蛋白結構組成和功能調節，COMP和FN1等6 個基因通過調節肌肉組織信號傳導和細胞黏附活動，維持肌肉正常生理功能。