柞蚕微粒子病潜在生物标志物的挖掘

孙 影,梁瑞业,张文娟,王 勇,3,姜义仁,3,*,秦 利,3,*

(1.沈阳工学院生命工程学院,辽宁抚顺 113122;2.沈阳农业大学生物科学技术学院,沈阳 110866; 3.辽宁省昆虫资源工程技术研究中心,沈阳 110866)

柞蚕Antheraeapernyi是我国重要的经济昆虫之一，近几年我国柞蚕放养面积基本稳定在80万hm2左右，年产柞蚕茧8万吨以上(封槐松和李建琴,2016)。柞蚕幼虫期生活在野外，各种病害严重影响柞蚕茧的产量及质量，往往造成巨大的经济损失。其中，柞蚕微粒子病是柞蚕生产中唯一的检疫性病害，该病害既可水平传播又可垂直传播，已成为制约柞蚕产业健康发展的瓶颈问题之一(秦利,2003)。目前，柞蚕生产上对于该病的防控仍一直沿用雌蛾显微镜检查的方法，该方法费时费力，且准确性受人为因素影响较大，不仅给柞蚕制种工作带来极大的负担，而且难以保证种卵质量。因此，开发柞蚕微粒子病的检测方法对于保障柞蚕产业健康发展具有重要的意义。

代谢组学是一种研究生物体内代谢的系统方法，主要以分子量在1 000以下的小分子代谢物为研究对象，通过比较不同状态生物样本中低分子量代谢物信息，反映生物体对病理生理刺激、遗传修饰和周围环境的反应，相比于转录组学、基因组学和蛋白质组学，代谢谱更接近于对疾病过程的反应，在很多疾病的早期诊断、监测及致病机理等方面具有广泛的应用潜力(Nicholsonetal.,1999,2002;Gowdaetal.,2008;Emwasetal.,2013;Zengetal.,2017)。利用代谢组学技术发掘新型生物标志物已是很多研究领域广泛采用的策略，包括癌症(Chengetal.,2011;Tumasetal.,2016)、糖尿病(Savolainenetal.,2017;Zengetal.,2017)、慢性病(Mastrangelo and Barbas,2017)及健康老龄化(Almeidaetal.,2021)等重要疾病研究。柞蚕微粒子病的病原物为柞蚕微孢子虫Nosemapernyi，微孢子虫侵染会造成柞蚕体内蛋白质在种类及含量上发生显著的变化，导致柞蚕生理或病理过程发生改变(姜义仁等,2012;臧敏等,2012)，推测在柞蚕微粒子病发生过程中柞蚕体内的代谢物也会发生相应的变化。本研究以健康及患柞蚕微粒子病雌成虫血淋巴为材料，采用代谢组学技术分析健康与患病柞蚕雌成虫血淋巴内代谢物差异，筛选差异代谢物，挖掘病害的潜在生物标志物，为开发柞蚕微粒子病的检测方法提供一种新的思路。

1 材料与方法

1.1 供试材料

供试柞蚕沈黄2号雌成虫由辽宁省昆虫资源工程技术研究中心提供。随机挑选刚羽化柞蚕雌成虫若干，冰上快速采集血淋巴，并对采血后雌成虫进行显微镜检查，镜检出微孢子虫的雌成虫(患微粒子病)的血淋巴标记为M组，镜检健康的雌成虫的血淋巴标记为H组，所有样品液氮冷冻后置于-80℃冰箱中备用。

1.2 代谢物提取

取25 μL冰上解冻后的柞蚕雌成虫血淋巴样品，加入250 μL 50%甲醇溶液，震动充分混匀，提取样品中的代谢物，常温静置10 min；提取液放置-20℃过夜，沉淀样品中的蛋白质，4 000 g离心20 min，转移上清液到96孔板；利用稀释液(异丙醇∶乙腈∶水=2∶1∶1，体积比)稀释代谢物提取液，样品稀释液作为代谢物上样样品，每组各取6个生物学重复。每个样品中等量取出10 μL稀释液混合作为质控(quality control,QC)样品。所有样品代谢物提取液上样前保存于-80℃冰箱。

1.3 液相色谱及质谱条件

利用超高效液相色谱仪(ExionLC,SCIEX)进行数据采集，色谱柱为ACQUITY UPLC HSS T3(100 mm×2.1 mm,1.8 μm,Waters)；采用高分辨率质谱仪TripleTOF5600plus飞行时间质谱(SCIEX，英国)进行代谢物采集，每个样本进行一次正离子模式采集和一次负离子模式采集，正离子模式的喷雾电压为5 000 V，负离子模式的喷雾电压为-4 500 V，仪器参数设置参照郑海英等(2019)方法。

1.4 代谢物注释和定量分析

1.4.1代谢物注释：利用Proteowizard的MSConver软件将原始文件转换成mzXML文件，将mzXML文件导入XCMS软件进行代谢物峰提取、保留时间矫正。利用CAMERA软件进行代谢物同位素、加和离子分析；利用metaX软件结合KEGG和HMDB数据库进行代谢物注释，质量误差为10 ppm；对质谱二级碎片数据与in-house代谢物标准品二级图谱库进行匹配、打分，提取相似性>80%的代谢物。

1.4.2代谢物定量分析：利用XCMS软件(Smithetal.,2006)提取所有样品的所有代谢物强度信息，采用metaX软件(Wenetal.,2017)对数据进行质控。去质量峰(QC样品中缺失超过50%或实际样品中缺失超过80%的离子)；利用KNN(k-nearest neighbor)方法进行缺失值补充；利用PQN(probabilistic quotient normalization)算法对所有样品进行数据归一化，再利用QC样品进行QC-RSC(QC-robust spline batch correction)校正；对校正后的数据进行过滤，去掉不稳定的代谢物，即在QC样品中相对标准差(relative standard deviation,RSD)>30%的离子。

1.5 差异表达代谢物筛选及其功能分析

使用metaX软件进行单变量和多变量分析获得健康及患微粒子病柞蚕雌成虫血淋巴样品间差异表达代谢物。采用单变量分析对每组6个生物学重复测代谢物强度的平均值计算其在健康及患微粒子病柞蚕雌成虫血淋巴样品间的比值(ratio)，即差异表达倍数(fold change)，利用T检验进行统计学分析，对检验结果作多重检验分析，主要利用BH(Benjamini-Hochberg)校正得到P值；结合多变量统计分析，运用偏最小二乘法判别分析(partial least squares discriminant analysis,PLS-DA)建立代谢物表达量和样品类别之间的关系模型，建立健康及患微粒子病柞蚕雌成虫血淋巴样品之间的PLS-DA模型，同时对该模型参数R2和Q2进行了置换检验，次数为200次，进而得到每个代谢物的变量重要性(variable importance,VIP)值；显著差异表达代谢物同时满足：ratio≥2或ratio≤1/2,P<0.05,VIP≥1。

2 结果

2.1 健康及患微粒子病柞蚕雌成虫血淋巴代谢物统计

为探讨柞蚕微粒子病雌成虫体内代谢物的变化，我们对柞蚕健康及患微粒子病雌成虫的血淋巴样品进行了代谢组学分析，获得样品的代谢物总离子流图(图1)。从图1可知，患病与健康雌成虫的血淋巴之间代谢物种类存在差异。在正离子模式下获得8 870个代谢物，其中注释代谢物5 390个，二级鉴定代谢物312个；在负离子模式下获得6 716个代谢物，其中注释代谢物3 848个，二级鉴定代谢物177个。

图1 健康与患微粒子病柞蚕雌成虫血淋巴样品代谢物质谱总离子流图Fig.1 Total ion chromatogram of metabolites from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometerA:正离子模式In positive ion model;B:负离子模式In negative ion model;H:健康Healthy;M:患病Sick;QC:质控Quality control.下同The same below.

2.2 健康及患微粒子病柞蚕雌成虫血淋巴代谢物鉴定

正、负离子模式下分别鉴定了2 629和1 983个一级鉴定代谢物，均被注释到40个KEGG二级代谢通路，两种模式鉴定的代谢物分别有2 109和1 606个被分类到新陈代谢(metabolism)中12个条目，134和107个被分类到人类疾病(human diseases)中11个条目，256和159个被分类到有机体系统(organismal systems)中9个条目，8和5个被分类到细胞过程(cellular processes)中3个条目，103和82个被分类到环境信息处理(environmental information processing)中的3个条目，19和24个被分类到遗传信息处理(genetic information processing)中的2个条目；在人类疾病二级通路中，正、负离子模式下代谢物数量最多的是癌症概况(cancers:overview)(分别为75和59个)，其次是抗肿瘤药物耐药性(drug resistance:antineoplastic)(分别为12和12个)和神经退行性疾病(neurodegenerative diseases)(分别为12和8个)，还包括药物依赖性(substance dependence)(分别为9和6个)、特定类型癌症(cancers:specific types)(分别为6和5个)、内分泌和代谢疾病(endocrine and metabolic diseases)(分别为6和6个)、细菌性传染病(infectious diseases:bacterial)(分别为4和3个)和传染性寄生虫病(infectious diseases:parasitic)(分别为6和4个)等(图2)。鉴定代谢物数量最多的前20个代谢通路如图3所示，正、负离子模式下前3位代谢通路依次为代谢通路(metabolic pathways)(分别为704和516个)、2-氧羧酸代谢(2-oxocarboxylic acid metabolism)(分别为75和72个)和氨基酸生物合成(biosynthesis of amino acids)(分别为69和68个)。正、负离子模式下，分别有28 622和20 683个一级鉴定代谢物被注释到HMDB数据库中的24个和22个功能分类中，前2位功能分类依次为脂和类脂分子(lipids and lipid-like molecules)(分别为16 530和11 223个)和苯丙素类和聚酮化合物(phenylpropanoids and polyketides)(分别为4 959和4 833个)(图4)。

图4 健康与患微粒子病柞蚕雌成虫血淋巴样品正离子(A)和负离子(B)模式下质谱检测的一级鉴定代谢物的HMDB注释Fig.4 HMDB super class annotation of metabolites primarily identified from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometer in positive ion (A)and negative ion (B)models

正离子模式下，312个二级鉴定代谢物中有293个被注释到HMDB数据库的11个super classes的42个条目，16个没有被注释(na)，3个为未知(unknown)，注释代谢物数量前3位的条目分别为甘油磷脂(glycerophospholipids)(67)、脂肪酰基化合物(fatty acyls)(62)和羧酸及其衍生物(carboxylic acids and derivatives)(34)(图5:A)；251个代谢物被注释到85个KEGG代谢通路中，其中排名前20位的通路如图6(A)所示，前3位通路分别为代谢通路(metabolic pathways)(42)、ABC转运蛋白(ABC transporters)(12)和蛋白质消化和吸收(protein digestion and absorption)(12)。负离子模式下，177个二级鉴定代谢物中有175个被注释到HMDB数据库的8个super classes的25个条目，2个没有被注释，前3位分别为甘油磷脂(glycerophospholipids)(71)、羧酸及其衍生物(carboxylic acids and derivatives)(27)和脂肪酰基化合物(fatty acyls)(16)(图5:B)；213个代谢物被注释到85个KEGG代谢通路中，代谢物数量前3位通路分别为代谢通路(metabolic pathways)(32)、癌症的中心碳代谢(central carbon metabolism in cancer)(9)和氨基酸生物合成(biosynthesis of amino acids)(9)(图6:B)。

图5 健康与患微粒子病柞蚕雌成虫血淋巴样品正离子(A)和负离子(B)模式下质谱检测的二级鉴定代谢物的HMDB注释Fig.5 HMDB annotation of metabolites secondarily identified from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometer in positive ion (A)and negative ion (B)models

图6 健康与患微粒子病柞蚕雌成虫血淋巴样品正离子(A)和负离子(B)模式下质谱检测的二级鉴定代谢物的前20个KEGG代谢通路注释Fig.6 Annotation of top 20 KEGG pathways of metabolites secondarily identified from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometer in positive ion (A)and negative ion (B)models

2.3 健康及患微粒子病柞蚕雌成虫血淋巴代谢物定量

正、负离子模式下，分别筛选出7 966个和 6 143个高质量代谢物归一化后用于统计分析，二级鉴定的高可信度代谢物在每个样品中的强度分布见图7。

图7 健康与患微粒子病柞蚕雌成虫血淋巴样品正离子(A)和负离子(B)模式下质谱检测的二级鉴定代谢物强度热图Fig.7 Heatmap of the intensity of metabolites secondarily identified from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometer in positive ion (A)and negative ion (B)models

2.4 健康及患微粒子病柞蚕雌成虫血淋巴差异表达代谢物的筛选

以健康组(H)为对照，筛选微粒子病组(M)差异表达的代谢物。结果表明，正离子模式下，筛选到差异表达代谢物472个，其中260个上调，212个下调；负离子模式下，差异表达代谢物301个，其中207个上调，94个下调。对两组各个样品的差异表达代谢物进行聚类，结果如图8所示，健康组(H)与微粒子病组(M)的差异表达代谢物彼此分开。

图8 健康和患微粒子病柞蚕雌成虫血淋巴样品正离子(A)和负离子(B)模式下质谱检测的差异表达代谢物表达强度热图Fig.8 Heatmap of the expression intensity of differentially expressed metabolites identified from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometer in positive ion (A)and negative ion (B)models

2.5 健康及患微粒子病柞蚕雌成虫血淋巴差异表达二级鉴定代谢物筛选

柞蚕健康雌蛾和患微粒子病雌蛾两组之间的PLS-DA模型结果表明，正离子模式下，建立的PLS-DA模型的R2为0.9949，Q2为0.6719(图9:A)；负离子模式下，建立的PLS-DA模型的R2为0.9976，Q2为0.7310(图9:B)，表明正、负离子模式下构建的PLS-DA模型较为可靠，差异表达代谢物在两组之间的表达量能够较好地将健康与患微粒子病柞蚕雌成虫血淋巴样品区分开。

图9 健康和患微粒子病柞蚕雌成虫血淋巴样品正离子(A)和负离子(B)模式下质谱检测的差异表达代谢物PLS-DA得分图Fig.9 PLS-DA score plot of differentially expressed metabolites identified from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometer in positive ion (A)and negative ion (B)models

筛选出差异表达的二级鉴定的代谢物(表1)中正离子模式下，筛选到差异表达代谢物12个，其中8个上调，4个下调；负离子模式下，筛选到差异表达代谢物9个，其中8个上调，1个下调。

3 讨论

本研究利用代谢组学技术分析了柞蚕健康及患微粒子病雌成虫血淋巴内代谢物的差异，了解了柞蚕雌成虫血淋巴内代谢物种类和柞蚕微粒子病对柞蚕雌成虫血淋巴内代谢物的影响，为深入开展柞蚕微粒子病的检测及致病机理等研究奠定了基础。本研究中，正离子模式下质谱检测获得8 870个代谢物，其中注释代谢物5 390个，筛选到微粒子病组相对于对照组健康的差异表达代谢物472个，其中260个上调，212个下调，这些差异表达代谢物中有二级鉴定代谢物12个(8个上调，4个下调)，分别为缬氨酸(valine)、苯并噻唑(benzothiazole),3-脱羟基肉碱(3-dehydroxycarnitine),1-甲基鸟嘌呤(1-methylguanine)、2-乙氧基萘(2-ethoxynaphthalene)、N6-乙酰基-L-赖氨酸(N6-acetyl-L-lysine)、生物素(biotin)、桑色素(morin)、噻吗洛尔(timolol)、酰基肉碱15∶0(acylcarnitine 15∶0)、酰基肉碱18∶4(acylcarnitine 18∶4)和异槲皮苷(isoquercitrin)；负离子模式下获得6 716个代谢物，其中注释代谢物3 848个，筛选到差异表达代谢物301个(207个上调，94个下调)，其中二级鉴定差异表达代谢物9个(8个上调，1个下调)，分别为二甲基丙二酸(dimethylmalonic acid)、戊二酸(glutaric acid)、2,5-二羟基苯甲酸(2,5-dihydroxybenzoic acid)、1,3-二乙酰基丙烷(1,3-diacetylpropane)、3-(4-羟基苯基)乳酸(DL-p-hydroxyphenyllactic acid)、泛酸(pantothenate)、荧光素(fluorescein)、飞燕草素-3-O-beta-吡喃葡萄糖苷(delphinidin-3-O-beta-glucopyranoside)和溶血磷酯酰肌醇16∶1(lysoPI 16∶1)(表1)。

表1 健康及患微粒子病柞蚕雌成虫血淋巴样品质谱二级鉴定的差异表达代谢物统计表Table 1 Statistical analysis of differentially expressed metabolites secondarily identified from the haemolymph samples of healthy and pebrine-infected female adults of Antheraea pernyi by mass spectrometer

目前，关于柞蚕代谢组方面的研究报道较少，本研究鉴定到这些差异代谢物可能与柞蚕微粒子病的病程相关，有待进一步研究。这21个二级鉴定差异表达代谢物分别注释到HMDB数据库中的8个super classes，其中包括苯型烃类(benzenoids)(2)、脂和类脂分子(lipids and lipid-like molecules)(3)、有机酸及其衍生物(organic acids and derivatives)(5)、有机化合物(organic compounds)(1)、有机氮化合物(organic nitrogen compounds)(1)、有机氧化合物(organic oxygen compounds)(1)、有机金属杂环化合物(organoheterocyclic compounds)(4)和苯丙素类和聚酮化合物(phenylpropanoids and polyketides)(4)(图5)。本研究发现，柞蚕微粒子病组的3-脱羟基肉碱显著上调表达，差异表达倍数达到了20.96，有意思的是在糖尿病(Zengetal.,2017)、红肉消费及其相关疾病(Romboutsetal.,2017)等研究中该代谢物同样也上调表达。有研究表明，3-脱羟基肉碱为肠道细菌分解代谢左旋肉碱(L-carnitine)的中间代谢物，推测肠道微生物可能参与到糖尿病进程(Kleber,2006;Koethetal.,2014)。柞蚕微粒子病会严重影响柞蚕幼虫的营养消化吸收功能(姜义仁等,2012)，推测柞蚕微孢子虫也可能通过影响柞蚕幼虫肠道微生物而参与到柞蚕微粒子病进程。这些差异代谢物在柞蚕微粒子病发生与发展过程中可能发挥重要作用，具有潜在的生物标志物功能。虽然本研究仅按照柞蚕制种生产的操作程序进行取样及初步分析了患柞蚕微粒子病及健康雌蛾血淋巴的差异代谢物，但发现代谢组学技术在柞蚕微粒子病研究中具有积极的应用价值及潜力。我们后续研究中将结合柞蚕微粒子病的分子检测手段来进一步提高检测灵敏性来确保取样的准确性，并计划以柞蚕幼虫为材料进行人工添食微孢子虫孢子，深入探讨侵染时期及孢子浓度等对柞蚕体内代谢物产生的影响，为进一步揭示柞蚕微粒子病的致病机理及开发病害诊断检测技术奠定基础。