杨苗，孙洁，沈富军，凌珊珊，吴虹林，张亮，侯蓉，黄炎，岳碧松，张修月*

(1.四川大学生命科学学院，生物资源与生态环境教育部重点实验室，成都610065；2.成都大熊猫繁育研究基地，四川省濒危动物保护生物学重点实验室，成都610081；3.中国大熊猫保护研究中心，大熊猫国家公园珍稀动物保护生物学国家林业和草原局重点实验室，四川都江堰623006)

大熊猫属于食肉目Carnivora，是我国珍稀濒危野生动物，幼年早期主要以母乳为食，成年却专性以竹子为食。大熊猫母乳营养丰富，含有大量的蛋白质、碳水化合物和脂肪(Nakamura.，2003)，而竹子属于低营养高纤维食物。不仅如此，大熊猫消化道短且无盲肠，食物通过时间短，依然保留着肉食性消化系统结构和遗传系统特征(邹兴淮等，1998)。因此，大熊猫生长过程中，营养利用与食物代谢及其变化可能与一般的草食动物和肉食动物均不同，具有其独特性。研究大熊猫生长过程中的营养利用及调控特征对其保护具有重要意义。

胰腺是动物重要的消化代谢器官，作为第二大消化腺，兼具内分泌和外分泌的双重功能。胰腺分泌胰液经胰导管输送至十二指肠的过程即为外分泌，胰液中含有的碱性碳酸氢盐及各种消化酶可中和胃酸，消化糖类、蛋白质及脂肪。胰腺的内分泌主要依靠内部的胰岛细胞产生胰岛素、胰高血糖素及生长激素释放抑制激素和胃泌素等(Leung，2010；Lorberbaum.，2020)。

转录组水平研究表明，大熊猫在幼年早期的胆固醇代谢、蛋白质消化吸收相关基因高表达，以满足其对营养的高需求，有利于出生后的快速生长。碳水化合物代谢、能量产生、氨基酸和蛋白质代谢等相关基因表达上调，表明成年大熊猫具有高的代谢水平(Ma.，2020)。但其具体代谢调控的分子机制有待进一步研究。microRNA(miRNA)作为一种非编码的短序列RNA，已被证实能够参与多种形式的生命活动，如疾病、免疫、细胞凋亡、代谢、生殖和发育等(Ambros，2004；Wienholds & Plasterk，2005)。胰腺组织相关miRNA的研究表明，一些miRNA可能参与胰腺的发生，包括miR-124(Baroukh.，2007)，miR-15a、miR-15b、miR-16、miR-195(Joglekar.，2007)，miR-376、miR-375(Kloosterman.，2007；Joglekar.，2009)和miR-7(Correa-Medina.，2009；Joglekar.，2009)等，这些miRNA的异常表达会影响胰腺的正常发育和生理功能(Dumortier & Obberghen，2012)。此外，miRNA可以调控营养物质的代谢过程，包括脂类代谢(Esau.，2006)、氨基酸代谢(Dang，2009；Sengupta.，2020)、糖代谢(Eichner.，2010；Kefas.，2010)等。因此，通过对大熊猫不同年龄阶段胰腺组织中差异miRNA的分析，可以了解在大熊猫从幼年到成年的食性转变过程中，miRNA对营养代谢的调控作用。

本研究探讨幼年早期与成年大熊猫胰腺组织中miRNA表达谱的差异，以了解miRNA对幼年早期和成年大熊猫营养代谢的调控差异。

1 实验材料和方法

1.1 样本采集

大熊猫胰腺组织样本由成都大熊猫繁育研究基地和中国大熊猫保护研究中心提供。从4只死亡个体中采集到了4份胰腺样本，其中，样本YR和LT为1日龄和6日龄，死亡原因为意外窒息；样本DN和CC为成年个体，为野外严重受伤，抢救无效死亡。所有样本的胰腺组织均未见病理改变。采样工作由专业兽医进行，项目研究经四川大学生命科学学院伦理委员会批准(批准号：20190506001)。

1.2 Small RNA提取和测序

对胰腺组织匀浆后，使用Trizol试剂按照说明提取RNA。使用Agilent 2100 Bioanalyzer检测样本的完整性和质量，RNA完整值(RIN)≥7.5，说明样本具有较好完整性，可用于后续small RNA建库测序。使用Small RNA Sample Pre Kit构建文库，将提取的总RNA过滤后得到 15～40个核苷酸的small RNA。用T4连接酶将small RNA的3’和5’端加上接头，PCR反转录扩增得到small RNA文库。扩增产物纯化检测合格后，由北京诺禾致源科技股份有限公司测序，测序平台为Illumina HiSeq2000。

1.3 Small RNA数据处理

测序得到的原始数据需经质控后才可用于后续分析，因此对原始数据进行以下处理：首先，删除质量值≤5的reads占总读段>50%的低质量序列，去除N(无法确定碱基信息)比例>10%的序列。其次，删除5’接头污染的序列和没有3’接头和插入片段的序列。修剪3’接头序列后，去除有poly A、poly T、poly C和poly G的序列。最后，去除过长或过短的序列，保留18～30 nt的序列用于后续分析。为避免非miRNA序列对结果的影响，使用BLASTN对非miRNA进行注释。为确保每个small RNA仅分配 1个注释，使用以下规则：rRNA>tRNA>snRNA>snoRNA>repeat sequences。删除非miRNA后，使用miRDeep2将剩余序列比对到大熊猫参考基因组，鉴定保守miRNA及预测新miRNA。

1.4 保守miRNA差异表达分析

使用miRDeep2中的quantifier.pl脚本计算每个样本miRNA的表达量，得到原始表达矩阵，导入 R软件的DEseq2程序包进行差异表达分析。DESeq2的estimateSizeFactors函数可对数据进行标准化处理，以消除测序文库的影响。使用标准化计数计算每个miRNA的表达量并采用Benjamini-Hochberg修正方法对值进行调整。以|log(fold change)|≥1且<0.05作为差异表达miRNA(differentially expressed miRNA，DE miRNA)筛选的阈值。

1.5 预测靶基因

使用miRTarBase数据库预测miRNA的靶基因，该数据库中的miRNA靶基因已通过实验验证(Huang.，2020)，可提高分析结果的准确性。

1.6 miRNA与mRNA的相互作用关系

为探究miRNA与mRNA之间的靶向调控关系，对已有mRNA样本的表达量进行差异表达分析。mRNA数据来自Ma等(2020)的研究，获得了该论文中与本研究一致的YR、LT与DN、CC 4只个体胰腺组织的原始表达量数据进行后续分析。将YR、LT与DN、CC的表达量数据同样采用DEseq2程序包进行差异表达分析，并认为符合|log(fold change)|≥1且<0.05筛选条件的基因是显著差异表达基因(differentially expressed gene，DEG)。miRNA主要通过使其靶基因沉默表达来起调控作用(Huntzinger & Izaurralde，2011)，因此，本研究主要讨论与其靶向基因表达有负相关的DE miRNA。

1.7 miRNA靶基因的富集分析

为了解miRNA的功能，分别对鉴定出的DE miRNA的靶向DEG进行富集分析。使用g:Profiler将其映射到基因本体论(GO)(http://www.geneontology.org/)和京都基因与基因组百科全书(KEGG)(http://www.genome.jp/kegg/)数据库，其中，<0.05的条目和通路被认为是显著富集的。

2 结果

2.1 Small RNA序列特征

测序获得4只大熊猫胰腺组织的small RNA数据，质控保留长度18～30 nt的序列后，共获得 46 357 895 条small RNA序列(图1)。

图1 幼年早期(YR、LT)和成年大熊猫(CC、DN)胰腺组织中small RNA分类注释及其所占比例统计

2.2 保守miRNA的识别和新miRNA的预测

将miRNA与miRBase数据库进行比对，共获得202个保守miRNA和8个新miRNA。其中，保守miRNA在幼年早期大熊猫胰腺组织中较为丰富，在成年大熊猫胰腺组织中较少。202个保守miRNA中，有36个在幼年早期大熊猫胰腺中特有表达，成年大熊猫胰腺组织未见表达(图2)。

图2 幼年早期(YR、LT)和成年大熊猫(CC、DN)胰腺组织中保守miRNA数目维恩图

2.3 幼年早期与成年大熊猫胰腺组织中miRNA的差异表达

miRNA的差异表达分析显示，幼年早期与成年大熊猫胰腺组织中miRNA表达存在较大差异。共获得77个DE miRNA，其中，在成年大熊猫胰腺组织中，39个上调，38个下调(图3)。

图3 DE miRNA的聚类热图

2.4 miRNA靶基因预测

通过miRTarBase数据库的预测，幼年早期个体中特有的36个miRNA共预测到3 488个有实验验证的靶基因。DE miRNA中，成年大熊猫胰腺组织中39个上调的miRNA预测到 5 836个有实验验证的靶基因，38个下调的预测到4 193个有实验验证的靶基因。

2.5 miRNA靶基因功能富集分析

幼年早期个体中特有的36个miRNA的靶基因GO富集分析显示，幼年早期个体特有的miRNA主要参与了物质代谢及合成等生物过程的调控，如，细胞代谢过程的调控(GO:0031323)、生物合成过程的调控(GO:0009889)和大分子代谢过程的调控(GO:0060255)等。

DE miRNA的靶基因富集分析显示，上调DE miRNA靶基因GO富集结果多与细胞、大分子、含氮化合物等代谢的正调控有关，KEGG富集结果多与信号通路和疾病通路等有关。下调DE miRNA靶基因GO富集结果与上调DE miRNA功能相反，主要参与细胞、大分子及含氮化合物等代谢的负调控，KEGG富集结果也多与信号通路和疾病通路相关。

2.6 miRNA与mRNA的相互作用关系的构建

mRNA差异表达分析显示，成年大熊猫胰腺组织中获得661个上调DEG和691个下调DEG。将39个上调DE miRNA对应的靶基因与691个下调DEG取交集处理，得到240个由上调DE miRNA调控的下调DEG(图4)。将38个下调DE miRNA对应的靶基因与661个上调DEG取交集处理，得到 138个由下调DE miRNA调控的上调DEG(图5)。

图4 上调DE miRNA调控的下调DEG

图5 下调DE miRNA调控的上调DEG

2.7 DE miRNA-DEG富集分析

为探究交集基因的生物学功能，对由DE miRNA调控的DEG进行了富集分析。KEGG富集结果显示，成年组中上调DE miRNA调控的下调DEG主要富集到蛋白质消化吸收等通路(图6：A)。成年组中下调DE miRNA调控的上调DEG主要富集到疾病、信号通路及氨基酸代谢等相关通路(图6：B)。

图6 DE miRNA调控的DEG的KEGG富集结果

2.8 消化代谢相关基因及调控基因表达的miRNA分析

胰腺是重要的消化代谢器官，为探究miRNA对不同年龄大熊猫消化代谢的调控作用，对参与糖代谢、脂质代谢以及蛋白质代谢等生物过程的相关基因及其调控miRNA进行分析。其中，编码电子传递链主要蛋白，在成年组中显著上调，是miR-1307的靶基因之一。与胆固醇代谢有关，在成年组中显著下调，是miR-1的靶基因之一。蛋白质消化代谢相关基因中，在成年组中显著下调，是miR-1的靶基因之一，这些miRNA与它们的靶基因表达均呈反向负调控趋势。

3 讨论

大熊猫的营养代谢一直是遗传学家和进化生物学家最感兴趣的问题。胰腺是动物重要的消化器官，在消化和营养代谢中起着至关重要的作用。分析幼年早期和成年大熊猫胰腺组织中miRNA的差异表达，有助于了解不同年龄阶段miRNA对大熊猫食物代谢的调控作用及其变化规律。

本文鉴别到幼年早期大熊猫胰腺组织中特有表达的36个保守miRNA，主要参与合成及代谢的调控，包括含氮化合物、大分子化合物、核酸代谢过程的调控，暗示幼年早期大熊猫大分子合成较弱；而在成年大熊猫体内，这些miRNA的沉默表达可能促进了蛋白质、脂质等生物大分子的合成过程。对DE miRNA的靶基因进行功能分析，发现成年组中上调DE miRNA的靶基因对大分子和含氮化合物代谢过程的调控是正调控作用，而成年组下调DE miRNA的靶基因对这些物质合成的生物学过程主要是负调控。这与mRNA分析结果相吻合，即成年大熊猫胰腺组织中参与这些大分子合成的相关基因表达上调(Ma.，2020)。

生物体在正常营养代谢过程中，会产生活性氧物质(Poyton.，2009)。幼年早期大熊猫胰腺组织中上调DE miRNA调控的下调DEG富集到了活性氧代谢过程的正调控通路上，说明相较于成年大熊猫，幼年大熊猫活性氧代谢的能力较弱。研究表明，当需求能量增加但摄入能量较少时，会激活体内脂肪的分解代谢，而造成细胞的活性氧代谢产物大量增加(张帆等，2020)。这一方面可能由于成年大熊猫代谢水平高(Ma.，2020)，产生的活性氧物质多，另一方面可能由于幼年早期大熊猫从母体乳汁中摄入的能量较多，因此体内活性氧化物质较少，而成年大熊猫以竹子为食，竹子中所含能量物质少，体内脂肪代谢途径可能被激活从而产生了较多的活性氧物质，提高活性氧物质的代谢反应以维持机体内正常的活性氧化物水平。本研究结果表明，miRNA参与了成年大熊猫活性氧代谢的调控作用。

幼年早期组中下调DE miRNA调控的上调DEG富集到了蛋白质消化吸收的通路上，而成年组中下调DE miRNA调控的上调DEG富集到了氨基酸代谢的相关通路上。mRNA分析结果也显示，成年大熊猫胰腺组织内上调DEG富集到了与氨基酸和蛋白质代谢有关的通路上，而在幼体大熊猫胰腺组织中的上调DEG富集到了与蛋白质消化吸收有关的通路上(Ma.，2020)，miRNA结果与mRNA分析结果一致。蛋白质是幼体生长发育期间所必需的重要营养物质(Wiedmeier.，2011)，蛋白质利用率下降会导致哺乳期幼崽生长受限(Aguinaga.，2011)。负调控蛋白质消化吸收相关基因的miRNA在幼年早期大熊猫胰腺组织中低表达，而在成年后高表达，可能由于乳汁中蛋白质含量高，幼年大熊猫对乳汁中蛋白质吸收充分，以保证出生后的快速生长。成年大熊猫调控氨基酸和蛋白质代谢相关基因的miRNA下调，可能是因为成年大熊猫代谢水平的全面提高(Ma.，2020)。

针对一些DE miRNA的功能分析发现，它们可以调节许多重要的代谢过程，将这些代谢过程分为脂质代谢、能量代谢以及蛋白质代谢。其中，miR-27a、miR-122和miR-370与脂质代谢有关。miR-27a可以促进脂肪分解，释放出更多的甘油及游离脂肪酸(Wang.，2011)，它在成年大熊猫中的表达上调表明了成年大熊猫胰腺组织中的脂质代谢较幼年早期个体有提高。miR-122可以参与脂类代谢，Esau等(2006)的研究表明，抑制 miR-122的表达对小鼠血浆的胆固醇和脂肪酸代谢均有促进作用。成年大熊猫体内miR-122的表达相较于幼年早期下降，也说明成年大熊猫体内的脂肪酸代谢较幼年早期旺盛。miR-370也在成年大熊猫体内表达下调。miR-370对miR-122的表达具有正向调控作用；此外，miR-370还与miR-122的调控作用类似，可通过调控基因的表达来调控脂质代谢过程(Iliopoulos.，2010)。这些miRNA的差异表达说明，成年大熊猫体内的脂质代谢水平高，与 mRNA的研究结果一致。Li等(2011)研究显示，包括miR-363、miR-382和miR-494等27个miRNA被预测为靶向参与氧化磷酸化的基因，且蛋白质组学实验显示miRNA与其靶基因呈反向表达。miR-494、miR-382和miR-363在成年大熊猫中的显著下调可能促进其能量代谢。此外，ATP合酶是氧化磷酸化过程中一个重要因子，而miR-127可以调控ATP合酶的亚基——ATP5B，miR-127靶向β-mRNA的3’UTR区域，并抑制ATP5B翻译(Willers.，2012)。在成年组中，miR-127的表达下调保证了ATP合酶的生成，促进了能量代谢。这一结果也与mRNA研究一致，即在成年大熊猫体内的能量代谢较幼年早期旺盛(Ma.，2020)。目前，有关miRNA对动物体内氨基酸代谢过程中的调控研究较少，在差异表达的miRNA中仅发现miR-1对氨基酸代谢有调控作用。有研究表明，仅摄取10 g必需氨基酸后3 h内，miR-1的表达会出现显著上调，从而加速氨基酸代谢与后期蛋白质的合成(Mcgregor.，2014)，miR-1在成年大熊猫胰腺组织中的表达上调可能与其相关代谢水平的提高有关。

综上，miRNA在大熊猫生长发育、营养代谢的动态变化中发挥了重要的调控作用，miRNA的差异表达模式与幼年早期大熊猫对高营养的需求和成年大熊猫胰腺组织中高水平的脂质代谢、氨基酸代谢和能量代谢相一致，这些结果与mRNA的分析一致。本研究丰富了大熊猫消化代谢的基础研究，为进一步了解大熊猫食性转变与营养利用提供了重要的参考。然而，由于研究对象及组织的特殊性，目前样本量较小，因此后续可增加样本量继续开展深入研究，进一步证实本研究的结果。