Author response: Competitive interactions between culturable bacteria are highly non-additive

Abstract

Article Figures and data Abstract Editor's evaluation eLife digest Introduction Results Discussion Materials and methods Appendix 1 Data availability References Decision letter Author response Article and author information Metrics Abstract Microorganisms are found in diverse communities whose structure and function are determined by interspecific interactions. Just as single species seldom exist in isolation, communities as a whole are also constantly challenged and affected by external species. Though much work has been done on characterizing how individual species affect each other through pairwise interactions, the joint effects of multiple species on a single (focal) species remain underexplored. As such, it is still unclear how single-species effects combine to a community-level effect on a species of interest. To explore this relationship, we assayed thousands of communities of two, three, and four bacterial species, measuring the effect of single, pairs of, and trios of 61 affecting species on six different focal species. We found that when multiple species each have a negative effect on a focal species, their joint effect is typically not given by the sum of the effects of individual affecting species. Rather, they are dominated by the strongest individual-species effect. Therefore, while joint effects of multiple species are often non-additive, they can still be derived from the effects of individual species, making it plausible to map complex interaction networks based on pairwise measurements. This finding is important for understanding the fate of species introduced into an occupied environment and is relevant for applications in medicine and agriculture, such as probiotics and biocontrol agents, as well as for ecological questions surrounding migrating and invasive species. Editor's evaluation This important study presents an interesting example of how complexities of communities may be reduced by showing that when partner species exert a negative effect on the focal species, the joint effects are generally not additive, but rather dominated by the strongest single effect. The evidence, enabled by thousands of measurements using nanodroplet-based microfluidics, is compelling, although the generality of the conclusion awaits further studies. This paper is of interest to microbial ecologists and synthetic biologists. https://doi.org/10.7554/eLife.83398.sa0 Decision letter Reviews on Sciety eLife's review process eLife digest Bacteria can be found almost everywhere on earth. Often, multiple species of bacteria live together in communities, which perform vital roles that affect everything from our health to the planet’s ecosystems. A single species within this community can sometimes be particularly important, for example if it is causing disease in its host or producing a vital nutrient. However, the other species within this community can influence the growth of this focus species, either by inhibiting or promoting it. It is challenging to predict how a certain species is going to fare within a bacterial community as it remains partly unclear how groups of bacteria affect each other. Some theory suggests that adding up or averaging the influences of all the bacteria in a community would be the best way to predict what will happen. To study these microorganism interactions, Baichman-Kass, Song and Friedman monitored thousands of bacterial communities, consisting of two to four different species. By using species that express fluorescent proteins, they were able to measure the abundance of the specific bacteria of interest in each of these communities. Baichman-Kass et al. found that in communities where all the species were only competing with or supporting the bacteria of interest, the individual affecting species with the strongest effect dominated the combined effect. This ‘strongest effect’ model offered accurate predictions for the joint effects of competitive communities, however predicting outcomes in supporting communities proved more difficult. This could indicate that the mechanisms for supporting other species are more intricate than the means of competition. The study of Baichman-Kass et al. brings us closer to understanding how the abundance of a given bacterium can be influenced through the actions of other bacterial species. Among other uses, it may be important in medicine, where it is desirable to reduce the amount of a bacteria that causes disease, or in agriculture where bacteria that protect plants from diseases and fungi, need to be boosted. Improving our ability to predict the outcome of introducing new species to an environment increases both the effectiveness and possible scope of such applications. Introduction Scarce are the environments on Earth not colonized by bacteria. In addition to naturally playing important roles from driving biogeochemical cycles at the ecosystem level (Cavicchioli et al., 2019; Arrigo, 2005; Falkowski et al., 2008) to supporting host health at the individual level (Berendsen et al., 2012; Manor et al., 2020; Gilbert et al., 2018), bacteria have also been harnessed for countless biotechnological applications across industries such as food preservation (Motarjemi, 2002), plant and animal health (Berendsen et al., 2012; Júnior et al., 2021), biocontrol of pathogens (Köhl et al., 2019), as well as decomposition of toxic compounds, and production of valuable ones (Varjani et al., 2017; Ro et al., 2006; Fang and Smith, 2016; Mainka et al., 2021). In natural environments, bacteria often form rich and complex communities, but understanding how these communities organize has proven difficult (Widder et al., 2016). Elucidating the rules that govern microbial ecology can both offer insight into larger ecological systems and allow us to better manipulate and design microbial communities for the desired functions. The structure of microbial communities is determined by the interactions between the involved species (Konopka et al., 2015; Barbier et al., 2018; Qian and Akçay, 2020). In recent years, much effort has been put into measuring pairwise interactions of different species from, and in, different environments (Foster and Bell, 2012; Vetsigian et al., 2011; Kehe et al., 2021). But it is still unclear to what extent the joint effects of multiple species on a focal species of interest (e.g., a pathogen) can be inferred from pairwise measurements. Such inference may be challenging due to the presence of indirect interactions: the affecting species may alter each other’s abundances (termed interaction chains) or may modify each other’s effect on the focal species (termed interaction modification, or higher-order interactions) (Sanchez, 2019; Wootton, 2002). Despite a strong theoretical foundation, empirical studies in recent years have shown conflicting results regarding the importance of higher order interactions and indirect effects (Levine et al., 2017). For some functions, such as degradation of complex molecules (Sanchez-Gorostiaga et al., 2019; Gralka et al., 2020) and antibiotic production (Tyc et al., 2014; Qi et al., 2021; Westhoff et al., 2021), clear evidence of such interactions has been found. Furthermore, in both empirical and theoretical studies, the presence or absence of an additional species has been shown to affect interactions, and even the outcome of invasion and coexistence in some systems (Mickalide and Kuehn, 2019; Chang et al., 2022; Hsu et al., 2019). Additional theoretical work has shown that commonly used ecological models (i.e., generalized Lotka-Voltera) do not properly capture microbial community interactions, partially due to the nature of these interactions (chemically mediated as opposed to predator-prey) (Momeni et al., 2017). However, other studies have shown that both structures of, and interactions within, larger communities can be accurately predicted from pairwise interactions alone using variations of said models (Friedman et al., 2017; Meroz et al., 2021; Guo and Boedicker, 2016; Venturelli et al., 2018). This being the case, how microbial interactions combine into the joint effect of multiple species on a single species of interest is still poorly understood. In our research, we used high-throughput nanodroplet-based microfluidics to measure over 14,000 bacterial communities composed of subsets of a library of 61 soil and leaf isolates of which six were fluorescently labeled (Figure 1). We quantified the effect of individual species and the joint effects of species pairs and trios on the growth of six focal bacterial species and found that the effects of multiple species are dominated by the strongest single-species effect, and specifically that negative effects combine non-additively. Figure 1 with 5 supplements see all Download asset Open asset Measuring effects of 61 affecting species, and their pairs and trios on 6 focal species. (A) A library of 61 soil and leaf-associated bacterial strains was used in this experiment. All strains are from four orders: Proteobacteria (orange), Firmicutes (blue), Bacteroidetes (purple), and Actinobacteria (green) (full list in Supplementary file 1a, Source data 1). Also, 6 of the 61 species were labeled with GFP and used as ‘focal’ species whose growth was tested in the presence of the other isolates (affecting species). These strains are labeled on the phylogenetic tree (Escherichia coli [EC], Ewingella americana [EA], Raoultella planticola [RP], Buttiauxella izardii [BI], Citrobacter freundii [CF], and Pantoea agglomerans [PA].) (B) Each focal species was grown in monoculture, with (between 18 and 52) single affecting species, and (between 153 and 1464) pairs of affecting species. Additionally, E. coli was grown with 3009 trios of affecting species. (C) Effects of pairs and trios were then predicted using the effects of single species and single species and pairs, respectively. Predictions were made using three different models: additive, mean, and strongest (detailed in ‘Results’ and ‘Materials and methods’). (D) Equation used for calculating the effect of an affecting species on the focal species. Results We conducted high-throughput assays involving 61 affecting species and 6 different focal species to understand the effects of single species, pairs, and trios on the growth of a given (focal) species. The 61 affecting species included soil and leaf isolates as well as lab strains representing 19 genera from four phyla: Proteobacteria (n = 14), Firmicutes (n = 2), Bacteroidetes (n = 2), and Actinobacteria (n = 1) (full list in Supplementary file 1, Source data 1). The focal species were a subset of six of these species (all proteobacteria) that were transformed to constitutively express a fluorescent protein: (Escherichia coli [EC], Ewingella americana [EA], Raoultella planticola [RP], Buttiauxella izardii [BI], Citrobacter freundii [CF], and Pantoea agglomerans [PA]) (see ‘Materials and methods’). Except for EC, which is a lab strain (E. coli K-12 substr. MG1655), all focals were isolated from soil samples (Kehe et al., 2021). First, we characterized each species phylogenetically by performing Sanger sequencing of their 16S ribosomal RNA gene and phenotypically by growing each species on each of 20 different carbon sources and 11 antibiotics. The species showed large variability in carbon utilization profiles with no species growing well on all carbon sources (Figure 1—figure supplement 3). There was also high variability in growth on antibiotics with 15 species showing little or no growth on any antibiotics, while 16 species showed resistance to at least seven antibiotics (Figure 1—figure supplement 4). We performed the interaction assays in the kChip microfluidics device (Kulesa et al., 2018; Kehe et al., 2019), allowing for extensive screening in parallel (see ‘Materials and methods,’ Figure 1—figure supplement 1). We measured the effects of 243 single species (18–52 for each of the six focal species), the joint effects of 5357 species pairs (between 153 and 1464 for each of the six focal species), and the joint effects of 3009 species trios (from a subset of 26 affecting species on one focal species). Cultures were normalized and mixed after pre-growth, such that the starting densities in the kChip were approximately 1:1 for all species in wells containing two droplets and two affecting species, but ratios varied in three droplet wells (see ‘Materials and methods,’ Figure 1—figure supplements 1 and 2). Interaction assays were carried out in minimal M9 media with 0.5% [w/v] glucose for 24 hr. The growth of the focal species was measured by fluorescence, and effects were calculated as the log ratio of growth in coculture to growth in monoculture (see ‘Materials and methods,’ Figure 1D). Positive and negative effects are defined as a net increase or decrease in growth compared to the monoculture respectively, while affecting species with no observable effect (see ‘Materials and methods’) were defined as neutral. Joint effects of species pairs tend to be stronger than those of individual affecting species We started our interaction assays by measuring the individual effects of single affecting species on each of the focal species (see ‘Materials and methods’). Individual effects covered a wide range (median = −0.15, interquartile range = 0.94) (Figure 2A), and positive effects (the focal species reaching a higher density in the presence of an affecting species than in monoculture) were common overall (32.9%, Figure 2B), in line with previous studies (Kehe et al., 2021). The distribution of effects varied based on the focal species, with E. coli and B. izardii showing the most negative (–2.83) and positive (+0.43) median effects, respectively (Figure 2D). Additionally, we found no affecting species that had strong effects across all focal species (Figure 2—figure supplement 1). Figure 2 with 1 supplement see all Download asset Open asset Pairs of affecting species have stronger effects than single species. (A) Distribution of the effects of single and pairs of affecting species on all focal species. Mann–Whitney–Wilcoxon test two-sided, p-value = 1e-9. Dots show individual effects, solid lines represent the median, boxes represent the interquartile range, and whiskers are expanded to include values no further than 1.5× interquartile range. (B, C) Distribution of qualitative effects of single and pairs of affecting species respectively on all focal species. (D) Distribution of the effect of single and pairs of affecting species for each focal species individually. Dots represent individual measurements, solid lines represent the median, boxes represent the interquartile range, and whiskers are expanded to include values no further than 1.5× interquartile range. Mann–Whitney–Wilcoxon test two-sided tests were performed for each focal species, and p-values are shown on the graph. The measured traits of individual species showed no consistent correlations with their effects on the focal species. In particular, the similarity of metabolic profile, resistance profile, or phylogeny between the focal and affecting species did not correlate strongly with the effect across focal species. Some traits showed little to no correlation for most focals (e.g., antibiotic resistance), while other traits were correlated with effect for a number focal species but not all (e.g., phylogenetic distance). Most of these correlations were not statistically significant (Figure 3—figure supplement 1A). After characterizing the individual effects of all single species, we assayed each pair of affecting species against the focal species. Overall, negative effects were significantly more prevalent in joint pair effects (77.1%) than in effects of individual species (60.5%) (p=1.3e-9, Fisher’s exact test) (Figure 2B and C). The median effect on each focal was more negative by 0.28 on average, though the difference was not significant in all cases; additionally, focals with mostly positive single-species interactions showed a small increase in median effect (Figure 2D). Despite this, the minimum and maximum effects for each focal species remained similar. As with the single affecting species, pairs’ joint effects did not correlate well with species traits, with similarity between the two affecting species, or with their similarity to the focal species (Figure 3—figure supplement 1B). These results indicate that it may be challenging to connect the effects of single and pairs of species on a focal strain to a specific trait of the involved strains using simple analysis. Negative effects combine non-additively and joint effects are dominated by the stronger single-species effect Next, we examined how the effects of individual species relate to their joint effect. In particular, we were interested in finding a model that describes the effects of pairs based on the data from single-species effects. Based on previous studies' success in predicting community structure from pairwise interactions (Friedman et al., 2017; Meroz et al., 2021; Guo and Boedicker, 2016; Venturelli et al., 2018), we posited that predicting how effects combine based solely on the effects of the single species should also be feasible. To do so, we considered three models: an additive effect model, a mean effect model, and a strongest effect model. The additive effect model proposes that the effects of each species on the focal are the same whether they act individually or jointly. Therefore, the combined effect will be equal to the sum of the effects of each species on their own. This is equivalent to additivity of effects between antibiotics, which is common in the drug combinations field (Bollenbach, 2015). The mean model represents a simple phenomenological model that assumes that the effects of different species will be diluted in the presence of a third species. By contrast, the strongest effect model posits that the species with the strongest effect dominates the effect of other affecting species, leaving the joint effect the maximum single effect, and not the sum or mean of single-species effects. When measured across all species and interaction types, we found that the model that best agrees with the measured effects is the strongest effect model (Figure 3B).Though supplementing the mean model with additional species information (i.e., carrying capacity) did improve the model accuracy, it was still less accurate than the strongest effect model (Figure 3—figure supplement 2). The accuracy of the models and identity of the best-fitting model varied across interaction types. The strongest effect model was the most accurate overall (nRMSE = 0.46, 0.32, and 0.16 for the additive, mean, and strongest models, correspondingly), and considerably more accurate when both species affected the focal negatively (nRMSE = 0.65, 0.25, and 0.16). The additive model was slightly more accurate when one effect was negative and the other positive (nRMSE = 0.14, 0.43, and 0.16). Overall, predictions when both effects were positive were less accurate, but here too the strongest model gave the most accurate predictions (nRMSE = 0.81, 0.78, and 0.69) (Figure 3C). Figure 3 with 6 supplements see all Download asset Open asset Strongest single-species effect offers the most accurate model for the combined effect of two species. (A) Graphical representation for each model. The additive model assumes that the effects of each species will accumulate, indicating they are acting independently, and are unaffected by one another. The mean model assumes the combined effect will be an average of the two single-species effects. The final model, strongest effect, assumes that whichever species had a stronger effect on its own will determine the joint effect when paired with an additional species. The y-axis represents the growth of the focal species in different conditions, and in these examples effects are negative. (B) Comparison of predicted effects and the experimental data, with their respective root mean squared error normalized to the interquartile range of the observed data (nRMSE). nRMSE values are calculated from 1000 bootstrapped datasets and represent the median and interquartile range in parentheses (see ‘Materials and methods’). Each dot represents the joint effect of a pair of affecting species on a focal species. Colors indicate the signs of the measured effects of the individual affecting species. (C) Similar to panel (B), but data is stratified by interaction signs of the individual affecting species. The distribution of errors further supported the strongest effect model (Figure 3B, Figure 3—figure supplement 3B): When both single-species effects were negative, the mean model was prone to underestimating the combined effect due to the reduction of the stronger effect by taking into account the weaker effect; while contrastingly, the additive model overestimated effects due to the addition of the weaker effect to the stronger effect, which was more accurate on its own. We saw the opposite trend when both single-species effects were positive, and no particular trend when there was one positive and negative effect. As with the effects themselves, model accuracy was not strongly correlated with any specific species trait (Figure 3—figure supplement 4). In regard to negative effects, support for the strongest model is also evident in how the difference in size of effect influences the model accuracy (Figure 3—figure supplement 3A). When effects are close to equal, the mean model is fairly accurate while the additive model does particularly poorly as these effects would be calculated as twice the strongest effect. Contrastly, when one effect is much stronger than the other, the additive model is accurate since the addition of the weak effect is negligible, whereas the mean model underestimates the joint effects by taking into account the weaker effect. The strongest effect model is also the most accurate for larger communities With this information in hand, we were interested to see whether the same rules held up for larger communities. To this end, we screened trios of a subset (i.e., 26) of the affecting species against a single focal species (E. coli) and found similar trends to all those seen for pairs of affecting species. Similar to what was observed in the move from single species to pairs, effects were stronger (in this case more negative effects) in trios than in the pairs (Figure 4C). Additionally, as with joint pairs’ effects, the strongest effect model was more accurate than the additive and mean models (nRMSE = 2.65, 1.23, and 0.63 for the additive, mean, and strongest models, correspondingly), which is consistent with the fact that the single-species effects in this subset were predominantly negative. Similar distributions of error were seen as in the pairs’ effects, but further exaggerated with the more extreme under and overestimation of the combined trios’ effects by the mean and additive models, respectively (Figure 4A). Figure 4 with 1 supplement see all Download asset Open asset The strongest effect model is also the most accurate for trios. (A, B) Correlation between three different models for how (A) single-species effects and (B) pairwise species effects combine into trio effects, and the experimental data. Root squared mean error normalized (nRMSE) to the interquartile range. nRMSE values are calculated from 1000 datasets and represent the median and interquartile range in parentheses (see ‘Materials and methods’). (C) Distribution of the effects of single, pairs, and trios of affecting species on E. coli. All Mann–Whitney–Wilcoxon two-sided tests were significant, p values are shown on plot. Dots show individual effects, solid lines represent the median, boxes represent the interquartile range, and whiskers are expanded to include values no further than 1.5× interquartile range. (D, E) Distribution of errors for each model based on (D) single-species data and (E) pairs data. Dots show individual effects, solid lines represent the median, boxes represent the interquartile range, and whiskers are expanded to include values no further than 1.5× interquartile range. We further explored the additive, mean, and strongest models in trios by basing the model on effects of the three pairs comprising each trio (i.e., joint effect of AB, AC, BC to predict effect of ABC), as opposed to only using single-effect data (i.e. ,effect of A, B, C on their own) (Figure 4A and B). The effects of single species and pairs were measured again independently in this experiment (see ‘Materials and methods,’ Figure 1—figure supplement 1). The move to pairs-based predictions improved the accuracy for both the mean and strongest model, while further pushing the additive predictions away from the observed effects (Figure 4B). These data suggest that even in the presence of additional species the strongest single-species effect still dominates the combined effect of a community. Discussion By measuring thousands of simplified microbial communities, we quantified the effects of single species, pairs, and trios on multiple focal species. The most accurate model, overall and specifically when both single-species effects were negative, was the strongest effect model. This is in stark contrast to models often used in antibiotic compound combinations, despite most effects being negative, where additivity is often the default model (Bollenbach, 2015). The additive model performed well for mixed effects (i.e., one negative and one positive), but only slightly better than the strongest model, and poorly when both species had effects of the same sign. When both single-species effects were positive, the strongest model was also the best, though the difference was less pronounced and all models performed worse for these interactions. This may be due to the small effect size seen with positive effects, as when we limited negative and mixed effects to a similar range of effects strength, their accuracy dropped to similar values (Figure 3—figure supplement 5). We posit that the difference in accuracy across species is affected mainly by the effect type dominating different focal species’ interactions, rather than by inherent species traits (Figure 3—figure supplement 6). We phenotypically and genetically profiled all species, but did not find strong correlations between the measured traits, or similarity in traits, to the effect on the focal species. Though positive effects were common, making up about one-third of the single-species effects, they became less common as the number of community members increased, making up only 16% of the effects of species pairs. Furthermore, we found similar trends in the larger communities of four species (three affecting species and one focal), both that effects combined in a non-additive manner, being dominated by the strongest single-species effect, and that effects became stronger in larger communities, which is consistent with previous studies (Cook et al., 2006; van Elsas et al., 2012 ; Jones et al., 2021; Piccardi et al., 2019; Gould et al., 2018; Palmer and Foster, 2022). The mechanistic basis underlying the joint effects of multiple species is still unclear. The additive model’s accuracy for mixed effects may indicate that negative and positive effects act independently. For negative effects, it is difficult to identify a single biological mechanism that could explain why the strongest effect model agreed best with our experimental data. Intuitively, we assumed this could be explained by resource competition (i.e., an affecting species that consumes resources quickly would negatively affect the focal species, as well as the other affecting species). However, this explanation is not consistent with the fact that the affecting species’ growth rate did not correlate well with their effect on some focal species (Figure 3—figure supplement 1A). Secondly, we thought effects could be saturating (either biologically or with regard to the detection limits in this experimental setup), but this would not explain why the model works for weaker effects. A hierarchical ranking, where each species affects all the species ranked below it could lead to the strongest affecting species affecting both the focal and the other affecting species, thus dominating their joint effects, but this does not coincide with the fact that we observed almost no single species or pairs with a strong effect across all focals (Figure 2—figure supplement 1). As we did not measure the abundance of all species in each community (only the focal), we cannot disentangle interaction modification (changes in per capita effect of specific species), from interaction chains (affecting the amount of an affecting species, and as such its effect on the focal), and further work is needed in order to pinpoint the exact mechanism(s) leading to the dominance of the strongest model for negative effects in our system. We also note that it is possible that the manner in which effects combine is affected by the mechanism of interaction; For instance, previous studies have shown that interference competition can combine additively, or even synergistically, results not seen in our work (Tyc et al., 2014; Westhoff et al., 2021). Understanding how microbial communities assemble and how large numbers of species interact is of both utmost importance and difficulty. Harnessing such information would open up a plethora of currently underutilized applications in food, medical, and agricultural industries. Specifically, understanding how the effects of multiple species on a single species combine is important for introduction of new species into a given environment. Our results suggest that when we want to affect a single