高通量数据库有哪些？细菌高通量数据荟萃分析的可行性研究

细菌(Bacteria)作为最重要的微生物类群,其分布涵盖地球主要生境，包括各类极端环境，参与各生境中的物质和能量循环,对于维持生态系统的功能和稳定具有极其重要的作用。

据估算，目前仅约2%的细菌被详细调查研究。

测序技术的兴起为细菌多样性研究带来了重要变革，第二、三代高通量测序的普及和应用，导致更多的细菌种类被发现，其更多的生态功能正不断被揭示。

随着NCBI、EBI等高通量数据库的开放，以及Muther、Qiime1、Qiime2等高通量数据分析软件的完善，通过收集已发表的细菌高通量数据，进行整理、再分析，可以打通不同研究的壁垒，扩大了其研究尺度，取得了一些全球性的研究成果。

然而，这种剪切得到的目标区段(可称“伪”区段)与“真”引物扩增得到的区段是否存在差异，这一核心问题在上述研究中均未被考虑过。

避免讨论该问题，将忽略不同引物所造成的测序影响，从而影响以上研究结果的准确性。

为此，本研究使用了78个包含V3_V4、V4以及V4_V5区段的16S rRNA高通量测序数据,并使用Qiime2平台比较了V3_V4和V4_V5串联区截取得到的“伪”V4区段与“真”V4区段(即原生V4扩增片段)所得到的细菌群落，来探究提取细菌统一扩增区段的荟萃分析方法的可行性和准确性，该研究能为今后类似分析方法提供重要的参考。

材料与方法

材料

使用Wear等测定的83个海水样品中的细菌16S rRNA数据，该数据已上传至NCBI,接收号为PRJNA412105。

该数据集包含四对不同引物的测序结果，选取了其中V3_V4、V4以及V4_V5区段的高通量数据进行分析，对应的测序引物信息如表1所示。

为排除较低测序对后续结果的影响，将包含任意测序深度低于5 000 reads的样品剔除，共得到234份高通量数据。

表1 扩增16S rRNA不同区段的引物信息

目标 区段

引物 名称

引物序列(5′-3′)

扩增 长度/bp

V3_V4

341F

CCTACGGGNGGCWGCAG

460

785R

GACTACHVGGGTATCTAATCC

V4

515F

GTGYCAGCMGCCGCGGTAA

292

806R

GGACTACNVGGGTWTCTAAT

V4_V5

515F

GTGYCAGCMGCCGCGGTAA

411

926R

CCGYCAATTYMTTTRAGTTT

方法

所有高通量分析使用Qiime2(2021.04)平台进行,主要步骤：第一，训练分类集(train feature classifier)。

使用Silva 138(2020.11)16S rRNA全长数据库作为参考数据库。

从全长数据库中截取目标区段，提高物种分类精度，主要使用RESCRIPt剔除原本重复序列和截取引物后目标序列的重复序列；使用feature-classifier extract-reads命令分别截取对应引物序列，设置chunk-size=1 000的运算通量来训练细菌分类数据集。

第二，合并双端测序结果。

使用Vsearch 2.17.1,按照-fastq_maxdiffs 30-fastq_minovlen 20-fastq_maxee 0.5的参数进行合并。

先合并序列能避免在下一步引物剪切后，各双端剩余序列过短，影响序列合并。

第三，剪切前段引物，以保留V4区段，该步骤主要对V3_V4区段使用515F引物进行剪切，保留后端序列，即为“伪”V4区段。使用cutadapt剪切对应引物。

第四，细菌群落分析。

该分析主要分为三部分：①单独分析各扩增区段，即对V3_V4区段、V4区段以及V4_V5分别截取至420 bp、260 bp以及360 bp长度后进行单独分析。

②单独分析V3_V4以及V4_V5区段截取后的V4区段，包括V3_V4区段剪切后的“伪”V4区段(剪切至260 bp),以及截取至260 bp的V4_V5区段(剩余序列即为V4区段)。

③将所有V4区段整合至一起进行分析，即将所有V4区段序列合并后，截取至260 bp进行分析。

序列使用DADA2进行嵌合体筛选和划分扩增子序列变异(Amplicon Sequence Variants, ASV,可简单理解为100%相似度进行聚类),ASV使用feature-classifier classify-sklearn命令，结合对应区段的训练分类集进行细菌ASV物种鉴定，其中置信度最低为0.8。

最后，将所有细菌矩阵reads数抽平至1 849,以避免不同测序深度对结果的影响。

数据分析

细菌群落多样性、群落组成等均在R语言平台分析完成,主要使用“vegan”程序包。

多样性差异使用配对的t检验进行，采用最小二乘法(OLS)拟合多样性和相对丰度之间的线性关系，使用Mantel test进行群落相似性分析。

结果与分析

不同扩增片段下细菌群落差异

同样测序深度下，不同扩增区段得到的细菌各纲相对丰度存在明显差异,得到的多样性也有显著不同。

其中，扩增V4区段得到的细菌丰富度(Richness)最高，V3_V4区段次之。

群落均匀度(Evenness)也有显著差异，扩增V3_V4区段的均匀度最高，V4区段次之。

三种扩增区段下，细菌丰富度相关性显著,其中V3_V4与V4_V5扩增下，细菌丰富度相关程度最高。

就均匀度来看，V4与V4_V5扩增下细菌群落均匀度相关性最高。

通过对群落相似性距离的相关性分析发现，三种扩增区段所构建的群落存在着一定程度的相似，其中扩增V4与V4_V5之间相似程度最高。

进一步对丰度最高的三种细菌纲之间的相关性分析发现，V4和V4_V5扩增得到的三种细菌纲相对丰度的相关性最高。

不同引物下细菌组成和多样性的差异

A:主要细菌纲；B:丰富度；C:均匀度。

表2 不同引物或剪切区段细菌多样性相关性

比较区段

目标区段

丰富度

均匀度

r

P

r

P

原始区段

V3_V4 vs V4

<0.001

V3_V4 vs V4_V5

<0.001

V4 vs V4_V5

<0.001

<0.001

原始区段与V4区段

V3_V4 vs V4(V3_V4)

<0.001

V4_V5 vs V4(V4_V5)

<0.001

<0.001

所有V4区段

V4(V3_V4) vs V4

V4(V3_V4) vs V4(V4_V5)

<0.001

<0.001

V4 vs V4(V4_V5)

<0.001

<0.001

注：V4后面括号代表原引物扩增区段，下同

不同引物和剪切区段细菌群落组成的相关性

表4 不同引物和剪切区段主要细菌纲丰度的相关性

比较内容

目标区段

Alphaproteobacteria

Bacteroidia

Gammaproteobacteria

r

P

r

P

r

P

原始区段

V3_V4 vs V4

<0.001

<0.001

<0.001

V3_V4 vs V4_V5

<0.001

<0.001

<0.001

V4 vs V4_V5

<0.001

<0.001

<0.001

原始与V4区段

V3_V4 vs V4 (V3_V4)

<0.001

<0.001

<0.001

V4_V5 vs V4 (V4_V5)

<0.001

<0.001

<0.001

所有V4区段

V4 (V3_V4) vs V4

<0.001

<0.001

<0.001

V4 (V3_V4) vs V4 (V4_V5)

<0.001

<0.001

<0.001

V4 vs V4 (V4_V5)

<0.001

<0.001

<0.001

同种引物下不同区段的细菌群落差异

通过对V3_V4以及V4_V5区段进行剪切，得到“伪”V4区段，分析发现，“伪”V4区段主要细菌纲的组成与原始区段差别不大,但物种多样性差异明显。

剪切后的“伪”V4区段丰富度高于剪切前区段,所得的均匀度低于剪切前区段。

剪切前后的区段所得的细菌多样性之间相关性显著，但V3_V4区段剪切前后均匀度相关程度相对较低。

剪切前后细菌群落之间具有极为显著的相关性，其中V4_V5区段剪切前后群落相似性更高,且主要的三种细菌纲相对丰度之间相关性也高于V3_V4剪切前后。

同一区段下不同引物细菌群落差异

V3_V4截取后的V4区段与V4_V5截取后的V4区段细菌群落丰富度相关性最高；V4区段与V4_V5截取后的V4区段均匀度相关性最高；而V3_V4截取后的V4区段与V4区段均匀度并无明显相关关系。

剪切至V4区段后细菌组成和多样性的差异

不同引物的V4区段细菌组成和多样性的差异

不同引物的V4区段细菌群落NMDS排序

椭圆圈代表对应颜色细菌群落的95%置信区间

16S rRNA基因的可变区具有属或种的特异性,通过单一可变区或可变区串联测序可以确定微生物物种组成。

然而，研究发现不同可变区或串联组合测序所得到的结果均相同。

本研究通过对不同扩增引物进行剪切，得到统一的V4区段数据分析结果来看，同种引物下，剪切前后的群落主要细菌纲相对丰度几乎相同。

该结果表明不同区段引物扩增导致的群落差异并非由不同可变区对细菌种鉴定差异引起，而是由引物特异性导致。

另外，Marcus等研究表明V3_V4、V4_V5这两个可变区串联组合对细菌的鉴定准确性较高，而本研究中V4区段引物扩增后细菌丰富度最高，表明即使V4区段引物515F/806R受到地球微生物组计划(EMP)的推荐,但该引物有可能会高估细菌多样性，并损失细菌物种鉴定的准确性。

本研究使用了相似度阈值在100%的ASV进行代表序列划分，该分析中单个碱基的差异将导致划分出新的ASV。

基于该分析方法，理论上序列长度越长，变异的ASV也会越多。

本研究V3_V4和V4_V5剪切至V4区段后，部分前段(V3_V4)或后段(V4_V5)的可变区序列被移除，理论上，序列变短将导致检测到的ASV数量下降，而实际结果与此相反，这与理论结果不相符。

这是因为DADA2会将丰度极低的物种当做测序错误而剔除,部分V3_V4前段和V4_V5后段涵盖不同碱基的低丰度序列被剪切后，剩余V4区段的低丰度序列将被被识别为“较高”丰度序列，从而纳入分析，最终提高了群落的丰度。

剪切后，群落更低的均匀度也能佐证部分低丰度的序列被纳入最终群落分析。

合并“伪”V4和“真”V4区段序列统一分析发现，不同引物的同段序列群落组成和多样性具有显著差异，该结果表明通过截取同一段序列进行荟萃分析的方法受原始引物的影响比较大，不能进行整合分析。

但考虑到群落组成的显著相关性，该方法对于群落组成分析具有一定的可行性。

总之，本研究发现不同引物对于群落结构和多样性的调查结果均会产生显著影响，采用序列截取的方法，得到同一区段的荟萃分析，也会引入引物差异，其结果值得重新分析考虑,部分真菌、原生生物等其他微生物高通量数据的荟萃分析结果是否会引入引物误差，还需要更多的研究来探讨。

Bar-On YM,Phillips R,Milo R.The biomass distribution on Earth[J].Proceedings of the National Academy of Sciences of the United States of America,2018,115(25):6506-6511.

Flemming HC,Wuertz S.Bacteria and archaea on Earth and their abundance in biofilms[J].Nature Reviews Micobiology,2019,17(4):247-260.

Bardgett RD,Putten WH,Van Der.Belowground biodiversity and ecosystem functioning[J].Nature,2014,515(515):505-511.