王鵬飛，丁燕玲，楊朝云，馬 應(yīng)，周小南，史遠(yuǎn)剛，康曉龍

（寧夏大學(xué)農(nóng)學(xué)院，銀川 750021）

飼料成本在家畜養(yǎng)殖成本中占比最高，提高飼料利用率、降低養(yǎng)殖成本為家畜規(guī)模化養(yǎng)殖及提高效益重要環(huán)節(jié)。剩余采食量（Residual feed intake，RFI）為家畜實際采食量與其維持生產(chǎn)（產(chǎn)蛋、產(chǎn)奶等）的預(yù)期采食量之差[1]，為評價動物飼料轉(zhuǎn)化率及代謝水平指標(biāo)[2]，本試驗開展家畜RFI表型研究，為高飼料效率表型的家畜選育提供理論基礎(chǔ)。

環(huán)狀RNA（circRNAs）為一種非編碼RNA，由線性RNA反向剪接產(chǎn)生[3]，與其他非編碼RNA相比，circRNAs具有封閉圓形結(jié)構(gòu)，更穩(wěn)定，不易降解[4-5]。circRNAs由內(nèi)含子和外顯子構(gòu)成，可調(diào)節(jié)多種生物學(xué)功能，如疾病、肌肉發(fā)育、脂肪代謝等[6]。circRNAs最重要功能為充當(dāng)miRNA的“海綿”，消除miRNAs對靶基因抑制[7]，如在沼澤水牛乳腺上皮細(xì)胞中，circ_015003和circ_014121的海綿miR-103過表達(dá)增加乳脂合成相關(guān)基因表達(dá)，導(dǎo)致脂滴形成及甘油三酯累積[8]；ssc_circ_0002807和ssc_circ_0009352為miRNA-433的海綿，miRNA-433能激活p38MAPK通路，而激活的p38MAPK促進(jìn)骨髓間充質(zhì)干細(xì)胞脂肪生成[9]。研究表明，miRNA參與調(diào)控動物機(jī)體采食及能量代謝，如miRNA-499和miRNA-15a可分別通過其靶基因TCF12和FOXO1介導(dǎo)AMPK-PGC-1α、TGF-β及PI3KAKT-mTOR等信號通路增強(qiáng)機(jī)體對營養(yǎng)物質(zhì)吸收利用，調(diào)控動物采食[10]。circRNAs可直接參與該生物學(xué)過程，或間接通過與其他分子（如miRNAs等）互作共同調(diào)控目標(biāo)性狀，表明circRNAs在調(diào)控表型變異過程中具有廣泛且獨特的生物學(xué)功能。

目前，circRNAs在畜禽主要表型性狀研究集中于脂肪沉積[11]、肌肉發(fā)育[12]、泌乳性能[13]等，但circRNAs是否參與調(diào)控牛RFI表型變異，或在牛營養(yǎng)物質(zhì)吸收和能量代謝過程中是否具有潛在生物學(xué)作用鮮有報道。下丘腦為調(diào)節(jié)食欲和能量平衡關(guān)鍵中樞，通過調(diào)節(jié)特定神經(jīng)遞質(zhì)/神經(jīng)肽表達(dá)對營養(yǎng)物質(zhì)/激素信號作出反應(yīng)，調(diào)控能量攝入和消耗[14]?；⌒魏耍ˋrcuate nucleus，ARC）被認(rèn)為是下丘腦控制食物攝入量主要中心，其兩個不同神經(jīng)元群體接收外周營養(yǎng)信號，通過不同食欲肽表達(dá)參與攝食調(diào)節(jié)，如促食類神經(jīng)肽有神經(jīng)肽Y（Neuropeptide Y，NPY）、刺鼠相關(guān)蛋白（Agoutirelated protein，AGRP）等。研究發(fā)現(xiàn)，活化蛋白激酶（AMP-activated protein kinase，AMPK）在下丘腦神經(jīng)元中扮演重要角色，作為身體能量穩(wěn)態(tài)的傳感器和調(diào)節(jié)器，將代謝狀態(tài)與神經(jīng)肽系統(tǒng)聯(lián)系起來，調(diào)節(jié)攝食[15]。鑒于上述背景，本研究采用高通量測序技術(shù)鑒定牛下丘腦中RFI表型相關(guān)差異circRNAs，并對靶基因富集分析，旨在揭示circRNAs在牛日糧攝入、能量代謝和消化吸收中潛在作用，為后續(xù)全面探究circRNAs在動物主要經(jīng)濟(jì)性狀表型變異中分子功能及調(diào)控作用提供參考依據(jù)。

1 材料與方法

1.1 研究對象

選取健康狀況良好初始體重相近（280±30）kg，月齡相近（16±1）mo高RFI（High RFI，HRFI）、低RFI（Low RFI，LRFI）秦川牛各5頭。試驗牛頸部放血屠宰開顱后，立即從下丘腦前段剪取1 cm×1 cm下丘腦組織，生理鹽水沖洗裝入凍存管，立即放入液氮保存。HRFI、LRFI牛均設(shè)5個重復(fù)，分別編 號 為HRFI_1、HRFI_2、HRFI_3、HRFI_4、HRFI_5和LRFI_1、LRFI_2、LRFI_3、LRFI_4、LRFI_5。依據(jù)標(biāo)準(zhǔn)方法計算RFI[16]。試驗牛飼養(yǎng)方式為單欄飼喂，每日兩次（8:00和15:00），采食1.5 h，自由飲水。

1.2 RNA提取及文庫構(gòu)建

總RNA提取采用Trizol法。對提取的RNA作濃度、純度（OD230，OD260和OD280）及完整性檢測，評估總RNA質(zhì)量。總RNA的OD260/OD280值≥1.8，完整性值在8.6以上，表明總RNA質(zhì)檢合格，用于文庫構(gòu)建。合格RNA樣品使用Covaris超聲波破碎儀隨機(jī)打斷成長度約350 bp片段，經(jīng)末端修復(fù)、加A尾、加測序接頭、純化、PCR擴(kuò)增等步驟完成文庫制備。使用Qubit2.0作初步定量，稀釋文庫至2 ng.μL-1，隨后使用Agilent 2100檢測文庫insert size，符合預(yù)期后，使用qRT-PCR方法準(zhǔn)確定量文庫有效濃度（文庫有效濃度高于3 ng.μL-1），保證文庫質(zhì)量。文庫質(zhì)檢合格后，用Illumina平臺作雙端測序，得到用于分析的原始測序片段（raw reads）。

1.3 差異表達(dá)circRNAs篩選

為確保數(shù)據(jù)分析可靠性，使用fastx_tool-kit軟件（v0.0.14）去除raw reads中帶接頭、低質(zhì)量序列及接頭污染等，得到clean reads。采用BWA中BWA-MEM算法將clean reads比對到參考基因組中，采用CIRI識別circRNAs，比對結(jié)果以SAM文件輸出后掃描該文件，找到PCC（Paired chiastic clipping）和PEM（Paired-endmapping）位點及剪接信號GT-AG，最后為確保識別circRNAs可靠性，將具有junction位點的序列重新比對。用TPM（Transcripts per million）對circRNAs表達(dá)量作歸一化處理，DEGseq分析樣本間circRNAs表達(dá)差異，以|Fold Change|≥1和P＜0.05作為表達(dá)標(biāo)準(zhǔn)確定差異表達(dá)的circRNAs。

1.4 差異表達(dá)circRNAs功能富集分析

差異表達(dá)circRNAs的GO及KEGG富集分析采用DAVID數(shù)據(jù)庫；GO富集分析包括3個層次，分為生物學(xué)過程（Biological process，BP）、細(xì)胞組分（Cellular component，CC）、分子功能（Molecular function，MF）。KEGG數(shù)據(jù)庫是將差異表達(dá)circRNAs所對應(yīng)基因作通路富集，確定來源基因主要參與的代謝途徑和信號通路。

1.5 差異表達(dá)circRNAs驗證

為驗證測序結(jié)果準(zhǔn)確性，以GAPDH為內(nèi)參基因，隨機(jī)選取6個差異表達(dá)circRNAs（novel_circ_0031230、novel_circ_0038739、novel_circ_0022168、novel_circ_0005798、novel_circ_0018798、novel_circ_0034790）作qRT-PCR驗證，引物由上海生工合成（見表1）。取RNA樣本反轉(zhuǎn)錄為cDNA，以此為模板作qRT-PCR檢測，具體操作方法依照TaKaRa試劑盒說明書。qRT-PCR反應(yīng)程序為95℃預(yù)變性30 sec，95℃變性10 sec，退火30 sec，4℃保存，共40個循環(huán)。采用公式2-ΔΔCt計算差異circRNAs相對表達(dá)量。

表1 基因及對應(yīng)引物Table 1 Gene and corresponding primer sequences

2 結(jié)果與分析

2.1 下丘腦組織circRNAs測序數(shù)據(jù)質(zhì)量分析

去除測序數(shù)據(jù)接頭序列、低質(zhì)量及錯配reads等，得到HRFI、LRFI組中clean reads（見表2）。由表2可知在所測10個個體中clean reads占總raw reads比例均在97%以上，表示測序數(shù)據(jù)質(zhì)量良好，可用于后續(xù)分析。

表2 數(shù)據(jù)過濾統(tǒng)計情況Table 2 Data filtering statistics

2.2 circRNA來源及長度統(tǒng)計

circRNAs可分為exon、intergenic、intron 3種類型，通過與參考基因組比對，發(fā)現(xiàn)circRNAs在exon中占比最高（約87%），intron中占比最少（約5%），結(jié)果見圖1A。證實circRNAs為閉合環(huán)狀RNA，大部分由外顯子構(gòu)成，故無多聚A尾巴；長度統(tǒng)計表明大部分circRNAs分布在100～500 nt（見圖1B）。

圖1 circRNAs來源（A）及長度（B）Fig.1 Source(A)and length(B)of circRNAs

2.3 差異circRNAs篩選

采用TPM對各樣本中已知和新circRNAs表達(dá)量作歸一化處理（見圖2A），通過差異倍數(shù)（|log2（Fold change）|＞1）和顯著水平（P＜0.05）兩個水平篩選兩組中差異表達(dá)circRNAs，差異基因分布情況以火山圖表示（見圖2B），從HRFI組中篩選得到736個與LRFI差異表達(dá)的circRNAs，其中表達(dá)上調(diào)有315個，下調(diào)有421個。差異顯著上調(diào)和下調(diào)前10個circRNAs相關(guān)信息見表3，其中novel_circ_0031230與bta-miR-103具有靶向關(guān)系。據(jù)報道，bta-miR-103及其靶基因（Dicer、PANK）可調(diào)節(jié)脂肪生成[17]等。

表3 顯著差異前10的circRNAsTable3 Top 10 significant differential expressed circRNAs

2.4 差異表達(dá)circRNAs靶基因GO富集分析

對736個差異表達(dá)circRNAs靶基因作GO分類（見圖3）。

由圖3可知，在BP中，差異表達(dá)circRNAs對應(yīng)基因主要涉及細(xì)胞器組織、細(xì)胞定位、蛋白質(zhì)定位及大分子定位等；在MF中主要涉及GTP酶結(jié)合、Ras GTP酶結(jié)合、酸性胺酸連接酶活性等過程；在CC中主要富集在細(xì)胞器、細(xì)胞內(nèi)膜、細(xì)胞質(zhì)及細(xì)胞核組成等；其中細(xì)胞內(nèi)細(xì)胞器、膜結(jié)合細(xì)胞器和細(xì)胞內(nèi)膜細(xì)胞器注釋高，說明與細(xì)胞器生成密切相關(guān)，而線粒體調(diào)控能量代謝，內(nèi)質(zhì)網(wǎng)促使蛋白質(zhì)加工及脂質(zhì)合成，核糖體將氨基酸合成蛋白質(zhì)，最終通過高爾基體加工合成多糖[18]，這對機(jī)體能量代謝發(fā)揮重要作用。以上結(jié)果表明，差異表達(dá)circRNAs可在不同RFI表型牛脂肪合成和能量代謝中起重要作用。

圖3 高、低RFI牛下丘腦差異circRNAs的GO功能富集分析Fig.3 GO enrichment analysis of hypothalamusdifferential circRNAs in bovine with high and low RFI

2.5 差異表達(dá)circRNAs靶基因的KEGG富集分析

使用KEGG分析探討差異表達(dá)circRNAs生物學(xué)作用。KEGG富集分析表明（見圖4），差異表達(dá)circRNAs靶基因主要富集在AMPK信號通路、mTOR信號通路及肌動蛋白細(xì)胞骨架形成通路中。AMPK信號通路在骨骼肌和肝臟中發(fā)揮重要作用，如丁酸通過增加肌肉和肝臟中AMP/ATP比率，直接激活A(yù)MPK信號通路[19]，該通路激活后可抑制糖原和蛋白質(zhì)合成，葡萄糖6-磷酸酶和磷酸烯醇式丙酮酸羧激酶基因表達(dá)下降，最終葡萄糖轉(zhuǎn)運和脂肪酸氧化增加[20]。這些信號通路的富集為后續(xù)開展circRNAs在牛食物攝入及能量代謝過程中分子功能提供重要信息。

圖4 高、低RFI牛下丘腦差異circRNAs的KEGG功能富集分析Fig.4 KEGG enrichment analysis of hypothalamus differential circRNAs in bovine with high and low RFI

2.6 差異表達(dá)circRNAs驗證

為驗證測序數(shù)據(jù)可靠性，對4個差異表達(dá)circRNAs（novel_circ_0031230、novel_circ_0038739、novel_circ_0022168、novel_circ_0005798、novel_circ_0018798、novel_circ_0034790）作qRT-PCR分析。結(jié)果表明（見圖5），在高RFI牛下丘腦組織中novel_circ_0031230、novel_circ_0038739、novel_circ_0018798表達(dá)上調(diào)，novel_circ_0022168、novel_circ_0005798、novel_circ_0034790表達(dá)下調(diào)，定量結(jié)果與測序結(jié)果趨勢一致，說明測序數(shù)據(jù)相對可靠，相關(guān)差異分子可用于后續(xù)功能驗證分析。

圖5 高、低RFI牛下丘腦6個差異表達(dá)circRNAs的qRT-PCR驗證Fig.5 qRT-PCR verification of six differentially expressed circRNAs in hypothalamus of brovine with high and low RFI

2.7 miRNA靶點預(yù)測及circRNA-miRNA-mRNA網(wǎng)絡(luò)分析

circRNAs富含miRNA結(jié)合位點，通過miRNA海綿效應(yīng)吸附miRNA，解除其對靶基因抑制作用。HRFI組315個上調(diào)和421個下調(diào)的差異表達(dá)circRNAs中共篩選到1 018個miRNAs靶點，其中在HRFI牛下丘腦差異表達(dá)上調(diào)最大的novel_circ_0037052靶向結(jié)合34個miRNAs，差異表達(dá)下調(diào)最大的novel_circ_0000923結(jié)合31個miRNAs靶點，提示上述miRNAs可能在牛RFI調(diào)控方面具有重要作用。通過篩選的miRNAs靶點制作circRNA-miRNA-mRNA靶向關(guān)系圖，圖6展示顯著差異表達(dá)的4個circRNAs及其靶向miRNAs，mRNA。結(jié)果顯示，差異circRNAs靶向結(jié)合的bta-miR-103、btamiR-33 a/b、bta-miR-499、bta-miR-15a等miRNAs調(diào)節(jié)動物脂肪沉積及能量代謝等相關(guān)mRNA。

圖6 差異表達(dá)circRNA-miRNA-mRNA調(diào)控網(wǎng)絡(luò)Fig.6 Differentially expressed of circRNA-miRNA-mRNA regulatory network

3 討論

本研究分析高、低RFI牛下丘腦組織中circRNAs，共得到736個差異表達(dá)的circRNAs，其中315個在HRFI牛中表達(dá)上調(diào)，421個表達(dá)下調(diào)。通過qRT-PCR驗證在HRFI牛下丘腦組織中novel_circ_0031230和novel_circ_0038739顯著上調(diào)；novel_circ_0022168和novel_circ_0005798顯著下調(diào)，與測序結(jié)果趨勢一致，說明本次測序數(shù)據(jù)真實可靠。對差異circRNAs的來源基因作GO和KEGG富集分析發(fā)現(xiàn)，差異circRNAs來源基因主要參與細(xì)胞質(zhì)及細(xì)胞器組成，能量代謝相關(guān)的AMPK信號通路、mTORx信號通路等生物學(xué)過程。

本研究預(yù)測circRNAs-miRNAs靶標(biāo)關(guān)系，得出bta-miR-103、bta-miR-33 a/b、bta-miR-499及bta-miR-15a分別為novel_circ_0031230、novel_circ_0020172、novel_circ_0002558、novel_circ_0019875及novel_circ_0020172的海綿。bta-miR-103屬于高度保守miRNA家族，已被證明bta-miR-103及其靶基因（如Dicer、PANK）可調(diào)節(jié)脂肪生成[21]。研究發(fā)現(xiàn)bta-miR-103可增加細(xì)胞內(nèi)乙酰輔酶A含量，并將乙酰輔酶A導(dǎo)入三羧酸循環(huán)，bta-miR-103也可通過與PANK蛋白共生作用調(diào)節(jié)細(xì)胞乙酰輔酶A和脂肪酸水平，影響機(jī)體能量代謝和脂肪沉積[22]。在山羊乳腺上皮細(xì)胞中，bta-miR-103被認(rèn)為是circ_015003和circ_014121的海綿，過表達(dá)bta-miR-103可增加乳脂合成相關(guān)基因表達(dá)，并抑制與脂肪分解和β-氧化相關(guān)基因（如瘦素、AMP活化蛋白激酶亞單位α），導(dǎo)致脂肪滴形成、甘油三酯積累增加[23]；Zhang等研究發(fā)現(xiàn)，circARF3可作為bta-miR-103的海綿而發(fā)揮作用，bta-miR-103抑制脂肪組織中有絲分裂，而bta-miR-103的海綿circARF3及bta-miR-103的靶基因TRAF3可促進(jìn)脂肪組織中有絲分裂，即circARF3通過阻斷btamiR-103效應(yīng)，導(dǎo)致TRAF3基因表達(dá)增加，調(diào)控小鼠脂肪合成[24]。Najafi-Shoushtari等發(fā)現(xiàn)miR-33 a/b為膽固醇控制的miRNA，其對膽固醇平衡的調(diào)節(jié)為通過與宿主基因協(xié)作，進(jìn)而抑制三磷酸腺苷結(jié)合盒轉(zhuǎn)運體A1（ATP binding cassette transporter A1，ABCAl）轉(zhuǎn)運蛋白表達(dá)[25-26]。當(dāng)膽固醇增加時，miR-33通過調(diào)節(jié)高密度脂蛋白消除血液中多余膽固醇[27]，改善動脈粥樣硬化；當(dāng)膽固醇減少時，miR-33下調(diào)ABCAl表達(dá)，減少膽固醇外流，升高細(xì)胞內(nèi)膽固醇水平。miR-499和miR-15a可分別通過其靶基因TCF12和FOXO1介導(dǎo)AMPKPGC-1α、TGF-β及PI3K-AKT-mTOR等信號通路調(diào)控動物采食，增強(qiáng)機(jī)體吸收利用，參與調(diào)控動物機(jī)體采食及能量代謝[28-29]。盡管上述miRNA與circRNA間共同作用參與眾多脂肪沉積調(diào)控通路，但其是否對牛食物攝入、能量代謝等生物學(xué)過程具有調(diào)控作用還需進(jìn)一步探討和驗證。

為深入探索circRNAs在能量代謝方面功能，對差異表達(dá)circRNAs的靶基因作GO及KEGG富集分析，結(jié)果表明靶基因主要富集在AMPK信號通路、mTORx信號通路及肌動蛋白細(xì)胞骨架形成通路中。AMPK信號通路在能量代謝中具有核心地位，該通路被激活后可促進(jìn)葡萄糖攝取和脂質(zhì)氧化，當(dāng)能量缺乏時產(chǎn)生能量，減緩能量消耗[30]；mTOR控制細(xì)胞代謝，包括氨基酸生物合成、葡萄糖穩(wěn)態(tài)和其他方面[31]，AMPK可通過mTOR通路調(diào)節(jié)細(xì)胞生長和蛋白質(zhì)翻譯過程，在能量充足條件下，AMPK無活性，但mTOR活躍，而在能量缺乏時，AMPK活性增加導(dǎo)致mTOR活性降低，機(jī)體通過加快攝取葡萄糖提供能量。在下丘腦中，AMPK和mTOR信號通路為食物攝入的決定因素，下丘腦AMPK的調(diào)節(jié)依賴于攝食狀態(tài)，攝入食物可降低下丘腦AMPK活性，而禁食導(dǎo)致下丘腦多個區(qū)域AMPK活性增加。下丘腦AMPK和mTOR活性改變可導(dǎo)致代謝紊亂，影響食物攝入量和體重及機(jī)體能量平衡[32]。

瘦素為一種脂肪細(xì)胞衍生的激素，可激活下丘腦mTOR信號通路活性，因而在減少食物攝入和長期控制體重方面發(fā)揮關(guān)鍵作用[33]；注射瘦素可抑制下丘腦AMPK活性，給予促食欲肽可增加AMPK活性[34]，如在自由采食的大鼠中，AMPK過度表達(dá)降低神經(jīng)肽Y的mRNA表達(dá)，從而控制攝食行為[35]。生物信息學(xué)分析顯示，circRNAs可能參與脂肪沉積及能量代謝等生物學(xué)過程，但關(guān)于circRNAs在脂肪沉積及能量代謝功能及circRNAmiRNA-mRNA互作關(guān)系還需進(jìn)一步在細(xì)胞水平進(jìn)行探討與驗證。

4 結(jié)論

本研究采用高通量測序技術(shù)獲得牛RFI極端表型個體下丘腦組織中circRNAs表達(dá)譜和差異表達(dá)信息。在樣本中檢測到豐富circRNAs，并統(tǒng)計其基因類型、長度分布及基因組特征。差異表達(dá)circRNAs靶基因主要參與細(xì)胞質(zhì)及細(xì)胞器組成，多富集于AMPK、mTOR等食物攝入和能量代謝相關(guān)信號通路。為后續(xù)探索circRNAs參與牛RFI的生物學(xué)過程和分子功能，研究circRNAs在哺乳動物食物攝入方面的調(diào)控機(jī)制提供數(shù)據(jù)支持。