# Integrating microbial abundance time series with fermentation dynamics of the rumen microbiome via mathematical modelling

Mohsen Davoudkhani, Francesco Rubino, Christopher J Creevey, Seppo Ahvenjärvi, Ali R Bayat, Ilma Tapio, Alejandro Belanche, Rafael Muñoz-Tamayo

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Mohsen Davoudkhani, Francesco Rubino, Christopher J Creevey, Seppo Ahvenjärvi, Ali R Bayat, et al.. Integrating microbial abundance time series with fermentation dynamics of the rumen microbiome via mathematical modelling. PLoS ONE, 2024, 19 (3), pp.e0298930. 10.1371/journal.pone.0298930. hal-04558514

HAL Id: hal-04558514

https://hal.inrae.fr/hal-04558514

Submitted on 25 Apr 2024

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
RESEARCH ARTICLE

# Integrating microbial abundance time series with fermentation dynamics of the rumen microbiome *via* mathematical modelling

Mohsen Davoudkhani<sup>1*</sup>, Francesco Rubino<sup>2*</sup>, Christopher J. Creevey<sup>2</sup>, Seppo Ahvenjärvi<sup>3</sup>, Ali R. Bayat<sup>3</sup>, Ilma Tapio<sup>4</sup>, Alejandro Belanche<sup>1,5</sup>, Rafael Muñoz-Tamayo<sup>1*</sup>

1 INRAE, AgroParisTech, UMR Modélisation Systémique Appliquée aux Ruminants, Université Paris-Saclay, Palaiseau, France, 2 Institute of Global Food Security, School of Biological Sciences, Queen's University Belfast, Northern Ireland, United Kingdom, 3 Animal Nutrition, Production Systems, Natural Resources Institute Finland (Luke), Jokioinen, Finland, 4 Genomics and Breeding, Production Systems, Natural Resources Institute Finland (Luke), Jokioinen, Finland, 5 Departamento de Producción Animal y Ciencia de los Alimentos, Universidad de Zaragoza, Zaragoza, Spain

© These authors contributed equally to this work.

* Rafael.munoz-tamayo@inrae.fr

<figure><img src="image_1.png" /></figure>

<figure><img src="image_2.png" /></figure>

## OPEN ACCESS

**Citation:** Davoudkhani M, Rubino F, Creevey CJ, Ahvenjärvi S, Bayat AR, Tapio I, et al. (2024) Integrating microbial abundance time series with fermentation dynamics of the rumen microbiome *via* mathematical modelling. PLoS ONE 19(3): e0298930. https://doi.org/10.1371/journal.pone.0298930

**Editor:** Adham A. Al-Sagheer, Zagazig University Faculty of Agriculture, EGYPT

**Received:** October 6, 2023

**Accepted:** February 2, 2024

**Published:** March 20, 2024

**Peer Review History:** PLOS recognizes the benefits of transparency in the peer review process; therefore, we enable the publication of all of the content of peer review and author responses alongside final, published articles. The editorial history of this article is available here: https://doi.org/10.1371/journal.pone.0298930

**Copyright:** © 2024 Davoudkhani et al. This is an open access article distributed under the terms of the [Creative Commons Attribution License](https://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** The scripts for the inference of the functional modules are available at https://github.com/frubino/cowpi and https://doi.org/10.1371/journal.pone.0298930

## Abstract

The rumen represents a dynamic microbial ecosystem where fermentation metabolites and microbial concentrations change over time in response to dietary changes. The integration of microbial genomic knowledge and dynamic modelling can enhance our system-level understanding of rumen ecosystem's function. However, such an integration between dynamic models and rumen microbiota data is lacking. The objective of this work was to integrate rumen microbiota time series determined by 16S rRNA gene amplicon sequencing into a dynamic modelling framework to link microbial data to the dynamics of the volatile fatty acids (VFA) production during fermentation. For that, we used the theory of state observers to develop a model that estimates the dynamics of VFA from the data of microbial functional proxies associated with the specific production of each VFA. We determined the microbial proxies using CowPi to infer the functional potential of the rumen microbiota and extrapolate their functional modules from KEGG (Kyoto Encyclopedia of Genes and Genomes). The approach was challenged using data from an *in vitro* RUSITEC experiment and from an *in vivo* experiment with four cows. The model performance was evaluated by the coefficient of variation of the root mean square error (CRMSE). For the *in vitro* case study, the mean CVRMSE were 9.8% for acetate, 14% for butyrate and 14.5% for propionate. For the *in vivo* case study, the mean CVRMSE were 16.4% for acetate, 15.8% for butyrate and 19.8% for propionate. The mean CVRMSE for the VFA molar fractions were 3.1% for acetate, 3.8% for butyrate and 8.9% for propionate. Ours results show the promising application of state observers integrated with microbiota time series data for predicting rumen microbial metabolism.
RESEARCH ARTICLE

# Integrating microbial abundance time series with fermentation dynamics of the rumen microbiome *via* mathematical modelling

Mohsen Davoudkhani<sup>1*</sup>, Francesco Rubino<sup>2*</sup>, Christopher J. Creevey<sup>2</sup>, Seppo Ahvenjärvi<sup>3</sup>, Ali R. Bayat<sup>3</sup>, Ilma Tapio<sup>4</sup>, Alejandro Belanche<sup>1,5</sup>, Rafael Muñoz-Tamayo<sup>1*</sup>

1 INRAE, AgroParisTech, UMR Modélisation Systémique Appliquée aux Ruminants, Université Paris-Saclay, Palaiseau, France, 2 Institute of Global Food Security, School of Biological Sciences, Queen's University Belfast, Northern Ireland, United Kingdom, 3 Animal Nutrition, Production Systems, Natural Resources Institute Finland (Luke), Jokioinen, Finland, 4 Genomics and Breeding, Production Systems, Natural Resources Institute Finland (Luke), Jokioinen, Finland, 5 Departamento de Producción Animal y Ciencia de los Alimentos, Universidad de Zaragoza, Zaragoza, Spain

© These authors contributed equally to this work.

* Rafael.munoz-tamayo@inrae.fr

<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /></figure>

OPEN ACCESS

**Citation:** Davoudkhani M, Rubino F, Creevey CJ, Ahvenjärvi S, Bayat AR, Tapio I, et al. (2024) Integrating microbial abundance time series with fermentation dynamics of the rumen microbiome *via* mathematical modelling. PLoS ONE 19(3): e0298930. https://doi.org/10.1371/journal.pone.0298930

**Editor:** Adham A. Al-Sagheer, Zagazig University Faculty of Agriculture, EGYPT

**Received:** October 6, 2023

**Accepted:** February 2, 2024

**Published:** March 20, 2024

**Peer Review History:** PLOS recognizes the benefits of transparency in the peer review process; therefore, we enable the publication of all of the content of peer review and author responses alongside final, published articles. The editorial history of this article is available here: https://doi.org/10.1371/journal.pone.0298930

**Copyright:** © 2024 Davoudkhani et al. This is an open access article distributed under the terms of the [Creative Commons Attribution License](https://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** The scripts for the inference of the functional modules are available at https://github.com/frubino/cowpi and https://doi.org/10.1371/journal.pone.0298930

## Abstract

The rumen represents a dynamic microbial ecosystem where fermentation metabolites and microbial concentrations change over time in response to dietary changes. The integration of microbial genomic knowledge and dynamic modelling can enhance our system-level understanding of rumen ecosystem's function. However, such an integration between dynamic models and rumen microbiota data is lacking. The objective of this work was to integrate rumen microbiota time series determined by 16S rRNA gene amplicon sequencing into a dynamic modelling framework to link microbial data to the dynamics of the volatile fatty acids (VFA) production during fermentation. For that, we used the theory of state observers to develop a model that estimates the dynamics of VFA from the data of microbial functional proxies associated with the specific production of each VFA. We determined the microbial proxies using CowPi to infer the functional potential of the rumen microbiota and extrapolate their functional modules from KEGG (Kyoto Encyclopedia of Genes and Genomes). The approach was challenged using data from an *in vitro* RUSITEC experiment and from an *in vivo* experiment with four cows. The model performance was evaluated by the coefficient of variation of the root mean square error (CRMSE). For the *in vitro* case study, the mean CVRMSE were 9.8% for acetate, 14% for butyrate and 14.5% for propionate. For the *in vivo* case study, the mean CVRMSE were 16.4% for acetate, 15.8% for butyrate and 19.8% for propionate. The mean CVRMSE for the VFA molar fractions were 3.1% for acetate, 3.8% for butyrate and 8.9% for propionate. Ours results show the promising application of state observers integrated with microbiota time series data for predicting rumen microbial metabolism.
org/10.5281/zenodo.8401851. The data and scripts of the implementation of the observer for each case study are available at https://doi.org/10.5281/zenodo.8386786 [48]. Sequencing data for the in vitro case study are accessible at the EBI Short Read Archive from the European Nucleotide Archive (accession number PRJEB20255). The sequencing data for the in vivo case study are accessible at NCBI SRA under the BioProject PRJNA1023082.

**Funding:** All authors receive funding from the MASTER project, an Innovation Action funded by the European Union's Horizon 2020 research and innovation programme under grant agreement No 818368. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

**Competing interests:** The authors have declared that no competing interests exist.

## Introduction

The function of the rumen microbiota affects animal production phenotypes, feed efficiency and methane emissions [1]. Our knowledge on the structure and function of the rumen microbiota has been greatly improved due to the progress on culture independent omic techniques [2, 3]. A powerful use of these techniques is the analysis of microbial time series data allowing to characterise the dynamics of the rumen microbial ecosystem. Applications include the study of microbial colonization of feed particles [4–6], microbial resilience in response to perturbations [7], impact of acidotic challenge on microbial function [8], activity of methanogens in response to the supplementation of a methane inhibitor [9] and evolution of gut microbiota through weaning transition [10, 11].

In parallel to the use of omics techniques to characterise rumen microbiota patterns, dynamic models have been developed to represent the rumen fermentation profile under *in vitro* [12, 13] and *in vivo* [14–16] conditions. Existing models of rumen function consider an aggregated representation of the rumen microbiota and its metabolic function. However, none of these models integrate microbial genomic knowledge and thus do not capitalize on the rich information that microbial genomic sequencing provides. Integration of dynamic modelling and microbial data has the potential to improve the understanding of the rumen ecosystem, to enhance predictive power of rumen models and to help the design of microbial manipulation strategies to improve rumen function [17]. Recently, some studies have applied the genome-scale metabolic approach to reconstruct metabolic networks of rumen microbes species [18–20] and to predict the metabolism of minimal rumen microbial consortium [21]. Another model approach consists in exploiting microbial time series data, with a variety of dedicated mathematical approaches [22] including the generalized Lotka-Volterra (gLV) model [23, 24]. In its standard form, the gLV approach determines interactions between microbes but it does not provide either mechanistic insights or predictions on the dynamics of microbial metabolism. An attempt to couple the gLV approach with a kinetic metabolic model that integrates both microbial interactions and metabolism in a nitrification reactor was developed in [25]. The model was useful to identify qualitative aspects of the system such as the coexistence of two competing bacteria. However, this type of modelling approach based on pairwise microbial interactions is limited to ecosystems with few species due to the high number of parameters that need to be estimated.

In this work, we aim at integrating microbial time series within a dynamic modelling framework *via* the application of state observers (also called software sensors). An observer is an algorithm that combines measurements and a mathematical model to estimate unmeasured variables (see, e.g., [26] for a review). Within the classes of state observers, asymptotic observers are of particular interest since they do not require knowledge on the kinetic functions representing the reaction rates. Asymptotic observers have been applied in simple microbial processes with few microbial species or in processes where the microbiota is represented in aggregated fashion by few microbial functional groups [27, 28]. Observers in combination with optimization routines have been used to assign the metabolic functions of species (Operational Taxonomic Units—OTUs) participating in a nitrification process [29, 30]. In these previous works, the functional assignment of species was translated into an optimization problem which assumes that one species participates only in one reaction. Such a hypothesis might hold for the nitrification process. However, the assumption of single microbial function does not seem to apply to the rumen microbiota, since microbes participate in various metabolic pathways simultaneously. To circumvent the single functional constraint, here we propose a novel approach that uses OTU data to derive microbial functional proxies for specific process of rumen metabolism such as volatile fatty acid (VFA) production. Our approach was tested
org/10.5281/zenodo.8401851. The data and scripts of the implementation of the observer for each case study are available at https://doi.org/10.5281/zenodo.8386786 [48]. Sequencing data for the in vitro case study are accessible at the EBI Short Read Archive from the European Nucleotide Archive (accession number PRJEB20255). The sequencing data for the in vivo case study are accessible at NCBI SRA under the BioProject PRJNA1023082.

**Funding:** All authors receive funding from the MASTER project, an Innovation Action funded by the European Union's Horizon 2020 research and innovation programme under grant agreement No 818368. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

**Competing interests:** The authors have declared that no competing interests exist.

## Introduction

The function of the rumen microbiota affects animal production phenotypes, feed efficiency and methane emissions [1]. Our knowledge on the structure and function of the rumen microbiota has been greatly improved due to the progress on culture independent omic techniques [2, 3]. A powerful use of these techniques is the analysis of microbial time series data allowing to characterise the dynamics of the rumen microbial ecosystem. Applications include the study of microbial colonization of feed particles [4–6], microbial resilience in response to perturbations [7], impact of acidotic challenge on microbial function [8], activity of methanogens in response to the supplementation of a methane inhibitor [9] and evolution of gut microbiota through weaning transition [10, 11].

In parallel to the use of omics techniques to characterise rumen microbiota patterns, dynamic models have been developed to represent the rumen fermentation profile under *in vitro* [12, 13] and *in vivo* [14–16] conditions. Existing models of rumen function consider an aggregated representation of the rumen microbiota and its metabolic function. However, none of these models integrate microbial genomic knowledge and thus do not capitalize on the rich information that microbial genomic sequencing provides. Integration of dynamic modelling and microbial data has the potential to improve the understanding of the rumen ecosystem, to enhance predictive power of rumen models and to help the design of microbial manipulation strategies to improve rumen function [17]. Recently, some studies have applied the genome-scale metabolic approach to reconstruct metabolic networks of rumen microbes species [18–20] and to predict the metabolism of minimal rumen microbial consortium [21]. Another model approach consists in exploiting microbial time series data, with a variety of dedicated mathematical approaches [22] including the generalized Lotka-Volterra (gLV) model [23, 24]. In its standard form, the gLV approach determines interactions between microbes but it does not provide either mechanistic insights or predictions on the dynamics of microbial metabolism. An attempt to couple the gLV approach with a kinetic metabolic model that integrates both microbial interactions and metabolism in a nitrification reactor was developed in [25]. The model was useful to identify qualitative aspects of the system such as the coexistence of two competing bacteria. However, this type of modelling approach based on pairwise microbial interactions is limited to ecosystems with few species due to the high number of parameters that need to be estimated.

In this work, we aim at integrating microbial time series within a dynamic modelling framework *via* the application of state observers (also called software sensors). An observer is an algorithm that combines measurements and a mathematical model to estimate unmeasured variables (see, e.g., [26] for a review). Within the classes of state observers, asymptotic observers are of particular interest since they do not require knowledge on the kinetic functions representing the reaction rates. Asymptotic observers have been applied in simple microbial processes with few microbial species or in processes where the microbiota is represented in aggregated fashion by few microbial functional groups [27, 28]. Observers in combination with optimization routines have been used to assign the metabolic functions of species (Operational Taxonomic Units—OTUs) participating in a nitrification process [29, 30]. In these previous works, the functional assignment of species was translated into an optimization problem which assumes that one species participates only in one reaction. Such a hypothesis might hold for the nitrification process. However, the assumption of single microbial function does not seem to apply to the rumen microbiota, since microbes participate in various metabolic pathways simultaneously. To circumvent the single functional constraint, here we propose a novel approach that uses OTU data to derive microbial functional proxies for specific process of rumen metabolism such as volatile fatty acid (VFA) production. Our approach was tested
firstly using published experimental data from an *in vitro* study [5]. We presented preliminary results of this approach at the 10th Workshop on Modelling Nutrient Digestion and Utilization in Farm Animals (MODNUT, 18th–21st September 2022, Sardinia, Italy). The corresponding abstract is published in the conference proceedings [31]. We performed further an experiment with four cows to assess the modelling approach under *in vivo* conditions. Following open science practices [32], the data and scripts are freely available.

## Material and methods

In this section, we provide a general overview of the experimental case studies. Next, we describe the modelling approach.

### *In vitro* case study

We used data from an *in vitro* experiment carried out in RUSITEC (Rumen Simulation Technique) systems [5] to study the feed microbial colonization dynamics of fresh ryegrass or rye-grass hay. Rumen inoculum was obtained from rumen fistulated cows, and the incubation ran for 18 consecutive days in 16 vessels. Fermentation vessels (800 mL effective volume) were inoculated with rumen fluid and incubated with the experimental diets (80:20 forage to concentrate ratio). One bag containing feed (11.25 g DM) was daily supplied to each vessel and incubated for 48h. Artificial saliva was continuously infused at a dilution rate of 3.35%/h (equivalent to 645 mL/d). After 14 days of adaptation to the experimental diets, the fermentation dynamics was monitored during 3 consecutive days. For each day, samples were taken at 2, 4, 8 and 24 h for determination of acetate, butyrate and propionate and for microbial characterization. Metabolite concentrations for each sampling time was reported as the mean value of the three samples. For microbial characterization, the samples were pooled per time point. The microbial RNA was extracted from the liquid and solid phases and the bacterial community structure was characterized by 16SrRNA (cDNA) Next Generation Sequencing (NGS) as previously described [5]. Raw sequence reads are accessible at the EBI Short Read Archive from the European Nucleotide Archive (accession number PRJEB20255).

### *In vivo* experiment

We conducted an experiment with four Nordic Red dairy cows, selected based on similar calving dates and equipped with rumen fistulas to provide dynamic data to assess our modelling approach *in vivo*. For 10 days of diet adaptation period, cows were offered total mixed ration consisting of grass silage (timothy-meadow fescue sward) preserved with formic acid based additive (AIV Ässä Na; 5 litres/tonne) provided at 45:55 forage to concentrate ratio on a dry matter basis. Concentrate mixture consisted of barley 210, oats 210, wheat 100, sugar beet pulp 220, rapeseed meal 230, and a mixture of minerals and vitamins 30 g/kg on as fed basis. Following the adaptation period, cows were located to metabolic chambers. Day 1 was dedicated to adaptation to the chamber conditions, while on d 2 and d 3 gas exchanges were measured for 48 h (for details see [33]). To determine circadian changes in rumen fermentation and rumen microbial community composition, on chamber d 4 and d 5 for 48 h period rumen liquid samples were collected every 3 h through ruminal fistula from the ventral site of the rumen using 500-mL bottle. Two sub-samples were taken for VFA and ammonia-N determination as described in [11] and were stored at -20°C for later analysis. Rumen liquid samples for microbial analysis were mixed to get homogenous distribution of microorganisms in the liquid, aliquoted into 2-mL sterile screw cap tubes, snap frozen on dry ice and stored at -80°C until DNA extraction. The total DNA was extracted from 500 µL of rumen liquid following protocol described in [34]. Libraries of the bacterial 16S ribosomal RNA (rRNA) V4 region were
firstly using published experimental data from an *in vitro* study [5]. We presented preliminary results of this approach at the 10th Workshop on Modelling Nutrient Digestion and Utilization in Farm Animals (MODNUT, 18th–21st September 2022, Sardinia, Italy). The corresponding abstract is published in the conference proceedings [31]. We performed further an experiment with four cows to assess the modelling approach under *in vivo* conditions. Following open science practices [32], the data and scripts are freely available.

## Material and methods

In this section, we provide a general overview of the experimental case studies. Next, we describe the modelling approach.

### *In vitro* case study

We used data from an *in vitro* experiment carried out in RUSITEC (Rumen Simulation Technique) systems [5] to study the feed microbial colonization dynamics of fresh ryegrass or rye-grass hay. Rumen inoculum was obtained from rumen fistulated cows, and the incubation ran for 18 consecutive days in 16 vessels. Fermentation vessels (800 mL effective volume) were inoculated with rumen fluid and incubated with the experimental diets (80:20 forage to concentrate ratio). One bag containing feed (11.25 g DM) was daily supplied to each vessel and incubated for 48h. Artificial saliva was continuously infused at a dilution rate of 3.35%/h (equivalent to 645 mL/d). After 14 days of adaptation to the experimental diets, the fermentation dynamics was monitored during 3 consecutive days. For each day, samples were taken at 2, 4, 8 and 24 h for determination of acetate, butyrate and propionate and for microbial characterization. Metabolite concentrations for each sampling time was reported as the mean value of the three samples. For microbial characterization, the samples were pooled per time point. The microbial RNA was extracted from the liquid and solid phases and the bacterial community structure was characterized by 16SrRNA (cDNA) Next Generation Sequencing (NGS) as previously described [5]. Raw sequence reads are accessible at the EBI Short Read Archive from the European Nucleotide Archive (accession number PRJEB20255).

### *In vivo* experiment

We conducted an experiment with four Nordic Red dairy cows, selected based on similar calving dates and equipped with rumen fistulas to provide dynamic data to assess our modelling approach *in vivo*. For 10 days of diet adaptation period, cows were offered total mixed ration consisting of grass silage (timothy-meadow fescue sward) preserved with formic acid based additive (AIV Ässä Na; 5 litres/tonne) provided at 45:55 forage to concentrate ratio on a dry matter basis. Concentrate mixture consisted of barley 210, oats 210, wheat 100, sugar beet pulp 220, rapeseed meal 230, and a mixture of minerals and vitamins 30 g/kg on as fed basis. Following the adaptation period, cows were located to metabolic chambers. Day 1 was dedicated to adaptation to the chamber conditions, while on d 2 and d 3 gas exchanges were measured for 48 h (for details see [33]). To determine circadian changes in rumen fermentation and rumen microbial community composition, on chamber d 4 and d 5 for 48 h period rumen liquid samples were collected every 3 h through ruminal fistula from the ventral site of the rumen using 500-mL bottle. Two sub-samples were taken for VFA and ammonia-N determination as described in [11] and were stored at -20°C for later analysis. Rumen liquid samples for microbial analysis were mixed to get homogenous distribution of microorganisms in the liquid, aliquoted into 2-mL sterile screw cap tubes, snap frozen on dry ice and stored at -80°C until DNA extraction. The total DNA was extracted from 500 µL of rumen liquid following protocol described in [34]. Libraries of the bacterial 16S ribosomal RNA (rRNA) V4 region were
prepared using 515F and 806R primers [35] and sequenced on Illumina MiSeq (Finnish Functional Genomics Centre, Turku) using the Paired-End approach and 2 x 250 bp chemistry. The sequencing data are accessible at NCBI SRA under the BioProject PRJNA1023082.

The cow experiments were conducted at the Luke Minkiö research dairy barn (Finland) and all experimental procedures were approved by the Project Authorisation Board (Regional Administrative Agency for Southern Finland, Hämeenlinna, Finland; ESAVI/34265/2019) in accordance with the guidelines established by the European Community Council Directive 86/609/EEC.

## Asymptotic observer and microbial functional proxies

Let us consider a metabolic process occurring in a continuous reactor where *n* microbial species *x*<sub>*i*</sub> grow at specific kinetic rates *r*<sub>*i*</sub> and produce the compound *s*. We consider that the reactor has a known dilution rate *D* and that the microbes leave the system at a lower rate than the dilution rate, which implies that the residence time of microbes is higher than the residence time of soluble compounds. The higher microbial residence time results from biological phenomena such as attachment to solid particles and formation of microbial aggregates. A simple way to model this phenomena is by adding a residence time factor *α* [36]. The dynamics of *s* and a microbial species *x*<sub>*i*</sub> follow:

$$ \frac{dx_i}{dt} = r_i - \alpha Dx_i $$

$$ \frac{ds}{dt} = \sum_{i=1}^{n} Y_{s,i} r_i - Ds $$

with $Y_{s,i}$ the production yield for the reaction $r_i$. The kinetic rates $r_i$ can be represented by a mathematical function like the Monod or Haldane equations. These kinetic functions are defined by specific parameters for each microbe $x_i$. The model equations [1, 2] assume that the microbe $x_i$ participates only in reaction $r_i$ which does not hold for many rumen microbes. We developed here an alternative approach adapted to the rumen microbiota that overpasses the need of assigning functionalities to microbial species. As observed in equation [2], the production of *s* results from the collective metabolic activity of the microbial consortium. In our approach, we assumed that the production of *s* can be related to a subunit of the specific metabolic pathways involved in the production of *s*. This subunit gathers the functional activity of the full consortium and can be used as a proxy of microbial activity of the rumen microbiome. The construction of the microbial proxy will be detailed later on. We will call $m_j$ the microbial proxy associated with the production of $s_j$. Here, our compounds $s_j$ are the major VFA from rumen fermentation: acetate ($s_{ac}$), butyrate ($s_{bu}$) and pro-pionate ($s_{pr}$).

For the compound $s_j$, the resulting model is defined by two equations

$$ \frac{dm_j}{dt} = \rho_j - \alpha Dm_j $$

$$ \frac{ds_j}{dt} = Y_j \rho_j - Ds_j $$

Where $\rho_j$ represents the reaction rate catalysed by the microbial proxy $m_j$. The next step of the approach consists of building an asymptotic observer that will enable us to estimate the dynamics of $s_j$ from measurements of $m_j$. For that, let us consider the following state transformation:

$$ z_j = s_j - Y_j m_j $$
prepared using 515F and 806R primers [35] and sequenced on Illumina MiSeq (Finnish Functional Genomics Centre, Turku) using the Paired-End approach and 2 x 250 bp chemistry. The sequencing data are accessible at NCBI SRA under the BioProject PRJNA1023082.

The cow experiments were conducted at the Luke Minkiö research dairy barn (Finland) and all experimental procedures were approved by the Project Authorisation Board (Regional Administrative Agency for Southern Finland, Hämeenlinna, Finland; ESAVI/34265/2019) in accordance with the guidelines established by the European Community Council Directive 86/609/EEC.

## Asymptotic observer and microbial functional proxies

Let us consider a metabolic process occurring in a continuous reactor where n microbial species $x_i$ grow at specific kinetic rates $r_i$ and produce the compound s. We consider that the reactor has a known dilution rate D and that the microbes leave the system at a lower rate than the dilution rate, which implies that the residence time of microbes is higher than the residence time of soluble compounds. The higher microbial residence time results from biological phenomena such as attachment to solid particles and formation of microbial aggregates. A simple way to model this phenomena is by adding a residence time factor $\alpha$ [36]. The dynamics of s and a microbial species $x_i$ follow:

$$ \frac{dx_i}{dt} = r_i - \alpha Dx_i $$

$$ \frac{ds}{dt} = \sum_{i=1}^{n} Y_{s,i} r_i - Ds $$

with $Y_{s,i}$ the production yield for the reaction $r_i$. The kinetic rates $r_i$ can be represented by a mathematical function like the Monod or Haldane equations. These kinetic functions are defined by specific parameters for each microbe $x_i$. The model equations [1, 2] assume that the microbe $x_i$ participates only in reaction $r_i$ which does not hold for many rumen microbes. We developed here an alternative approach adapted to the rumen microbiota that overpasses the need of assigning functionalities to microbial species. As observed in equation [2], the production of s results from the collective metabolic activity of the microbial consortium. In our approach, we assumed that the production of s can be related to a subunit of the specific metabolic pathways involved in the production of s. This subunit gathers the functional activity of the full consortium and can be used as a proxy of microbial activity of the rumen microbiome. The construction of the microbial proxy will be detailed later on. We will call $m_j$ the microbial proxy associated with the production of $s_j$. Here, our compounds $s_j$ are the major VFA from rumen fermentation: acetate ($s_{ac}$), butyrate ($s_{bu}$) and pro-pionate ($s_{pr}$).

For the compound $s_j$, the resulting model is defined by two equations

$$ \frac{dm_j}{dt} = \rho_j - \alpha Dm_j $$

$$ \frac{ds_j}{dt} = Y_j \rho_j - Ds_j $$

Where $\rho_j$ represents the reaction rate catalysed by the microbial proxy $m_j$. The next step of the approach consists of building an asymptotic observer that will enable us to estimate the dynamics of $s_j$ from measurements of $m_j$. For that, let us consider the following state transformation:

$$ z_j = s_j - Y_j m_j $$
By deriving in time [5], we get

$$ \frac{dz_j}{dt} = -D[z_j + (1 - \alpha)Y_j m_j] $$

Let us denote $\hat{z}_j$ an on-line estimate of $z_j$. The dynamics of $\hat{z}_j$ follows

$$ \frac{d\hat{z}_j}{dt} = -D[\hat{z}_j + (1 - \alpha)Y_j m_j] $$

To simulate $\hat{z}_j$, we require dynamic data of $m_j$ and the parameters $\alpha$, $Y_j$ to be known. Under these conditions, we can use the estimate $\hat{z}_j$ and the dynamic data $m_j$ to provide an estimate $\hat{s}_j$ at each sampling time:

$$ \hat{s}_j = \hat{z}_j + Y_j m_j $$

The asymptotic observer is the conjunction of the measurements $m_j$ with equations [7, 8]. As mentioned before, the approach does not require an explicit definition of the kinetic rate function $\rho_j$. We applied this development for acetate, butyrate and propionate, which requires measurements of the microbial proxies $m_{ac}$, $m_{bu}$, $m_{pr}$ and the yields $Y_{ac}$, $Y_{bu}$, $Y_{pr}$.

For both *in vitro* and *in vivo* case studies, the specific microbial proxy associated to the production of each VFA was derived from rumen microbial time series data determined by 16S rRNA gene amplicon sequencing. Raw data were analyzed using QIIME2 (version 2021.8) [37] using the script (run-samples) in https://github.com/frubino/cowpi and https://doi.org/10.5281/zenodo.8401851. Data were imported and denoised using the trim length option set at 400 and clustered using vsearch [38], in particular the 'cluster-features-de-novo' command with a percent identity set to 0.99. Classification was performed using the 'feature-classifier classify-consensus-blast' that uses BLAST [39] to classify the representative sequences, using default options and the Silva database (version 138) [40]. Finally, the OTU table representing microbial abundance time series was exported to be used in subsequent analysis steps. The resulting OTU table and representative sequences were then analyzed using CowPI [41] and modified using the information available at https://github.com/frubino/cowpi. The approach used here includes scaling data using DESeq2 [42] and infers modules from KEGG [43] instead of pathways, using the information in CowPI data. The pipeline for the inference of the functional modules is available at https://github.com/frubino/cowpi. KEGG modules resolve more specific processes and allow for greater precision when performing metabolic analysis, and like pathways, genes are their building blocks, so CowPI data can be used. For our analysis, modules allowed us to focus on more detailed aspects of metabolism for our model that pathways would not permit. Consequently, we used information about modules instead of pathways. However, there was no module in KEGG to present the butyrate metabolism. Therefore, an additional module for the metabolism of butyrate was introduced using information from [44]. These modules are the microbial proxies $m_j$. The abundances of the microbial proxies are relative measurements. However, since the data are scaled, we assumed that the absolute concentrations of microbial proxies are proportional to their relative abundances. As mentioned above, the asymptotic observer requires $\alpha$ and the yield factors to be known. The factor $\alpha$ was set to 0.85 according to [45]. We estimated the yield factors via the maximum likelihood (ML) approach as implemented in the Matlab IDEAS toolbox [46], which is freely available at http://genome.jouy.inra.fr/logiciels/IDEAS. For the optimization step, the Nelder-Mead Simplex method [47] implemented in the fminseach function was used. The model performance was
By deriving in time [5], we get

$$ \frac{dz_j}{dt} = -D[z_j + (1 - \alpha)Y_j m_j] $$

Let us denote $\hat{z}_j$ an on-line estimate of $z_j$. The dynamics of $\hat{z}_j$ follows

$$ \frac{d\hat{z}_j}{dt} = -D[\hat{z}_j + (1 - \alpha)Y_j m_j] $$

To simulate $\hat{z}_j$, we require dynamic data of $m_j$ and the parameters $\alpha$, $Y_j$ to be known. Under these conditions, we can use the estimate $\hat{z}_j$ and the dynamic data $m_j$ to provide an estimate $\hat{s}_j$ at each sampling time:

$$ \hat{s}_j = \hat{z}_j + Y_j m_j $$

The asymptotic observer is the conjunction of the measurements $m_j$ with equations [7, 8]. As mentioned before, the approach does not require an explicit definition of the kinetic rate function $\rho_j$. We applied this development for acetate, butyrate and propionate, which requires measurements of the microbial proxies $m_{ac}$, $m_{bu}$, $m_{pr}$ and the yields $Y_{ac}$, $Y_{bu}$, $Y_{pr}$.

For both *in vitro* and *in vivo* case studies, the specific microbial proxy associated to the production of each VFA was derived from rumen microbial time series data determined by 16S rRNA gene amplicon sequencing. Raw data were analyzed using QIIME2 (version 2021.8) [37] using the script (run-samples) in https://github.com/frubino/cowpi and https://doi.org/10.5281/zenodo.8401851. Data were imported and denoised using the trim length option set at 400 and clustered using vsearch [38], in particular the 'cluster-features-de-novo' command with a percent identity set to 0.99. Classification was performed using the 'feature-classifier classify-consensus-blast' that uses BLAST [39] to classify the representative sequences, using default options and the Silva database (version 138) [40]. Finally, the OTU table representing microbial abundance time series was exported to be used in subsequent analysis steps. The resulting OTU table and representative sequences were then analyzed using CowPI [41] and modified using the information available at https://github.com/frubino/cowpi. The approach used here includes scaling data using DESeq2 [42] and infers modules from KEGG [43] instead of pathways, using the information in CowPI data. The pipeline for the inference of the functional modules is available at https://github.com/frubino/cowpi. KEGG modules resolve more specific processes and allow for greater precision when performing metabolic analysis, and like pathways, genes are their building blocks, so CowPI data can be used. For our analysis, modules allowed us to focus on more detailed aspects of metabolism for our model that pathways would not permit. Consequently, we used information about modules instead of pathways. However, there was no module in KEGG to present the butyrate metabolism. Therefore, an additional module for the metabolism of butyrate was introduced using information from [44]. These modules are the microbial proxies $m_j$. The abundances of the microbial proxies are relative measurements. However, since the data are scaled, we assumed that the absolute concentrations of microbial proxies are proportional to their relative abundances. As mentioned above, the asymptotic observer requires $\alpha$ and the yield factors to be known. The factor $\alpha$ was set to 0.85 according to [45]. We estimated the yield factors via the maximum likelihood (ML) approach as implemented in the Matlab IDEAS toolbox [46], which is freely available at http://genome.jouy.inra.fr/logiciels/IDEAS. For the optimization step, the Nelder-Mead Simplex method [47] implemented in the fminseach function was used. The model performance was
assessed by computing the coefficient of variation of the root mean squared error (CVRMSE). The Matlab scripts are available at [48] https://doi.org/10.5281/zenodo.8386786.

## Results

The scripts for the inference of the functional microbial modules are available at https://github.com/frubino/cowpi and https://doi.org/10.5281/zenodo.8401851. The abundances of each module for each case study are in the Tables modules.Rusitec.xls and modules.Cows in [48]. A total of 308 modules were identified. The microbial proxies for VFA production are named as M00579 ($m_{ac}$), M99999 ($m_{bu}$), M00013 ($m_{pr}$). The implementation of the observer for each case study is available at [48].

### In vitro case study

Table 1 shows the estimated abundance of the microbial proxies of VFA production. The values are the sum of the abundances from the liquid and solid phase. Fig 1 shows the dynamics of VFA (acetate, butyrate, and propionate) compared to the estimated VFA concentrations by the observer for the *in vitro* experiments with grass and hay using the rumen inoculum from two cows. Table 2 shows the model evaluation in terms of the CVRMSE. The model performance is satisfactory with average CVRMSE for acetate, butyrate and propionate of 9.8%, 14% and 14.5%. Table 3 shows the estimated yield factors for each diet and inoculum. The average values of the yields of acetate, butyrate and propionate are $0.13 \times 10^{-6}$, $0.39 \times 10^{-7}$ and $0.49 \times 10^{-6}$.

### In vivo case study

Table 4 shows the estimated abundance of the microbial proxies of VFA production from the liquid phase. Fig 2 shows the dynamics of VFA compared to the estimated VFA concentrations by the observer.

Table 5 shows the CVRMSE for the VFA concentrations. The model performance is adequate with average CVRMSE for acetate, butyrate and propionate of 16.4%, 15.8% and 19.8%. As shown in Table 6, the model performance is very satisfactory with regard to prediction of VFA molar proportions. The CVRMSE for the molar proportions of acetate, butyrate and propionate are $0.13 \times 10^{-6}$, $0.39 \times 10^{-7}$ and $0.49 \times 10^{-6}$.

<table><caption>**Table 1. Microbial abundance time series of functional proxies of VFA production for the *in vitro* case study.** Values are the sum of the abundances from the liquid and solid phase.</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th colspan="4">Sampling time (h)</th></tr><tr><th>2</th><th>4</th><th>8</th><th>24</th></tr></thead><tbody><tr><td rowspan="3">Rumen inoculum 1-Grass</td><td>m<sub>ac</sub></td><td>33139</td><td>37601</td><td>60680</td><td>70338</td></tr><tr><td>m<sub>bu</sub></td><td>52574</td><td>57719</td><td>116469</td><td>109282</td></tr><tr><td>m<sub>pr</sub></td><td>3503</td><td>3669</td><td>8181</td><td>8971</td></tr><tr><td rowspan="3">Rumen inoculum 1-Hay</td><td>m<sub>ac</sub></td><td>14487</td><td>42303</td><td>36663</td><td>82646</td></tr><tr><td>m<sub>bu</sub></td><td>18600</td><td>61808</td><td>49934</td><td>118876</td></tr><tr><td>m<sub>pr</sub></td><td>640</td><td>4323</td><td>3673</td><td>13470</td></tr><tr><td rowspan="3">Rumen inoculum 2-Grass</td><td>m<sub>ac</sub></td><td>87093</td><td>91845</td><td>105360</td><td>83174</td></tr><tr><td>m<sub>bu</sub></td><td>184833</td><td>206660</td><td>227768</td><td>235760</td></tr><tr><td>m<sub>pr</sub></td><td>23114</td><td>27944</td><td>24257</td><td>36697</td></tr><tr><td rowspan="3">Rumen inoculum 2-Hay</td><td>m<sub>ac</sub></td><td>111166</td><td>106780</td><td>101168</td><td>102789</td></tr><tr><td>m<sub>bu</sub></td><td>240822</td><td>240218</td><td>194157</td><td>218833</td></tr><tr><td>m<sub>pr</sub></td><td>31287</td><td>31758</td><td>22945</td><td>22005</td></tr></tbody></table>

$m_{ac}$, $m_{bu}$, $m_{pr}$: abundances of microbial functional proxies for acetate, butyrate and propionate production.

https://doi.org/10.1371/journal.pone.0298930.t001
assessed by computing the coefficient of variation of the root mean squared error (CVRMSE). The Matlab scripts are available at [48] https://doi.org/10.5281/zenodo.8386786.

## Results

The scripts for the inference of the functional microbial modules are available at https://github.com/frubino/cowpi and https://doi.org/10.5281/zenodo.8401851. The abundances of each module for each case study are in the Tables modules.Rusitec.xls and modules.Cows in [48]. A total of 308 modules were identified. The microbial proxies for VFA production are named as M00579 ($m_{ac}$), M99999 ($m_{bu}$), M00013 ($m_{pr}$). The implementation of the observer for each case study is available at [48].

### In vitro case study

Table 1 shows the estimated abundance of the microbial proxies of VFA production. The values are the sum of the abundances from the liquid and solid phase. Fig 1 shows the dynamics of VFA (acetate, butyrate, and propionate) compared to the estimated VFA concentrations by the observer for the *in vitro* experiments with grass and hay using the rumen inoculum from two cows. Table 2 shows the model evaluation in terms of the CVRMSE. The model performance is satisfactory with average CVRMSE for acetate, butyrate and propionate of 9.8%, 14% and 14.5%. Table 3 shows the estimated yield factors for each diet and inoculum. The average values of the yields of acetate, butyrate and propionate are $0.13 \times 10^{-6}$, $0.39 \times 10^{-7}$ and $0.49 \times 10^{-6}$.

### In vivo case study

Table 4 shows the estimated abundance of the microbial proxies of VFA production from the liquid phase. Fig 2 shows the dynamics of VFA compared to the estimated VFA concentrations by the observer.

Table 5 shows the CVRMSE for the VFA concentrations. The model performance is adequate with average CVRMSE for acetate, butyrate and propionate of 16.4%, 15.8% and 19.8%. As shown in Table 6, the model performance is very satisfactory with regard to prediction of VFA molar proportions. The CVRMSE for the molar proportions of acetate, butyrate and

<table><caption>**Table 1. Microbial abundance time series of functional proxies of VFA production for the *in vitro* case study.** Values are the sum of the abundances from the liquid and solid phase.</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th colspan="4">Sampling time (h)</th></tr><tr><th>2</th><th>4</th><th>8</th><th>24</th></tr></thead><tbody><tr><td rowspan="3">Rumen inoculum 1-Grass</td><td>m<sub>ac</sub></td><td>33139</td><td>37601</td><td>60680</td><td>70338</td></tr><tr><td>m<sub>bu</sub></td><td>52574</td><td>57719</td><td>116469</td><td>109282</td></tr><tr><td>m<sub>pr</sub></td><td>3503</td><td>3669</td><td>8181</td><td>8971</td></tr><tr><td rowspan="3">Rumen inoculum 1-Hay</td><td>m<sub>ac</sub></td><td>14487</td><td>42303</td><td>36663</td><td>82646</td></tr><tr><td>m<sub>bu</sub></td><td>18600</td><td>61808</td><td>49934</td><td>118876</td></tr><tr><td>m<sub>pr</sub></td><td>640</td><td>4323</td><td>3673</td><td>13470</td></tr><tr><td rowspan="3">Rumen inoculum 2-Grass</td><td>m<sub>ac</sub></td><td>87093</td><td>91845</td><td>105360</td><td>83174</td></tr><tr><td>m<sub>bu</sub></td><td>184833</td><td>206660</td><td>227768</td><td>235760</td></tr><tr><td>m<sub>pr</sub></td><td>23114</td><td>27944</td><td>24257</td><td>36697</td></tr><tr><td rowspan="3">Rumen inoculum 2-Hay</td><td>m<sub>ac</sub></td><td>111166</td><td>106780</td><td>101168</td><td>102789</td></tr><tr><td>m<sub>bu</sub></td><td>240822</td><td>240218</td><td>194157</td><td>218833</td></tr><tr><td>m<sub>pr</sub></td><td>31287</td><td>31758</td><td>22945</td><td>22005</td></tr></tbody></table>

$m_{ac}$, $m_{bu}$, $m_{pr}$: abundances of microbial functional proxies for acetate, butyrate and propionate production.

https://doi.org/10.1371/journal.pone.0298930.t001
<figure><img src="image_3.png" /><figcaption>Fig 1. Experimental data of VFA concentrations from a RUSITEC experiment with two feeds (grass and hay) and rumen fluid from two cows (○: cow 1; ◊: cow 2) are compared against the estimated concentrations by the state observer (red solid line for cow 1, dashed blue line for cow 2).</figcaption></figure>

https://doi.org/10.1371/journal.pone.0298930.g001

propionate are 3.1%, 3.8% and 8.9%. Table 7 shows the estimated yield factors for each cow. The average values of the yields of acetate, butyrate and propionate are 1.43 x10⁻⁶, 1.62x10⁻⁷ and 2.43x10⁻⁶.

Discussion

The objective of this work was to integrate rumen microbiota time series determined by 16S rRNA gene amplicon sequencing into a mathematical model linking microbial data to the dynamics of the volatile fatty acids (VFA) production during rumen fermentation. This objective followed the rational that microbial data can be used to enhance predictive capabilities of rumen fermentation models. Our model development provided satisfactory results for estimating the dynamics of VFA from microbial data. The *in vivo* study showed that yield factors for the specific VFA production are similar between the cows. It will be useful to test our approach with animals fed at different diets including diets supplemented with additives impacting rumen fermentation patterns.

Given the limited number of animals used in our study, it is difficult to provide a fair comparison with existing rumen models. Furthermore, few rumen modelling studies have compared model predictions against dynamic data of VFA concentrations, with the exception of

<table><caption>Table 2. Coefficient of variation of the root mean square error (%) for the observer applied to *in vitro* data of VFA concentrations.</caption><thead><tr><th></th><th>Acetate</th><th>Butyrate</th><th>Propionate</th></tr></thead><tbody><tr><td>Rumen inoculum 1- Grass</td><td>7.1</td><td>11.8</td><td>13.2</td></tr><tr><td>Rumen inoculum 2- Grass</td><td>4.6</td><td>7.8</td><td>1.3</td></tr><tr><td>Rumen inoculum 1- Hay</td><td>16.9</td><td>20</td><td>25.6</td></tr><tr><td>Rumen inoculum 2- Hay</td><td>10.5</td><td>16.4</td><td>18</td></tr><tr><td>Mean ± s.d.</td><td>9.8 ± 5.3</td><td>14 ± 5.3</td><td>14.5 ± 10.2</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0298930.t002
<figure><img src="image_5.png" /><figcaption>Fig 1. Experimental data of VFA concentrations from a RUSITEC experiment with two feeds (grass and hay) and rumen fluid from two cows (○: cow 1; ◊: cow 2) are compared against the estimated concentrations by the state observer (red solid line for cow 1, dashed blue line for cow 2).</figcaption></figure>

https://doi.org/10.1371/journal.pone.0298930.g001

propionate are 3.1%, 3.8% and 8.9%. Table 7 shows the estimated yield factors for each cow. The average values of the yields of acetate, butyrate and propionate are 1.43 x10⁻⁶, 1.62x10⁻⁷ and 2.43x10⁻⁶.

Discussion

The objective of this work was to integrate rumen microbiota time series determined by 16S rRNA gene amplicon sequencing into a mathematical model linking microbial data to the dynamics of the volatile fatty acids (VFA) production during rumen fermentation. This objective followed the rational that microbial data can be used to enhance predictive capabilities of rumen fermentation models. Our model development provided satisfactory results for estimating the dynamics of VFA from microbial data. The *in vivo* study showed that yield factors for the specific VFA production are similar between the cows. It will be useful to test our approach with animals fed at different diets including diets supplemented with additives impacting rumen fermentation patterns.

Given the limited number of animals used in our study, it is difficult to provide a fair comparison with existing rumen models. Furthermore, few rumen modelling studies have compared model predictions against dynamic data of VFA concentrations, with the exception of

<table><caption>Table 2. Coefficient of variation of the root mean square error (%) for the observer applied to *in vitro* data of VFA concentrations.</caption><thead><tr><th></th><th>Acetate</th><th>Butyrate</th><th>Propionate</th></tr></thead><tbody><tr><td>Rumen inoculum 1- Grass</td><td>7.1</td><td>11.8</td><td>13.2</td></tr><tr><td>Rumen inoculum 2- Grass</td><td>4.6</td><td>7.8</td><td>1.3</td></tr><tr><td>Rumen inoculum 1- Hay</td><td>16.9</td><td>20</td><td>25.6</td></tr><tr><td>Rumen inoculum 2- Hay</td><td>10.5</td><td>16.4</td><td>18</td></tr><tr><td>Mean ± s.d.</td><td>9.8 ± 5.3</td><td>14 ± 5.3</td><td>14.5 ± 10.2</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0298930.t002
<table><caption>Table 3. Estimated yield parameters for the *in vitro* data.</caption><thead><tr><th></th><th>Y<sub>ac</sub></th><th>Y<sub>bu</sub></th><th>Y<sub>pr</sub></th></tr></thead><tbody><tr><td>Rumen inoculum 1- Grass</td><td>0.19x10<sup>-6</sup></td><td>0.77x10<sup>-7</sup></td><td>0.93x10<sup>-6</sup></td></tr><tr><td>Rumen inoculum 2-Grass</td><td>0.17x10<sup>-6</sup></td><td>0.26x10<sup>-7</sup></td><td>0.14x10<sup>-6</sup></td></tr><tr><td>Rumen inoculum 1-Hay</td><td>0.07x10<sup>-6</sup></td><td>0.30x10<sup>-7</sup></td><td>0.29x10<sup>-6</sup></td></tr><tr><td>Rumen inoculum 2-Hay</td><td>0.11x10<sup>-6</sup></td><td>0.22x10<sup>-7</sup></td><td>0.23x10<sup>-6</sup></td></tr><tr><td>Mean ± s.d.</td><td>0.13x10<sup>-6</sup> ± 5.3x10<sup>-8</sup></td><td>0.39x10<sup>-7</sup> ± 2.54x10<sup>-8</sup></td><td>0.40x10<sup>-6</sup> ± 3.61x10<sup>-7</sup></td></tr></tbody></table>

Y<sub>ac</sub>, Y<sub>bu</sub>, Y<sub>pr</sub>: yields of acetate, butyrate and propionate production (mol VFA / microbial functional proxy abundance).

https://doi.org/10.1371/journal.pone.0298930.t003

the work of [16]. Under this context, however, it should be mentioned that the CVRMSE of our modelling approach are lower than those reported in [16], which strengthens our model development. From the modelling perspective, the structure of the model developed here is very simple compared to the structure of the existing rumen models. This parsimonious property is of usefulness for the development of *in silico* tools to predict rumen function from microbial data within a Precision Livestock Farming context. It should be said, however, that our study is a theoretical work that demonstrates the proof of concept that rumen microbial data can be used to predict rumen fermentation variables. However, the application of our modelling approach in farm conditions is currently unfeasible due to the need of rumen sampling and the costs of microbial analysis. To overcome the obstacle of rumen accessibility, it will be useful to test our approach using buccal microbial samples as proxies of the rumen microbiota [49, 50]. Buccal sampling might be a solution for getting quick samples and if it is coupled with easy sequencing instrument, it might become possible to use this approach to estimate rumen VFA using our modelling approach.

The fundamental aspect of our approach is the definition of the microbial proxies for each VFA. The microbial proxies here developed followed the principle that the rumen ecosystem operates as supra-organism provided with all the metabolic capabilities of its species. This approach does not account for species connectivity which is a relevant aspect in gut microbial ecosystems [51]. It is indeed interesting that the aggregated microbial proxies provided an

<table><caption>Table 4. Microbial abundance time series of functional proxies of VFA production for the *in vivo* case study. Values are from the liquid phase.</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th colspan="16">Sampling time (h)</th></tr><tr><th>7</th><th>10</th><th>13</th><th>16</th><th>19</th><th>22</th><th>25</th><th>28</th><th>31</th><th>34</th><th>37</th><th>40</th><th>43</th><th>46</th><th>49</th><th>52</th></tr></thead><tbody><tr><td rowspan="3">Cow 1</td><td>m<sub>ac</sub></td><td>47708</td><td>45126</td><td>55602</td><td>57647</td><td>64097</td><td>56031</td><td>54965</td><td>65061</td><td>54497</td><td>69391</td><td>57204</td><td>67378</td><td>63766</td><td>67489</td><td>57861</td><td>56384</td></tr><tr><td>m<sub>bu</sub></td><td>76650</td><td>72740</td><td>92161</td><td>95845</td><td>108538</td><td>89023</td><td>90047</td><td>108791</td><td>90950</td><td>117048</td><td>94074</td><td>115549</td><td>109147</td><td>112391</td><td>96751</td><td>97451</td></tr><tr><td>m<sub>pr</sub></td><td>7229</td><td>6963</td><td>10052</td><td>9955</td><td>11753</td><td>8104</td><td>8651</td><td>11608</td><td>9396</td><td>12755</td><td>10278</td><td>13197</td><td>11730</td><td>11249</td><td>10067</td><td>11329</td></tr><tr><td rowspan="3">Cow 2</td><td>m<sub>ac</sub></td><td>49745</td><td>48073</td><td>52536</td><td>49837</td><td>59831</td><td>58632</td><td>58376</td><td>40248</td><td>66988</td><td>69821</td><td>67276</td><td>54960</td><td>64366</td><td>50587</td><td>58103</td><td>54381</td></tr><tr><td>m<sub>bu</sub></td><td>80804</td><td>80864</td><td>87082</td><td>84439</td><td>98867</td><td>102169</td><td>97994</td><td>67514</td><td>113688</td><td>118291</td><td>114353</td><td>94216</td><td>110819</td><td>87848</td><td>100618</td><td>96853</td></tr><tr><td>m<sub>pr</sub></td><td>8064</td><td>8626</td><td>9009</td><td>8991</td><td>9808</td><td>11502</td><td>9774</td><td>6934</td><td>12201</td><td>12928</td><td>12633</td><td>10080</td><td>12292</td><td>9732</td><td>11186</td><td>10745</td></tr><tr><td rowspan="3">Cow 3</td><td>m<sub>ac</sub></td><td>42144</td><td>48535</td><td>50274</td><td>55855</td><td>45496</td><td>58116</td><td>49365</td><td>54277</td><td>52852</td><td>66731</td><td>60735</td><td>53184</td><td>54231</td><td>65121</td><td>49727</td><td>51256</td></tr><tr><td>m<sub>bu</sub></td><td>73426</td><td>84226</td><td>85419</td><td>97337</td><td>78949</td><td>102467</td><td>80924</td><td>93929</td><td>89409</td><td>111399</td><td>102816</td><td>90607</td><td>92020</td><td>110914</td><td>84625</td><td>91281</td></tr><tr><td>m<sub>pr</sub></td><td>8149</td><td>9554</td><td>9264</td><td>11066</td><td>8440</td><td>11366</td><td>7839</td><td>10324</td><td>9596</td><td>11580</td><td>11227</td><td>10341</td><td>9335</td><td>11322</td><td>9163</td><td>10974</td></tr><tr><td rowspan="3">Cow 4</td><td>m<sub>ac</sub></td><td>46298</td><td>12944</td><td>58665</td><td>60244</td><td>54556</td><td>55437</td><td>53430</td><td>54539</td><td>60154</td><td>64751</td><td>56813</td><td>56088</td><td>58200</td><td>53288</td><td>52679</td><td>57640</td></tr><tr><td>m<sub>bu</sub></td><td>77495</td><td>24482</td><td>99615</td><td>100443</td><td>90413</td><td>95868</td><td>86833</td><td>90640</td><td>103506</td><td>114801</td><td>93832</td><td>93850</td><td>101479</td><td>89818</td><td>91936</td><td>98357</td></tr><tr><td>m<sub>pr</sub></td><td>8128</td><td>2847</td><td>10811</td><td>10838</td><td>9539</td><td>10574</td><td>8927</td><td>9399</td><td>12079</td><td>13578</td><td>10265</td><td>9756</td><td>11216</td><td>9414</td><td>10771</td><td>11167</td></tr></tbody></table>

m<sub>ac</sub>, m<sub>bu</sub>, m<sub>pr</sub>: abundances of microbial functional proxies for acetate, butyrate and propionate production.

https://doi.org/10.1371/journal.pone.0298930.t004
<table><caption>Table 3. Estimated yield parameters for the *in vitro* data.</caption><thead><tr><th></th><th>Y<sub>ac</sub></th><th>Y<sub>bu</sub></th><th>Y<sub>pr</sub></th></tr></thead><tbody><tr><td>Rumen inoculum 1- Grass</td><td>0.19x10<sup>-6</sup></td><td>0.77x10<sup>-7</sup></td><td>0.93x10<sup>-6</sup></td></tr><tr><td>Rumen inoculum 2-Grass</td><td>0.17x10<sup>-6</sup></td><td>0.26x10<sup>-7</sup></td><td>0.14x10<sup>-6</sup></td></tr><tr><td>Rumen inoculum 1-Hay</td><td>0.07x10<sup>-6</sup></td><td>0.30x10<sup>-7</sup></td><td>0.29x10<sup>-6</sup></td></tr><tr><td>Rumen inoculum 2-Hay</td><td>0.11x10<sup>-6</sup></td><td>0.22x10<sup>-7</sup></td><td>0.23x10<sup>-6</sup></td></tr><tr><td>Mean ± s.d.</td><td>0.13x10<sup>-6</sup> ± 5.3x10<sup>-8</sup></td><td>0.39x10<sup>-7</sup> ± 2.54x10<sup>-8</sup></td><td>0.40x10<sup>-6</sup> ± 3.61x10<sup>-7</sup></td></tr></tbody></table>

Y<sub>ac</sub>, Y<sub>bu</sub>, Y<sub>pr</sub>: yields of acetate, butyrate and propionate production (mol VFA / microbial functional proxy abundance).

https://doi.org/10.1371/journal.pone.0298930.t003

the work of [16]. Under this context, however, it should be mentioned that the CVRMSE of our modelling approach are lower than those reported in [16], which strengthens our model development. From the modelling perspective, the structure of the model developed here is very simple compared to the structure of the existing rumen models. This parsimonious property is of usefulness for the development of *in silico* tools to predict rumen function from microbial data within a Precision Livestock Farming context. It should be said, however, that our study is a theoretical work that demonstrates the proof of concept that rumen microbial data can be used to predict rumen fermentation variables. However, the application of our modelling approach in farm conditions is currently unfeasible due to the need of rumen sampling and the costs of microbial analysis. To overcome the obstacle of rumen accessibility, it will be useful to test our approach using buccal microbial samples as proxies of the rumen microbiota [49, 50]. Buccal sampling might be a solution for getting quick samples and if it is coupled with easy sequencing instrument, it might become possible to use this approach to estimate rumen VFA using our modelling approach.

The fundamental aspect of our approach is the definition of the microbial proxies for each VFA. The microbial proxies here developed followed the principle that the rumen ecosystem operates as supra-organism provided with all the metabolic capabilities of its species. This approach does not account for species connectivity which is a relevant aspect in gut microbial ecosystems [51]. It is indeed interesting that the aggregated microbial proxies provided an

<table><caption>Table 4. Microbial abundance time series of functional proxies of VFA production for the *in vivo* case study. Values are from the liquid phase.</caption><thead><tr><th></th><th></th><th colspan="16">Sampling time (h)</th></tr><tr><th></th><th></th><th>7</th><th>10</th><th>13</th><th>16</th><th>19</th><th>22</th><th>25</th><th>28</th><th>31</th><th>34</th><th>37</th><th>40</th><th>43</th><th>46</th><th>49</th><th>52</th></tr></thead><tbody><tr><td rowspan="3">Cow 1</td><td>m<sub>ac</sub></td><td>47708</td><td>45126</td><td>55602</td><td>57647</td><td>64097</td><td>56031</td><td>54965</td><td>65061</td><td>54497</td><td>69391</td><td>57204</td><td>67378</td><td>63766</td><td>67489</td><td>57861</td><td>56384</td></tr><tr><td>m<sub>bu</sub></td><td>76650</td><td>72740</td><td>92161</td><td>95845</td><td>108538</td><td>89023</td><td>90047</td><td>108791</td><td>90950</td><td>117048</td><td>94074</td><td>115549</td><td>109147</td><td>112391</td><td>96751</td><td>97451</td></tr><tr><td>m<sub>pr</sub></td><td>7229</td><td>6963</td><td>10052</td><td>9955</td><td>11753</td><td>8104</td><td>8651</td><td>11608</td><td>9396</td><td>12755</td><td>10278</td><td>13197</td><td>11730</td><td>11249</td><td>10067</td><td>11329</td></tr><tr><td rowspan="3">Cow 2</td><td>m<sub>ac</sub></td><td>49745</td><td>48073</td><td>52536</td><td>49837</td><td>59831</td><td>58632</td><td>58376</td><td>40248</td><td>66988</td><td>69821</td><td>67276</td><td>54960</td><td>64366</td><td>50587</td><td>58103</td><td>54381</td></tr><tr><td>m<sub>bu</sub></td><td>80804</td><td>80864</td><td>87082</td><td>84439</td><td>98867</td><td>102169</td><td>97994</td><td>67514</td><td>113688</td><td>118291</td><td>114353</td><td>94216</td><td>110819</td><td>87848</td><td>100618</td><td>96853</td></tr><tr><td>m<sub>pr</sub></td><td>8064</td><td>8626</td><td>9009</td><td>8991</td><td>9808</td><td>11502</td><td>9774</td><td>6934</td><td>12201</td><td>12928</td><td>12633</td><td>10080</td><td>12292</td><td>9732</td><td>11186</td><td>10745</td></tr><tr><td rowspan="3">Cow 3</td><td>m<sub>ac</sub></td><td>42144</td><td>48535</td><td>50274</td><td>55855</td><td>45496</td><td>58116</td><td>49365</td><td>54277</td><td>52852</td><td>66731</td><td>60735</td><td>53184</td><td>54231</td><td>65121</td><td>49727</td><td>51256</td></tr><tr><td>m<sub>bu</sub></td><td>73426</td><td>84226</td><td>85419</td><td>97337</td><td>78949</td><td>102467</td><td>80924</td><td>93929</td><td>89409</td><td>111399</td><td>102816</td><td>90607</td><td>92020</td><td>110914</td><td>84625</td><td>91281</td></tr><tr><td>m<sub>pr</sub></td><td>8149</td><td>9554</td><td>9264</td><td>11066</td><td>8440</td><td>11366</td><td>7839</td><td>10324</td><td>9596</td><td>11580</td><td>11227</td><td>10341</td><td>9335</td><td>11322</td><td>9163</td><td>10974</td></tr><tr><td rowspan="3">Cow 4</td><td>m<sub>ac</sub></td><td>46298</td><td>12944</td><td>58665</td><td>60244</td><td>54556</td><td>55437</td><td>53430</td><td>54539</td><td>60154</td><td>64751</td><td>56813</td><td>56088</td><td>58200</td><td>53288</td><td>52679</td><td>57640</td></tr><tr><td>m<sub>bu</sub></td><td>77495</td><td>24482</td><td>99615</td><td>100443</td><td>90413</td><td>95868</td><td>86833</td><td>90640</td><td>103506</td><td>114801</td><td>93832</td><td>93850</td><td>101479</td><td>89818</td><td>91936</td><td>98357</td></tr><tr><td>m<sub>pr</sub></td><td>8128</td><td>2847</td><td>10811</td><td>10838</td><td>9539</td><td>10574</td><td>8927</td><td>9399</td><td>12079</td><td>13578</td><td>10265</td><td>9756</td><td>11216</td><td>9414</td><td>10771</td><td>11167</td></tr></tbody></table>

m<sub>ac</sub>, m<sub>bu</sub>, m<sub>pr</sub>: abundances of microbial functional proxies for acetate, butyrate and propionate production.

https://doi.org/10.1371/journal.pone.0298930.t004
<figure><img src="image_4.png" /><figcaption>Fig 2. Experimental data of VFA concentrations (○) from an *in vivo* experiment with four cows are compared against the estimated concentrations by the state observer (blue solid line).</figcaption></figure>

https://doi.org/10.1371/journal.pone.0298930.g002

<table><caption>Table 5. Coefficient of variation of the root mean square error (%) for the observer applied to *in vivo* data of VFA concentrations.</caption><thead><tr><th></th><th>Acetate</th><th>Butyrate</th><th>Propionate</th></tr></thead><tbody><tr><td>Cow 1</td><td>10.1</td><td>10.2</td><td>14.7</td></tr><tr><td>Cow 2</td><td>15.4</td><td>14.1</td><td>17.1</td></tr><tr><td>Cow 3</td><td>16.5</td><td>18.4</td><td>25.9</td></tr><tr><td>Cow 4</td><td>23.8</td><td>20.5</td><td>21.3</td></tr><tr><td>Mean± s.d.</td><td>16.4 ± 5.6</td><td>15.8 ± 4.6</td><td>19.8 ± 4.9</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0298930.t005
<figure><img src="image_6.png" /><figcaption>Fig 2. Experimental data of VFA concentrations (○) from an *in vivo* experiment with four cows are compared against the estimated concentrations by the state observer (blue solid line).</figcaption></figure>

https://doi.org/10.1371/journal.pone.0298930.g002

<table><caption>Table 5. Coefficient of variation of the root mean square error (%) for the observer applied to *in vivo* data of VFA concentrations.</caption><thead><tr><th></th><th>Acetate</th><th>Butyrate</th><th>Propionate</th></tr></thead><tbody><tr><td>Cow 1</td><td>10.1</td><td>10.2</td><td>14.7</td></tr><tr><td>Cow 2</td><td>15.4</td><td>14.1</td><td>17.1</td></tr><tr><td>Cow 3</td><td>16.5</td><td>18.4</td><td>25.9</td></tr><tr><td>Cow 4</td><td>23.8</td><td>20.5</td><td>21.3</td></tr><tr><td>Mean± s.d.</td><td>16.4 ± 5.6</td><td>15.8 ± 4.6</td><td>19.8 ± 4.9</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0298930.t005
<table><caption>**Table 6.** Coefficient of variation of the root mean square error (%) for the observer applied to *in vivo* data of VFA molar proportions.</caption><thead><tr><th></th><th>Acetate</th><th>Butyrate</th><th>Propionate</th></tr></thead><tbody><tr><td>Cow 1</td><td>2.9</td><td>3.4</td><td>8.3</td></tr><tr><td>Cow 2</td><td>1.8</td><td>2.8</td><td>5.7</td></tr><tr><td>Cow 3</td><td>3.6</td><td>4.6</td><td>10.9</td></tr><tr><td>Cow 4</td><td>4.2</td><td>4.5</td><td>10.5</td></tr><tr><td>Mean ± s.d.</td><td>3.1 ± 1.0</td><td>3.8 ± 0.9</td><td>8.9 ± 2.4</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0298930.t006

<table><caption>**Table 7.** Estimated yield parameters for the *in vivo* data.</caption><thead><tr><th></th><th>Y<sub>ac</sub></th><th>Y<sub>bu</sub></th><th>Y<sub>pr</sub></th></tr></thead><tbody><tr><td>Cow 1</td><td>1.38x10<sup>-6</sup></td><td>1.73x10<sup>-7</sup></td><td>2.37x10<sup>-6</sup></td></tr><tr><td>Cow 2</td><td>1.52x10<sup>-6</sup></td><td>1.55x10<sup>-7</sup></td><td>2.16x10<sup>-6</sup></td></tr><tr><td>Cow 3</td><td>1.43x10<sup>-6</sup></td><td>1.50x10<sup>-7</sup></td><td>2.64x10<sup>-6</sup></td></tr><tr><td>Cow 4</td><td>1.39x10<sup>-6</sup></td><td>1.83x10<sup>-7</sup></td><td>2.54x10<sup>-6</sup></td></tr><tr><td>Mean ± s.d.</td><td>1.43x10<sup>-6</sup> ± 6.28x10<sup>-8</sup></td><td>1.62x10<sup>-7</sup> ± 1.56x10<sup>-8</sup></td><td>2.43x10<sup>-6</sup> ± 2.10x10<sup>-7</sup></td></tr></tbody></table>

Y<sub>ac</sub>, Y<sub>bu</sub>, Y<sub>pr</sub>: yields of acetate, butyrate and propionate production (mol VFA / microbial functional proxy abundance).

https://doi.org/10.1371/journal.pone.0298930.t007

adequate representation of VFA dynamics. Although the results are promising, future research is needed to refine the definition of these proxies, by taking into account quantitative tran scriptomic information to identify the actual metabolic activity.

Our study has taken a step in the direction of capitalizing on microbial data to predict rumen function. In this first step, we used asymptotic state observers due to its simple struc ture and the advantage of not requiring information on the mathematical functions describing the metabolic kinetic rates. However, asymptotic observers have the limitation that the conver gence rate cannot be tuned because it is determined by operational conditions (e.g., ruminal fractional passage rate). This limitation can be overcome by other observers. However, more sophisticated technicalities are needed for the design and implementation of such observers [26]. Our approach was applied to estimate VFA dynamics. We think that the approach can be extended to estimate other metabolites including methane production by the rumen methano- genic archaea. Such an application will be of great interesting for the monitoring of methane emissions and for the evaluation of methane inhibition strategies.

## Conclusions

Existing mechanistic models of rumen fermentation consider an aggregated representation of the rumen microbiota and its metabolic function. However, none of these models integrate microbial genomic knowledge and thus do not capitalize on the rich information of microbial genomic sequencing. In this work, we integrated microbial time series of the rumen micro biota determined by 16S rDNA into a dynamic model of the rumen microbiome. Our results showed the estimated VFA concentrations converge towards the real VFA concentration dynamics demonstrating the promising potential of our approach.

## Acknowledgments

Rafael Muñoz-Tamayo thanks Pablo Ugalde-Salas (Inria, France) and Jerôme Harmand (INRAE, France) for their helpful explanations on the use of state observers for microbial
<table><caption>**Table 6.** Coefficient of variation of the root mean square error (%) for the observer applied to *in vivo* data of VFA molar proportions.</caption><thead><tr><th></th><th>Acetate</th><th>Butyrate</th><th>Propionate</th></tr></thead><tbody><tr><td>Cow 1</td><td>2.9</td><td>3.4</td><td>8.3</td></tr><tr><td>Cow 2</td><td>1.8</td><td>2.8</td><td>5.7</td></tr><tr><td>Cow 3</td><td>3.6</td><td>4.6</td><td>10.9</td></tr><tr><td>Cow 4</td><td>4.2</td><td>4.5</td><td>10.5</td></tr><tr><td>Mean ± s.d.</td><td>3.1 ± 1.0</td><td>3.8 ± 0.9</td><td>8.9 ± 2.4</td></tr></tbody></table>

https://doi.org/10.1371/journal.pone.0298930.t006

<table><caption>**Table 7.** Estimated yield parameters for the *in vivo* data.</caption><thead><tr><th></th><th>Y<sub>ac</sub></th><th>Y<sub>bu</sub></th><th>Y<sub>pr</sub></th></tr></thead><tbody><tr><td>Cow 1</td><td>1.38x10<sup>-6</sup></td><td>1.73x10<sup>-7</sup></td><td>2.37x10<sup>-6</sup></td></tr><tr><td>Cow 2</td><td>1.52x10<sup>-6</sup></td><td>1.55x10<sup>-7</sup></td><td>2.16x10<sup>-6</sup></td></tr><tr><td>Cow 3</td><td>1.43x10<sup>-6</sup></td><td>1.50x10<sup>-7</sup></td><td>2.64x10<sup>-6</sup></td></tr><tr><td>Cow 4</td><td>1.39x10<sup>-6</sup></td><td>1.83x10<sup>-7</sup></td><td>2.54x10<sup>-6</sup></td></tr><tr><td>Mean ± s.d.</td><td>1.43x10<sup>-6</sup> ± 6.28x10<sup>-8</sup></td><td>1.62x10<sup>-7</sup> ± 1.56x10<sup>-8</sup></td><td>2.43x10<sup>-6</sup> ± 2.10x10<sup>-7</sup></td></tr></tbody></table>

Y<sub>ac</sub>, Y<sub>bu</sub>, Y<sub>pr</sub>: yields of acetate, butyrate and propionate production (mol VFA / microbial functional proxy abundance).

https://doi.org/10.1371/journal.pone.0298930.t007

adequate representation of VFA dynamics. Although the results are promising, future research is needed to refine the definition of these proxies, by taking into account quantitative tran scriptomic information to identify the actual metabolic activity.

Our study has taken a step in the direction of capitalizing on microbial data to predict rumen function. In this first step, we used asymptotic state observers due to its simple struc ture and the advantage of not requiring information on the mathematical functions describing the metabolic kinetic rates. However, asymptotic observers have the limitation that the conver gence rate cannot be tuned because it is determined by operational conditions (e.g., ruminal fractional passage rate). This limitation can be overcome by other observers. However, more sophisticated technicalities are needed for the design and implementation of such observers [26]. Our approach was applied to estimate VFA dynamics. We think that the approach can be extended to estimate other metabolites including methane production by the rumen methano- genic archaea. Such an application will be of great interesting for the monitoring of methane emissions and for the evaluation of methane inhibition strategies.

## Conclusions

Existing mechanistic models of rumen fermentation consider an aggregated representation of the rumen microbiota and its metabolic function. However, none of these models integrate microbial genomic knowledge and thus do not capitalize on the rich information of microbial genomic sequencing. In this work, we integrated microbial time series of the rumen micro biota determined by 16S rDNA into a dynamic model of the rumen microbiome. Our results showed the estimated VFA concentrations converge towards the real VFA concentration dynamics demonstrating the promising potential of our approach.

## Acknowledgments

Rafael Muñoz-Tamayo thanks Pablo Ugalde-Salas (Inria, France) and Jerôme Harmand (INRAE, France) for their helpful explanations on the use of state observers for microbial
ecosystems. The Finnish Functional Genomics Centre supported by the University of Turku, Åbo Akademi University, and Biocenter Finland is acknowledged for sequencing. The authors wish to acknowledge CSC-IT Center for Science, Finland, for computational resources.

## Author Contributions

**Conceptualization:** Mohsen Davoudkhani, Francesco Rubino, Christopher J. Creevey, Seppo Ahvenjärvi, Ali R. Bayat, Ilma Tapio, Alejandro Belanche, Rafael Muñoz-Tamayo.

**Data curation:** Francesco Rubino, Ilma Tapio, Alejandro Belanche.

**Formal analysis:** Mohsen Davoudkhani, Francesco Rubino, Rafael Muñoz-Tamayo.

**Investigation:** Seppo Ahvenjärvi, Ali R. Bayat, Ilma Tapio, Alejandro Belanche.

**Methodology:** Mohsen Davoudkhani, Francesco Rubino, Christopher J. Creevey, Rafael Muñoz-Tamayo.

**Project administration:** Rafael Muñoz-Tamayo.

**Software:** Mohsen Davoudkhani, Francesco Rubino, Rafael Muñoz-Tamayo.

**Supervision:** Christopher J. Creevey, Rafael Muñoz-Tamayo.

**Writing – original draft:** Mohsen Davoudkhani, Francesco Rubino, Rafael Muñoz-Tamayo.

**Writing – review & editing:** Christopher J. Creevey, Seppo Ahvenjärvi, Ali R. Bayat, Ilma Tapio, Alejandro Belanche.

## References

1. Wallace JR, Sasson G, Garnsworthy PC, Tapio I, Gregson E, Bani P, et al. A heritable subset of the core rumen microbiome dictates dairy cow productivity and emissions. Sci Adv. 2019; 5: eaav8391. https://doi.org/10.1126/sciadv.aav8391 PMID: 31281883

2. Huws SA, Creevey CJ, Oyama LB, Mizrahi I, Denman SE, Popova M, et al. Addressing global ruminant agricultural challenges through understanding the rumen microbiome: past, present, and future. Front Microbiol. 2018; 9: 2161. https://doi.org/10.3389/fmicb.2018.02161 PMID: 30319557

3. Gruninger RJ, Ribeiro GO, Cameron A, McAllister TA. Invited review: Application of meta-omics to understand the dynamic nature of the rumen microbiome and how it responds to diet in ruminants. Animal. 2019; 13: 1843–1854. https://doi.org/10.1017/S1751731119000752 PMID: 31062682

4. Edwards JE, Kingston-smith AH, Jimenez HR, Huws SA, Skøt KP. Dynamics of initial colonization of nonconserved perennial ryegrass by anaerobic fungi in the bovine rumen Dynamics of initial colonization of nonconserved perennial ryegrass by anaerobic fungi in the bovine rumen. 2008. https://doi.org/10.1111/j.1574-6941.2008.00563.x PMID: 18673390

5. Belanche A, Newbold CJ, Lin W, Rees Stevens P, Kingston-Smith AH. A systems biology approach reveals differences in the dynamics of colonization and degradation of grass vs. hay by rumen microbes with minor effects of vitamin E supplementation. Front Microbiol. 2017;8. https://doi.org/10.3389/fmicb.2017.01456 PMID: 28824585

6. Huws SA, Edwards JE, Lin W, Rubino F, Alston M, Swarbreck D, et al. Microbiomes attached to fresh perennial ryegrass are temporally resilient and adapt to changing ecological niches. Microbiome. 2021; 9: 1–17. https://doi.org/10.1186/S40168-021-01087-W/FIGURES/8

7. Li RW, Wu S, Baldwin VI RL, Li W, Li C. Perturbation Dynamics of the Rumen Microbiota in Response to Exogenous Butyrate. PLoS One. 2012; 7: e29392. https://doi.org/10.1371/journal.pone.0029392 PMID: 22253719

8. Petri RM, Pourazad P, Khiaosa-ard R, Klevenhusen F, Metzler-Zebeli BU, Zebeli Q. Temporal dynamics of in-situ fiber-adherent bacterial community under ruminal acidotic conditions determined by 16S rRNA gene profiling. PLoS One. 2017; 12: e0182271. https://doi.org/10.1371/journal.pone.0182271 PMID: 28763489

9. Pitta DW, Melgar A, Hristov AN, Indugu N, Narayan KS, Pappalardo C, et al. Temporal changes in total and metabolically active ruminal methanogens in dairy cows supplemented with 3-nitrooxypropanol. J Dairy Sci. 2021; 104: 8721–8735. https://doi.org/10.3168/jds.2020-19862 PMID: 34024597
ecosystems. The Finnish Functional Genomics Centre supported by the University of Turku, Åbo Akademi University, and Biocenter Finland is acknowledged for sequencing. The authors wish to acknowledge CSC-IT Center for Science, Finland, for computational resources.

## Author Contributions

**Conceptualization:** Mohsen Davoudkhani, Francesco Rubino, Christopher J. Creevey, Seppo Ahvenjärvi, Ali R. Bayat, Ilma Tapio, Alejandro Belanche, Rafael Muñoz-Tamayo.

**Data curation:** Francesco Rubino, Ilma Tapio, Alejandro Belanche.

**Formal analysis:** Mohsen Davoudkhani, Francesco Rubino, Rafael Muñoz-Tamayo.

**Investigation:** Seppo Ahvenjärvi, Ali R. Bayat, Ilma Tapio, Alejandro Belanche.

**Methodology:** Mohsen Davoudkhani, Francesco Rubino, Christopher J. Creevey, Rafael Muñoz-Tamayo.

**Project administration:** Rafael Muñoz-Tamayo.

**Software:** Mohsen Davoudkhani, Francesco Rubino, Rafael Muñoz-Tamayo.

**Supervision:** Christopher J. Creevey, Rafael Muñoz-Tamayo.

**Writing – original draft:** Mohsen Davoudkhani, Francesco Rubino, Rafael Muñoz-Tamayo.

**Writing – review & editing:** Christopher J. Creevey, Seppo Ahvenjärvi, Ali R. Bayat, Ilma Tapio, Alejandro Belanche.

## References

1. Wallace JR, Sasson G, Garnsworthy PC, Tapio I, Gregson E, Bani P, et al. A heritable subset of the core rumen microbiome dictates dairy cow productivity and emissions. Sci Adv. 2019; 5: eaav8391. https://doi.org/10.1126/sciadv.aav8391 PMID: 31281883

2. Huws SA, Creevey CJ, Oyama LB, Mizrahi I, Denman SE, Popova M, et al. Addressing global ruminant agricultural challenges through understanding the rumen microbiome: past, present, and future. Front Microbiol. 2018; 9: 2161. https://doi.org/10.3389/fmicb.2018.02161 PMID: 30319557

3. Gruninger RJ, Ribeiro GO, Cameron A, McAllister TA. Invited review: Application of meta-omics to understand the dynamic nature of the rumen microbiome and how it responds to diet in ruminants. Animal. 2019; 13: 1843–1854. https://doi.org/10.1017/S1751731119000752 PMID: 31062682

4. Edwards JE, Kingston-smith AH, Jimenez HR, Huws SA, Skøt KP. Dynamics of initial colonization of nonconserved perennial ryegrass by anaerobic fungi in the bovine rumen Dynamics of initial colonization of nonconserved perennial ryegrass by anaerobic fungi in the bovine rumen. 2008. https://doi.org/10.1111/j.1574-6941.2008.00563.x PMID: 18673390

5. Belanche A, Newbold CJ, Lin W, Rees Stevens P, Kingston-Smith AH. A systems biology approach reveals differences in the dynamics of colonization and degradation of grass vs. hay by rumen microbes with minor effects of vitamin E supplementation. Front Microbiol. 2017;8. https://doi.org/10.3389/fmicb.2017.01456 PMID: 28824585

6. Huws SA, Edwards JE, Lin W, Rubino F, Alston M, Swarbreck D, et al. Microbiomes attached to fresh perennial ryegrass are temporally resilient and adapt to changing ecological niches. Microbiome. 2021; 9: 1–17. https://doi.org/10.1186/S40168-021-01087-W/FIGURES/8

7. Li RW, Wu S, Baldwin VI RL, Li W, Li C. Perturbation Dynamics of the Rumen Microbiota in Response to Exogenous Butyrate. PLoS One. 2012; 7: e29392. https://doi.org/10.1371/journal.pone.0029392 PMID: 22253719

8. Petri RM, Pourazad P, Khiaosa-ard R, Klevenhusen F, Metzler-Zebeli BU, Zebeli Q. Temporal dynamics of in-situ fiber-adherent bacterial community under ruminal acidotic conditions determined by 16S rRNA gene profiling. PLoS One. 2017; 12: e0182271. https://doi.org/10.1371/journal.pone.0182271 PMID: 28763489

9. Pitta DW, Melgar A, Hristov AN, Indugu N, Narayan KS, Pappalardo C, et al. Temporal changes in total and metabolically active ruminal methanogens in dairy cows supplemented with 3-nitrooxypropanol. J Dairy Sci. 2021; 104: 8721–8735. https://doi.org/10.3168/jds.2020-19862 PMID: 34024597
10. Hennessy ML, Indugu N, Vecchiarelli B, Bender J, Pappalardo C, Leibstein M, et al. Temporal changes in the fecal bacterial community in Holstein dairy calves from birth through the transition to a solid diet. PLoS One. 2020; 15: e0238882. https://doi.org/10.1371/journal.pone.0238882 PMID: 32898158

11. Huuki H, Ahvenjärvi S, Lidauer P, Popova M, Vilkki J, Vanhatalo A, et al. Fresh Rumen Liquid Inoculant Enhances the Rumen Microbial Community Establishment in Pre-weaned Dairy Calves. Front Microbiol. 2022; 12: 758395. https://doi.org/10.3389/fmicb.2021.758395 PMID: 35095788

12. Muñoz-Tamayo R, Giger-Reverdin S, Sauvant D. Mechanistic modelling of in vitro fermentation and methane production by rumen microbiota. Anim Feed Sci Technol. 2016; 220: 1–21. https://doi.org/10.1016/j.anifeedsci.2016.07.005

13. Muñoz-Tamayo R, Chagas JC, Ramin M, Krizsan SJ. Modelling the impact of the macroalgae Asparagopsis taxiformis on rumen microbial fermentation and methane production. Peer Community J. 2021; 1: e7. https://doi.org/10.24072/PCJOURNAL.11

14. Huhtanen P, Ramin M, Udén P. Nordic dairy cow model Karoline in predicting methane emissions: 1. Model description and sensitivity analysis. Livest Sci. 2015; 178: 71–80. https://doi.org/10.1016/j.livsci.2015.05.009

15. Gregorini P, Beukes P, Waghorn G, Pacheco D, Hanigan M. Development of an improved representation of rumen digesta outflow in a mechanistic and dynamic model of a dairy cow, Molly. Ecol Modell. 2015; 313: 293–306. https://doi.org/10.1016/j.ecolmodel.2015.06.042

16. van Lingen HJ, Fadel JG, Moraes LE, Bannink A, Dijkstra J. Bayesian mechanistic modeling of thermodynamically controlled volatile fatty acid, hydrogen and methane production in the bovine rumen. J Theor Biol. 2019; 480: 150–165. https://doi.org/10.1016/j.jtbi.2019.08.008 PMID: 31401059

17. Muñoz-Tamayo R, Davoudkhani M, Fakih I, Robles-Rodriguez CE, Rubino F, Creevey CJ, et al. Review: Towards the next-generation models of the rumen microbiome for enhancing predictive power and guiding sustainable production strategies. animal. 2023; 17: 100984. https://doi.org/10.1016/j.animal.2023.100984 PMID: 37821326

18. Pereira B, Miguel J, Vilaça P, Soares S, Rocha I, Carneiro S. Reconstruction of a genome-scale metabolic model for Actinobacillus succinogenes 130Z. BMC Syst Biol. 2018; 12: 61. https://doi.org/10.1186/s12918-018-0585-7 PMID: 29843739

19. Lee NR, Lee CH, Lee DY, Park JB. Genome-scale metabolic network reconstruction and in silico analysis of hexanoic acid producing Megasphaera elsdenii. Microorg 2020, Vol 8, Page 539. 2020; 8: 539. https://doi.org/10.3390/MICROORGANISMS8040539 PMID: 32283671

20. Fakih I, Got J, Robles-Rodriguez CE, Siegel A, Forano E, Muñoz-Tamayo R. Dynamic genome-based metabolic modeling of the predominant cellulolytic rumen bacterium Fibrobacter succinogenes S85. mSystems. 2023; 8: e01027-22. https://doi.org/10.1128/msystems.01027-22 PMID: 37289026

21. Islam MM, Fernando SC, Saha R. Metabolic modeling elucidates the transactions in the rumen microbiome and the shifts upon virome interactions. Front Microbiol. 2019; 10: 2412. https://doi.org/10.3389/fmicb.2019.02412 PMID: 31866953

22. Faust K, Lahti L, Gonze D, de Vos WM, Raes J. Metagenomics meets time series analysis: unraveling microbial community dynamics. Curr Opin Microbiol. 2015; 25: 56–66. https://doi.org/10.1016/j.mib.2015.04.004 PMID: 26005845

23. Stein RR, Bucci V, Toussaint NC, Buffie CG, Rätsch G, Pamer EG, et al. Ecological Modeling from Time-Series Inference: Insight into Dynamics and Stability of Intestinal Microbiota. PLoS Comput Biol. 2013; 9: e31–36. https://doi.org/10.1371/journal.pcbi.1003388 PMID: 24348232

24. Gonze D, Coyte KZ, Lahti L, Faust K. Microbial communities as dynamical systems. Current Opinion in Microbiology. Elsevier Ltd; 2018. pp. 41–49. https://doi.org/10.1016/j.mib.2018.07.004 PMID: 30041083

25. Dumont M, Godon JJ, Harmand J. Species coexistence in nitrifying chemostats: A model of microbial interactions. Processes. 2016; 4: pr4040051. https://doi.org/10.3390/pr4040051

26. Dochain D. State and parameter estimation in chemical and biochemical processes: a tutorial. J Process Control. 2003; 13: 801–818. https://doi.org/10.1016/S0959-1524(03)00026-X

27. Selişteanu D, Tebbani S, Roman M, Petre E, Georgeanu V. Microbial production of enzymes: Nonlinear state and kinetic reaction rates estimation. Biochem Eng J. 2014; 91: 23–36. https://doi.org/10.1016/J.BEJ.2014.07.010

28. Aceves-Lara CA, Latrille E, Steyer JP. Optimal control of hydrogen production in a continuous anaerobic fermentation bioreactor. Int J Hydrogen Energy. 2010; 35: 10710–10718. https://doi.org/10.1016/J.IJHYDENE.2010.02.110

29. Ugalde-Salas P, Harmand J, Quemener ED Le. Asymptotic observers and integer programming for functional classification of a microbial community in a chemostat. 2019 18th European Control
10. Hennessy ML, Indugu N, Vecchiarelli B, Bender J, Pappalardo C, Leibstein M, et al. Temporal changes in the fecal bacterial community in Holstein dairy calves from birth through the transition to a solid diet. PLoS One. 2020; 15: e0238882. https://doi.org/10.1371/journal.pone.0238882 PMID: 32898158

11. Huuki H, Ahvenjärvi S, Lidauer P, Popova M, Vilkki J, Vanhatalo A, et al. Fresh Rumen Liquid Inoculant Enhances the Rumen Microbial Community Establishment in Pre-weaned Dairy Calves. Front Microbiol. 2022; 12: 758395. https://doi.org/10.3389/fmicb.2021.758395 PMID: 35095788

12. Muñoz-Tamayo R, Giger-Reverdin S, Sauvant D. Mechanistic modelling of in vitro fermentation and methane production by rumen microbiota. Anim Feed Sci Technol. 2016; 220: 1–21. https://doi.org/10.1016/j.anifeedsci.2016.07.005

13. Muñoz-Tamayo R, Chagas JC, Ramin M, Krizsan SJ. Modelling the impact of the macroalgae Asparagopsis taxiformis on rumen microbial fermentation and methane production. Peer Community J. 2021; 1: e7. https://doi.org/10.24072/PCJOURNAL.11

14. Huhtanen P, Ramin M, Udén P. Nordic dairy cow model Karoline in predicting methane emissions: 1. Model description and sensitivity analysis. Livest Sci. 2015; 178: 71–80. https://doi.org/10.1016/j.livsci.2015.05.009

15. Gregorini P, Beukes P, Waghorn G, Pacheco D, Hanigan M. Development of an improved representation of rumen digesta outflow in a mechanistic and dynamic model of a dairy cow, Molly. Ecol Modell. 2015; 313: 293–306. https://doi.org/10.1016/j.ecolmodel.2015.06.042

16. van Lingen HJ, Fadel JG, Moraes LE, Bannink A, Dijkstra J. Bayesian mechanistic modeling of thermodynamically controlled volatile fatty acid, hydrogen and methane production in the bovine rumen. J Theor Biol. 2019; 480: 150–165. https://doi.org/10.1016/j.jtbi.2019.08.008 PMID: 31401059

17. Muñoz-Tamayo R, Davoudkhani M, Fakih I, Robles-Rodriguez CE, Rubino F, Creevey CJ, et al. Review: Towards the next-generation models of the rumen microbiome for enhancing predictive power and guiding sustainable production strategies. animal. 2023; 17: 100984. https://doi.org/10.1016/j.animal.2023.100984 PMID: 37821326

18. Pereira B, Miguel J, Vilaça P, Soares S, Rocha I, Carneiro S. Reconstruction of a genome-scale metabolic model for Actinobacillus succinogenes 130Z. BMC Syst Biol. 2018; 12: 61. https://doi.org/10.1186/s12918-018-0585-7 PMID: 29843739

19. Lee NR, Lee CH, Lee DY, Park JB. Genome-scale metabolic network reconstruction and in silico analysis of hexanoic acid producing Megasphaera elsdenii. Microorg 2020, Vol 8, Page 539. 2020; 8: 539. https://doi.org/10.3390/MICROORGANISMS8040539 PMID: 32283671

20. Fakih I, Got J, Robles-Rodriguez CE, Siegel A, Forano E, Muñoz-Tamayo R. Dynamic genome-based metabolic modeling of the predominant cellulolytic rumen bacterium Fibrobacter succinogenes S85. mSystems. 2023; 8: e01027-22. https://doi.org/10.1128/msystems.01027-22 PMID: 37289026

21. Islam MM, Fernando SC, Saha R. Metabolic modeling elucidates the transactions in the rumen microbiome and the shifts upon virome interactions. Front Microbiol. 2019; 10: 2412. https://doi.org/10.3389/fmicb.2019.02412 PMID: 31866953

22. Faust K, Lahti L, Gonze D, de Vos WM, Raes J. Metagenomics meets time series analysis: unraveling microbial community dynamics. Curr Opin Microbiol. 2015; 25: 56–66. https://doi.org/10.1016/j.mib.2015.04.004 PMID: 26005845

23. Stein RR, Bucci V, Toussaint NC, Buffie CG, Rätsch G, Pamer EG, et al. Ecological Modeling from Time-Series Inference: Insight into Dynamics and Stability of Intestinal Microbiota. PLoS Comput Biol. 2013; 9: e31–36. https://doi.org/10.1371/journal.pcbi.1003388 PMID: 24348232

24. Gonze D, Coyte KZ, Lahti L, Faust K. Microbial communities as dynamical systems. Current Opinion in Microbiology. Elsevier Ltd; 2018. pp. 41–49. https://doi.org/10.1016/j.mib.2018.07.004 PMID: 30041083

25. Dumont M, Godon JJ, Harmand J. Species coexistence in nitrifying chemostats: A model of microbial interactions. Processes. 2016; 4: pr4040051. https://doi.org/10.3390/pr4040051

26. Dochain D. State and parameter estimation in chemical and biochemical processes: a tutorial. J Process Control. 2003; 13: 801–818. https://doi.org/10.1016/S0959-1524(03)00026-X

27. Selişteanu D, Tebbani S, Roman M, Petre E, Georgeanu V. Microbial production of enzymes: Nonlinear state and kinetic reaction rates estimation. Biochem Eng J. 2014; 91: 23–36. https://doi.org/10.1016/J.BEJ.2014.07.010

28. Aceves-Lara CA, Latrille E, Steyer JP. Optimal control of hydrogen production in a continuous anaerobic fermentation bioreactor. Int J Hydrogen Energy. 2010; 35: 10710–10718. https://doi.org/10.1016/J.IJHYDENE.2010.02.110

29. Ugalde-Salas P, Harmand J, Quemener ED Le. Asymptotic observers and integer programming for functional classification of a microbial community in a chemostat. 2019 18th European Control
Conference, ECC 2019. Institute of Electrical and Electronics Engineers Inc.; 2019. pp. 1665–1670.
https://doi.org/10.23919/ECC.2019.8795854

30. Dumont M, Harmand J, Rapaport A, Godon J. Towards functional molecular fingerprints. Environ Microbiol. 2009; 11: 1717–1727. https://doi.org/10.1111/j.1462-2920.2009.01898.x PMID: 19453611

31. Davoudkhani M, Rubino F, Creevey CJ, Belanche A, Muñoz-Tamayo R. Integration of microbial time series into a mechanistic model of the rumen microbiome under the RUSITEC condition. Anim—Sci Proc. 2022; 13: 572–573. https://doi.org/10.1016/J.ANSCIP.2022.07.443

32. Muñoz-Tamayo R, Nielsen BL, Gagaoua M, Gondret F, Krause ET, Morgavi DP, et al. Seven steps to enhance open science practices in animal science. PNAS Nexus. 2022; 1: pgac106. https://doi.org/10.1093/pnasexus/pgac106 PMID: 36741429

33. Bayat AR, Vilkki J, Razzaghi A, Leskinen H, Kettunen H, Khurana R, et al. Evaluating the effects of high-oil rapeseed cake or natural additives on methane emissions and performance of dairy cows. J Dairy Sci. 2022; 105: 1211–1224. https://doi.org/10.3168/jds.2021-20537 PMID: 34799103

34. Rius AG, Kittelmann S, Macdonald KA, Waghorn GC, Janssen PH, Sikkema E. Nitrogen metabolism and rumen microbial enumeration in lactating cows with divergent residual feed intake fed high-digestibility pasture. J Dairy Sci. 2012; 95: 5024–5034. https://doi.org/10.3168/jds.2012-5392 PMID: 22916906

35. Caporaso JG, Lauber CL, Walters WA, Berg-Lyons D, Lozupone CA, Turnbaugh PJ, et al. Global patterns of 16S rRNA diversity at a depth of millions of sequences per sample. Proc Natl Acad Sci U S A. 2011; 108: 4516–4522. https://doi.org/10.1073/pnas.1000080107 PMID: 20534432

36. Bernard O, Hadj-Sadok Z, Dochain D, Genovesi A, Steyer JP. Dynamical model development and parameter identification for an anaerobic wastewater treatment process. Biotechnol Bioeng. 2001; 75: 424–438. https://doi.org/10.1002/bit.10036 PMID: 11668442

37. Bolyen E, Rideout JR, Dillon MR, Bokulich NA, Abnet CC, Al-Ghalith GA, et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat Biotechnol 2019 378. 2019; 37: 852–857. https://doi.org/10.1038/s41587-019-0209-9 PMID: 31341288

38. Rognes T, Flouri T, Nichols B, Quince C, Mahé F. VSEARCH: a versatile open source tool for metagenomics. 2016; 1–22. https://doi.org/10.7717/peerj.2584 PMID: 27781170

39. Camacho C, Coulouris G, Avagyan V, Ma N, Papadopoulos J, Bealer K, et al. BLAST +: architecture and applications. 2009; 9: 1–9. https://doi.org/10.1186/1471-2105-10-421 PMID: 20003500

40. Quast C, Pruesse E, Yilmaz P, Gerken J, Schweer T, Glo FO, et al. The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. 2013; 41: 590–596. https://doi.org/10.1093/nar/gks1219 PMID: 23193283

41. Wilkinson TJ, Huws SA, Edwards JE, Kingston-Smith AH, Siu-Ting K, Hughes M, et al. CowPI: A rumen microbiome focussed version of the PICRUSt functional inference software. Front Microbiol. 2018; 9: 1095. https://doi.org/10.3389/fmicb.2018.01095 PMID: 29887853

42. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15. https://doi.org/10.1186/s13059-014-0550-8 PMID: 25516281

43. Furumichi M, Sato Y, Ishiguro-watanabe M. KEGG: integrating viruses and cellular organisms. 2021; 49: 545–551. https://doi.org/10.1093/nar/gkaa970 PMID: 33125081

44. Hackmann TJ, Firkins JL. Electron transport phosphorylation in rumen butyrvibrios: unprecedented ATP yield for glucose fermentation to butyrate. 2015; 6: 1–11. https://doi.org/10.3389/fmicb.2015.00622 PMID: 26157432

45. Mills JA, Dijkstra J, Bannink A, Cammell SB, Kebreab E, France J. A mechanistic model of whole-tract digestion and methanogenesis in the lactating dairy cow: model development, evaluation, and application The online version of this article, along with updated information and services, is located on the World Wide W. 2001; 1584–1597.

46. Muñoz-Tamayo R, Laroche B, Leclerc M, Walter E. IDEAS: A parameter identification toolbox with symbolic analysis of uncertainty and its application to biological modelling. IFAC Proceedings Volumes. 2009. pp. 1271–1276. https://doi.org/10.3182/20090706-3-FR-2004.0211

47. Lagarias JC, Reeds JA, Wright MH, Wright PE. Convergence Properties of the Nelder-Mead Simplex Method in Low Dimensions. SIAM J Optim. 1998; 9: 112–147.

48. Davoudkhani M, Rubino F, Creevey CJ, Ahvenjärvi S, Bayat AR, Tapio I, et al. Implementation of state observers linking rumen microbial abundance time series with fermentation products. 2023. https://doi.org/10.5281/ZENODO.8386786

49. Kittelmann S, Kirk MR, Jonker A, McCulloch A, Janssen PH. Buccal swabbing as a noninvasive method to determine bacterial, archaeological, and eukaryotic microbial community structures in the rumen. Appl Environ Microbiol. 2015; 81: 7470–7483. https://doi.org/10.1128/AEM.02385-15 PMID: 26276109
Conference, ECC 2019. Institute of Electrical and Electronics Engineers Inc.; 2019. pp. 1665–1670.
https://doi.org/10.23919/ECC.2019.8795854

30. Dumont M, Harmand J, Rapaport A, Godon J. Towards functional molecular fingerprints. Environ Microbiol. 2009; 11: 1717–1727. https://doi.org/10.1111/j.1462-2920.2009.01898.x PMID: 19453611

31. Davoudkhani M, Rubino F, Creevey CJ, Belanche A, Muñoz-Tamayo R. Integration of microbial time series into a mechanistic model of the rumen microbiome under the RUSITEC condition. Anim—Sci Proc. 2022; 13: 572–573. https://doi.org/10.1016/J.ANSCIP.2022.07.443

32. Muñoz-Tamayo R, Nielsen BL, Gagaoua M, Gondret F, Krause ET, Morgavi DP, et al. Seven steps to enhance open science practices in animal science. PNAS Nexus. 2022; 1: pgac106. https://doi.org/10.1093/pnasexus/pgac106 PMID: 36741429

33. Bayat AR, Vilkki J, Razzaghi A, Leskinen H, Kettunen H, Khurana R, et al. Evaluating the effects of high-oil rapeseed cake or natural additives on methane emissions and performance of dairy cows. J Dairy Sci. 2022; 105: 1211–1224. https://doi.org/10.3168/jds.2021-20537 PMID: 34799103

34. Rius AG, Kittelmann S, Macdonald KA, Waghorn GC, Janssen PH, Sikkema E. Nitrogen metabolism and rumen microbial enumeration in lactating cows with divergent residual feed intake fed high-digestibility pasture. J Dairy Sci. 2012; 95: 5024–5034. https://doi.org/10.3168/jds.2012-5392 PMID: 22916906

35. Caporaso JG, Lauber CL, Walters WA, Berg-Lyons D, Lozupone CA, Turnbaugh PJ, et al. Global patterns of 16S rRNA diversity at a depth of millions of sequences per sample. Proc Natl Acad Sci U S A. 2011; 108: 4516–4522. https://doi.org/10.1073/pnas.1000080107 PMID: 20534432

36. Bernard O, Hadj-Sadok Z, Dochain D, Genovesi A, Steyer JP. Dynamical model development and parameter identification for an anaerobic wastewater treatment process. Biotechnol Bioeng. 2001; 75: 424–438. https://doi.org/10.1002/bit.10036 PMID: 11668442

37. Bolyen E, Rideout JR, Dillon MR, Bokulich NA, Abnet CC, Al-Ghalith GA, et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat Biotechnol 2019 378. 2019; 37: 852–857. https://doi.org/10.1038/s41587-019-0209-9 PMID: 31341288

38. Rognes T, Flouri T, Nichols B, Quince C, Mahé F. VSEARCH: a versatile open source tool for metagenomics. 2016; 1–22. https://doi.org/10.7717/peerj.2584 PMID: 27781170

39. Camacho C, Coulouris G, Avagyan V, Ma N, Papadopoulos J, Bealer K, et al. BLAST +: architecture and applications. 2009; 9: 1–9. https://doi.org/10.1186/1471-2105-10-421 PMID: 20003500

40. Quast C, Pruesse E, Yilmaz P, Gerken J, Schweer T, Glo FO, et al. The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. 2013; 41: 590–596. https://doi.org/10.1093/nar/gks1219 PMID: 23193283

41. Wilkinson TJ, Huws SA, Edwards JE, Kingston-Smith AH, Siu-Ting K, Hughes M, et al. CowPI: A rumen microbiome focussed version of the PICRUSt functional inference software. Front Microbiol. 2018; 9: 1095. https://doi.org/10.3389/fmicb.2018.01095 PMID: 29887853

42. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15. https://doi.org/10.1186/s13059-014-0550-8 PMID: 25516281

43. Furumichi M, Sato Y, Ishiguro-watanabe M. KEGG: integrating viruses and cellular organisms. 2021; 49: 545–551. https://doi.org/10.1093/nar/gkaa970 PMID: 33125081

44. Hackmann TJ, Firkins JL. Electron transport phosphorylation in rumen butyrvibrios: unprecedented ATP yield for glucose fermentation to butyrate. 2015; 6: 1–11. https://doi.org/10.3389/fmicb.2015.00622 PMID: 26157432

45. Mills JA, Dijkstra J, Bannink A, Cammell SB, Kebreab E, France J. A mechanistic model of whole-tract digestion and methanogenesis in the lactating dairy cow: model development, evaluation, and application The online version of this article, along with updated information and services, is located on the World Wide W. 2001; 1584–1597.

46. Muñoz-Tamayo R, Laroche B, Leclerc M, Walter E. IDEAS: A parameter identification toolbox with symbolic analysis of uncertainty and its application to biological modelling. IFAC Proceedings Volumes. 2009. pp. 1271–1276. https://doi.org/10.3182/20090706-3-FR-2004.0211

47. Lagarias JC, Reeds JA, Wright MH, Wright PE. Convergence Properties of the Nelder-Mead Simplex Method in Low Dimensions. SIAM J Optim. 1998; 9: 112–147.

48. Davoudkhani M, Rubino F, Creevey CJ, Ahvenjärvi S, Bayat AR, Tapio I, et al. Implementation of state observers linking rumen microbial abundance time series with fermentation products. 2023. https://doi.org/10.5281/ZENODO.8386786

49. Kittelmann S, Kirk MR, Jonker A, McCulloch A, Janssen PH. Buccal swabbing as a noninvasive method to determine bacterial, archaeological, and eukaryotic microbial community structures in the rumen. Appl Environ Microbiol. 2015; 81: 7470–7483. https://doi.org/10.1128/AEM.02385-15 PMID: 26276109
50. Tapio I, Shingfield KJ, McKain N, Bonin A, Fischer D, Bayat AR, et al. Oral samples as non-invasive proxies for assessing the composition of the rumen microbial community. PLoS One. 2016; 11: e0151220. https://doi.org/10.1371/journal.pone.0151220 PMID: 26986467

51. Walker AW, Duncan SH, Louis P, Flint HJ. Phylogeny, culturing, and metagenomics of the human gut microbiota. Trends Microbiol. 2014; 22: 267-74. https://doi.org/10.1016/j.tim.2014.03.001 PMID: 24698744
50. Tapio I, Shingfield KJ, McKain N, Bonin A, Fischer D, Bayat AR, et al. Oral samples as non-invasive proxies for assessing the composition of the rumen microbial community. PLoS One. 2016; 11: e0151220. https://doi.org/10.1371/journal.pone.0151220 PMID: 26986467

51. Walker AW, Duncan SH, Louis P, Flint HJ. Phylogeny, culturing, and metagenomics of the human gut microbiota. Trends Microbiol. 2014; 22: 267-74. https://doi.org/10.1016/j.tim.2014.03.001 PMID: 24698744