通过饮食中添加灭活的微生物细胞和马尾藻属（Sargassum sp.）免疫刺激剂来研究其对虾免疫调节的转录组学分析

《Comparative Immunology Reports》：Transcriptomic analysis of shrimp immune modulation by dietary inactivated microbial cells and Sargassum sp. immunostimulants

　　M.A. Amatul-Samahah | N.F.M Ikhsan | M.N.A Amal | M.Y. Ina-Salwany
马来西亚普特拉大学农业学院水产养殖系，43400，Serdang，Selangor，马来西亚

摘要
近年来，关于虾免疫系统的研究引起了广泛关注，因为更好地理解虾免疫系统中的事件有助于改善虾的整体健康状况，从而提高虾的产量。本研究以Litopenaeus vannamei为对象，采用转录组分析方法来研究基因表达和免疫反应。实验中使用含有Sargassum sp.的饲料基灭活微生物细胞对虾进行免疫刺激。接受免疫刺激的虾在免疫后没有表现出任何物理变化，对其肝胰腺的组织学检查显示，各组中的细胞结构均正常。通过RNA测序和质量评估，共获得84,545个转录本和47,763个独特基因。利用多个数据库（NR、NT、KO、SwissProt、PFAM、GO、KOG）进行的功能注释发现，63.36%的独特基因与已知免疫相关基因有显著匹配。差异表达分析确定处理组与对照组之间有2,868个上调基因和4,721个下调基因。关键的免疫相关基因，包括参与病原体反应、抗菌肽生成和多种免疫信号通路的基因（如Toll、JAK-STAT、P13K-AKT、Wnt），显示出显著的差异表达。本研究为虾免疫反应的分子机制提供了宝贵的见解，强调了功能性基因组学在通过靶向免疫刺激剂增强虾健康和产量方面的潜力。

引言
由于虾养殖的重要经济价值以及疾病爆发带来的挑战，特别是在Litopenaeus vannamei这样的物种中，免疫刺激已成为一个关键的研究领域。传统上使用抗生素和化学消毒剂进行疾病控制的方法由于对细菌耐药性和环境影响的担忧而逐渐被淘汰[1]。采用来自海藻、螺旋藻和植物提取物等天然来源的免疫刺激剂作为一种有前景的替代方法，可以增强虾的非特异性免疫反应，从而提高水产养殖中的疾病抵抗力和总体产量[2]。这种海藻富含多糖、多酚和抗氧化剂等生物活性化合物，已知它们能够增强水生动物的免疫反应。假设将Sargassum sp.与灭活细菌结合使用可以通过利用这两种成分的免疫刺激特性来协同增强虾的免疫力[2]。这些天然免疫刺激剂可以通过口服、浸泡或注射等多种方式给药，并已被证明可以增加虾的总血细胞计数、噬菌作用活性和酚氧化酶活性，为虾养殖中的疾病控制提供了一种更可持续和环保的方法[3]。

转录组分析在理解虾的免疫反应中起着关键作用，尤其是在由Vibrio parahaemolyticus引起的急性肝胰腺坏死病（AHPND）的情况下[[4]、[5]、[6]、[7]、[8]]。研究表明，使用Pandanus tectorius叶提取物等天然免疫刺激剂可以增强白脚虾（L. vannamei）对V. parahaemolyticus等病原体的免疫反应，提高存活率并上调Hsp70和crustin等免疫相关基因的表达[9]。通过分析虾不同组织样本中的RNA转录本，研究人员可以了解哪些基因正在活跃表达以及它们如何与特定特征相关[10,11]。这也是了解该物种免疫反应分子机制的有趣方法，尽管在这方面相关研究较少。在虾的研究中，表达序列标签（ESTs）分析有助于将知识与来自其他生物（如蛋白酶抑制剂）的已知免疫功能基因联系起来，这些基因可以对虾的免疫刺激作出反应[12,13]。大规模的EST和高通量基因表达研究将增加制定关于候选免疫功能基因在虾中作用的有效假设的可能性。此外，分析还发现，大多数来自虾的EST与任何已知序列没有显著相似性[14,15]。此外，在L. vannamei中鉴定出的免疫相关基因揭示了虾免疫和病原体反应的分子机制，突显了理解宿主代谢在对抗AHPND等感染中的重要性[16]。近年来，我们在甲壳类动物免疫系统的关键方面取得了重大进展，包括吞噬细胞的作用、prophenoloxidase级联反应、黑色素化系统和抗菌肽[17]、[18]、[19]、[20]]。虽然一些保守的免疫效应途径（如黑色素化和抗菌肽生成）在生化水平上已经得到了较好理解，但大多数甲壳类动物（如虾）的免疫反应背后的分子事件仍不清楚。因此，在本研究中，我们希望了解使用含有益生元Sargassum sp.的饲料基灭活微生物免疫刺激剂对虾免疫系统在转录水平上的影响。此外，我们还希望阐明处理样品中涉及免疫相关基因的通路。

材料与方法
本研究的实验严格遵循马来西亚普特拉大学（UPM）动物护理与使用委员会（IACUC）制定的动物使用规定。本研究的虾处理和实验方案得到了UPM伦理委员会的批准（20190416151439AUP/101）。

**含有Sargassum sp.的饲料基灭活微生物细胞的制备**
根据先前描述的方法[21]，使用Vibrio parahaemolyticus C4B菌株制备了灭活微生物细胞。实验在马来西亚普特拉大学农业学院水产养殖系的生物技术水生实验室进行。选择V. parahaemolyticus C4B菌株作为制备灭活微生物免疫刺激剂的菌株，是因为该菌株与AHPND疾病有关并且已有相关研究[22]。V. parahaemolyticus C4B菌株在含有1.5% NaCl（重量/体积）的盐溶液胰蛋白酶大豆肉汤（STSB）中培养24至36小时。细菌细胞通过13,500 g离心机收集，沉降后的细胞重新悬浮在无菌PBS中，并在70°C下热灭活30分钟[23]。通过将热处理后的细菌悬浮液接种在盐溶液胰蛋白酶大豆琼脂平板上进行培养确认其完全灭活。经过48小时培养后未观察到细菌生长。使用的灭活细菌细胞密度为10^9 CFU/mL，对应的OD值为575 nm处的0.85，用于稀释或浓缩细菌悬浮液。同时，根据Farhan等人的方法（2020），将实验饲料（Cargill，马来西亚商业颗粒：45%）与2%的Sargassum sp.混合。

**实验饲料的制备**
实验饲料通过将灭活微生物细胞悬浮液稀释至每公斤饲料10^5 CFU的最终浓度来制备。该悬浮液与0.1%的瓜尔胶（Merck，德国）混合后均匀涂抹在预先添加了2% Sargassum sp.的饲料颗粒上。本实验使用的对照饲料为Cargill，马来西亚商业颗粒（No.2，45%）。

**虾和饲料基免疫刺激剂的处理**
本研究共使用了500只健康的P. vannamei后幼体（PL15），大小均匀（0.24 ± 0.01 g）。这些虾由马来西亚Selangor Darul Ehsan的Oasis Long Diann Bio-Tech Sdn Bhd公司提供。实验在马来西亚普特拉大学生物科学研究所（IBS）的孵化场单元进行。到达后，虾被浸泡在装有充气海水的500升玻璃纤维容器中30分钟，然后分成几个容器进行适应程序[24]。虾被饲养在有恒定空气供应的封闭循环容器中，每天喂食两次商业虾饲料（Cargill，马来西亚）。实验容器中的水温和为28°C。海水通过Resun空气泵LP100（深圳Xing Risheng工业有限公司，中国）进行充气，每个操作容器和水族箱都连接了空气石。海水的盐度保持在20至26 ppt之间。

**免疫刺激剂的施用**
免疫刺激剂施用持续四周（每个处理组n=40只虾PL）。免疫刺激剂通过口服方式给予虾幼体（喂食量为虾体重的5%），每天喂食两次，一次在清晨，另一次在傍晚[23]。对照组喂食Cargill，马来西亚商业颗粒（No.2，45%），喂食量相同（表1）。一致的喂食时间表有助于减少同类相食现象。每天在第二次喂食后一小时更换大约10%的容器水量。使用虹吸管和小网手动清除剩余的沉没饲料。容器中的水质参数保持在适合虾养殖的标准范围内。

**表1. 用于转录组分析的虾样本处理组**

| 处理组 | 施用方法 | 剂量 | |
|----------------|------------------------|-------------------------------------------|
| 商业饲料 | 每公斤饲料1×10^5 CFU + 含2% Sargassum sp.的商业饲料 | |
| 口服 | 每公斤饲料1×10^5 CFU | |

**组织病理学**
进行了组织学分析，以研究病理变化。取样包括虾的全身部分，如头部、尾部、肝胰腺、肌肉组织和外骨骼。所有样本均用10%的福尔马林（v/v）固定。组织学分析在马来西亚普特拉大学兽医学院实验室诊断系（UPM）的兽医诊断实验室按照常规组织制备和石蜡包埋技术进行[25]。石蜡块切成4 μm厚的切片，并用苏木精和伊红（H&E）染色。切片在Axioskop 2光学显微镜（Carl Zeiss Jena GmbH，德国）下观察，并配备AxioCamERc5s显微镜相机（Carl Zeiss Jena GmbH，德国）。显微镜还连接了ZEN 2.3 Lite软件平台，以便控制相机并将czi格式图像导出为jpeg格式。

**RNA提取和质量评估**
来自对照组和接受免疫刺激处理的虾样本储存在RNAlater溶液中（Thermofisher Scientific，美国）。使用Trizol RNA分离（Invitrogen，美国）协议从虾样本中提取总RNA[26]。在RNA纯化过程中进行了柱上DNase处理以去除残留的gDNA。使用1%琼脂糖凝胶测量RNA降解和污染情况。RNA纯度使用NanoPhotometer?分光光度计（IMPLEN，CA，美国）确定。使用Agilent Bioanalyzer 2100系统的RNA Nano 6000 Assay Kit（Agilent Technologies，CA，美国）评估RNA完整性数（RIN）和定量。

**文库构建和测序**
每个样本使用1 μg的RNA作为输入材料进行RNA样本制备。使用NEBNext? Ultra? RNA Library Prep Kit（NEB，美国）根据制造商的建议生成测序文库，并为每个样本添加索引码。首先使用poly-T寡核苷酸磁珠从总RNA中纯化mRNA。在NEBNext First Strand Synthesis Reaction Buffer（5X）中，在高温下使用二价阳离子进行片段化。然后使用随机六聚体引物和M-MuLV反转录酶（RNase H-）合成第一链cDNA。随后使用DNA聚合酶I和RNase H合成第二链cDNA。剩余的突出端通过外切核酸酶/聚合酶活性转化为平末端。在DNA片段的3’端加上腺苷酸化后，使用NEBNext Adaptor和发夹环结构进行连接，为杂交做准备。为了选择长度约250～300 bp的cDNA片段，使用AMPure XP系统（Beckman Coulter，Beverly，美国）进行纯化。然后使用3 μl USER Enzyme（NEB，美国）在37°C下处理大小选择后的适配子连接的cDNA 15分钟，再在95°C下处理5分钟，之后进行PCR。最后使用Phusion High-Fidelity DNA聚合酶、通用PCR引物和Index（X）引物进行PCR。PCR产物使用AMPure XP系统纯化，并在Agilent Bioanalyzer 2100系统上评估文库质量。

**测序的质量控制**
高通量测序的数据文件通过CASAVA 1.6进行碱基识别（base calling）转换为测序读段。原始数据存储在FASTQ（fq）格式文件中，其中包含读段序列和相应的碱基质量信息。通过分析错误率、GC含量分布和数据过滤来控制数据质量。测序错误率和碱基质量因测序仪、试剂残留物和不同的样本类型而异。对于RNA-seq技术，错误率随着测序读取量的增加而增加，这是因为测序试剂的消耗[27]。GC含量分布分析用于检测潜在的AT/GC分离情况，这会影响后续的基因表达定量。原始数据（原始读取）以fastq格式首先通过内部脚本进行处理。在这一步中，通过移除包含接头序列、poly-N序列和低质量读取的数据，获得干净的数据（干净读取）。同时，计算了干净数据的Q20、Q30、GC含量和序列重复水平。给定碱基的测序质量分数Q由以下公式定义：Q = -10log10(e)，其中e是碱基调用错误的估计概率。质量分数20（Q20）表示1%的错误率，相应的调用准确率为99%。类似地，Q30表示99.9%的调用准确率。所有后续分析都基于高质量的数据。

**unigene的从头组装和功能注释**

干净读取使用Trinity进行从头组装，Trinity是一个专业的转录组组装工具，包含三个组件：Inchworm、Chrysalis和Butterfly。Inchworm通过搜索k-mer图中的路径将读取片段组装成线性contigs。Chrysalis如果这些contigs共享k-1-mers，则将它们合并并构建单独的de Bruijn图。Butterfly通过修剪虚假边、压缩线性路径并将其与读取片段对照来处理这些图，以生成每个剪接形式和旁系转录本的线性序列[28]。聚类使用Corset 1.0进行，它根据共享读取和表达模式对contigs进行聚类，选择最长的转录本作为unigenes，并过滤掉映射读取少于十个的contigs[29]。BUSCO v5.4.4通过分析剪接结果并根据单拷贝同源物评估剪接的质量和准确性来评估转录本的完整性[30]。

**注释分析**

使用了七个数据库进行注释分析：NR（非冗余蛋白质数据库）、NT（核苷酸序列数据库）、Pfam、COG/KOG、Swiss-Prot、KEGG和GO。NR和NT分别是NCBI的蛋白质和核苷酸序列的正式数据库，包括来自GenBank和SwissProt等仓库的信息。对于NCBI NT数据库，使用的软件和参数是NCBI BLAST（基本局部比对搜索工具）2.13.0+（https://blast.ncbi.nlm.nih.gov/Blast.cgi），e值阈值为1e-5。对于NCBI NR、SwissProt和KOG，使用DIAMOND v0.8.31[31]，切线e值阈值为1e-5。对于PFAM（蛋白质结构域的预测），使用HMMER版本3.0（http://hmmer.org/）[32]，e值阈值为0.01。对于GO（基于NR和Pfam的蛋白质注释结果），使用Blast2GO（https://www.blast2go.com/）[33]，e值阈值为1e-6。对于KEGG，使用KAAS（KEGG自动注释服务器）（https://www.genome.jp/kegg/kaas/）[34]，e值阈值为1e-5。

**差异表达基因（DEGs）的鉴定**

差异基因表达分析的输入数据来自基因表达水平分析的读取计数。差异基因表达分析包含三个步骤：首先是读取计数标准化，其次是模型依赖的p值估计，第三步是基于多重假设检验的假发现率（FDR）值估计。FDR定义为被拒绝的假设中错误拒绝的预期比例。在不同情况下应用不同的软件和参数集。使用DESeq R包（1.10.1）（https://bioconductor.orgpackages/release/bioc/html/DESeq2.html）[35]对两组进行了差异表达分析[35]。DESeq提供了基于负二项分布的模型来确定数字基因表达数据中的差异表达。使用Benjamini-Hochberg方法[35,36]调整得到的P值，以控制FDR。DESeq发现的调整后P值< 0.05的基因被归类为差异表达基因。

**GO富集分析和KEGG通路富集分析**

进行GO富集分析以确定一组基因或蛋白质中GO术语的统计显著过表达。GO术语（http://geneontology.org/docs/go-enrichment-analysis/）[38]中padj（调整后的p值）<0.05的被认为是显著富集的。三种不同的分类代表GO术语的三个基本分类：生物过程（BP）、细胞组分（CC）和分子功能（MF）。GO分类将基因与其GO功能相关联以获得基因注释信息。GO富集检查了一组基因中具有相似GO功能的基因的富集情况，用于分析具有相似功能的基因。使用DAG（有向无环图）直观显示差异表达基因的富集GO术语及其层次结构。KEGG（京都基因组与基因组百科全书）通路富集分析（https://www.genome.jp/kegg/pathway.html）研究了可能参与某些生物功能的多个基因之间的相互作用[37]。通路富集分析识别了与整个基因组背景相比显著富集的代谢通路或信号转导通路。使用KOBAS软件测试差异表达基因在KEGG通路中的统计富集[38]。

**虾中免疫相关基因和信号通路的分析**

基于差异基因表达分析，鉴定出免疫相关基因并确定了参与的信号通路。log2(倍数变化)>1的免疫相关基因被认为是差异表达的，并且基于padj<0.05（具有生物学重复）在统计上是显著的。还确定了参与免疫相关基因和信号通路的通路。

**统计分析**

通过单因素方差分析（ANOVA）分析了实验组和对照组之间均值值的显著性倍数变化。使用IBM SPSS统计软件版本22（SPSS Inc., Chicago, IL, USA）进行了Duncan的新多重范围测试进行数据分析，P < 0.05的差异被认为是显著的。

**基于饲料的免疫刺激剂对肝胰腺整体和组织形态的影响**

免疫后四周，所有虾都没有显示出任何物理变化，看起来正常，没有明显的病变或物理变化。没有观察到与AHPND感染一致的体征变化或感染症状，如嗜睡、游泳行为不规律、饲料摄入量减少、身体颜色苍白和厌食。根据组织形态学分析，肝胰腺细胞和结构正常，未因免疫刺激剂而发生改变（图1）。未观察到AHPND感染虾所预期的组织病理学特征，如上皮管脱落、血细胞浸润、细菌定植、核巨大化和白色粪便综合征。

**Illumina测序、质量评估和从头组装**

为了进一步确认免疫刺激后的免疫反应表达，对免疫虾和对照虾的双份样本的总RNA进行了测序。在没有参考基因组的情况下，使用Trinity从头组装干净读取以获得转录组，将所有样本的读取片段组合在一起以获得参考基因组。共获得了84,545个转录本和47,763个unigene（表2）。

**表2. 不同长度区间内的转录本和unigene数量概述**

| 转录本长度区间 | 300-500bp | 500-1kbp | 1k-2kbp | >2kbp |
| --- | --- | --- | --- | --- |
| 转录本数量 | 27,906 | 20,894 | 16,881 | 18,864 |
| Unigene数量 | 20,101 | 12,816 | 7,830 | 7,016 |
| 47,763 |

根据读取统计摘要，Library G1_A产生了38,465,352个干净读取，Library G1_B产生了35,324,298个干净读取，Library G6_A产生了31,049,285个干净读取，Library G6_B产生了30,313,061个干净读取（表3）。对于Q20（Q Phred值大于20的碱基百分比）(具有Q Phred值> 20的碱基数量) / (总碱基数量) *100)，每个库的所有干净读取的值都超过97%。同时，对于Q30（Q Phred值大于30的碱基百分比) (具有Q Phred值> 30的碱基数量) / (总碱基数量) *100)，每个库的所有干净读取的值都超过93%（表6.3）。较高的Q分数表示较低的错误概率，而较低的Q分数可能导致大量读取片段不可用，从而增加假阳性变异的调用，导致结论不准确。作为下一代测序（NGS）质量的基准，Q30的质量分数表示错误的碱基调用概率仅为1/1000。根据表3，Library G1_A和Library G1_B的干净读取的Guanine-Cytosine（GC）含量接近46%。Library G6_A和Library G6_B的干净读取的GC含量接近45%。上述结果表明序列数据质量良好，适合进一步分析。

**表3. 转录组分析的测序数据摘要**

| 库ID | 原始读取 | 干净读取 | 干净碱基 | 错误率（%） | Q20 | Q30 | GC（%） |
| --- | --- | --- | --- | --- | --- | --- |
| G1_A | 38,994 | 771 | 11.73 | 8,465 | 352 | 11.50 | 0.03 | 97.46 | 93.15 | 46.62 |
| G1_B | 35,785 | 588 | 10.73 | 5,324 | 298 | 10.60 | 0.03 | 97.43 | 93.08 | 46.42 |
| G6_A | 31,425 | 033 | 9.43 | 1,049 | 285 | 9.30 | 0.03 | 97.25 | 92.66 | 45.47 |
| G6_B | 30,742 | 265 | 9.23 | 0,313 | 061 | 9.10 | 0.03 | 97.59 | 93.41 | 44.92 |

**功能注释和分类**

在这项研究中，我们使用七个数据库成功注释了47,763个unigene（63.36%）：非冗余蛋白质数据库（NR）、核苷酸序列数据库（NT）、基因本体（GO）、真核生物同源群（KOG）、KEGG同源群（KO）、SwissProt和蛋白质家族（PFAM）（表4）。根据注释结果，22,930个unigene（48.00%）在NR数据库中有相似匹配，而21,524个unigene（45.06%）在核苷酸序列（NT）数据库中有显著相似性。unigene分别在KO、SwissProt、PFAM、GO和KOG数据库中独立注释。在KO、SwissProt、PFAM、GO和KOG数据库中成功注释的基因百分比分别为20.03%、30.99%、39.27%和17.09%。图2展示了unigene映射结果的Venn图，重叠区域显示了能够映射到重叠数据库的总unigene数量。Venn图展示了五个数据库（NR、NT、KOG、GO和Pfam）中注释序列的分布和重叠情况。大量序列（5810）在所有数据库中共享，表明这些序列在序列相似性、功能注释和保守域方面具有高度一致性，这表明这些序列的预测生物学角色具有很高的可信度。此外，NT（2975）和NR（2169）包含大量独特注释，反映了通过相似性搜索可识别的序列，但在功能或基于域的数据库中缺乏进一步分类。

**表4. 在每个功能数据库中成功注释的基因百分比**

| 功能数据库 | 注释的unigene数量 | 百分比（%） |
| --- | --- | --- |
| NR | 22,930 | 48.00 |
| NT | 21,524 | 45.06 |
| KO | 9,567 | 20.03 |
| SwissProt | 14,803 | 30.99 |
| PFAM | 18,760 | 39.27 |
| GO | 18,757 | 39.27 |
| KOG | 8,164 | 17.09 |
| 所有数据库 | 45,339 | 99.49 |

**下载：下载高分辨率图像（156KB）**
**下载：下载全尺寸图像**

**图2. 使用Venn图映射unigene。重叠区域显示成功映射到重叠数据库的unigene数量。结果仅显示非冗余蛋白质数据库（NR）、真核生物同源群（KOG）、基因本体（GO）和核苷酸序列数据库（NT）以及蛋白质家族（PFAM）。**

数据库之间的部分重叠，如NR和NT之间或Pfam、GO和KOG之间，表示注释水平中等的序列。这些序列可能有可识别的域或部分功能信息，但不是在所有数据库中都得到完全表征。总体而言，基于序列的数据库（NR和NT）提供了更广泛的覆盖范围，而功能数据库（GO、KOG、Pfam）提供了更具体的但有限的注释。因此，结合多个数据库可以提高功能解释的总体可靠性和深度。图3显示了来自非冗余蛋白质数据库（NR）的unigene与所研究物种P. vannamei之间的百分比相似性分布（%）。根据NR数据库中的unigene注释，unigene与P. vannamei的相似性显著（76.4%）。

**图3. 饼图显示了从NR数据库注释的unigene的物种分布相似性（%）**

**GO、KOG和KEGG对转录组序列的分类**

基因本体注释有效地将基因分类为三个GO领域：细胞组分（CC）、生物过程（BP）和分子功能（MF）（图4）。生物过程基因大约有26个数量，但主要涉及细胞过程（10,865个基因）、代谢过程（9,013个基因）、生物调控（4,162个基因）、生物过程的调节（3,818个基因）、定位（3,750个基因）、对刺激的反应（2,813个基因）和信号传导（1,802个基因）。166个基因被归类为免疫系统过程（图4）。与细胞成分相关的基因分类共有5种，包括细胞解剖实体（9,667个基因）、细胞内成分（4,924个基因）、含蛋白复合物（3,650个基因）、病毒颗粒（687个基因）以及病毒颗粒的组成部分（687个基因）。分子功能方面，基因被分为12个类别，涉及结合功能（9,041个基因）、催化活性（7,027个基因）、转运活性（2,115个基因）和结构分子活性（1,027个基因）（图4）。下载：下载高分辨率图片（327KB）下载：下载全尺寸图片

图4. 组装后的单基因的基因本体分类（GO）。结果被分为三个主要类别：细胞成分、分子功能和生物过程。右侧的y轴表示每个类别中的基因数量。左侧的y轴表示特定类别在该主要类别中的基因百分比。

根据这项研究，共有8,164个单基因（占17.09%）被注释到KOG中，并被分类为26个KOG功能簇。其中，通用功能预测簇的单基因数量最多（1,213个），信号传导机制簇的单基因数量位居第二（1,055个），其次是翻译、核糖体结构和生物发生簇（923个基因）（图5）。接下来是翻译后修饰、蛋白质周转、伴侣蛋白簇（922个基因）、细胞内运输、分泌和囊泡运输簇（485个基因）、转录簇（448个基因）以及RNA加工和修饰簇（347个基因）。此外，KOG还注释了大约44个归类为防御机制簇的单基因（图5）。

下载：下载高分辨率图片（449KB）下载：下载全尺寸图片

图5. 根据euKaryotic Ortholog Group（KOG）对组装后的单基因进行了功能分类。X轴表示26个KOG群的名称；Y轴表示该组中注释的基因占所有注释基因的百分比。

通过使用KO注释，我们根据基因参与KEGG代谢途径的情况将基因分为32组（图6）。在这项研究中，涉及的途径包括信号传导（1,229个基因）、翻译途径（932个基因）、运输和分解代谢途径（830个基因）、内分泌系统途径（682个基因）以及折叠、排序和降解途径（512个基因）。研究发现，来自生物系统的免疫系统基因（512个基因）也被识别出来。与KEGG注释相关的免疫系统途径在表5中呈现。这包括B细胞受体信号通路、抗原处理和呈递、趋化因子信号通路以及补体和凝血通路。因此，从这项分析中可以清楚地看到保守的免疫系统通路。

下载：下载高分辨率图片（527KB）下载：下载全尺寸图片

图6. 组装后的单基因的KEGG途径功能分类。KEGG途径被总结为五个主要类别：A. 细胞过程；B. 环境信息处理；C. 遗传信息处理；D. 代谢；E. 生物系统。y轴表示KEGG代谢途径的名称。x轴表示该途径中注释的基因占所有注释基因的百分比。

表5. 在P. vannamei中发现的与免疫系统相关的KEGG途径代码。

通过使用KO注释，我们根据基因参与KEGG代谢途径的情况将基因分为32组（图6）。在这项研究中，涉及的途径包括信号传导（1,229个基因）、翻译途径（932个基因）、运输和分解代谢途径（830个基因）、内分泌系统途径（682个基因）以及折叠、排序和降解途径（512个基因）。研究发现，来自生物系统的免疫系统基因（512个基因）也被识别出来。与KEGG注释相关的免疫系统途径在表5中呈现。这包括B细胞受体信号通路、抗原处理和呈递、趋化因子信号通路以及补体和凝血通路。因此，从这项分析中可以清楚地看到保守的免疫系统通路。

下载：下载高分辨率图片（527KB）下载：下载全尺寸图片

图6. 组装后的单基因的KEGG途径功能分类。KEGG途径被总结为五个主要类别：A. 细胞过程；B. 环境信息处理；C. 遗传信息处理；D. 代谢；E. 生物系统。y轴表示KEGG代谢途径的名称。x轴表示该途径中注释的基因占所有注释基因的百分比。

通过使用KO注释，我们根据基因参与KEGG代谢途径的情况将基因分为32组（图6）。在这项研究中，涉及的途径包括信号传导（1,229个基因）、翻译途径（932个基因）、运输和分解代谢途径（830个基因）、内分泌系统途径（682个基因）以及折叠、排序和降解途径（512个基因）。研究发现，来自生物系统的免疫系统基因（512个基因）也被识别出来。与KEGG注释相关的免疫系统途径在表5中呈现。

图6. 组装后的单基因的KEGG途径功能分类。KEGG途径被总结为五个主要类别：A. 细胞过程；B. 环境信息处理；C. 遗传信息处理；D. 代谢；E. 生物系统。y轴表示KEGG代谢途径的名称。x轴表示该途径中注释的基因占所有注释基因的百分比。

使用KO注释，我们根据基因参与KEGG代谢途径的情况将基因分为32组（图6）。在这项研究中，涉及的途径包括信号传导（1,229个基因）、翻译途径（932个基因）、运输和分解代谢途径（830个基因）、内分泌系统途径（682个基因）以及折叠、排序和降解途径（512个基因）。研究发现，来自生物系统的免疫系统基因（512个基因）也被识别出来。与KEGG注释相关的免疫系统途径在表5中呈现。这包括B细胞受体信号通路、抗原处理和呈递、趋化因子信号通路以及补体和凝血通路。因此，从这项分析中可以清楚地看到保守的免疫系统通路。

图7. 表达基因的Venn图。每个圆圈中的数字总和表示该组中表达的基因总数，重叠部分表示不同组之间共表达的基因。图8显示了差异表达基因的结果。在比较组6（处理组）和组1（对照组）时，有2,868个基因上调，4,721个基因下调，条件是对数变化值（log2FoldChange）大于1且padj小于0.05（具有生物学重复性）（图8）。

从基因表达水平分析中获得的读数被用作生成差异表达基因的输入数据。火山图显示了组1（对照组）和组6（处理组）中差异表达基因的总体分布（图7）。具有统计显著性的上调差异基因用红色点表示，绿色点表示下调差异基因，蓝色点表示无差异基因。聚类分析了组1和组6样本中的基因聚类模式（图8）。聚类分析应用于所有差异表达基因的公共基因，并计算了它们的聚类结果。FPKM聚类分析使用log10(FPKM+1)值进行聚类。红色表示高表达水平的基因，绿色表示低表达水平的基因。颜色范围从红色到绿色表示log10(FPKM+1)值从小到大。

下载：下载高分辨率图片（81KB）下载：下载全尺寸图片

图7. 表达基因的Venn图。每个圆圈中的数字总和表示该组中表达的基因总数，重叠部分表示不同组之间共表达的基因。使用Fpkm>0.3作为标准。

下载：下载高分辨率图片（572KB）下载：下载全尺寸图片

图8. 火山图。x轴显示不同样本之间基因表达的倍数变化，y轴显示差异的统计显著性。统计显著的差异用红色点表示。每个圆圈中的数字总和表示该组中表达的基因总数，重叠部分表示不同组之间共表达的基因。

图9. 组装后的单基因的GO分类（GO）。结果被分为三个主要类别：细胞成分、分子功能和生物过程。右侧的y轴表示每个类别中的基因数量。左侧的y轴表示特定类别在该主要类别中的基因百分比。

表6. 差异表达基因（DEGs）中显著富集的GO术语。对于细胞成分（CC），在细胞外区域、核糖体和外部包封结构下聚集的单基因表现出显著差异表达。对于分子功能（MF），在结构分子活性、核糖体的结构组分和酶调节因子下聚集的单基因也表现出显著差异表达。对于生物过程，核糖体生物发生和翻译下聚集的单基因表现出显著差异表达。上调和下调的DEGs数量在表6中列出。图10展示了与GO富集相关的上调基因的GO富集点图，图11展示了与GO富集相关的下调基因的GO富集点图。其中最上调的基因包括跨膜转运活性、氧化还原酶活性、细胞外区域和碳水化合物代谢过程（图10）。而最下调的基因包括跨膜转运活性、氧化还原酶活性、细胞外区域和碳水化合物代谢过程、结构分子活性、核糖体生物发生、翻译和细胞外区域（图11）。

下载：下载高分辨率图片（176KB）下载：下载全尺寸图片

图9. 组装后的单基因的基因本体分类（GO）。结果被分为三个主要类别：细胞成分、分子功能和生物过程。右侧的y轴表示每个类别中的基因数量。左侧的y轴表示特定类别在该主要类别中的基因百分比。

表6. 差异表达基因（DEGs）中显著富集的GO术语。

ID 描述 pvalue padj 计数
GO:0005 结构分子活性 1.46E-3 4467 372
GO:0005 细胞外区域 4.55E-1 350 120
GO:0003 核糖体的结构组分 9.70E-06 433 167
GO:0005 酶调节因子活性 1.75E-05 87 151
GO:0030 核糖体生物发生 6.99E-05 231
GO:0030 核糖体 0.0001 110 23
GO:0006 翻译 0.0003 90 223
GO:0003 外部包封结构 0.0008 80 137

下载：下载高分辨率图片（379KB）下载：下载全尺寸图片

图10. 与GO富集相关的上调基因的GO富集点图。

下载：下载高分辨率图片（346KB）下载：下载全尺寸图片

图11. 与GO富集相关的下调基因的GO富集点图。

对差异表达基因（DEGs）进行了KEGG通路富集分析，以进一步了解其生物学功能。校正后的P值小于0.05的KEGG通路被认为是显著富集的。在接种微生物免疫刺激剂后，共有371个DEGs的KEGG通路被诱导富集。根据组1和组6之间的DEGs结果，有六条显著富集的通路，其中最显著富集的通路包括细胞色素P450的外来物质代谢（ID:KO00980）、化学致癌作用（ID:KO05204）、药物代谢-其他酶（ID:KO00983）、药物代谢-细胞色素（KO00982）和流体剪切应力与动脉粥样硬化（ID:KO05418）（表7）。表7展示了组1和组6中检测到的前20个富集的KEGG通路。这些上位富集的KEGG通路涉及代谢通路、核糖体、微生物在多样化环境中的代谢、凋亡、氨基糖和核苷糖代谢、肝细胞癌和化学致癌作用。

表7. 在P. vannamei中接种免疫刺激剂后差异表达基因（DEGs）的KEGG通路富集列表。

ID 描述 pvalue padj 计数
KO00980 细胞色素P450的外来物质代谢 1.66E-05 152
KO05204 化学致癌作用 2.24E-05 152
KO00983 药物代谢-其他酶 9.84E-05 135
KO00982 药物代谢-细胞色素 0.0001 32
KO05418 流体剪切应力与动脉粥样硬化 0.0001 146

下载：下载高分辨率图片（290KB）下载：下载全尺寸图片

图12. KEGG通路富集散点图。y轴表示通路名称，x轴表示富集因子。点大小表示不同基因的数量，颜色表示q值。

在这项研究中，我们还识别了差异表达的免疫相关基因，如病原体诱导的交互反应基因、病原体相关分子模式基因、抗菌肽、proPhenoloxidases系统基因、抗氧化基因、吞噬相关基因、凝固系统基因、氧化应激基因、蛋白水解活性基因、DSCAM（唐氏综合征细胞粘附分子异构体）以及其他免疫相关基因。表8列出了（处理组）和组1（对照组）之间免疫相关DEGs及其注释身份。log2(Fold Change)>1值的基因被认为是差异表达的，并且基于padj<0.05（具有生物学重复性）具有统计显著性。对于病原体诱导的交互反应基因，几丁质酶（几丁质酶1前体，部分）和粘蛋白（粘蛋白-5AC类似物）在处理组中上调。而对于病原体相关分子模式基因，C型凝集素、von Willebrand因子A域含有蛋白7、几丁质结合蛋白和beta-甘露糖苷酶在处理组中上调。Penaeidin、crustin和抗脂多糖因子也是高表达的抗菌基因（表8）。proPhenoloxidase激活酶III、超氧化物歧化酶和谷胱甘肽过氧化物酶也是其他高表达的基因，涉及proPhenoloxidases系统和抗氧化基因。参与先天免疫途径的基因，如吞噬相关基因、凝固系统基因、氧化应激基因和蛋白水解活性基因也在处理组中上调。其他与免疫相关通路相关的基因，如Wnt16、Wnt2、Wntless蛋白和toll蛋白以及toll样受体在处理组中也显著过表达（表8）。

表8. 在P. vannamei中接种免疫刺激剂后差异表达基因（DEGs）的免疫相关DEGs列表及其注释身份、注释和log2（倍数变化）。

类别/基因ID（簇） 注释的同源性 log2（倍数变化）
病原体诱导的交互反应基因 438 3.0
几丁质酶 4.61 94
几丁质酶1前体 4.61 105
粘蛋白-5AC类似物 3.81 94
病原体相关分子模式基因 219 1.0
C型凝集素 5.15 65
C型凝集素5 4.36 30
von Willebrand因子A域含有蛋白7 4.46 24
几丁质结合蛋白 2.72 90
beta-甘露糖苷酶 2.72 15

下载：下载高分辨率图片（346KB）下载：下载全尺寸图片

图10. 与GO富集相关的上调基因的GO富集点图。

下载：下载高分辨率图片（346KB）下载：下载全尺寸图片

图11. 与GO富集相关的下调基因的GO富集点图。

对差异表达基因（DEGs）进行了KEGG通路富集分析，以进一步了解其生物学功能。校正后的P值小于0.05的KEGG通路被认为是显著富集的。接种微生物免疫刺激剂后，共有371个DEGs的KEGG通路被富集。根据组1和组6之间的DEGs结果，有六条显著富集的通路，其中最显著过表达的通路包括细胞色素P450的外来物质代谢（ID:KO00980）、化学致癌作用（ID:KO05204）、药物代谢-其他酶（ID:KO00983）、药物代谢-细胞色素（KO00982）和流体剪切应力与动脉粥样硬化（ID:KO05418）（表7）。表7展示了组1和组6中检测到的前20个富集的KEGG通路。

表7. 在P. vannamei中接种免疫刺激剂后差异表达基因（DEGs）的KEGG通路富集列表。

ID 描述 pvalue padj 计数
KO00980 细胞色素P450的外来物质代谢 1.66E-05 152
KO05204 化学致癌作用 2.24E-05 152
KO00983 药物代谢-其他酶 9.84E-05 135
KO00982 药物代谢-细胞色素 0.0001 32
KO05418 流体剪切应力与动脉粥样硬化 0.0001 146在本研究中，已免疫的L. vannamei在免疫过程中的整体免疫反应模式表明，病原体相关分子模式（C型凝集素）激活了先天免疫，并被宿主细胞上的模式识别受体（Toll样受体）所识别[39]。这导致信号通路（Wnt16, Wnt2）的激活以及抗菌肽（抗脂多糖因子、penaeidin、crustin）的表达升高，进一步靶向并杀死入侵的病原体[40]。免疫虾体内黏蛋白、几丁质酶和几丁质脱乙酰酶基因表达的增加表明了几丁质的相关作用[41]。这种相互作用对于应对细菌致病性非常重要，因为像V. parahaemolyticus这样的细菌能够利用菌毛（例如IV型菌毛）与几丁质结合，在宿主表面（通常是外骨骼或肠道腔）定植[42]。免疫还激活了转录活动，触发抗菌反应以消除入侵的病原体。抗菌反应包括proPO系统的激活、AMPs的释放和吞噬作用的激活[43]。proPO基因表达也上调，以支持和维持proPO系统的级联反应。在甲壳类动物中，prophenoloxidase（proPO）系统是重要的先天免疫防御机制，它涉及一系列丝氨酸蛋白酶将无活性的proPO转化为活性的phenoloxidase（PO），从而引发下游的免疫反应，如Toll通路激活、免疫基因合成和黑色素形成[44]。抗氧化基因表达的上调表明了一线抗氧化防御的参与[45]。超氧化物歧化酶、谷胱甘肽过氧化物酶和Down综合征细胞粘附分子（DSCAM）基因表达的上调表明血小板或类似机制以及免疫记忆的激活[43]。此外，在免疫虾中还发现Down综合征细胞粘附分子（DSCAM）的表达上调，这是一种高度可变的蛋白质，可以在虾中作为病原体特异性识别分子发挥作用[46]。该分子的上调表明本研究开发的失活微生物刺激物能够触发虾的替代性适应性免疫，并进一步增强了对AHPND的免疫反应[46]。这些与免疫相关的信号通路的激活通过同样上调的AMPs得到验证，其中包括抗脂多糖因子（ALF）、penaeidin和crustin[46]。上调的基因（Wnt/β-连环蛋白通路和Toll样受体）的检测表明相应的TLRs、IMD、JAK-STAT和细胞质感应通路的激活[40]。目前已知的对疾病抵抗重要的虾免疫信号通路包括Janus激酶-信号转导子和转录激活因子（JAK-STAT）通路、免疫缺陷（IMD）通路、TLRs通路、RNA干扰（RNAi）通路和P38丝裂原活化蛋白激酶（MAPK）通路[40]。因此，这项研究能够识别出使用含益生元的基于饲料的失活微生物免疫刺激剂进行免疫后的不同类型和时间的免疫反应顺序，为虾免疫学知识的发展提供了宝贵的信息。近年来，我们在甲壳类动物免疫系统的几个重要方面获得了实质性见解，包括吞噬细胞的作用、prophenoloxidase级联反应、黑色素形成系统和抗菌肽[47]。虽然一些保守的免疫效应通路（如黑色素形成和抗菌肽的产生）在生化水平上已经得到了相当好的理解，但大多数甲壳类动物（如虾）的免疫反应背后的分子事件仍然未知。我们对甲壳类动物免疫了解的最大空白之一是其对免疫刺激反应的分子基础。这些发现强调了综合管理虾健康的重要性，结合转录组分析的见解和免疫刺激剂的应用来增强虾的免疫系统。通过加深对L. vannamei免疫反应的理解，我们可以开发出更有效的策略来减轻AHPND等疾病的影响，从而支持虾养殖业的可持续发展。这些发现共同促进了虾养殖中抗病策略的发展，强调了虾免疫系统、转录组反应和免疫刺激之间的复杂关系，以提高L. vannamei养殖产品的健康和抵抗力。然而，这项研究承认当前工作的一个局限性是缺乏qRT-PCR验证。未来的研究将重点对选定的免疫相关基因进行针对性的qRT-PCR验证，以进一步证实转录组发现并加强基因表达分析的可靠性。

结论
使用基于饲料的失活微生物细胞和Sargassum sp作为免疫刺激剂对L. vannamei进行免疫后，根据KEGG通路数据库中显示的不同表达基因（DEGs），发现了几种显著的免疫相关基因变化。观察到的免疫调节主要涉及先天免疫信号通路，特别是Toll样受体（TLR）通路、免疫缺陷（IMD）通路和JAK–STAT通路，这得到了模式识别受体和关键信号基因上调的支持。此外，prophenoloxidase（proPO）系统和抗菌肽反应（如penaeidin、crustin和ALF）也被强烈激活，表明这些通路在介导L. vannamei中观察到的免疫刺激效果方面起着核心作用。此外，我们的研究还能够确定使用失活微生物免疫刺激剂后免疫反应的各种类型和时间顺序，显著促进了人们对虾免疫学理解的深入。这些发现有助于通过强调虾免疫系统、转录组反应和免疫刺激之间的复杂关系来开发虾养殖中的抗病策略，从而提高L. vannamei养殖产品的健康和抵抗力。

未引用的参考文献[[48], [49], [50]]

CRediT作者贡献声明：
Md. Ali Amatul-Samahah：概念化、虾养殖、方法学、软件、数据管理、撰写——初稿。
Natrah Fatin Mohd Ikhsan：方法学、研究、监督。
Mohamad Noor Amal Azmai：方法学、研究、监督。
Md. Yasin Ina-Salwany：概念化、可视化、方法学、研究、监督。