李 嬡,張秀秀,黃萬龍,解領(lǐng)麗,苗向陽
(中國農(nóng)業(yè)科學(xué)院北京畜牧獸醫(yī)研究所,動物營養(yǎng)學(xué)國家重點(diǎn)實(shí)驗室,北京 100193)
豬肉是人民膳食結(jié)構(gòu)的主要組成部分,據(jù)美國農(nóng)業(yè)部經(jīng)濟(jì)學(xué)家和經(jīng)合組織-糧農(nóng)組織預(yù)測,與2014年相比,到2023年全球豬肉產(chǎn)量將會至少增加10.5%[1]。隨著人們需求量的增加,對豬肉肉質(zhì)的要求也越來越高,主要包括肉色、嫩度、風(fēng)味和多汁性等。肌內(nèi)脂肪是影響肉質(zhì)嫩度、風(fēng)味以及多汁性的重要因素[2-3],有研究發(fā)現(xiàn),某些生物活性物質(zhì)如共軛亞油酸在提高肌內(nèi)脂肪含量的同時,卻降低或維持了皮下脂肪的沉積量[4],由此推測,肌內(nèi)脂肪可能具有不同于皮下脂肪的生物學(xué)功能以及分子調(diào)控機(jī)制。因此,了解豬肉肌內(nèi)脂肪沉積的分子機(jī)制有助于改善豬肉品質(zhì)。
環(huán)狀RNA(circular RNA, circRNA)是一類由特殊的選擇性剪切產(chǎn)生的環(huán)形內(nèi)源性分子,呈共價閉合的環(huán)形結(jié)構(gòu),沒有5′和3′極性,不能夠被核糖核酸酶降解。在哺乳動物細(xì)胞中,circRNA分布廣泛,含量豐富,且具有穩(wěn)定性、保守性和組織特異性。已有研究發(fā)現(xiàn),circRNA與阿爾茨海默病[5]、大腸癌和卵巢癌[6]、肝細(xì)胞癌[7]、心血管疾病[8]、中樞神經(jīng)紊亂[9]等疾病相關(guān),除此之外,還與細(xì)胞增殖[10]、卵子發(fā)生與受精卵的形成[11]有關(guān)。在真核細(xì)胞中,circRNA可以作為miRNA海綿來調(diào)控基因表達(dá)[12-13]。研究發(fā)現(xiàn),miRNA和多種轉(zhuǎn)錄因子參與調(diào)控與脂肪細(xì)胞分化和脂質(zhì)代謝相關(guān)基因的表達(dá)[14-18]。但是,關(guān)于circRNA是否參與調(diào)控脂肪沉積和脂代謝還未見報道。
萊蕪豬肉質(zhì)優(yōu)良,其肌內(nèi)脂肪含量高達(dá)10.32%[19];而大白豬屬于典型的瘦肉型豬種,其肌內(nèi)脂肪含量為(1.51±0.19)%[20],因此本研究利用RNA-Seq技術(shù)對萊蕪豬和大白豬肌內(nèi)脂肪中circRNA進(jìn)行鑒定分析,以找出調(diào)控豬脂肪沉積的circRNA,進(jìn)一步了解circRNA調(diào)控脂肪沉積和脂代謝作用的機(jī)制,為豬肉品質(zhì)的提高奠定理論基礎(chǔ)。
本試驗以脂質(zhì)沉積存在顯著差異的大白豬和地方品種萊蕪豬為材料,均飼養(yǎng)于萊蕪市大千農(nóng)牧有限公司,在相同飼養(yǎng)條件和環(huán)境下生長育肥,參照滿足當(dāng)前營養(yǎng)需要標(biāo)準(zhǔn)(National Research Council, NRC, 1998)飼喂日糧,選擇處于屠宰期(約180日齡),健康無病、體況良好,種內(nèi)個體體重相近的雌性大白豬(約100 kg)和萊蕪豬(約35 kg)各3頭,采集其背最長肌肌內(nèi)脂肪組織,參照之前報道的肌內(nèi)脂肪剝離方法[21],用消毒的眼科鑷和手術(shù)刀,小心地分離出肌束間的脂肪,以避免肌肉纖維的污染。并且為減少RNA的降解,該過程在冰上進(jìn)行。試驗設(shè)置大白豬肌內(nèi)脂肪組織和萊蕪豬肌內(nèi)脂肪組織兩組樣本,每組樣本3個重復(fù),分別為D-JN-1、D-JN-2、D-JN-3和L-JN-1、L-JN-2、L-JN-3。
1.2.1 RNA 的提取及質(zhì)控 取等量低溫保存的脂肪組織樣本,根據(jù)說明書,利用mirVanaTM RNA提取試劑盒(#AM1561,Ambion, USA)提取每個脂肪組織樣本的總RNA,將分離的總RNA樣本置于-80 ℃保存。用NanoDrop 2000 (Thermo Scientific, USA)分光光度計測定OD260 nm/OD280 nm值以及RNA樣品濃度,并控制在1.9~2.1之間。Bioanalyzer 2100(Agilent, Santa Clara, CA)用于評估總RNA質(zhì)量,檢測RNA完整性(RNA integrity number, RIN),并控制RIN≥7且 28S/18S≥0.7。
1.2.2 cDNA文庫構(gòu)建和 RNA測序 根據(jù)Illumina?TruSeqTMRNA樣品制備流程,構(gòu)建 6個cDNA 文庫(包括 D-JN-1、D-JN-2、D-JN-3、L-JN-1、L-JN-2 和L-JN-3)。構(gòu)建cDNA文庫的步驟主要包括 Oligo(dT)磁珠分離富集純化 mRNA,酶促RNA片段化,cDNA 的合成,測序接頭銜接及 PCR 擴(kuò)增等。應(yīng)用Agilent DNA 1000 試劑盒和Agilent 2100生物分析儀(Agilent technologise, Santa Clara, CA)測定 cDNA 文庫大小和純度。應(yīng)用 ABI stepOnePlus 實(shí)時熒光定量PCR系統(tǒng)對cDNA文庫有效濃度進(jìn)行準(zhǔn)確定量(文庫有效濃度>2 nmol·L-1),以保證文庫質(zhì)量。文庫質(zhì)檢合格后,利用IlluminaHiSeqTM2500,采用雙端測序,對每個cDNA文庫進(jìn)行測序,得到原始測序數(shù)據(jù)raw reads。
1.2.3 circRNA測序分析 使用 fastx_toolkit[22]軟件 (v0.0.14)對 raw reads 進(jìn)行質(zhì)控,去低質(zhì)量序列、去接頭污染、去除rRNA等得到 clean reads。采用BWA中的BWA-MEM算法[23]將clean reads比對到參考基因組上,用軟件CIRI[24]對circRNA進(jìn)行識別,對比對結(jié)果SAM文件進(jìn)行掃描,尋找PCC(paired chiastic clipping)和PEM(paired-end mapping)位點(diǎn),以及GT-AG剪接信號,最后將具有junction位點(diǎn)的序列用動態(tài)規(guī)劃算法進(jìn)行重新比對,來確保識別circRNA的可靠性。采用DEseq2[25]進(jìn)行circRNA差異表達(dá)分析,比較處理組與參考組,并選取|Fold Change| ≥1.5和P<0.05 作為差異表達(dá)閾值,鑒定差異表達(dá)的circRNA。
1.2.4 差異表達(dá)circRNA的靶基因預(yù)測 運(yùn)用miRanda和Targetscan對差異表達(dá)的circRNA進(jìn)行靶標(biāo)關(guān)系預(yù)測[26-27],由于circRNA的miRNA結(jié)合位點(diǎn)較多,我們在差異表達(dá)circRNAs中選出了能較多結(jié)合miRNA的circRNAs之后,從Ensemebel數(shù)據(jù)庫中提取了豬所有基因的3′UTR序列,在Linux系統(tǒng)下利用miRanda[26](http://www.microrna.org/microrna)、RNAhybrid[28](http://bibiserv.techfak.uni-bielefeld.de/rnahybrid)、PITA[29](http://genie.weizmann.ac.il/pubs/mir07/mir07_data.html)3種軟件對這些miRNAs進(jìn)行靶基因預(yù)測,并取三者預(yù)測結(jié)果的交集用于后續(xù)進(jìn)一步分析。
1.2.5 circRNA-miRNA-mRNA互作網(wǎng)絡(luò)圖 結(jié)合本實(shí)驗室前期所做的萊蕪豬與大白豬肌內(nèi)脂肪組織中差異表達(dá)mRNA的鑒定與分析工作[30],篩選出與差異表達(dá)circRNA相關(guān)的差異表達(dá)基因,利用Cytoscape 軟件對circRNA-miRNA-mRNA構(gòu)建互作網(wǎng)絡(luò)圖[31]。
1.2.6 靶基因的功能富集分析 應(yīng)用 ClueGO[32]軟件對與差異表達(dá)circRNA相關(guān)的差異表達(dá)基因進(jìn)行富集分析?;虮倔w論(gene ontology, GO, http://www. geneontology.org/)是基因功能國際分類標(biāo)準(zhǔn),由分子功能(molecular function)、生物過程(biological process)和細(xì)胞組分(cellular component)組成。通路富集分析能確定差異表達(dá)基因參與的主要代謝途徑和信號通路,KEGG (kyoto encyclopedia of genes and genomes)數(shù)據(jù)庫[33](http://www.genome.jp/kegg)作為相關(guān)的主要公共數(shù)據(jù)庫,是進(jìn)行生物體代謝分析、代謝網(wǎng)絡(luò)研究的強(qiáng)有力工具?;诔瑤缀畏植?,計算差異表達(dá)基因顯著富集的 GO 條目和信號通路,利用 Benjamini-Hochberg 法對P_value進(jìn)行校正,當(dāng) Q_value (correctedP_value)≤0.05時,則顯著富集。
1.2.7 實(shí)時熒光定量 PCR驗證 qRT-PCR方法用于驗證基因的表達(dá)水平。取RNA樣品反轉(zhuǎn)錄合成cDNA模板。利用QuantiFast?SYBR?Green PCR Kit (Qiagen, Germany)和LightCycler?480 ⅡReal-time PCR Instrument (Roche, Swiss)進(jìn)行 qRT-PCR 分析。對4個circRNAs進(jìn)行qRT-PCR的驗證,4個circRNAs對應(yīng)的引物序列見表1。其反應(yīng)體系:1 μL cDNA,0.5 μL正向引物,0.5 μL反向引物,SYBR mix 10 μL,ddH2O 10 μL。反應(yīng)條件:94 ℃預(yù)變性10 min;然后45個循環(huán):94 ℃變性15 s、60 ℃ 60 s、72 ℃ 10 min。以甘油醛-3-磷酸脫氫酶(glyceraldehyde-3-phosphate dehydrogenase,GAPDH)基因為內(nèi)參,應(yīng)用 2-△△Ct法計算樣本間circRNA 的相對表達(dá)量,用t-檢驗對相對表達(dá)量進(jìn)行統(tǒng)計分析,數(shù)據(jù)表示為“平均數(shù)±標(biāo)準(zhǔn)差(Mean±SD)”,P<0.05表示差異顯著。
表1基因以及對應(yīng)的引物
Table1Geneandthecorrespondingprimersequences
基因Gene引物序列(5'→3')Primer sequenceGAPDH(內(nèi)參)F:ACCCAGAAGACTGTGGATGGR:TGAGCTCAGGGATGACCTTGssc_circ_0002807F:TCTCCTCGCAGTCATAGCCTR:TCACCTGGTGGATGTTCCCTssc_circ_0005382F:TGTAGTTCCCGGAGCCACTTR:TGATGTGTCTATGCCAAGCAGAssc_circ_0009172F:GGCCACTTTGCAGAAGAACATR:GGCTGGAGGCAAGGACATTssc_circ_0004824F:TGCTCCAGGATATCTTCAGGGAR:TGTGAAGTGGTCCCGTCTTG
通過環(huán)狀RNA測序分析,共檢測得到9 649個circRNAs,包括CLASSIC、ALTER_EXON、INTRON、OVERLAP_EXON、ANTISENSE和INTERGENIC 6種類型,所占比例如圖1所示。大白豬的3個肌內(nèi)脂肪組織分別得到88 890 352、89 083 366、91 659 768 條reads,比對到總讀段數(shù)的reads為88 358 935、88 674 923、91 184 863 條;萊蕪豬的3個肌內(nèi)脂肪組織分別得到90 180 707、90 392 480、91 093 142條reads,比對到總讀段數(shù)的reads為89 444 077、89 871 730、90 524 685 條,兩品種豬比對率均達(dá)到95%以上(表2)。這些結(jié)果表明,測序數(shù)據(jù)質(zhì)量可靠,為下游數(shù)據(jù)分析的可靠性奠定了堅實(shí)的基礎(chǔ)。
以|Fold Change| ≥1.5和P<0.05作為差異表達(dá)閾值,篩選得到181個差異表達(dá)circRNAs,其中102個上調(diào)circRNAs(LvsD),79個下調(diào)circRNAs(LvsD)。由于ssc_circ_0002807(下調(diào))、ssc_circ_0009352(上調(diào))有較多miRNA結(jié)合位點(diǎn),因此對其進(jìn)行后續(xù)分析。
圖1 各類型circRNAs占總circRNAs的比例Fig.1 The proportion of various circRNAs in total circRNAs
表26個樣本比對參考基因組結(jié)果
Table2Theresultof6samplesaftermappingtothereferencegenome
讀段類型Library樣品名稱 Sample nameD_JN_1D_JN_2D_JN_3L_JN_1L_JN_2L_JN_3總讀段數(shù)Total reads88 890 35289 083 36691 659 76890 180 70290 392 48091 093 142比對到基因組的總讀段數(shù)(比例/%) Mapped reads88 358 935(99.40)88 674 923(99.54)91 184 863(99.48)89 444 077(99.18)89 871 730(99.42)90 524 685(99.38)比對到基因組多個位置的讀段數(shù)(比例/%) MultiMap reads16 922 841(19.04)18 449 165(20.71)17 038 956(18.59)19 311 577(21.41)17 812 841(19.71)17 408 170(19.11)
對具有較多miRNA結(jié)合位點(diǎn)的差異表達(dá)circRNAs:ssc_circ_0002807和ssc_circ_0009352進(jìn)行miRNA結(jié)合位點(diǎn)的預(yù)測(圖2),每四行作為一個靶向關(guān)系展示。其中自由能越小,打分值越高,該靶向關(guān)系越可靠。ssc-miR-136、ssc-miR-4331、ssc-miR-4339和ssc-miR-874可以結(jié)合ssc_circ_0002807,ssc-miR-4331和ssc-miR-9822-3p可以結(jié)合ssc_circ_0009352,對ssc_circ_0002807、ssc_circ_0009352及其靶基因構(gòu)建了互作網(wǎng)絡(luò)圖(圖3)。
為了深入探究這些靶基因的功能,運(yùn)用ClueGO進(jìn)行GO和KEGG富集分析。在生物學(xué)過程,發(fā)現(xiàn)ssc_circ_0009352的靶基因主要富集于蛋白激酶C活性、蛋白去磷酸化、膽固醇平衡、磷脂平衡、甘油酯平衡和脂質(zhì)代謝平衡等;在分子生物功能,富集于蛋白磷酸酶、脂肪酸結(jié)合、磷脂酰肌醇磷酸結(jié)合(圖4A)。而在通路富集中也富集到與脂代謝相關(guān)的通路,如過氧化物酶體、脂肪酸降解(圖4B)。對于ssc_circ_0002807,在生物學(xué)過程,靶基因中差異表達(dá)基因富集于JAK-STAT級聯(lián)、蛋白激酶D信號調(diào)控、脂蛋白脂質(zhì)氧化、氧化應(yīng)激反應(yīng)、膽固醇平衡、脂蛋白修飾及氧化、磷脂平衡、脂質(zhì)代謝過程等。在分子功能,富集于JAK通道信號轉(zhuǎn)導(dǎo)適配器活動、過氧化物酶活性、脂肪酸結(jié)合、MAP激酶活性等(圖5A)。通路富集分析顯示,Ras信號通路、細(xì)胞因子受體相互作用、血管內(nèi)皮細(xì)胞生長因子信號通路、二型糖尿病、AMPK信號通路、不飽和脂肪酸的生物合成等通路顯著富集(圖5B)。
選擇4個差異表達(dá)circRNAs進(jìn)行驗證, 結(jié)果顯示,在萊蕪肌內(nèi)脂肪組織中ssc_circ_0002807(chr14:52135828..52247318:-)、ssc_circ_0005382(chr1:216338622..216447657:-)顯著下調(diào);ssc_circ_0009172(chr6:148414520..148422837:+)、ssc_circ_0004824(chr1:120698326..120698912:-)顯著上調(diào),與測序結(jié)果完全一致。
第1行每個字段分別表示miRNA名稱、circRNA id、打分值、最小自由能、miRNA與circRNA基因結(jié)合的起始位點(diǎn)、miRNA與circRNA基因的結(jié)合終止位點(diǎn)、circRNA基因與miRNA結(jié)合的起始位點(diǎn)、circRNA基因與miRNA結(jié)合的終止位點(diǎn);第2行(sRNA)表示反向互補(bǔ)的miRNA序列;第3行表示miRNA每個堿基和circRNA的互補(bǔ)配對情況,豎線表示AT或GC匹配,冒號表示GT匹配;第4行(Ref)表示miRNA結(jié)合的circRNA區(qū)域The first line represent miRNA names, circRNA id, score, minimum free energy, starting sites of miRNA combined with circRNA, ending sites of miRNA combined with circRNA, starting sites of circRNA combined with miRNA, ending sites of circRNA combined with miRNA, respectively; The second line (sRNA) represent the reverse complementary miRNA sequence; The third line represent complementary base-pairing of miRNA and circRNA: the vertical indicate matching of A and T or G and C, the colon represent matching of G and T; The fourth line (Ref) represent the circRNA region combined with miRNA圖2 circRNA與miRNA靶標(biāo)關(guān)系示意圖Fig.2 Schematic diagram of the relationship between circRNA and miRNA
箭頭表示circRNAs,圓形表示miRNAs,菱形表示miRNA靶基因The arrows represent circRNAs, the circles represent miRNAs, the rhomboids represent the target genes of miRNAs圖3 生物信息學(xué)預(yù)測的circRNA-miRNA-gene互作網(wǎng)絡(luò)圖Fig.3 Predicted circRNA-miRNA-gene network for two selected circRNAs by bioinformatics
圖4 ssc_circ_0009352靶基因GO和KEGG pathway富集分析Fig.4 GO(A)and KEGG pathway (B) enrichment analysis of the ssc_circ_0009352 targeted genes
圖5 ssc_circ_0002807靶基因GO和KEGG pathway富集分析Fig.5 GO(A)and KEGG pathway (B) enrichment analysis of the ssc_circ_0002807 targeted genes
不同豬種間,*. P<0.05;**. P<0.01;***. P<0.001*. P<0.05;**. P<0.01;***. P<0.001 between different pig breeds圖6 大白和萊蕪肌內(nèi)脂肪組織間4個差異表達(dá)circRNAs的qRT-PCR驗證Fig.6 qRT-PCR validation for the expression of 4 selected differentially expressed circRNAs in intramuscular adipose tissue between Large White and Laiwu pigs
本研究對大白豬和萊蕪豬肌內(nèi)脂肪組織中的circRNA進(jìn)行鑒定和分析,共鑒定到9 649個circRNAs,有181個circRNAs差異表達(dá),其中102個circRNAs在萊蕪豬肌內(nèi)脂肪組織中上調(diào)表達(dá),79個下調(diào)。通過qRT-PCR,驗證了ssc_circ_0002807、ssc_circ_0005382、ssc_circ_0009172和ssc_circ_0004824,與測序結(jié)果一致。circRNA是一類穩(wěn)定的非編碼RNA,廣泛存在于不同的細(xì)胞中,在基因的表達(dá)中發(fā)揮了關(guān)鍵的調(diào)節(jié)作用[34]。circRNA具有高度穩(wěn)定性,多樣性,豐富的miRNA反應(yīng)元件(miRNA Response Elements, MRE),序列特異性和物種間保守性[35],這些特征表明,circRNA在轉(zhuǎn)錄及轉(zhuǎn)錄后水平發(fā)揮重要作用,可能參與調(diào)控脂質(zhì)沉積,為代謝相關(guān)性疾病提供潛在的預(yù)防和治療靶標(biāo)。
最近研究表明,在許多真核細(xì)胞中circRNA可以作為miRNA海綿來調(diào)控基因表達(dá)。Hansen等[12-13]報道,ciRs7在人腦組織中豐富表達(dá),它含有多個串聯(lián)的miR-7結(jié)合位點(diǎn),可以作為內(nèi)源性的miRNA海綿吸附miR-7,在正常生理條件下ciRs7-miR-7-miR-7的靶基因處于動態(tài)平衡,當(dāng)這種穩(wěn)態(tài)被打破時則可能誘發(fā)疾病。本研究進(jìn)行了circRNA-miRNA靶標(biāo)關(guān)系預(yù)測。這些miRNA中,miR-136作為ssc_circ_0002807的靶向,與ASCL1、PPARα、C/EBPα呈負(fù)相關(guān)。而在脂肪組織,PPARα和PPARγ可以誘導(dǎo)ASCL1基因,從而酯化長鏈脂肪酸[36-38]。miR-4331作為ssc_circ_0002807和ssc_circ_0009352的靶向,能夠激活p38MAPK通路[39],激活的p38MAPK可以促進(jìn)骨髓間充質(zhì)干細(xì)胞的脂肪形成[40]。由此推測,ssc_circ_0002807可能參與成脂分化的調(diào)控。還有報道指出在感染HIV患者的脂肪組織中,一些miRNA異常表達(dá),其中miR-874上調(diào)表達(dá),這些miRNA的表達(dá)量變化可能導(dǎo)致與HIV相關(guān)的脂肪代謝障礙[41],而miR-874是ssc_circ_0002807的靶向,由此推測,ssc_circ_0002807可能參與脂質(zhì)代謝的調(diào)控。
競爭性內(nèi)源RNA(competing endogenous RNAs, ceRNA)是通過miRNA反應(yīng)元件競爭結(jié)合共同的miRNA來實(shí)現(xiàn)相互調(diào)控作用的轉(zhuǎn)錄產(chǎn)物,該機(jī)制也就是所謂的“ceRNA”假說。根據(jù)這一理論,circRNA和mRNA可以結(jié)合相同的miRNA[42],為了深入探索circRNA在脂肪沉積和代謝方面的功能,我們構(gòu)建了circRNA-miRNA-gene互作網(wǎng)絡(luò)圖(圖3)。通過GO和KEGG富集分析發(fā)現(xiàn),ssc_circ_0009352和ssc_circ_0002807的靶基因富集于與脂類代謝相關(guān)的生物學(xué)過程,如氧化應(yīng)激反應(yīng)、膽固醇平衡、磷脂平衡、脂質(zhì)代謝過程等。由此推斷,ssc_circ_0009352、ssc_circ_0002807可能在脂肪代謝中發(fā)揮作用。另外,ssc_circ_0002807的靶基因通路中還富集到二型糖尿病、胰島素信號通路和非酒精性脂肪肝等。脂肪組織能夠儲存大量的甘油三酯和脂溶性物質(zhì)以及機(jī)體代謝產(chǎn)生的多余能量,如果脂肪組織吸收和儲存循環(huán)脂質(zhì)的能力失調(diào),則可能會導(dǎo)致其它非脂肪組織中脂質(zhì)的積累,導(dǎo)致肥胖和能量代謝平衡失調(diào),引發(fā)肥胖相關(guān)性疾病的產(chǎn)生,如二型糖尿病、動脈粥樣硬化、胰島素抵抗、非酒精性脂肪肝炎、心血管疾病等[43-44]。生物信息學(xué)分析顯示,circRNA可能參與調(diào)控脂肪沉積和脂質(zhì)代謝,盡管如此,關(guān)于circRNA在脂肪沉積和脂質(zhì)代謝的功能以及circRNA-miRNA-gene的互作關(guān)系還需要進(jìn)一步探索和驗證。
本研究鑒定分析了大白和萊蕪豬肌內(nèi)脂肪組織中circRNA的表達(dá)情況,發(fā)現(xiàn)181個差異表達(dá)circRNAs。對具有較多miRNA結(jié)合位點(diǎn)的ssc_circ_0002807和ssc_circ_0009352靶基因進(jìn)行了GO和KEGG富集分析,發(fā)現(xiàn)靶基因富集于膽固醇平衡和磷脂平衡等與脂代謝相關(guān)的生物學(xué)過程。并通過qRT-PCR驗證了ssc_circ_0002807、ssc_circ_0005382、ssc_circ_0009172和ssc_circ_0004824,與測序結(jié)果一致,推測ssc_circ_0002807、ssc_circ_0009352可能參與調(diào)控脂肪沉積和脂質(zhì)代謝。