奶牛與肉牛不同生長發(fā)育階段肝臟轉(zhuǎn)錄組比較分析

周 樂 石彩霞 常晨城 王 宇 滿 達(dá) 高 靜 張文廣*

(1.內(nèi)蒙古農(nóng)業(yè)大學(xué) 動物科學(xué)學(xué)院，呼和浩特 010018；2.內(nèi)蒙古自治區(qū)農(nóng)業(yè)基因組大數(shù)據(jù)工程研究中心，呼和浩特 010018；3.內(nèi)蒙古農(nóng)業(yè)大學(xué) 獸醫(yī)學(xué)院，呼和浩特 010018;4.內(nèi)蒙古農(nóng)業(yè)大學(xué) 計(jì)算機(jī)與信息工程學(xué)院，呼和浩特 010018)

奶牛與肉牛的早期發(fā)育對奶牛群和肉牛群的后續(xù)發(fā)展很重要，因?yàn)樯L良好的奶?；蛉馀？赡軙绠a(chǎn)犢并開始產(chǎn)奶或產(chǎn)肉，這不僅會增加利潤也會提高生長效應(yīng)。對于牛場來說，奶牛或肉牛的身體大小會被檢測生長和跟蹤每只牛在整個(gè)育成期期間發(fā)育進(jìn)展，以更好的進(jìn)行育種選擇和管理。其中大多數(shù)肉牛育種過程中會以初生重和斷奶重等作為選擇標(biāo)準(zhǔn)，但已有研究發(fā)現(xiàn)與生長性狀相關(guān)的基因在不同年齡階段、品種及環(huán)境中存在差異[1-3]。為提高奶牛與肉牛生產(chǎn)的經(jīng)濟(jì)效應(yīng)和完善奶牛與肉牛的生長性狀，挖掘與牛生長性狀相關(guān)的基因，研究動物生長發(fā)育和生產(chǎn)力的基本作用機(jī)制以及探索生長發(fā)育調(diào)控機(jī)制對促進(jìn)奶牛與肉牛分子育種具有很重要的意義。近年的研究主要是簡單的基于奶?；蛉馀Ｆ贩N進(jìn)行與生長相關(guān)基因的挖掘和探索。Sadkowski等[4]研究以6～12月齡的波蘭黑白公牛骨骼肌為研究樣本，發(fā)現(xiàn)有53個(gè)基因在6～12月齡之間的表達(dá)方式可能表明它們在出生后會參與到牛骨骼肌的生長，發(fā)育和代謝變化。Lu等[5]研究將2個(gè)肉牛品種的垂體組織進(jìn)行轉(zhuǎn)錄組分析，發(fā)現(xiàn)有8個(gè)基因在垂體中通過調(diào)節(jié)生長相關(guān)激素分泌來調(diào)節(jié)肉牛生長的候選基因。Zhang等[6]研究基于轉(zhuǎn)錄組測序?qū)?個(gè)肉牛品種的背闊肌進(jìn)行比較，發(fā)現(xiàn)他們之間的生長差異主要是代謝途徑中的脂肪分解和糖代謝。還有一些與肌內(nèi)脂肪生成等有關(guān)研究，都是以奶?；蛉馀Ｆ贩N進(jìn)行比較分析，發(fā)現(xiàn)與牛生長相關(guān)的基因或信號通路[7-8]。但除此之外，應(yīng)該基于奶牛與肉牛品種之間轉(zhuǎn)錄組的差異比較，進(jìn)一步探索發(fā)育機(jī)制和挖掘與生長性狀相關(guān)聯(lián)的基因及信號通路。

肝臟是哺乳動物體內(nèi)最大的消化腺和代謝器官，在機(jī)體的消化、代謝等許多生理活動中存在不可缺少的重要功能：參與多種重要生物化合物的代謝和合成、生長信號通路內(nèi)分泌控制、宿主對入侵病原體的防御以及外源物質(zhì)的解毒[9-11]。由于肝臟在圍產(chǎn)期和育肥期時(shí)調(diào)節(jié)代謝的關(guān)鍵作用，很多研究都是探究肝臟對成年奶?；蛉馀５牟町惡妥饔脵C(jī)制。這些研究表明，在荷斯坦奶牛產(chǎn)后早期時(shí)脂質(zhì)代謝、脂肪酸代謝和線粒體起重要作用；能量不平衡的奶牛和肉牛更難維持葡萄糖穩(wěn)態(tài)，造成肝臟脂質(zhì)代謝紊亂；肉牛在兩種成熟時(shí)期時(shí)可能造成肝臟對飼料效率的不同[12-14]。綜上所述，通過對不同品種及用途的奶?；蛉馀８闻K轉(zhuǎn)錄組進(jìn)行比較，不僅可以研究奶牛與肉牛肝臟在生長發(fā)育中的基因表達(dá)變化差異，也可以揭示不同品種牛肝臟在生長發(fā)育階段中的變化趨勢為后期飼養(yǎng)管理和選育提供重要依據(jù)。

因此，本研究以荷斯坦奶牛與海福特肉牛分別作為奶牛和肉牛的代表，通過對6、9和12月齡3個(gè)不同生長發(fā)育階段的肝臟轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行生物信息學(xué)分析，獲得差異表達(dá)基因及功能通路。通過不同生長階段的海福特肉牛與荷斯坦奶牛肝臟轉(zhuǎn)錄組數(shù)據(jù)進(jìn)行生物信息學(xué)技術(shù)的分析，篩選在奶牛與肉牛不同生長階段的肝臟差異表達(dá)基因，旨在挖掘影響奶牛與肉牛在生長發(fā)育階段中產(chǎn)生差異的重要基因與潛在的信號通路，并探究奶牛與肉牛肝臟組織在不同生長發(fā)育階段的表達(dá)變化趨勢，有助于揭示肉牛與奶牛肝臟在生長發(fā)育階段中有關(guān)代謝途徑的變化與差異。

1 材料與方法

1.1 數(shù)據(jù)來源

牛肝臟轉(zhuǎn)錄組數(shù)據(jù)來源于NCBI 的高通量測序數(shù)據(jù)存儲庫SRA (Sequence read archive)，獲取碼PRJNA312148，共獲得12 個(gè)肝臟轉(zhuǎn)錄組數(shù)據(jù)，如表1所示，分別是3個(gè)不同發(fā)育階段的海福特肉牛(Hereford，Here，n=6)和荷斯坦奶牛(Holstein，Hols，n=6)，其中 6、9和12月齡是指不同發(fā)育階段，每一月齡中有2個(gè)重復(fù)樣本，6月齡海福特肉牛的樣本定義為Here_6MA和Here_6MB，6月齡荷斯坦奶牛的樣本定義為Hols_6MA和Hols_6MB，9月齡和12月齡以此類推。

表1 海福特肉牛與荷斯坦奶牛肝臟轉(zhuǎn)錄組數(shù)據(jù)下載信息Table 1 Download information of liver transcriptome data for Hereford and Holstein cows

1.2 海福特肉牛與荷斯坦奶牛肝臟基因表達(dá)量的計(jì)算

利用SRAtoolkit 軟件的fast-dump參數(shù)將高通量測序數(shù)據(jù)轉(zhuǎn)變?yōu)閒astq文件，利用FastQC軟件對原始測序數(shù)據(jù)進(jìn)行定量分析，將低質(zhì)量以及未知堿基數(shù)目過多的reads過濾掉，獲得高質(zhì)量的clean data。

利用Hisat2軟件將clean data與牛參考基因組(ARS-UCD1.2)比對，利用Stringtie 軟件對轉(zhuǎn)錄本組裝，之后采用R軟件的Ballgown包計(jì)算基因在每個(gè)樣本中的轉(zhuǎn)錄本表達(dá)量。使用FPKM(Fragment per kilobase of per million mapped reads)值來衡量基因的表達(dá)水平。

1.3 海福特肉牛與荷斯坦奶牛肝臟差異表達(dá)基因(DEGs)分析

牛肝臟DEGs數(shù)據(jù)集來自6頭不同發(fā)育階段6(6月)、9(9月)和12月齡(12月)Here以及6頭不同發(fā)育階段6(6月)、9(9月)和12月齡(12月)Hols,通過不同品種相同發(fā)育月齡為標(biāo)準(zhǔn)，分為Here_6 M vs Hols_6 M；Here_9 M vs Hols_9 M；Here_12 M vs Hols_12 M這3組。使用DEseq2軟件對以上三組進(jìn)行DEGs分析，在差異基因檢測過程中，將|log2(Fold Change)|≥2且P<0.05作為篩選標(biāo)準(zhǔn)，最后為了對比分析不同組別中DEGs的表達(dá)狀況，需繪制火山圖判斷DEGs的分布情況。

1.4 海福特肉牛與荷斯坦奶牛短時(shí)序表達(dá)分析

STEM(Short time-series expression miner)是通過針對短時(shí)間序列的基因表達(dá)數(shù)據(jù)(少于或等于8個(gè)時(shí)間點(diǎn))進(jìn)行聚類和可視化分析，其每個(gè)模式圖代表每組基因的時(shí)序性軌跡。得到的聚類模塊中有不同顏色的是顯著性表達(dá)模式(Bonferroni矯正后的P<0.05)，顏色相同的表達(dá)趨勢一樣。為了呈現(xiàn)更好的聚類效果，選擇每組之間差異明顯的cluster。

1.5 海福特肉牛與荷斯坦奶牛差異基因GO和KEGG功能富集

GO (Gene Ontology)富集分析是國際標(biāo)準(zhǔn)化的基因功能分類體系。KEGG(Kyoto encyclopedia of genes and genomes)是與Pathway相關(guān)的主要公共數(shù)據(jù)庫。本研究主要利用網(wǎng)站DAVID Functional Annotation Bioinformatics Microarray Analysis v6.8 (https:∥david.ncifcrf.gov/tools.jsp)進(jìn)行功能富集分析，以P-value<0.05且FDR<0.1作為顯著富集的GO條目和Pathway通路的篩選標(biāo)準(zhǔn)。

2 結(jié)果與分析

2.1 海福特肉牛與荷斯坦奶牛轉(zhuǎn)錄組數(shù)據(jù)比對分析

利用Hisat2軟件將所有樣品的fastq文件分別與牛參考基因組ARS-UCD1.2進(jìn)行比對。平均每個(gè)個(gè)體將近86%左右的reads比對到?；蚪M。每個(gè)樣品的比對到參考基因組的信息如表2所示。

表2 海福特肉牛與荷斯坦奶牛RNA-Seq數(shù)據(jù)對比分析Table 2 RNA-Seq data mapping results statistics in Hereford and Holstein cattle

海福特肉牛與荷斯坦奶牛不同生長階段的轉(zhuǎn)錄組數(shù)據(jù)顯示，同一品種牛的不同生長階段中表達(dá)基因的數(shù)量不同，將每組重復(fù)樣本中表達(dá)基因數(shù)總和除去重復(fù)基因的數(shù)目，最終在海福特肉牛6月、9月和12月中分別鑒定到14 224、13 784和15 037個(gè)表達(dá)基因，在荷斯坦奶牛6月、9月和12月中分別鑒定到13 829、14 542和14 364個(gè)基因表達(dá)。在圖1中顯示了海福特肉牛和荷斯坦奶牛每個(gè)樣品中表達(dá)基因的數(shù)量，以FPKM表示(FPKM>0.5)，其中FPKM低于0.5的基因在后續(xù)分析中刪除。

由圖2可見(樣品命名詳情見材料與方法中的1.1)，可以通過主成分分析(PCA)探索所有牛肝臟樣品之間的全局關(guān)系。第一個(gè)主成分(PC1)解釋了基因表達(dá)的大部分差異都是通過品種而分離成兩大部分。第二主成分(PC2)是通過發(fā)育階段(6月、9月和12月)分離樣品，但分離的結(jié)果不是很明顯。結(jié)果發(fā)現(xiàn)肝臟組織中基因表達(dá)的主要差異可能是通過不同品種造成的，而在不同發(fā)育階段上不明顯。

圖1 海福特肉牛(a)和荷斯坦奶牛(b)基因表達(dá)量分布圖Fig.1 Distribution of gene expression of Hereford group (a) and Holstein group (b)

圖2 荷斯坦奶牛和海福特肉牛樣本PCA圖Fig.2 PCA group of Holstein and Hereford cattle

2.2 海福特肉牛與荷斯坦奶牛不同發(fā)育時(shí)期肝臟差異表達(dá)基因分析

基于PCA分析發(fā)現(xiàn)肝臟組織可能因品種不同而分離，但在不同生長階段中分離不明顯。因此，從海福特肉牛與荷斯坦奶牛在相同發(fā)育月齡中進(jìn)行分析，通過DEseq2軟件包在肝臟組織中篩選出Here_6 M vs Hols_6 M、Here_9 M vs Hols_9 M 和Here_12 M vs Hols_12 M這3組中的DEGs，在火山圖(圖3)中可以很直觀的觀察到這三組DEGs的上下調(diào)情況，x軸為log2FC值，y軸為-lgP值。在不同發(fā)育階段中的海福特肉牛和荷斯坦奶牛中顯著差異表達(dá)的基因用綠點(diǎn)(下調(diào))和紅點(diǎn)(上調(diào))表示。經(jīng)過差異基因的分析得到如下結(jié)果：1)海福特肉牛與荷斯坦奶牛在6月齡的DEGs為699個(gè)，其中373個(gè)DEGs在肉牛中表達(dá)上調(diào)，326個(gè)DEGs在奶牛中表達(dá)上調(diào)，最大log2FC值達(dá)到9.47(圖3(a))；2) 海福特肉牛與荷斯坦奶牛在9月齡的DEGs為486個(gè)，其中210個(gè)DEGs在肉牛中表達(dá)上調(diào)，276個(gè)DEGs在奶牛中表達(dá)上調(diào)，最大log2FC值達(dá)到23.38(圖3(b))；3) 海福特肉牛與荷斯坦奶牛在12月齡的DEGs為769個(gè)，其中405個(gè)DEGs在肉牛中表達(dá)上調(diào)，364個(gè)DEGs在奶牛中表達(dá)上調(diào)，最大log2FC值達(dá)到22.65(圖3(c))。根據(jù)FC值越大差異越顯著為標(biāo)準(zhǔn)，表3展示了這3組中差異顯著高表達(dá)前20的基因、P值和FC值。

圖3 海福特肉牛與荷斯坦奶牛6月(a)、9月(b)和12月(c)齡差異表達(dá)基因火山圖Fig.3 Volcano map of differentially expressed genes of Hereford and Holstein cattle in six (a),nine (b) and twelve (c) months of age

2.3 海福特肉牛與荷斯坦奶牛不同發(fā)育時(shí)期差異表達(dá)基因功能富集分析

2.3.1GO富集分析

將相同發(fā)育月齡不同品種牛的DEGs進(jìn)行GO功能注釋，表4分別列出奶牛與肉牛在6月、9月和12月這3組中富集到最顯著的條目。結(jié)果如下：1)海福特肉牛與荷斯坦奶牛在6月齡中的699個(gè)DEGs富集后得到4個(gè)GO條目，屬于生物學(xué)過程(Biological process，BP)；2)海福特肉牛與荷斯坦奶牛在9月齡中的486個(gè)DEGs富集后得到11個(gè)GO條目，其中1個(gè)屬于生物學(xué)過程、10個(gè)屬于細(xì)胞組分(Cellular component，CC)；3)海福特肉牛與荷斯坦奶牛在12月齡中的769個(gè)DEGs富集后得到5個(gè)GO條目，其中1個(gè)屬于生物學(xué)過程、3個(gè)屬于細(xì)胞組分、1個(gè)屬于分子功能。

表3 海福特肉牛與荷斯坦奶牛6、9和12月齡Top20差異表達(dá)基因Table 3 Top 20 differentially expressed genes of Hereford and Holstein cattle in six,nine and twelve months of age

2.3.2KEGG富集分析

對海福特肉牛與荷斯坦奶牛肝臟組織在相同發(fā)育月齡中的DEGs進(jìn)行KEGG信號通路分析，結(jié)果如下：1) 海福特肉牛與荷斯坦奶牛在6月齡中的699個(gè)DEGs富集后得到17個(gè)信號通路，依據(jù)P-value排序，前10個(gè)顯著富集的信號通路分別為類固醇生物合成、代謝途徑、生物合成的抗生素、萜類化合物生物合成、藥物代謝-細(xì)胞色素P450、細(xì)胞色素P450對外源性藥物的代謝作用、化學(xué)致癌作用、甾類激素生物合成、藥物代謝-其他酶和戊糖和葡萄糖醛酸的相互轉(zhuǎn)換過程(圖4(a))；2) 海福特肉牛與荷斯坦奶牛在9月齡中的486個(gè)DEGs富集后得到8個(gè)信號通路，依據(jù)PValue排序，顯著富集的信號通路分別為吞噬體、類風(fēng)濕性關(guān)節(jié)炎、細(xì)胞周期、造血細(xì)胞譜系、PPAR信號通路、I型糖尿病、肺結(jié)核和金黃色葡萄球菌感染(圖4(b))；3) 海福特肉牛與荷斯坦奶牛在12月齡中的769個(gè)DEGs富集后得到6個(gè)信號通路，依據(jù)P-value排序，顯著富集的信號通路分別為代謝途徑、生物合成的抗生素、PPAR信號通路、藥物代謝-細(xì)胞色素P450、類固醇生物合成和脂肪酸代謝過程(圖4(c))。

表4 6、9和12月齡功能富集分析GO條目統(tǒng)計(jì)圖Table 4 Statistical table of go items for function enrichment analysis in six，nine and twelve months of age

圖4 海福特肉牛與荷斯坦奶牛KEGG信號通路氣泡圖Fig.4 Bubble map of KEGG pathway for Hereford and Holstein cattle

通過比較海福特肉牛與荷斯坦奶牛參與脂質(zhì)代謝相關(guān)通路基因表達(dá)情況，發(fā)現(xiàn)細(xì)胞色素P450信號通路和PPAR信號通路2個(gè)通路中差異表達(dá)基因的表達(dá)量存在差異，并且不同生長階段表達(dá)也不同。例如細(xì)胞色素P450信號通路包含的16個(gè)基因中，GSTM1、GSTM2和GSTM3這3個(gè)基因的表達(dá)量在6月齡中海福特肉牛遠(yuǎn)遠(yuǎn)高于荷斯坦奶牛，但在12月齡中荷斯坦奶牛卻遠(yuǎn)遠(yuǎn)高于海福特肉牛(表5)；PPAR信號通路包含的18個(gè)基因在海福特肉牛中表達(dá)均高于荷斯坦奶牛，但其中ACOX1、APOA1和ACADM這3個(gè)基因表達(dá)差異最顯著，均在海福特肉牛生長階段中表達(dá)量最高(表6)；還發(fā)現(xiàn)與天然免疫相關(guān)的通路中所包含的8個(gè)基因中BOLA-DQB、BLA-DQB和BOLA-DQA5這3個(gè)基因在海福特肉牛中的表達(dá)均高于荷斯坦奶牛(表7)。

表5 海福特肉牛與荷斯坦奶牛中細(xì)胞色素P450信號通路DEGs的表達(dá)量Table 5 The expression of DEGs in cytochrome P450 signaling pathway for Hereford and Holstein cattle

表6 海福特肉牛與荷斯坦奶牛中PPAR信號通路DEGs的表達(dá)量Table 6 The expression of DEGs in PPAR signaling pathway for Hereford and Holstein cattle

表7 海福特肉牛與荷斯坦奶牛中抗原呈遞信號通路DEGs的表達(dá)量Table 7 The expression of DEGs in Antigen presentation for Hereford and Holstein cattle

2.4 奶牛與肉牛時(shí)序性動態(tài)表達(dá)模塊分析

通過GO和KEGG功能富集分析結(jié)果發(fā)現(xiàn)海福特肉牛與荷斯坦奶牛肝臟組織在不同發(fā)育階段中基因表達(dá)不同，基于此，對DEGs進(jìn)行了短時(shí)序表達(dá)分析，以了解在動態(tài)時(shí)序范圍內(nèi)基因表達(dá)趨勢的變化。對海福特肉牛和荷斯坦奶牛6、9和12月齡共1 666個(gè)DEGs進(jìn)行STEM分析。結(jié)果顯示，海福特肉牛中777個(gè)DEGs呈現(xiàn)出2種表達(dá)模式，其中6號模塊中579個(gè)DEGs在9～12月齡中的表達(dá)趨勢發(fā)生了變化，3號模塊中198個(gè)DEGs在6～9月齡和9～12月齡中的表達(dá)趨勢都發(fā)生了變化；而在荷斯坦奶牛中969個(gè)DEGs呈現(xiàn)出4種表達(dá)模式，其中7號模塊中467個(gè)DEGs和3號模塊中222個(gè)DEGs在6～9月齡和9～12月齡中的表達(dá)趨勢都發(fā)生了變化，6號模塊中222個(gè)DEGs在9～12月齡中的表達(dá)趨勢發(fā)生了變化，5號模塊115個(gè)DEGs在6～9月齡中的表達(dá)趨勢發(fā)生了變化(圖5)。

圖5 海福特肉牛與荷斯坦奶牛時(shí)序性表達(dá)基因的序列模塊Fig.5 Time-series modules of genes for Hereford and Holstein cattle

為了揭示奶牛與肉牛之間時(shí)序模塊中基因的表達(dá)趨勢和生長發(fā)育的關(guān)系，接下來重點(diǎn)選取兩組表達(dá)趨勢相同的6號模塊。通過KEGG富集，條形圖的長度代表通路中的-lgP值，值越大表示越顯著。結(jié)果表明，隨著時(shí)序動態(tài)變化6號模塊的表達(dá)趨勢為先升高到平穩(wěn)，再降低，海福特肉牛肝臟中氨基酸生物合成、類固醇激素生物合成、乙醛酸和二羧酸代謝、藥物代謝-細(xì)胞色素P450和脂肪酸的新陳代謝過程表達(dá)在6～9月齡中處于正常表達(dá)狀態(tài)，在9～12月齡中表達(dá)下降。荷斯坦奶牛肝臟中細(xì)胞周期、p53信號通路、卵母細(xì)胞減數(shù)分裂過程表達(dá)在6～9月齡中處于正常表達(dá)狀態(tài)，在9～12月齡中表達(dá)下降。由此可見在奶牛與肉牛肝臟之間的代謝過程存在差異(圖6)。

3 討 論

3.1 奶牛與肉牛不同發(fā)育時(shí)期肝臟差異表達(dá)基因分析

數(shù)據(jù)標(biāo)簽表示在此通路中所包含的基因數(shù)。The data label indicates the number of genes contained in this pathway.圖6 海福特肉牛與荷斯坦奶牛時(shí)序性表達(dá)基因的KEGG富集Fig.6 KEGG pathway of STEM genes for Hereford and Holstein cattle

通過差異分析、富集分析以及短時(shí)序模塊的聯(lián)合分析，本研究獲得了海福特肉牛與荷斯坦奶牛不同品種在相同發(fā)育階段的差異以及表達(dá)趨勢的不同。其中，發(fā)現(xiàn)6和12月齡中海福特肉牛與荷斯坦奶牛的DEGs表達(dá)個(gè)數(shù)較明顯，分別是699個(gè)差異基因和769個(gè)差異基因，而在9月齡中差異相比于其他2個(gè)月，差異不明顯，差異基因數(shù)量為486個(gè)。就差異分析結(jié)果來說，海福特肉牛與荷斯坦奶牛的肝臟組織在同一發(fā)育階段中確實(shí)存在差異。此外，通過KEGG富集分析，還發(fā)現(xiàn)PPAR信號通路和細(xì)胞色素P450這2個(gè)通路在生長階段中奶牛與肉牛脂質(zhì)代謝上存在差異的信號通路，還有1條與天然免疫有關(guān)的信號通路。

細(xì)胞色素P450是代表著一個(gè)很大且可自身氧化的亞鐵血紅素蛋白家族，肝臟中的細(xì)胞色素P450參與維持脂質(zhì)穩(wěn)態(tài)(膽固醇、維生素D、氧化甾醇和膽汁酸代謝)以及內(nèi)源性化合物(即膽汁酸)和外源物質(zhì)(藥物)的解毒過程[16]。當(dāng)部分細(xì)胞色素P450中的基因和酶表現(xiàn)出差異表達(dá)時(shí)，通常取決于牛的年齡和性別[17]。這與本文的結(jié)果皆證明在不同牛品種中細(xì)胞色素P450的差異表達(dá)較明顯。本研究中細(xì)胞色素P450信號通路中所包含的DEGs有16個(gè)，GSTM1、GSTM2和GSTM3這3個(gè)基因的表達(dá)量在6月齡中海福特肉牛遠(yuǎn)遠(yuǎn)高于荷斯坦奶牛，但在12月齡中荷斯坦奶牛卻遠(yuǎn)遠(yuǎn)高于海福特肉牛。GSTM1、GSTM2和GSTM3屬于谷胱甘肽轉(zhuǎn)氨酶(GST)基因家族。GST基因家族是解毒酶的關(guān)鍵基因家族之一，主要是基于谷胱甘肽和各種物質(zhì)的親電中心進(jìn)行結(jié)合，進(jìn)而導(dǎo)致毒性損失并產(chǎn)生親水能力更強(qiáng)的物質(zhì)，在不催化的情況下會隔離致癌物以及調(diào)控與解毒有關(guān)的信號通路[18]。

PPAR信號通路是一條與脂肪酸代謝、甾醇代謝以及脂肪分化密切相關(guān)的關(guān)鍵通路。大量的研究證實(shí)，在肝臟組織中與脂質(zhì)代謝和脂肪生成相關(guān)的關(guān)鍵基因會參與到PPAR信號通路中[19-20]。本文結(jié)果與上述相似，并且還闡明了PPAR信號通路在參與奶牛和肉牛的肝臟代謝過程中發(fā)揮不同的功能。本研究中PPAR信號通路所包含的18個(gè)基因在海福特肉牛中表達(dá)均高于荷斯坦奶牛，但其中ACOX1、APOA1和ACADM這3個(gè)基因表達(dá)差異最顯著，均在海福特肉牛生長階段中表達(dá)量最高。其中有研究發(fā)現(xiàn)，?；o酶A氧化酶1(ACOX1)不僅可通過調(diào)節(jié)過氧化物酶體脂肪酸β-氧化促進(jìn)牛前脂肪細(xì)胞的脂肪生成，還是通過標(biāo)記輔助選擇進(jìn)行肉質(zhì)鑒定的重要候選基因[21-22]。

除了細(xì)胞色素P450信號通路和PPAR信號通路外，本研究還發(fā)現(xiàn)了與免疫相關(guān)的信號通路，包含了8個(gè)基因，其中BoLA-DQB，BLA-DQB，BoLA-DQA5這3個(gè)基因在海福特肉牛中的表達(dá)均高于荷斯坦奶牛，并且已被研究報(bào)道是牛組織相容性復(fù)合體(MHC)Ⅱ類基因家族成員，代表了免疫應(yīng)答、抗病性和適應(yīng)性的關(guān)鍵候選基因[23-24]。例如牛淋巴細(xì)胞抗原(BoLA)經(jīng)過進(jìn)化在中國耗牛品種中產(chǎn)生了應(yīng)對各種病原體的免疫反應(yīng)能力[25]。結(jié)果表明海福特肉牛與荷斯坦奶牛在生長階段中與脂肪生成、肉品質(zhì)、免疫、抗病和適應(yīng)性中存在差異，表明奶牛與肉牛這兩品種在生長過程中肉品質(zhì)性狀、脂肪生成、免疫系統(tǒng)及適應(yīng)性可能存在一定的差異。

3.2 奶牛與肉牛時(shí)序性動態(tài)表達(dá)模塊分析

STEM支持對具有相同時(shí)間表達(dá)模式的基因集進(jìn)行基因本體(GO)注釋，為有效地對典型時(shí)間表達(dá)模式進(jìn)行統(tǒng)計(jì)上嚴(yán)格的生物學(xué)解釋提供了一種方法[26]。使用STEM軟件研究奶牛與肉牛肝臟在不同發(fā)育過程中基因表達(dá)如何隨時(shí)間變化而變化，最后確定了1 666個(gè)DEGs的表達(dá)趨勢圖，獲得了它們在奶牛與肉牛不同發(fā)育階段的動態(tài)表達(dá)模式。其中在表達(dá)模式相同的6號模塊中，海福特肉牛有579個(gè)DEGs，而荷斯坦奶牛只有222個(gè)DEGs，這表明在6、9和12月齡這3個(gè)階段中奶牛與肉牛DEGs的個(gè)數(shù)在表達(dá)趨勢中存在差異，相對于奶牛來說，肉牛表達(dá)的基因數(shù)量更多。因此，選取奶牛和肉牛的肝臟在不同發(fā)育階段中相同的6號模塊進(jìn)行探索，發(fā)現(xiàn)奶牛與肉牛之間時(shí)序動態(tài)表達(dá)模式存在差異。結(jié)果表明肉牛組中5個(gè)信號通路(包括40個(gè)DEGs)與氨基酸合成和脂質(zhì)代謝有關(guān)，說明隨著肉牛生長進(jìn)入成熟期時(shí)，氨基酸合成代謝和脂質(zhì)代謝相關(guān)基因的表達(dá)水平呈現(xiàn)出下調(diào)的趨勢。氨基酸合成代謝可調(diào)節(jié)影響生物體健康、生長發(fā)育的代謝途徑，其中谷氨酸是肝臟尿素再生、谷胱甘肽生成和新蛋白合成的中心底物[27]，其作用在肉牛育成期肌肉生長階段尤為關(guān)鍵。同品種牛谷氨酸合成代謝和蛋白質(zhì)含量的變化在生長階段和成熟期相比，成熟期的脂肪生成明顯高于育成階段[28]，這說明氨基酸合成代謝在成熟期脂肪增生階段的作用效果比育成階段的作用效果更顯著。脂質(zhì)代謝是影響肉牛生長發(fā)育速度的重要生物學(xué)功能。這結(jié)果與Mukiibi等[29]研究結(jié)果一致，在夏洛萊牛中發(fā)現(xiàn)19個(gè)與脂質(zhì)代謝有關(guān)的基因在快速生長階段中表達(dá)下調(diào)以及在安格斯牛中發(fā)現(xiàn)8個(gè)與脂肪生成有關(guān)的基因在月齡較大的肝臟組織中表達(dá)也下調(diào)。12月齡的肉牛已接近體成熟，身體發(fā)育也接近成熟并且有快速增重傾向，這可能是與脂質(zhì)代謝有關(guān)的DEGs在該階段表達(dá)下調(diào)的原因。這也與Mitchell[30]綜述中總結(jié)到快速生長階段中脂質(zhì)合成和脂肪沉積過程減弱，進(jìn)而產(chǎn)生更多的肉這一結(jié)論更貼合。綜上，氨基酸合成和脂質(zhì)代謝在肉牛生長發(fā)育階段中發(fā)揮著重要的作用，其表達(dá)變化可能是在早期發(fā)育階段中影響肌肉生長及脂肪沉積而上升，但在后期發(fā)育階段中為不影響育肥期快速增重而下降。

奶牛組中結(jié)果顯示3個(gè)信號通路(包括22個(gè)DEGs)與細(xì)胞周期、p53信號通路、卵母細(xì)胞減數(shù)分裂有關(guān)，其隨著時(shí)序表達(dá)的變化而變化，說明隨著奶牛生長進(jìn)入成熟期時(shí)，細(xì)胞周期、p53信號通路相關(guān)基因的表達(dá)水平呈現(xiàn)出下調(diào)的趨勢。細(xì)胞周期是細(xì)胞生命活動的基本過程，其通過調(diào)控G1、S、G2和M個(gè)時(shí)期中的細(xì)胞周期調(diào)控因子完成周期[31]。當(dāng)機(jī)體應(yīng)對外界刺激時(shí)細(xì)胞周期過程中會出現(xiàn)細(xì)胞周期G1阻滯現(xiàn)象產(chǎn)生保護(hù)作用，有助于維持基因穩(wěn)定以及修復(fù)受損的DNA，從而減少突變的發(fā)生，避免將突變基因遺傳到子代細(xì)胞[32]。安志芳等[33]研究已證明當(dāng)受到外界刺激時(shí)，機(jī)體通過抑制細(xì)胞周期下游基因的表達(dá)阻滯細(xì)胞周期G1期，從而提供充足時(shí)間進(jìn)行DNA修復(fù)，保證DNA復(fù)制的準(zhǔn)確性。奶牛在發(fā)育期向成熟期過渡時(shí)機(jī)體處于敏感時(shí)期，這可能會導(dǎo)致細(xì)胞周期下游基因的表達(dá)降低，從而影響細(xì)胞周期表達(dá)下調(diào)。P53信號通路主要功能是調(diào)節(jié)細(xì)胞周期和細(xì)胞凋亡，進(jìn)而促進(jìn)DNA修復(fù)[34]。在接近成熟期時(shí)表達(dá)不斷下調(diào)，與細(xì)胞周期表達(dá)下降結(jié)果保持一致，說明這一階段中奶牛的生理狀態(tài)較不穩(wěn)定，體內(nèi)很多系統(tǒng)正在建立過程中?？偟膩碚f，這些發(fā)現(xiàn)表明奶牛在生長發(fā)育階段中的狀態(tài)較不穩(wěn)定，因此在育成期時(shí)應(yīng)提供良好的飼養(yǎng)管理和營養(yǎng)以確保高效且穩(wěn)定的生長發(fā)育。

4 結(jié) 論

本研究集通過荷斯坦奶牛和海福特肉牛肝臟轉(zhuǎn)錄組數(shù)據(jù)在6、9和12月齡這3個(gè)生長階段進(jìn)行生物信息學(xué)技術(shù)分析，結(jié)果表明：奶牛與肉牛肝臟在不同生長階段中存在差異，主要表現(xiàn)在信號通路中，例如PPAR信號通路、細(xì)胞色素P450信號通路和免疫反應(yīng)中存在明顯差異，其次還表現(xiàn)在基因動態(tài)表達(dá)變化也存在差異，這進(jìn)一步為奶牛與肉牛肝臟生長發(fā)育的遺傳機(jī)制提供了依據(jù)。