原文链接：Emergent coexistence in multispecies microbial communities | Science

Editor’s summary

The question of how organisms coexist in communities is a pivotal issue for ecologists. To answer this question in a natural ecosystem would require the imposing task of isolating and competing all coexisting members. To render the question experimentally tractable, Chang et al. isolated organisms from stable synthetic bacterial communities and competed all possible combination of pairs of organisms to test their ability to live together. Competitive exclusion occurred in a majority of pairs, whereas a minority of pairs coexisted. Therefore, species coexistence is in part the result of networks of interactions and is an emergent property of community assembly, pointing to the importance of sustaining biodiversity. —Caroline Ash

编辑总结

生物如何在群落中共存是生态学家面临的一个关键问题。要在自然生态系统中回答这个问题，需要分离所有共存成员并对其进行竞争，这是一项艰巨的任务。为了在实验中解决这个问题，Chang 等人从稳定的合成细菌群落中分离出生物，并对所有可能的'生物对'组合进行竞争，以测试它们的共生能力。结果发现，大多数'生物对'发生了竞争性排斥，而少数生物对共存。因此，物种共存在一定程度上是相互作用网络的结果，也是群落组合的一个新兴属性，这说明了维持生物多样性的重要性。---卡罗琳 阿什

Abstract

Understanding the mechanisms that maintain microbial biodiversity is a critical aspiration in ecology. Past work on microbial coexistence has largely focused on species pairs, but it is unclear whether pairwise coexistence in isolation is required for coexistence in a multispecies community. To address this question, we conducted hundreds of pairwise competition experiments among the stably coexisting members of 12 different enrichment communities in vitro. To determine the outcomes of these experiments, we developed an automated image analysis pipeline to quantify species abundances. We found that competitive exclusion was the most common outcome, and it was strongly hierarchical and transitive. Because many species that coexist within a stable multispecies community fail to coexist in pairwise co-culture under identical conditions, we concluded that multispecies coexistence is an emergent phenomenon. This work highlights the importance of community context for understanding the origins of coexistence in complex ecosystems.

摘要

了解维持微生物生物多样性的机制是生态学的一个重要愿望。过去有关微生物共存的研究主要集中在'物种对'上，但目前还不清楚多物种群落中的共存是否需要孤立的成对共存。为了解决这个问题，我们在体外对 12 个不同富集群落中稳定共存的成员进行了数百次配对竞争实验。为了确定这些实验的结果，我们开发了一个自动图像分析管道来量化物种丰度。我们发现，竞争排斥是最常见的结果，而且具有很强的层次性和传递性。由于许多在稳定的多物种群落中共存的物种无法在相同条件下进行成对共培养，因此我们得出结论，多物种共存是一种涌现现象。这项研究强调了群落背景对于理解复杂生态系统中共存起源的重要性。

图 1. 富集微生物群落使我们能够检验物种共存的复杂性。(A）在我们的研究中测试的物种共存的两种假设。(B）为了区分这两种假说，我们使用了一个经验系统，该系统由之前在连续生长和稀释周期下组装的富集体外细菌群落构建而成（9）。在插图 I中，我们展示了一个代表性群落的完整组装动态，显示了每个生长期（转移）结束时每个 ESV 的频率。我们只显示频率大于 2% 的 ESV，每个 ESV 用不同的颜色表示。我们选择了 12 个具有代表性的群落，在第 12 次转移时，其 ESV 的丰富度介于 N = 5 和 N = 13 之间（插图 II），并分离了大多数群落成员（彩色条），平均覆盖了 89.4% 的丰富度。灰条代表我们未能分离的 ESV（见材料和方法）。原始数据来自先前的研究（9, 34, 35）。(C) 频率相关动力学预测了经验观察到的平衡频率。经验平衡频率（横轴）被量化为群落组装过程最后四次转移（第九至第十二次转移）中 ESV 的平均频率。为了确定预测的平衡频率 x*（x 轴），我们首先量化了每次转移中每种 ESV 的入侵适应度 F_i = log (x_i/x_i-1)，然后将 Fi 与 ESV 频率进行回归。回归结果显示，95/99 个 ESV 在各自群落的平衡点附近出现了负斜率（图 S3），表明这些 ESV 受到了频率依赖性负选择的影响。在这些情况下，我们将平衡频率 x* 估计为回归线的 x-截距（图 S3 和 S4）。(D) 插图 I 中群落入侵适应性分析的两个例子，显示了负的频率依赖性选择。黄线表示通过最小二乘回归确定的线性拟合结果（上图和下图分别为 N = 11，R² = 0.92 和 N = 11，R² = 0.70）。x 截距用于估计平衡频率 x*，以垂直虚线表示。

图 2. 多物种共存是群落的新兴属性。(A）为了确定孤立的物种配对是共存还是相互竞争，我们以三种不同的初始频率培养每对物种配对。配对物种在与其原生群落相同的培养条件下连续培养八次。所有 12 个富集群落的配对竞争结果见 (B)，群落按每个群落中菌株的数量排序，从最小的（3 个类群）到最大的（10 个类群）。每个条形图上方的数字分别表示 ESV 的数量、分离菌株的数量和测试配对的数量。请注意，有些群落中缺少配对，因为这些配对要么在共培养中没有任何菌落，要么分类模型准确率较低。(C) 144 对共培养的竞争结果。平均频率和 95% 的置信区间由泊松采样确定（N = 1000；见材料和方法）。为清晰起见，我们绘制了所有情况下在时间点 8（T8）平均频率较低的分离物的频率。在共存配对中，最终转移的平均平衡频率用水平虚线表示，95% 的置信区间（根据泊松采样计算，N = 1000）用其周围的阴影区域表示。每个插入网格表示从初始时间点（T0）到最终时间点（T8）的频率变化。背景颜色代表竞争结果，线条颜色表示三个初始频率。为了确定每个实验中 T0 和 T8 之间频率的显著变化，我们使用了 Wilcoxon-Mann-Whitney 检验，N = 2000，显著性阈值为 P <0.05（见材料和方法）。

图 3. 在稳定共存的群落中，物种对之间普遍存在竞争等级。(A) 利用图 2 所示的实验数据，将 12 个群落中的所有分离物根据它们在成对竞争中排除的其他分离物的数量从上到下排序。网络中的灰色节点代表每个分离物。红色箭头从配对竞争中获胜的分离物指向失败的分离物。蓝线连接共存的分离物。(B) N = 77、P = 0.25 的二项分布显示了如果我们随机交换 (A) 中的共存和排斥链接，表现出反向性的配对的预期数量。红色圆圈表示实验确定的破坏反式性的三元组数量（在我们的例子中为零）。

图 S1. 从土壤和植物样本中分离的微生物接种体的组成和多样性。图 S1 是我们之前研究成果（9）的翻版。请注意，我们未能获得接种体 11 的扩增子文库。在本图中，我们按照参考文献（9）中描述的稀释方案，从读数池中随机抽取 N=10000 个读数，创建了一个样本群落。 (A) 叠加条形图显示了 12 个初始接种体中 11 个接种体（在 x 轴上标注为 1-10 和 12）在目分类水平上的 ESV 群落组成。为清晰起见，我们显示了相对丰度超过 1%的目，并将稀有目归入 "其他"（浅灰色条）。(B) 11 个接种物样本的稀释度曲线。y 轴表示固定采样大小（x 轴）的 100 个随机子样本的平均稀有丰富度，采样读数范围为 10-10000 个。11 个样本中唯一 ESV 读数的数量范围为 110-1290 个。

图 S2. 葡萄糖 M9 盐培养基中 26 个群落的时间动态。用于生成此图的数据集最初由参考文献（9）测序，ESV 表由参考文献（34, 35）编辑。堆叠条形图显示了 26 个群落的时间动态，如参考文献（9）图 S2 和 S7 所示。这些群落包括来自同一接种体 2 的 8 个群落（C2R1、 C2R2、C2R3、C2R4、C2R5、C2R6、C2R7、C2R8）和来自同一接种体6 的另外 8 个群落（C6R1、C6R2、C6R3、C6R4、C6R5、C6R6、C6R7、C6R8）、 以及来自不同接种物的另外 10 个群落（C1R4、C3R4、C4R4、C5R4、C7R4、C8R4、C9R4、C10R4、C11R4、C12R4）。请注意，本次研究包括 C1R4、C2R6、C2R8 和 C8R4 四个群落。请注意，C5R4 第 9 次转移以及 C4R4 和 C9R4 第 11 次转移的时间点缺失 (9)。

图 S3. 26 个稳定群落中 ESV 的负的频率依赖性选择。我们显示了在平衡动态（转移 9-12）附近发现的 99 种 ESV 的相对丰度（x 轴）与入侵适存度（y 轴）的关系。请注意，C5R4 第 9 次转移以及 C4R4 和 C9R4 第 11 次转移的时间点缺失 (9)。这些点代表第 i 次转移时的 ESV 丰度 xi 及其入侵适合度 F=log(xi/xi-1)，其中 i = 2, 3,..., 12，log 为自然对数。蓝色垂直虚线表示平衡频率，即 xi 的平均值，其中 i = 9、10、11、12。在 99 个 ESV 中，有 95 个呈现负斜率（在 52 个标有红线的案例中，负相关具有统计学意义，Spearman's ρ： P<0.05，而不显著的相关性为灰色）。本文数据的最初排序见 (9)，ESV 表的编制见 (34)。

图 S4. 在 26 个稳定群落的组装过程中瞬时观察到的 46 个 ESV 的负频率依赖性。为了与图 S3 进行比较，我们在此显示了在最后 4 个时间点（i=9、10、11、12）中都不存在的 ESV 的相对丰度（x 轴）与入侵适配度（y 轴）的关系，这些 ESV 至少有 5 个入侵适配度数据点（即在两个连续时间点的 5 个实例中都存在），因此有 N=46 个 ESV。在 46 个 ESV 中，有 10 个呈现显著负相关（红线）（Spearman's ρ：P<0.05）。

图 S5. 在组装过程中瞬时出现的 ESV 具有负的平均适合度和/或接近于 0 的平衡频率。(A) 对于图 S3 （即 99 个稳定的 ESV）、图 S4 （即 46 个瞬时 ESV，至少有 5 个适合度数据点）所示26 个群落中的每个 ESV，及至少有 3 个适合度数据点的 64 个瞬时 ESV，我们绘制了其从组装时间序列（x 轴）中获得的平均经验适合度值与其线性模型预测的平衡频率（入侵适合度=0 的丰度）。颜色表示分离物是存在于最终的稳定群落中（红色），还是仅在组装过程中短暂存在（灰色）。我们只显示具有负趋势的 ESV（N=175，其中 95 个在最终群落中稳定共存，80 个在组装时间序列中瞬时出现，但在最终稳定群落中不存在 (方法））。这里显示的瞬时 ESV 至少有三个数据点。为清晰起见，我们省略了误差条，误差条显示在面板中。请注意，当回归的 x-截距为负时，预测频率值可能为负；这实际上意味着预测的平衡频率为零。(B) 我们显示了面板（A）中显示的每个适合度值的平均值 ± S.E.M。

图 S6. 来自 12 个代表性群落的 65 个细菌分离物的菌落形态。这些图像是从 65 个单培养平板的原始图像中剪切出来的。请注意，群落 C11R1 的配对竞争实验分两批（B2 和 C）进行，这里我们选择了 C 批的图像。为了使菌落大小具有可比性，我们将焦点菌落置于中心位置，并剪切出一个边长为 200 像素的固定方形窗口，从而剪切出菌落图像。请注意，我们后来删除了三个因配准质量低而与任何 ESV 都不匹配的分离物（图 S7），包括 C2R6 分离物 4 以及 C11R2 分离物 9 和 11（红色方框）。

图 S7. 12 个代表性群落中 65 个分离物（全 16S rRNA，桑格测序）与 102 个 ESV 之间的序列比对。热图显示了 65 个分离物（x 轴）的所有全长 Sanger 序列与 12 个群落中每个群落在转移 12（y 轴）时在最终稳定群落中发现的 102 个 ESV 之间的碱基对错配情况。颜色强度代表 0（最深的蓝色）、1、2、3、4 的错配，以及所有错配等于或超过 5 的排列（最浅的蓝色）。比对使用 R 包 Biostrings 函数 pairwiseAlignment 进行，比对类型设置为局部比对。为了确定 Sanger 序列和 ESV 之间的一对一匹配，我们首先使用了一个通用过滤器，以去除错配 > 4 和一致长度 <200 的低质量比对。根据这些标准，有三个分离物与任何 ESV 都不匹配，包括 C2R6 分离物 4 以及 C11R2 分离物 9 和 11。然后，我们给每个 Sanger 序列分配了一个最佳配对的 ESV 序列，该序列的配对具有最小的碱基对错配（黄色矩形）。请注意，一个 ESV 可能匹配多个菌株。在 44 个匹配的 ESV 中，32 个 ESV 与一个菌株唯一匹配，8 个 ESV 与两个菌株匹配，2 个 ESV 与 3 个菌株匹配，2 个 ESV 与 4 个菌株匹配。

图 S8. 图像处理和菌落特征计算的自动流水线。流程图显示了在该流程中处理一张图像的示例。主要步骤（1）-（8）已在 "材料与方法 "一节中说明。

图 S9. 图像分析中的菌落形态示例。(A) 这里显示的是原始图像（左图）、经过阈值处理的二值图像（中图）以及进一步应用第一个对象过滤器后的图像（右图）。分割后，平板边缘的物体可能会被误认为是物体（红色）。在分水岭步骤后，使用对物体形状有更严格标准的第二滤波器将这些假阳性物体去除。(B) 当菌落密度非常高时，就会出现菌落团块。一个菌落块的形状和强度属性可能与其他同种菌落块截然不同，这可能会在分类过程中误导模型训练。我们在分水岭之后也会使用第二个对象过滤器来去除这种情况，该过滤器对对象的形状有更严格的标准。图中所示为原始图像（左）、分水岭前图像（中）和分水岭后图像（右）。(C) 在某些情况下，一个菌落可能会出现在靠近平板边缘的位置。这种情况下，菌落会被阴影覆盖，从而在很大程度上改变了像素的分布。我们通过滤波器去除这些菌落，滤波器考虑了像素强度变化的七个指标中的异常值。几个菌落形态示例见 (D-H)。最常见的情况是菌落具有（D）均匀的像素强度或（E）径向强度模式的渐变。除了这两种情况外，有些菌落还具有独特的形态，例如（F）浅色外圈、（G）同心圆或（H）深色星形中心。这些特征是区分两种形态的关键证据。为了捕捉这些特征，我们计算了横断面上的像素强度变化。

图S10 使用随机森林进行共培养分类的示例。该图举例说明了153×3=459共培养分类结果之一。159对无菌落的配对中有6对无菌落被排除在外。最上面一行显示了两种单一培养物和共培养物的绿色通道、扣除背景图像，其中分离物1在这种情况下代表C1R2菌株1，分离物2代表C1R2毒株3，当分离物1以95%混合，分离物2中以5%混合时，即为共培养物。下面的面板显示了随机森林的输出：（A） 在重复交叉验证期间以八个值测试超参数MTRY所报告的准确性。（B） 40个菌落特征根据它们的重要性进行排序。（C） -（D）对于共同培养图像中分割的每个对象，来自随机森林的预测概率。（E） 散点图使用了物体的前两个重要特征。浅红色点表示来自分离物1（C1R2菌株1）单培养图像的菌落，而浅绿色点表示分离物2（C1R2毒株3）单培养。暗红色和绿色点分别表示来自共培养的菌落的预测形态类型。（F） 使用顶部两个主成分的类似散点图。

图 S11. 随机森林模型对共培养图像的准确率。模型准确率是指经过训练的模型在验证数据集中正确预测对象的分数。本文中报告的每种球菌的准确率值都是通过计算 30 个验证集的平均值得出的。我们对 159 对物种进行了配对竞争，其中有 6 对没有形成菌落。这样就有 153 对物种，即 153×3 = 459 个共同培养物。(A) 所有 459 张共培养图像的准确率。(B) 从分析中剔除了模型准确度低于 0.9 的 9×3=27 对共培养。

图 S12. 机器和人类结果之间的 CFU 结果。在成对竞争实验的 159 对图像中，我们使用机器学习方法得到了 144 对图像的结果，其中有 6 对图像没有显示菌落，9 对图像的随机森林模型准确率较低。同样，我们通过人工统计图像菌落获得了 128 对的结果，其中有 6 对没有显示菌落，另外 25 对很难用肉眼区分。两组结果的交集为 127 对，因此这里显示的事件总数为 127×3=381。(A) 通过图像分割得出的 CFU 计数与人类计数的比较（R2=0.85；RMSE=17.67；N=381）。(B) 随机森林预测的 CFU 频率与人工结果的对比（R2=0.87；RMSE=0.17；N=381）。红色虚线表示 y=x。

图 S13. 12 个群落中高丰度菌株之间的配对竞争。与图 2B 类似，这里我们只显示了 ESV 在原始群落中丰度大于 5%的菌株的竞争结果总结。排斥仍然是最常见的结果（竞争性排斥：26/84=31%；竞争性排斥路径：35/84=41.7%；总计：61/84=72.6%），而共存则不太常见（稳定共存：11/84=13.1%；竞争性排斥路径：35/84=41.7%；总计：61/84=72.6%）：11/84=13.1%；共存但无互不侵犯证据：6/84=7.1%，总计：17/84=20.2%），6 个案例尚无定论（6/84=7.14%）。

图 S14. 62 个菌株的竞争力等级与丰度等级之间的相关性。我们通过与其匹配的 ESV 的丰度来确定菌株在群落中的丰度（方法）。(A) 由于多个 ESV 与多个分离株匹配（方法），因此与同一 ESV 匹配的菌株丰度相同。为了在对菌株丰度进行排序时考虑到这一因素，我们采用了随机排序的方法。丰度排名时，为了考虑到这一点，我们在菌株丰度相同的情况下随机分配等级。这是用 R 函数 rank 进行的，tie.method='random'。根据这些排名，我们进行了双尾斯皮尔曼相关性分析，以确定它们之间的关系。然后，我们重复这一过程 1000 次。这里我们绘制了 Spearman's ρ 的分布图，其范围在 0.32 到 0.52 之间，表明菌株丰度与其竞争排名之间存在正相关关系。我们展示了（B）全部 62 个菌株（Spearman's ρ=0.437，P<0.001，N=62）、（C）属于 10 个小型群落的 43 个菌株子集（Spearman's ρ=0.395，P=0.009，N=43）和（D）属于两个大型群落的 19 个菌株子集（Spearman's ρ=0.0248，P=0.920，N=19）的 1000 次重复中的一个实例。红色虚线表示 y=x，圆圈的大小表示重叠点的数量。

图 S15. 竞争等级与菌株在葡萄糖培养基中的生长速率相关。(A) 如方法中所描述和总结的，我们在本研究使用的相同 M9 葡萄糖培养基中培养了各菌株的细菌培养物。除一种菌株外，其他所有菌株都能在葡萄糖培养基中生长。(B）我们绘制了竞争等级与生长速度的对比图，并计算了两者之间的相关性（Pearson's r = -0.314，P = 0.0129，N = 62）。红色虚线表示线性拟合，斜率 = -3.177，截距 = 5.389。

图 S16. 呼吸发酵菌的竞争排名高于呼吸菌。我们使用细菌分离物的全长 16S rRNA 基因进行了分类（方法）。按照参考文献（28）的方法，我们根据菌株的分类学科系为其分配了发酵特征--呼吸发酵菌或呼吸菌。例如，属于肠杆菌科和气单胞菌科的菌种被称为呼吸发酵菌，而我们的菌株中的呼吸菌包括属于假单胞菌科、摩拉菌科、黄单胞菌科、天牛科和科摩单胞菌科的菌种。然后，我们绘制了本研究中使用的 62 株菌株的竞争等级图。呼吸发酵菌株（平均值=2.54）的竞争等级高于呼吸发酵菌株（平均值=4.72）（Wilcoxon Mann-Whitney 检验，P < 0.001，N = 62）。

图 S17. 竞争实验的结果与序列相似性无关。x 轴上的系统发育距离是通过它们全长 16S rRNA 基因中核苷酸的差异数量来衡量的，简而言之，我们首先根据每个碱基的质量得分来修剪 16S 序列，使用的是 Richard Mott 修改后的修剪算法，该算法改编自 Biopython.SeqIO (47)。其次，我们对两个分离物的修剪序列进行了比对，并计算了两个修剪序列之间的碱基对错配。我们使用 EMBOSS 针算法（48）进行了比对和碱基对比较（方法）。这里我们绘制了 N=132 个碱基对的错配值分布图，其结果是竞争性排斥或共存（Wilcoxon-Mann-Whitney 检验，P=0.13）。

表 S1. 物体特征。从图像处理管道中检测出 40 个菌落的物体特征，并在随机森林分类中加以考虑。

表 S2。群落概述和内部标签。从数百个稳定群落中选出的 12 个代表性群落如下。为便于正文阅读，我们将这些群落编号为 1 至 12。在整个实验、扫描仪图像以及分析脚本中，我们对这些群落使用了内部标签。详见 "培养物和图像的命名规则 "一节（方法）。测试的物种对是指主图中显示的物种对，在去除不适用的物种对之后，这些物种对要么在共培养物中没有群落，要么随机森林模型准确率较低。

表 S3。用于单培养和菌落计数的分离物和图像文件名。时间 "一栏代表拍摄单培养图像的转移时间（方法）。