《Frontiers in Microbiology》：Mosquito viromes across England and Wales reveal hidden arbovirus signals and limited ecological structuring
引言
近年来，欧洲虫媒病毒活动加剧，以往零星的检测已被某些地区的持续传播所取代，并且历史上未曾出现的病毒也开始出现在新的区域。乌苏图病毒和西尼罗河病毒现已在中欧和南欧部分地区稳定存在，其中乌苏图病毒在英国于2020年在鸟类和蚊子中被检测到，成为第一个地方性蚊媒病毒；西尼罗河病毒则于2023年首次在英国蚊子中被检测到，这反映了英国与欧洲更广泛的虫媒病毒趋势日益一致。尽管乌苏图病毒和西尼罗河病毒目前被视为主要的蚊媒威胁，但其他虫媒病毒，包括甲病毒和正布尼亚病毒等，也已在欧洲蚊群中报告，凸显了在欧洲大陆传播的虫媒病毒范围广泛，这可能对英国构成未来的出现风险。
尽管有这些检测，英国尚未进行大规模的蚊媒相关病毒多样性调查。与大多数地区一样，当前的监测仍然集中在少数已确定的威胁上，主要通过靶向PCR或脊椎动物血清学进行。为了弥补这一差距，高通量测序提供了一个强大且无偏见的替代方案，能够同时检测已知的虫媒病毒、高度分化的分类群以及与蚊子先前没有关联的病毒。来自欧洲、亚洲和美洲的最新宏基因组研究证实，蚊群蕴藏着异常丰富的病毒群落，包括昆虫特异性病毒和宿主范围或致病潜力不确定的新谱系。
在这些检测中，昆虫特异性病毒因其已被证明能够调节虫媒病毒的复制和传播而受到越来越多的关注。例如，几种昆虫特异性黄病毒已被证明可以降低寨卡病毒、西尼罗河病毒和登革热病毒在伊蚊和库蚊中的传播或复制。提出的机制包括超级感染排斥，即密切相关的病毒竞争相似的复制生态位和细胞因子，或通过激活宿主抗病毒途径产生更广泛的免疫启动效应。尽管其在野生种群中的作用尚不确定，但昆虫特异性病毒正越来越多地被研究作为生物防治策略的候选者。
先前的研究曾利用宏基因组测序调查了英国两个动物园的蚊子病毒组，分析了超过4000只库蚊pipiens复合群和环纹库蚊，鉴定出26种病毒，包括首次在英国报告了两种具有潜在虫媒病毒能力的新型正布尼亚病毒。然而，这项基于动物园的研究地理范围有限，限制了对整个英国景观中病毒流行率和多样性生态驱动因素的推断。
本研究在此基础上，首次对英格兰和威尔士的库蚊物种进行了全面的病毒组调查，分析了从93个地点收集的蚊子。本研究的目标是：描述本地库蚊种群相关病毒的多样性和系统发育关系；检查病毒组组成的空间和生态模式；识别可能影响媒介能力的候选昆虫特异性病毒；检测可能与动物或公共卫生相关的病毒。通过这样做，首次提供了英国蚊媒相关病毒多样性的国家级评估，突出了该地区蚊子种群中存在的病毒多样性，并为未来的监测策略提供了信息。
方法
蚊子采集与合并
蚊子采集于2023年7月，对应于英国库蚊成虫活动的高峰期，并允许在可比较的季节条件下同步采样所有200个地点。采集使用BG-PRO?陷阱，辅以BG-Lure?和BG-CO2发生器，同时使用BG-GAT?产卵陷阱，每个地点运行72小时。蚊子储存于-80°C直至处理。
属于库蚊pipiens复合群和洪流库蚊的标本通过形态学鉴定并经PCR确认。每个地点将1到10个个体合并形成一个池，具体取决于地点产量。如果从一个地点收集到超过10只蚊子，则制备多个重复样本。使用珠磨仪在含有2毫米二氧化硅珠的100 μL蛋白酶K缓冲液中匀浆全蚊。其中50 μL留作物种鉴定，50 μL用于合并。合并体积在需要时用1×PBS调整至500 μL。只有雌蚊被用于病毒组测序。
总共为测序生成了代表93个地点的151个池（总计948个个体）。包含一个仅含PBS的样本作为阴性对照。对于阳性对照，一只库蚊pipiens molestus雌蚊喂食含有乌苏图病毒的血液餐，最终浓度为4.0 × 107pfu/mL，对应于每只蚊子约4.0 × 104pfu的估计剂量（假设摄入约1 μL血液）。
核酸提取和病毒RNA富集
合并的匀浆液在4°C下以16,000 × g离心5分钟，300 μL澄清上清液通过0.45 μm无菌旋转过滤器过滤。如果发生堵塞，将剩余材料转移至新的旋转柱，直至所有上清液处理完毕。
过滤后的匀浆液用2单位TURBO DNase处理以去除宿主和细菌DNA。RNA使用RNAClean xp磁珠按照制造商说明进行纯化。使用NEBNext rRNA去除试剂盒，辅以靶向库蚊保守rRNA区域的自定义探针进行核糖体RNA去除。去除遵循制造商方案，并添加了蚊子特异性探针。
使用安捷伦5300 Fragment Analyzer评估RNA质量和片段大小分布，并使用Qubit? RNA HS（高灵敏度）检测试剂盒测定浓度。进行逆转录和序列非依赖性单引物扩增以富集病毒RNA，遵循所述的改良方案。
使用NEBNext Ultra II FS DNA Library Prep Kit for Illumina制备文库，包括片段化、末端修复、适配器连接和索引添加。使用AMPure XP磁珠进行清理。在测序前，使用Qubit? 1X dsDNA高灵敏度检测试剂盒对文库进行定量，并使用安捷伦5300 Fragment Analyzer验证片段分布。
Illumina测序
所有文库在Illumina NovaSeq X Plus平台的两个通道上使用25B化学试剂进行测序，产生150 bp双端读长，生成了33.32亿条读长。
读长处理和组装
使用Cutadapt版本4.5从原始FASTQ文件中修剪Illumina接头和SISPA引物序列。进一步修剪读长以去除最低窗口质量分数为20的低质量碱基。然后去除短于15 bp的读长，并使用FastQC v0.12.1评估测序质量。使用默认参数进行从头组装，仅保留长于1,000个核苷酸的重叠群用于进一步分析。
初始病毒信号检测
使用同源性和基于特征的方法组合鉴定推定的病毒序列。首先，使用BLAST+ v2.15.0将重叠群与Virus-Host DB病毒蛋白数据库进行比较，使用1 × 10?5的e值截止值和每个高得分片段最小查询覆盖度30%，作为故意宽松的第一轮过滤器，然后在后续进行更严格的下游验证和整理。同时，运行VirSorter2 v2.2.3以默认参数检测RNA病毒。通过任一方法被识别为病毒的任何重叠群都被带入后续步骤。
重叠群延伸和病毒基因检测
为了提高连续性，使用基于重叠群重叠的重组处理候选病毒重叠群。使用Prodigal v2.6.3以“meta”模式从这些延伸的重叠群预测蛋白质编码基因，并使用HMMsearch将所得蛋白质序列针对RVDB-prot v29.0进行筛选。保留含有与病毒家族具有显著相似性（e值 ≤ 1 × 10?5）的蛋白质的重叠群。
完整性评估和过滤
使用ViralQC评估病毒重叠群的基因组完整性。保留估计完整性至少为50%的重叠群。对于未被ViralQC评分的重叠群，使用救援流程进行评估，其中预测的蛋白质使用MMseqs2 v14.7e284针对自定义的ICTV衍生NR蛋白质数据库进行查询，并使用最低共同祖先方法分配分类学。根据MMseqs2分配估计完整性，并应用相同的≥50%阈值。整合两种方法的结果以产生高置信度的病毒重叠群集合。
去重复和基因组过滤
为了减少冗余，使用dRep v3.4.0对高置信度重叠群进行去重复，使用最小重叠群长度1,000 bp，主要聚类阈值90%平均核苷酸一致性，次要阈值95%平均核苷酸一致性。对去重复集合重新进行Prodigal分析，仅保留含有至少一个完整开放阅读框（partial = 00标志）的基因组。
临时分类学注释和验证
通过使用BLASTx针对Virus-Host DB进行临时注释，e值截止值为1 × 10?5，最小查询覆盖度为30%，每个重叠群最多保留五个命中。这些分配用于指导系统发育定位。为了验证组装质量，使用bwa-mem2 v2.2.1将读长映射回保留的基因组，并在IGV v2.12.3中检查覆盖度。在下游分析之前，修剪读长支持不一致的末端区域。
系统发育分析和分类学分配
对于树重建，仅使用含有完整标记基因的去重复重叠群。选择RNA病毒的RNA依赖性RNA聚合酶和DNA病毒的复制相关蛋白。使用NCBI ORFfinder预测开放阅读框。每个标记ORF使用BLASTp与NCBI nr数据库进行比较，并检索顶级命中以及根据国际病毒分类委员会参考物种策划的代表性序列。
使用MAFFT v7.525为每个病毒科或目生成多序列比对，使用--maxiterate 1000 --globalpair选项以最大化比对准确性。使用trimAl v1.5去除比对不佳的位置，使用间隙阈值0.75和块大小10。然后使用IQ-TREE2 v2.3.4重建最大似然系统发育树，使用1,000次超快自举重复评估分支支持。所得树在FigTree v1.4.4中手动重新根以优化可解释性。使用ggtree包v3.17.1在RStudio v4.3.2中可视化树。
对于分类学分配，物种级别的调用遵循ICTV指定的科特异性划分标准，使用推荐的氨基酸一致性阈值（如果可用）。对于尚未定义正式ICTV物种划分阈值的给定物种，与其最接近的参考序列在RdRp（或等效的标志性蛋白）中共享<90%氨基酸一致性的谱系被视为不同的临时物种，与ICTV对谱系给出的最常见截止值一致。
病毒丰度估计和可视化
遵循所述方法量化病毒丰度。使用bwa-mem2 v2.2.1将读长映射回最终去重复的病毒重叠群，并使用CoverM v0.7.0估计重叠群水平的丰度。如果重叠群长度的至少50%被映射读长覆盖，则认为该重叠群存在于池中。对每个池的病毒重叠群读长计数求和以生成丰度矩阵。随后使用pheatmap包在Rstudio中可视化该矩阵。
病毒群落的生态和地理分布
为了探索空间和生态模式，使用上述映射读长和覆盖度标准定义地点水平的病毒存在-缺失。对于每个合并文库，通过对分配给该病毒的所有重叠群的读长求和获得病毒水平的读长计数。对于生态分析，来自同一采样位置的池被视为生物学重复，并通过对来自同一地点的所有池的病毒读长计数求和汇总到地点水平；如果该地点水平的总计数大于零（即，在该地点的至少一个合并库中检测到），则认为病毒存在于该地点。这些地点水平的计数（跨池的总和）支撑了所有后续基于丰度的分析，然后将所得的地点水平存在-缺失矩阵映射到监测位置，使用R中的散点饼图包在地图上叠加基底图来可视化群落组成。将分配给包含公认虫媒病毒物种的属的病毒定义为推定虫媒病毒，并分别绘制以突出其跨地点的分布。
计算病毒的相对丰度以评估跨环境和区域梯度的变化。病毒读长计数在科水平汇总，在每个地点内标准化，并在组间平均。对于生态比较，每个陷阱地点根据英国生态与水文中心土地覆盖图分配一个土地类型类别，使用陷阱位置周围1 km × 1 km正方形内的主导土地覆盖。原始调查设计在21个UKCEH土地覆盖类别中分配地点；为了分析，这些被汇总为城市和农村类别，“农村”包括所有非城市/郊区类别。区域梯度使用英格兰和威尔士的一级国际领土级别区域定义，遵循英国国家统计局的分类。使用ggplot2 v3.5.1进行可视化。
为了测试生态关联，使用地点水平的存在-缺失矩阵比较城市和农村地点之间的检测频率。对每个病毒物种和科独立执行Fisher精确检验，并应用错误发现率校正。同时，进行地点水平的差异丰度测试，以评估特定的病毒分类群是否在城市与农村地点中富集。读长计数跨重复汇总，按物种或科聚合，并使用DESeq2（v1.36）进行分析。为了减少稀有分类群的影响以及由少数样本中高度偏斜的读长分布驱动的虚假富集，在差异丰度测试之前应用了流行度过滤器。具体来说，我们仅保留在城市和农村组中至少约20%的地点（每组≥9个地点）存在的科，这是在基于DESeq2的微生物组和宏基因组分析中常用的流行度阈值，以限制极端稀有、零膨胀特征对模型的影响。仅对通过此过滤器的科检查物种水平模式，以帮助识别科水平信号的潜在贡献者。此外，量化了合并工作量（每个地点的池数）以评估每个地点合并库数量的变化是否会偏倚地点水平比较；因此测试了合并工作量是否随土地类型（Kruskal-Wallis检验）和纬度/经度（Spearman等级相关）存在系统性变化，未检测到系统性差异。
Alpha和Beta多样性分析
为了检查地点内和地点间的病毒多样性，将病毒读长计数稀释到对应于非零文库大小第5百分位的共同深度，稀释重复1,000次并保留平均多样性值。
使用观察到的丰富度（不同病毒分类群的数量）和香农多样性量化Alpha多样性，在vegan v2.6-4中计算。使用Wilcoxon秩和检验比较跨土地类型（农村和城市）和蚊子物种的多样性值，并使用Benjamini-Hochberg程序调整p值以进行多重比较。首先使用Spearman等级相关评估Alpha多样性与地理坐标（纬度和经度）之间的关联，随后使用最小二乘回归评估线性趋势的强度。
使用根据相对丰度矩阵计算的Bray-Curtis相异性评估Beta多样性。通过主坐标分析和非度量多维标度在vegan中进行排序，并在ggplot2中可视化排序图。使用PERMANOVA（9,999次排列）测试土地类型、纬度和经度对病毒群落组成的影响，重点关注库蚊pipiens以允许跨土地类型的平衡比较。使用基于质心的距离评估多元离散度的同质性。
结果
分类学广度和系统发育定位
病毒序列处理流程在151个文库中产生了253个去重复的重叠群，每个重叠群包含至少一个完整的ORF。其中，为41个不同的分类群恢复了完整的标志基因（RNA病毒的RdRp或DNA病毒的Rep），涵盖RNA病毒和一个DNA病毒。在英格兰和威尔士采样的93个地点中，有86个地点检测到至少一种这些病毒。这些包括负链RNA病毒（10个）、正链RNA病毒（22个）、双链RNA病毒（8个）和单链DNA病毒（1个）。系统发育重建确认了大多数谱系在已识别的病毒科中的位置，包括Iflaviridae（5个）、Solemoviridae（3个）、Tymoviridae（3个）、Peribunyaviridae（2个）、Partitiviridae（2个）、Rhabdoviridae（2个）、Sedoreoviridae（2个）、Orthomyxoviridae（2个）、Xinmoviridae（2个）、Amalgaviridae（1个）、Chrysoviridae（1个）、Chuviridae（1个）、Dicistroviridae（1个）、Draupnirviridae（1个）、Mesoniviridae（1个）、Nodaviridae（1个）。此外，几个序列聚集在已建立的ICTV指定病毒科之外，包括4个Negev样病毒、库蚊布尼亚病毒2（目：Hareavirales）、Daeseongdong样病毒2、两个Ghabrivirales物种、两个Tolivirales物种和一个Tymovirales物种。
总共有11种病毒符合ICTV的新物种标准，其RNA依赖性RNA聚合酶氨基酸与其最接近的已知亲属的一致性在31%到84%之间。除Ghabrivirales sp. 1和Ghabrivirales sp. 2在RdRp中共享95%氨基酸一致性，因此根据ICTV划分标准被视为单个临时物种外，所有分类群基于去重复和系统发育标准都是不同的。
分类学亮点
有12种病毒在本次全国调查和先前的基于动物园的调查中均被检测到。其余的检测代表了先前数据集中未观察到的分类群。在RNA病毒中，Picornavirales（5个Iflaviridae和1个Dicistroviridae）和Mononegavirales（4个分类群）的成员很突出；系统发育分析将所有这些都是昆虫特异性分支中。鉴定出两个来自Quaranjavirus属的分类群（武汉蚊子病毒4和武汉蚊子病毒6）和四个Negev样病毒，与已建立的蚊媒分支分组。除了昆虫相关分类群之外，还存在几种通常与植物或真菌相关的谱系，包括Solemoviridae、Chrysoviridae、Ghabrivirales、Partitiviridae和Amalgaviridae的成员。一种partitivirus与先前在英国报道的库蚊pipiens betapartitivirus 2匹配，而另一种代表了一种新的deltapartitivirus。
虫媒病毒相关检测
除了昆虫特异性谱系外，还检测到几种与公认虫媒病毒密切相关的病毒。Hedwig病毒（科：Peribunyaviridae；属：Gryffinivirus）分布最广，在英格兰南部和威尔士的10个地点观察到：斯旺西、纽波特、布里斯托尔、伊普斯维奇、梅尔顿莫布雷、纽马基特、斯利姆布里奇、上斯托克以及两个伦敦地点（哈尔斯登和沃尔沃思）。在这些检测中的六个中，恢复了完整的RdRp ORF，每个都与先前报道的Hedwig病毒序列显示>98%的氨基酸一致性。有些与德国分离株聚类最密切，其他则与瑞典或法国报告的病毒聚类，表明与多个欧洲谱系关系密切。
Umatilla病毒（科：Reoviridae；属：Orbivirus）在三个地点发现，包括斯利姆布里奇、纽波特和普利茅斯。这两个地点的代表性序列与2019年监测期间在德国捕获的鸟类中获得的其他序列聚类。
Atherstone病毒（科：Peribunyaviridae；属：Orthobunyavirus）仅限于斯温顿和剑桥的两个地点，其RdRp基因与先前研究中报告的病毒显示近乎相同的氨基酸一致性，并与最近从法国2015年检测发布的部分序列密切相关。
对于这些病毒中的每一种，所有预期的基因组片段都被恢复（Umatilla病毒的片段3除外），并且共存在于单个池中，确认了近乎完整基因组的组装而不是部分检测。
病毒分布模式
病毒检测范围从广泛分布到高度受限的分类群。只有三种病毒，Daeseongdong病毒2、武汉蚊子病毒4和Alphamesonivirus fluvideense，分布广泛，每种在超过三分之一的所有地点以及所有10个ITL区域被检测到。十六种病毒显示中间分布，出现在6-25个地点，包括Chrysoviridae sp.（25个地点）、库蚊Negev样病毒1（16个地点）和Marma病毒（18个地点）。相比之下，大多数分类群（41个中的22个）是受限的，在四个或更少的地点发现，其中八个仅观察到一次。虽然大多数单例先前未在库蚊中报道过（例如，Amalgaviridae sp.、库蚊Negev样病毒3、库蚊pipiens Tymo样病毒1和库蚊pipiens Tymovirales sp.），但其他如库蚊pipiens betapartitivirus 2、Almendravirus Chester和Valmbacken virus已在早期研究中记录，支持它们可能的蚊媒关联，尽管此处流行度低。
在跨越86个地点的312个分类群检测中，六个科（Mesoniviridae、Chrysoviridae、Orthomyxoviridae、Iflaviridae、Xinmoviridae和Solemoviridae）占所有记录的半数以上（59%）。另外23%的检测属于“未分类”类别，反映了无法置于已建立科中的病毒。这些未指定检测中的大多数反映了Daeseongdong病毒2，它分布广泛，在所有10个ITL区域的70个地点被检测
