(C) PLOS One This story was originally published by PLOS One and is unaltered. . . . . . . . . . . Quantifying microbial interactions based on compositional data using an iterative approach for solving generalized Lotka-Volterra equations [1] ['Yue Huang', 'Department Of Quantitative', 'Computational Biology', 'University Of Southern California', 'Los Angeles', 'California', 'United States Of America', 'Tianqi Tang', 'Department Of Neurology', 'University Of California'] Date: 2026-02 Understanding microbial interactions is fundamental for exploring population dynamics, particularly in microbial communities where interactions affect stability and host health. Generalized Lotka-Volterra (gLV) models have been widely used to investigate system dynamics but depend on absolute abundance data, which are often unavailable in microbiome studies. To address this limitation, we introduce an iterative Lotka-Volterra (iLV) model, a novel framework tailored for compositional data that leverages relative abundances and iterative refinements for parameter estimation. The iLV model features two key innovations: an adaptation of the gLV framework to compositional constraints and an iterative optimization strategy combining linear approximations with nonlinear refinements to enhance parameter estimation accuracy. Using simulations and real-world datasets, we demonstrate that iLV surpasses existing methodologies, such as the compositional LV (cLV) and the generalized LV (gLV) model, in recovering interaction coefficients and predicting species trajectories under varying noise levels and temporal resolutions. Applications to the lynx-hare predator-prey, Stylonychia pustula-P. caudatum mixed culture, and cheese microbial systems revealed consistency between predicted and observed relative abundances showcasing its accuracy and robustness. In summary, the iLV model bridges theoretical gLV models and practical compositional data analysis, offering a robust framework to infer microbial interactions and predict community dynamics using relative abundance data, with significant potential for advancing microbial research. Microbes and animals often live in communities where species interact in complex ways. Understanding these interactions is essential for predicting population dynamics and how ecosystems respond to environmental changes. A widely used mathematical tool for modeling these interactions is the generalized Lotka-Volterra model. However, it requires information on absolute population sizes, which is rarely available in microbiome studies. We developed the iterative Lotka-Volterra (iLV) model to overcome this challenge. iLV adapts the classical framework to work directly with relative abundance data—the kind most often collected in microbiome studies. Through an iterative optimization process, iLV improves the accuracy of inferred species interactions and population trajectories. We validated our approach using simulated data, along with real-world systems including the snowshoe hare–Canadian lynx system, a microbial co-culture experiment, and a cheese microbial community. In each case, iLV provided more accurate results than existing methods. Our work offers a practical framework for studying species interactions using relative abundance data. Funding: This research was supported in part by the US National Science Foundation (EF-2125142 to F.S.) and the US National Institute of Diabetes and Digestive and Kidney Diseases (R01DK142026 to X.D.). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. Copyright: © 2025 Huang et al. This is an open access article distributed under the terms of the Creative Commons Attribution License , which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. These previous advancements have highlighted the potential of compositional data analysis, but they also have some major limitations, such as moderate accuracy due to their nature of linear approximations. For cLV, the coefficients in the gLV model cannot be fully recovered [ 13 ], making it harder to interpret the results. Gaps remain in integrating theoretical and empirical approaches for compositional data. To bridge the gap, we proposed a novel method (called iterative Lotka-Volterra or iLV) of compositional data analysis based on the traditional gLV model, emphasizing the compositional constraints and their implications for model parameterization. There are two major innovative components of iLV, including defining the classical generalized Lotka-Volterra model with relative abundances and inter-species sum of absolute abundances, and iterative linear approximations followed by non-linear optimizations. By leveraging computational simulations, we compared the performance of iLV with other existing methods including generalized Lotka–Volterra (gLV) and compositional Lotka-Volterra (cLV). We also applied iLV in three real-world datasets to estimate the coefficients of gLV model and accurately recover community relative abundance trajectories. However, the widespread use of sequencing technologies in microbial studies introduces a unique challenge: most data are compositional, providing relative rather than absolute abundance profiles [ 10 , 11 ]. This poses significant limitations for traditional gLV models, which rely on absolute abundance inputs for accurate parameterization and prediction [ 12 , 13 ]. To address this problem, compositional adaptations of gLV models have been developed, such as the compositional Lotka-Volterra (cLV) framework, which maps dynamics onto a constrained simplex, accommodating the summation constraint of relative abundances [ 13 ]. Quantifying microbial interactions is essential for understanding the dynamics of microbial communities. This task has gained increased attention in microbiome research, as microbial interactions directly impact community stability and host health [ 1 – 3 ]. The Lotka-Volterra (LV) modeling framework was originally developed to describe predator-prey systems [ 4 ], based on which the generalized Lotka-Volterra (gLV) was proposed to model systems of more than two species. gLV has been widely adopted to infer these interactions due to its versatility in modeling nonlinear dynamics across diverse ecosystems, from microbial communities to macroscopic populations [ 5 – 7 ]. The generalized Lotka-Volterra (gLV) model remains a key tool for describing microbial interactions and forecasting community trajectories. Applications range from predicting microbiome responses to antibiotics and dietary shifts to understanding microbial coexistence in ecological contexts [ 8 , 9 ]. Results Performance and stability of different non-linear optimization methods When applying the iLV algorithm to parameter estimation, the choice of non-linear optimization methods in Subroutine 2 can influence both accuracy and numerical stability, particularly when the problem is ill-conditioned. Fig 1 presents an example where instability occurred for different non-linear optimization methods. We simulated data points using the following parameter setting that reflects periodic oscillations among three species: , , , , , , , , , , , , , with a time length of 20 and a time interval of 0.4. We executed the iLV Algorithm 20 times using the simulated data and estimated the parameters with three different non-linear optimization methods: leastsq(), least_squares(method = ‘lm’), and least_squares(method=’trf’). was set to 200. As shown in Fig 1, all three methods exhibited instabilities potentially due to rounding errors and the ill-conditioned nature of this dataset, and leastsq() achieved the lowest trajectory RMSE of . Therefore, to mitigate the effects of numerical instabilities and the variable performance of optimization methods across different datasets, we compared the trajectory RMSE returned by least_squares (method = “trf” or method = “lm”) and leastsq() in Subroutine 2, and retained the results with the lowest RMSE. Additionally, during benchmarking, we repeated the algorithm 20 times and reported the parameter that yielded the lowest RMSE among these 20 runs. PPT PowerPoint slide PNG larger image TIFF original image Download: Fig 1. Trajectory Root Mean Square Error (RMSE) distributions of three optimization methods: (A) SQ refers to leastsq(), (B) LM refers to least_squares(method = ‘lm’), and (C) TRF refers to least_squares(method=’trf’). Data points were simulated using the following parameter setting that reflects periodic oscillations among three species: , , , , , , , , , , , , , with a time length of 20 and a time interval of 0.4. All three methods exhibited instabilities for this dataset, while leastsq() achieved the lowest trajectory RMSE of . The mean RMSEs across 20 runs for leastsq(), least_squares(method = ’lm’), and least_squares(method = ’trf’) were 0.0687, 0.135, and 0.103, respectively, while their corresponding median RMSEs were 0.0386, 0.0958, and 0.0945. https://doi.org/10.1371/journal.pcbi.1013691.g001 Two subroutines of the iLV Algorithm jointly improve parameter estimation An accurate initial guess of parameters is crucial for the success of optimization functions like leastsq() or least_squares() in finding an optimal local minimum of the cost function. Subroutine 1 (the iterative subroutine) of the iLV Algorithm provides an iterative approach to generate such initial guesses effectively. To illustrate the impact of the iterative subroutine, we generated simulated data using the following parameter settings: , , , , , , , , , , , , , time range = 20, time interval = 1. was set to 200. Fig 2A shows the trajectory RMSE as a function of the number of iterations in Subroutine 1 (the iterative subroutine). The figure indicates that Subroutine 1 iteratively refined the initial guess, with the optimal guess achieved in the 13th iteration out of 100. Fig 2B and 2C illustrate how Subroutine 1 substantially improved the fit between predicted and observed relative abundances. Without applying Subroutine 1, using the parameter estimations derived from the gLV linear approximation as the starting point of Subroutine 2 (the least square estimation subroutine) resulted in poor optimization performance, with an RMSE of 0.122 for leastsq() (Fig 2C). In contrast, using the parameters from the iteration with the lowest trajectory RMSE within the first 100 iterations significantly improved optimization performance, yielding an RMSE of with leastsq() (Fig 2B). Under this parameter setting, least_squares(method = ‘lm’ or ‘trf’) did not perform as well as leastsq(). Notably, the iterative subroutine guarantees non-increasing trajectory RMSE values, as the first iteration corresponds to the result without applying the iterative subroutine, and only the iteration with the lowest trajectory RMSE (e.g., out of the first 100 iterations) is selected as the starting point of the least square estimation subroutine. PPT PowerPoint slide PNG larger image TIFF original image Download: Fig 2. Two subroutines of the iLV Algorithm jointly improve parameter estimation. Panel A depicts the trajectory RMSE values across 100 iterations of the iterative subroutine in the iLV Algorithm, illustrating its ability to iteratively identify parameter estimations with the lowest RMSE within the first 100 runs. The last three panels compare the predicted relative abundances (solid lines) to observed values (symbols) for each species, with both Subroutines (Panel B), without Subroutine 1 (Panel C), and without Subroutine 2 (Panel D), respectively, demonstrating significantly improved alignment after joint effects of both Subroutines in the iLV Algorithm. https://doi.org/10.1371/journal.pcbi.1013691.g002 Subroutine 2 helps to improve parameter estimation as leastsq() and least_squares() functions find a local minimum of the cost function F(x), near the initial guess returned by Subroutine 1. As shown in Fig 2D, the trajectory RMSE was 0.169 when only Subroutine 1 was applied without Subroutine 2. A cheese microbial community Mounier et al. previously quantified cell counts of five microbial groups within a cheese microbial community over a 21-day period [15]. We converted these absolute counts to relative abundances for modeling purposes. These groups include D. hansenii, Y. lipolytica, G. candidum, Leucobacter sp., and a bacterial group composed of Arthrobacter arilaitensis, Hafnia alvei, Corynebacterium casei, Brevibacterium aurantiacum, and Staphylococcus xylosus. The system plateaued around day 10, and we modeled the dynamic phase between days 0 and 10. We applied our iLV model to this dataset with set to 0.001 (see Table 7), and the trajectory RMSE is 0.149. For comparison, we also applied the gLV_relative and gLV_absolute models, which yielded trajectory RMSEs of 0.327 and 0.296, respectively. Significant interactions among the microbial species, as predicted by the iLV model, are shown in Fig 5D. Some of these interactions are supported by previous studies [15]. For example, Y. lipolytica has been shown to reduce the viability of D. hansenii during the stationary phase [15]. However, not all predicted interactions align with prior findings. For instance, Y. lipolytica was reported to inhibit the mycelial growth of G. candidum, rather than promoting it. This inhibition disrupted the typical mold-like mycelial structure of G. candidum, leading to the formation of spaghetti-like morphologies [15]. It is important to note, however, that Mounier et al. investigated pairwise interactions through two-species co-culture experiments [15], which may not fully capture the complexity of interactions within a multispecies community. PPT PowerPoint slide PNG larger image TIFF original image Download: Fig 5. Predicted abundance trajectory and interactions of the cheese microbial community. We set according to Table 7. Panels A, B, and C show the trajectories generated by iLV, gLV_relative, and gLV_absolute models, respectively, with observed data overlaid for comparison. iLV aligns closely with observed trajectories (RMSE = 0.149), whereas gLV_relative (RMSE = 0.327) and gLV_absolute (RMSE = 0.296) display higher errors. Panel D shows the predicted interaction network. The green arrows mean promotion, and the red arrows mean inhibition. The thickness of the arrows is proportional to the effect size of promotion or inhibition, indicated by the magnitude of b ij . Only interactions that are at least 50 times larger than the minimum predicted interaction were plotted. Dh is D. hansenii, Yl is Y. lipolytica, Gc is G. candidum, Ls is Leucobacter sp., and C is a bacterial group composed of Arthrobacter arilaitensis, Hafnia alvei, Corynebacterium casei, Brevibacterium aurantiacum, and Staphylococcus xylosus. https://doi.org/10.1371/journal.pcbi.1013691.g005 [END] --- [1] Url: https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1013691 Published and (C) by PLOS One Content appears here under this condition or license: Creative Commons - Attribution BY 4.0. via Magical.Fish Gopher News Feeds: gopher://magical.fish/1/feeds/news/plosone/