奶牛白细胞氧代谢VEGF通路和新基因散热调节

刘松松，王 剑

(西南大学 动物科技学院，重庆 400703)

热应激引起奶牛生理和免疫反应，两者相辅相成。高温散热导致低氧环境应激生理反应，而白细胞是免疫系统的重要免疫细胞，高温胁迫奶牛血液白细胞增加而免疫力下降[1]，受基因调控同时表达相关蛋白[2-3]。蔡亚非等[4]发现37.5℃时奶牛体内白细胞数量达到最低值。另有研究报道，非热应激奶牛血液中白细胞数量与热应激组差异极显著[5-6]。MIN等[7]用代谢组学技术分析夏季热应激奶牛免疫调控。

高温散热所致的低氧环境引起的奶牛氧化应激生理反应。氧化应激是一种毒性物质作用[8]，是体内抗氧化剂之间的不平衡[9]，当活性氧产生的速度超过其被抗氧化中和速度时就会产生氧化应激[10]。抗氧化防御作用减少所致的氧化应激引起生物大分子的损伤和正常代谢及生理反应紊乱[11]。而抗氧化剂能清除细胞内和细胞外的超氧化物，并抑制细胞膜的脂质过氧化[12]。 BERNABUCCI等[13]认为夏季热应激会导致奶牛能量负平衡而诱发氧化应激，同时，奶牛血浆抗氧化水平受温度和湿度影响，具体表现为自由基合成升高[14]。因此，用热环境引起娟珊牛氧化应激从而筛选信号通路相关基因可以作为潜在生物标记。

转录组学是继人类基因组项目后一种技术手段，实现筛选基因的差异表达、基因功能注释、基因功能富集和筛选出相关信号通路[15]。目前，转录组分析技术已经广泛应用到动物生理调控、免疫功能和品种选育，如，猪[16]、鸡[17]、羊[18]、秦川牛[19]和鱼[20]等。有的学者对奶热应激牛乳腺和肝脏组织进行mRNA-seq分析[21]。

前期研究发现耐热性好的奶牛品种有皮肤毛细血管面积大的形态特征和热转运快、热散失多的生理反应[22-23]，然而其基因调节机制仍不清楚。假设认为在不同热应激环境下奶牛白细胞转录水平差异基因调控散热信号通路从而实现高温散热。为了解决这个问题，本研究通过转录组学对不同热环境娟珊牛白细胞热基因和信号通路进行筛选并且对基因进行q-PCR验证。

1 材料与方法

1.1 实验动物和样品采集选用健康耐热品娟姗牛为实验动物，n=10,(340±40) kg，(3.5±0.3)年，在重庆光大奶牛场采集血液样品，外界环境温度为25，30和35℃，从尾静脉用真空抗凝采血管收集血液样品4℃保存运回实验室，2 mL血液样品分离出白细胞(5×106)备用。

1.2 RNA分离，cDNA文库构建和测序采用TRIzol®试剂从白细胞样品中分离出总RNA，通过NanoDrop ND-2000(NanoDrop Technologies)美国德克萨斯州威尔明顿市(Wilmington,DE,USA)检测RNA的浓度、纯度和完整性。用Smart-Seq2方法[24]进行扩增：加入反应buffer、反转录酶、含公共序列的Oligo-dT引物和TSO引物，反应得到第1条cDNA 链,不转管直接加入含共同序列ISPCR引物和PCR扩增试剂，反应得到长度(1～2) kb的第2条cDNA 链，形成稳定的双链结构。对合成的cDNA进行末端修复、磷酸化和添加‘A’碱基，用于连接Y字形的接头。对连接接头产物进行纯化和片段分选，用分选产物进行PCR扩增，纯化得到最终的文库。最后用Illumina HiSeq Xten(2×150 bp读取长度)对配对末端RNA-seq测序文库进行测序[25]。

1.3 转录组数据分析llumina HisSeqTM2500测序产生的原始图像数据经 Base calling转化为序列数据(Raw Reads)，结果以FASTQ文件格式存储，包含Reads的序列以及碱基的测序质量。RawReads经公共序列污染及酶切位点、Adapter污染、含N过多(>5%)低质量(Q≤19)以及比对到rRNA 库上的Reads过滤后，获得纯净Reads(Clean Reads)。采用 TopHatv2.0.12软件[22]将测序过滤后所得的Clean Reads与黄牛基因组(http://www.ensembl.org/Bos_taurus/Info/Index)进行比对分析，获得与参考基因相匹配的Reads(Unique reads)[26]。采用RPKM法(每百万Read中来自于某基因每千碱基长度的Reads数)计算基因表达量[27]，参照无生物学重复样品的DEGseq法基因差异表达分析[25]，以log2 Ratio≥1和Q<0.05为阈值，筛选出不同发育阶段间的差异表达基因(Differentially expressed genes,DEGs)将DEGs与参考基因比较，得出DEGs显著富集的基因本体论(GO,Gene Ontology)功能条目，并进一步挖掘出与DEGs显著相关的生物学功能。将DEGs向GO数据库 (http://www.geneontology.org/)各个条目进行映射，计算其数目，用 BENJAMINI等[28]方法对P-value进行校正后，以Q-value<0.05为阈值定义DEGs中显著富集的GO条目。通过与京都基因与基因组百科全书(Kyoto encyclopedia of genes and genomes,KEGG)数据库(http://wego.genomics.org.cn)进行比对，对基因涉及的信号通路或代谢途径(Pathway)进行分析。

1.4 q-PCR验证差异表达基因选择热调控基因与差异表达基因相关的PRDM1、SAMSN1、BHLHE40、NR4A2、DUSP6、EDN1、MYL9、PTGS2、IL-8、IL1A和TGFB2等11个DEGs，采用Primer5.0软件进行基因的定量引物设计(表1)。Q-PCR反应体系为15 μL：上下游引物0.6 μL，cDNA模板1.2 μL，TB Green Premix Ex TaqⅡ (TliRNaseH Plus) 12.5 μL，ddH2O 5.1 μL。反应程序：95℃预变性3 min；95℃ 10 s,60℃ 20 s，35个循环。以β-actin作为内参，采用2-△△Ct计算基因的差异倍数[29]。每个样品做3次重复。

表1 验证RNA-Seq的q-PCR基因及其引物

2 结果

2.1 差异基因表达试验筛选到1 550个热调控基因、212个DEGs和16个热调控DEGs(图1c,1d)。3个温度梯度(25℃对30℃、25℃对35℃、30℃对35℃)差异基因分析如下：差异表达基因数量分别为170，65，30(表2)；其中上调基因为144，58，2个(图1a)，下调基因为26，7，28个(FC>2，P<0.01)(图1b)。根据|log2倍变化|值≥1且P<0.05，使用DESeq2获得16个DEGs(包括13个上调的DEGs和3个下调的DEGs)(图1d、表3)。

图1 娟珊牛白细胞差异基因韦恩分析和热调控差异基因热图 a.不同温度梯度的娟珊牛白细胞转录组差异上调表达基因韦恩图；b.不同温度梯度的娟珊牛白细胞转录组差异下调表达基因韦恩图；c.不同温度下娟珊牛差异表达基因与热调控基因韦恩图;d.16个差异表达基因的热图

表2 不同温度娟珊牛的转录组热调控上、下调差异表达基因数量 个

表3 不同温度娟珊牛的转录组热调控差异表达基因统计

2.2 差异基因GO和KEGG功能注释本试验筛选出16个热调控DEGs(表3)，并对16个基因进行GO和KEGG功能注释分析。GO分析表明(图2)，3个温度的差异基因归类注释涉及生物过程(BP)、细胞组分(CC)和分子功能(MF)3大类44个二级条目。

图2 不同温度下娟珊牛中白细胞热调控差异基因的GO功能分析

KEGG分析表明(图3)，3个温度的差异基因归类注释涉及代谢、环境信息处理、细胞过程和有机体系4大类13个二级条目，其中11个基因注释在信号转导(signal transduction)、5个基因注释在内分泌(endocrine system)和5个基因注释在免疫系统 (immune system)。

2.3 差异表达基因的GO和KEGG功能富集分析对这16条热调控DEGs GO和KEGG功能富集分析。GO功能富集表明(图4a),这些基因富集在BP、CC和MP 3大类1 004条二级条目，其中主要富集于磷代谢进程调控、氧水平反应、运动、细胞成分正调控、氧水平下降反应、缺氧反应、细胞运动正调控、发育过程调节和细胞因子活性等。

从GO富集可以看出主要富集于关于磷代谢和氧有关的通路。KEGG富集分析表明(4b),这些基因涉及到87条通路，主要富集于AGE-RAGE、VEGF、TNF和HIF-1信号通路等。

2.4 热调控差异表达基因相关性分析对16条热调控DEGs进行相关性分析，发现只有14条基因互相具有相关性。尤其是ENSBTAG00000005339与其他6个基因具有相关性，而且该条基因还参与血管内皮生长因子(VEGF)信号通路；PRDM1与其他5个基因具有相关性；SAMSN1和DUSP6与其他4个基因具有相关性；NSBTAG00000010303、EDN1和IL1与其他3个基因具有相关性；IL8、PTGS2、NR4A2和MYL8和其他2个基因具有相关性(图5，6)。

2.5 q-PCR验证RNA-Seq为了验证本研究转录组测序结果的可靠性，利用q-PCR检测DEGs。与对照组相比，试验组中PRDM1、SAMSN1、BHLHE40、NR4A2、DUSP6、EDN1、MYL9、PTGS2和IL-8的mRNA表达增加；而IL1A和TGFB2的 mRNA表达下降。RNA-Seq结果中，PRDM1、SAMSN1、BHLHE40、NR4A2、DUSP6、EDN1、MYL9、PTGS2和IL-8上调DEGs；而IL1A和TGFB2下调DEGs，与q-PCR检测结果一致(图7)。

图3 娟珊牛白细胞热调控差异上调基因KEGG功能分析

图4 GO和KEGG富集娟珊牛白细胞中差异表达基因途径 a.GO富集;b.KEGG富集

图5 热调控差异基因共表达网络图

图6 VEGF信号通路

图7 mRNA-Sep和q-PCR上调和下调基因

3 讨论

目前，有很多的关于高温对奶牛生理、生化和免疫反应影响的研究。BERNABUCCI等[17]发现氧化应激是奶牛高温生理调节机制之一。 有学者报道热应激增加细胞内外过氧化物而损伤细胞从而诱导血浆氧化应激[13]。同时，奶牛血浆中的氧化剂受到温度湿度影响[14]。而WANG等[22-23]报道的散热性好的奶牛皮肤毛细血管面积大提供夏季高温奶牛低氧应激皮肤形态证据。 闵力等[7]利用mRNA-seq分析热应激奶牛代谢和生理调控。本试验目前主要用mRNA-seq分析研究奶牛轻微、中度和严重热应激氧化应激信号通路和基因调控机制。筛选出Jersey白细胞在不同温度刺激下mRNA转录水平有212个差异基因、1 550和16个热调控差异基因。3个温度梯度(25℃对30℃、25℃对35℃、30℃对35℃)差异基因(DEGs)分析如下：DEGs分别为170，65，30个差异基因；上调基因为144，58，2个；下调基因为26，7，28个。 本研究发现严重热应激上调差异基因数量少，但下调基因数量多。可能在极端热环境中奶牛除了散热，其他生理反应在基因转录水平处于抑制状态。

GO分析表明，3个温度的差异基因功能注释涉及BP、CC和MF 3大类44个二级条目。这些基因富集在BP、CC和MP 3大类1 004条二级条目，其中主要富集在和氧代谢有关的氧水平反应、氧水平下降反应、缺氧反应和磷代谢进程的调控等。3个功能注释直接和氧有关，而且磷代谢又需要氧的参与。GO富集分析体现了不同温度下娟珊牛对抗热应激过程在mRNA转录水平上存在一定的差异，而氧代谢是娟珊牛重要的一个调节功能。这些功能充分证明耐热品种娟珊牛适应高温与血氧代谢密切相关。热应激时奶牛机体处于一个能量负平衡状态，氧调控就显得尤为重要。

目前很多研究报道热应激诱导奶牛细胞凋亡，如，乳腺上皮细胞[30]、外周血液单核细胞和卵巢颗粒细胞[31]。有学者报道凋亡细胞通过TNF信号通路诱导凋亡[32-33]。MILLIN等[34]报道，轻度热休克刺激TNF相关凋亡，诱导配体介导白血病T淋巴细胞和原生细胞凋亡。本研究结果显示，在轻微、中度和严重热应激环境中TNF是奶牛最高激活调节者，可作为潜在的热应激反应标记分子。

本研究显示KEGG通路富集于VEGF、HIF-1和TNF信号通路等。VEGF是一种血管内皮细胞的特异性有丝分裂原，在体外可促进血管内皮细胞的生长，在体内可诱导血管增生[35-36]。尤其是在低氧环境下，VEGF与内皮细胞膜上VEGF受体结合，引起受体的自身磷酸化，从而激活有丝分裂原活化蛋白激酶(MAPK)，实现VEGF的有丝分裂原特性，诱导内皮细胞增生[37]。VEGF可增加血管通透性，其特点是作用迅速、持续时间短[38]。VEGF在低氧环境下可以改变细胞外基质，诱导血浆蛋白溶酶原激活物和血浆溶酶原激活物抑制剂-1、基质胶原酶和诱导组织因子等在内皮细胞表达，激发V3因子从内皮细胞中释放出来，从而改变细胞外基质，使其更易于血管生长[39]。因此VEGF信号通路增强血管的通透性使物质更利于氧等交换和运输，在缺氧条件下可以改变细胞外基质引起受体自身的磷酸化。

本试验筛选到16条热调控DEGs，发现有14条基因互相具有相关性，这14条基因有12条为上调表达。筛选出新基因ENSBTAG00000005339、ENSBTAG00000047561、ENSBTAG00000010303和ENSBTAG00000009987等。而ENSBTAG00000-005339与其他6个基因具有相关性，并且参与VEGF信号通路调节。本研究利用 RNA-Seq技术对娟珊牛血液中白细胞在热环境不同温度下转录组进行分析，获得了众多差异基因和有关通路的富集。研究发现高热环境耐热品种娟珊牛通过mRNA水平差异基因调控血氧运输和代谢实现高温调控，可能是DUSP6、SAMSN1、ENSBTAG00-000005339和PRDM1等基因调控VEGF通路增加氧代谢从而实现蒸发散热。