Stochastic SIR model predicts the evolution of COVID-19 epidemics from public health and wastewater data in small and medium-sized municipalities: A one year study.

The level of unpredictability of the COVID-19 pandemics poses a challenge to effectively model its dynamic evolution. In this study we incorporate the inherent stochasticity of the SARS-CoV-2 virus spread by reinterpreting the classical compartmental models of infectious diseases (SIR type) as chemical reaction systems modeled via the Chemical Master Equation and solved by Monte Carlo Methods. Our model predicts the evolution of the pandemics at the level of municipalities, incorporating for the first time (i) a variable infection rate to capture the effect of mitigation policies on the dynamic evolution of the pandemics (ii) SIR-with-jumps taking into account the possibility of multiple infections from a single infected person and (iii) data of viral load quantified by RT-qPCR from samples taken from Wastewater Treatment Plants. The model has been successfully employed for the prediction of the COVID-19 pandemics evolution in small and medium size municipalities of Galicia (Northwest of Spain).

SELANSI: a toolbox for simulation of stochastic gene regulatory networks.

Motivation: Gene regulation is inherently stochastic. In many applications concerning Systems and Synthetic Biology such as the reverse engineering and the de novo design of genetic circuits, stochastic effects (yet potentially crucial) are often neglected due to the high computational cost of stochastic simulations. With advances in these fields there is an increasing need of tools providing accurate approximations of the stochastic dynamics of gene regulatory networks (GRNs) with reduced computational effort. Results: This work presents SELANSI (SEmi-LAgrangian SImulation of GRNs), a software toolbox for the simulation of stochastic multidimensional gene regulatory networks. SELANSI exploits intrinsic structural properties of gene regulatory networks to accurately approximate the corresponding Chemical Master Equation with a partial integral differential equation that is solved by a semi-lagrangian method with high efficiency. Networks under consideration might involve multiple genes with self and cross regulations, in which genes can be regulated by different transcription factors. Moreover, the validity of the method is not restricted to a particular type of kinetics. The tool offers total flexibility regarding network topology, kinetics and parameterization, as well as simulation options. Availability and implementation: SELANSI runs under the MATLAB environment, and is available under GPLv3 license at https://sites.google.com/view/selansi. Contact: antonio@iim.csic.es.

Toward predictive food process models: A protocol for parameter estimation.

Mathematical models, in particular, physics-based models, are essential tools to food product and process design, optimization and control. The success of mathematical models relies on their predictive capabilities. However, describing physical, chemical and biological changes in food processing requires the values of some, typically unknown, parameters. Therefore, parameter estimation from experimental data is critical to achieving desired model predictive properties. This work takes a new look into the parameter estimation (or identification) problem in food process modeling. First, we examine common pitfalls such as lack of identifiability and multimodality. Second, we present the theoretical background of a parameter identification protocol intended to deal with those challenges. And, to finish, we illustrate the performance of the proposed protocol with an example related to the thermal processing of packaged foods.

Stochastic modeling and numerical simulation of gene regulatory networks with protein bursting.

Gene expression is inherently stochastic. Advanced single-cell microscopy techniques together with mathematical models for single gene expression led to important insights in elucidating the sources of intrinsic noise in prokaryotic and eukaryotic cells. In addition to the finite size effects due to low copy numbers, translational bursting is a dominant source of stochasticity in cell scenarios involving few short lived mRNA transcripts with high translational efficiency (as is typically the case for prokaryotes), causing protein synthesis to occur in random bursts. In the context of gene regulation cascades, the Chemical Master Equation (CME) governing gene expression has in general no closed form solution, and the accurate stochastic simulation of the dynamics of complex gene regulatory networks is a major computational challenge. The CME associated to a single gene self regulatory motif has been previously approximated by a one dimensional time dependent partial integral differential equation (PIDE). However, to the best of our knowledge, multidimensional versions for such PIDE have not been developed yet. Here we propose a multidimensional PIDE model for regulatory networks involving multiple genes with self and cross regulations (in which genes can be regulated by different transcription factors) derived as the continuous counterpart of a CME with jump process. The model offers a reliable description of systems with translational bursting. In order to provide an efficient numerical solution, we develop a semilagrangian method to discretize the differential part of the PIDE, combined with a composed trapezoidal quadrature formula to approximate the integral term. We apply the model and numerical method to study sustained stochastic oscillations and the development of competence, a particular case of transient differentiation attained by certain bacterial cells under stress conditions. We found that the resulting probability distributions are distinguishable from those characteristic of other transient differentiation processes. In this way, they can be employed as markers or signatures that identify such phenomena from bacterial population experimental data, for instance. The computational efficiency of the semilagrangian method makes it suitable for purposes like model identification and parameter estimation from experimental data or, in combination with optimization routines, the design of gene regulatory networks under molecular noise.

Modeling bacterial population growth from stochastic single-cell dynamics.

A few bacterial cells may be sufficient to produce a food-borne illness outbreak, provided that they are capable of adapting and proliferating on a food matrix. This is why any quantitative health risk assessment policy must incorporate methods to accurately predict the growth of bacterial populations from a small number of pathogens. In this aim, mathematical models have become a powerful tool. Unfortunately, at low cell concentrations, standard deterministic models fail to predict the fate of the population, essentially because the heterogeneity between individuals becomes relevant. In this work, a stochastic differential equation (SDE) model is proposed to describe variability within single-cell growth and division and to simulate population growth from a given initial number of individuals. We provide evidence of the model ability to explain the observed distributions of times to division, including the lag time produced by the adaptation to the environment, by comparing model predictions with experiments from the literature for Escherichia coli, Listeria innocua, and Salmonella enterica. The model is shown to accurately predict experimental growth population dynamics for both small and large microbial populations. The use of stochastic models for the estimation of parameters to successfully fit experimental data is a particularly challenging problem. For instance, if Monte Carlo methods are employed to model the required distributions of times to division, the parameter estimation problem can become numerically intractable. We overcame this limitation by converting the stochastic description to a partial differential equation (backward Kolmogorov) instead, which relates to the distribution of division times. Contrary to previous stochastic formulations based on random parameters, the present model is capable of explaining the variability observed in populations that result from the growth of a small number of initial cells as well as the lack of it compared to populations initiated by a larger number of individuals, where the random effects become negligible.

CRNreals: a toolbox for distinguishability and identifiability analysis of biochemical reaction networks.

UNLABELLED: Chemical reaction network theory is widely used in modeling and analyzing complex biochemical systems such as metabolic networks and cell signalling pathways. Being able to produce all the biologically and chemically important qualitative dynamical features, chemical reaction networks (CRNs) have attracted significant attention in the systems biology community. It is well-known that the reliable inference of CRN models generally requires thorough identifiability and distinguishability analysis together with carefully selected prior modeling assumptions. Here, we present a software toolbox CRNreals that supports the distinguishability and identifiability analysis of CRN models using recently published optimization-based procedures. AVAILABILITY AND IMPLEMENTATION: The CRNreals toolbox and the associated documentation are available at http://www.iim.csic.es/~gingproc/CRNreals/. The toolbox runs under the popular MATLAB computational environment and supports several free and commercial linear programming and mixed integer linear programming solvers.

CO2 budget of cultured mussels metabolism in the highly productive Northwest Iberian upwelling system.

Assessing the carbon footprint of marine bivalve aquaculture demands an accurate estimation of the CO2 release associated to capital goods and aquaculture operations but also to the metabolic CO2 budget of the cultured species. Nowadays, there are discrepancies on the processes to include in that budget, how to estimate them, and which scale should be applied, from individual to ecosystem. Site-specific environmental conditions and culture methods also affect significantly the estimates. Here, we have gathered environmental, biochemical and metabolic data from published scientific articles, reports and existing databases to present the metabolic CO2 budget for mussel aquaculture in the coastal inlets of the Northwest Iberian upwelling. We analyse the contribution of mussel flesh and shell production jointly and separately. At the individual scale, the shell CO2 budget is estimated from CO2 removal by shell matrix protein synthesis and CO2 release during calcification and respiration to support shell maintenance. Organic carbon in mussel flesh and CO2 released by respiration to support flesh maintenance contribute to the flesh CO2 budget. Only calcification and respiration processes are considered when estimating the metabolic carbon footprint of individual mussels because organic carbon in mussel flesh and shell returns to the atmosphere as CO2 in a relatively short period. While the metabolic carbon footprint associated to mussel shell remains constant at 365 kg CO2 per ton of shell, it varies from 92 to 578 kg CO2 per ton of mussel flesh. This large variability depends on mussel seeding time and harvesting size, due to the differential seasonal growth patterns of flesh and shell. Inclusion of the CO2 potentially immobilised in mussel faeces buried in the sediments would lead to a reduction of the metabolic carbon footprint estimates by up to 6 % compared with the individual estimates.

Wastewater and marine bioindicators surveillance to anticipate COVID-19 prevalence and to explore SARS-CoV-2 diversity by next generation sequencing: One-year study.

This study presents the results of SARS-CoV-2 surveillance in sewage water of 11 municipalities and marine bioindicators in Galicia (NW of Spain) from May 2020 to May 2021. An integrated pipeline was developed including sampling, pre-treatment and biomarker quantification, RNA detection, SARS-CoV-2 sequencing, mechanistic mathematical modeling and forecasting. The viral load in the inlet stream to the wastewater treatment plants (WWTP) was used to detect new outbreaks of COVID-19, and the data of viral load in the wastewater in combination with data provided by the health system was used to predict the evolution of the pandemic in the municipalities under study within a time horizon of 7 days. Moreover, the study shows that the viral load was eliminated from the treated sewage water in the WWTP, mainly in the biological reactors and the disinfection system. As a result, we detected a minor impact of the virus in the marine environment through the analysis of seawater, marine sediments and, wild and aquacultured mussels in the final discharge point of the WWTP.

Corrigendum: Stochastic Individual-Based Modeling of Bacterial Growth and Division Using Flow Cytometry.

[This corrects the article on p. 2626 in vol. 8, PMID: 29354110.].

Stochastic Individual-Based Modeling of Bacterial Growth and Division Using Flow Cytometry.

A realistic description of the variability in bacterial growth and division is critical to produce reliable predictions of safety risks along the food chain. Individual-based modeling of bacteria provides the theoretical framework to deal with this variability, but it requires information about the individual behavior of bacteria inside populations. In this work, we overcome this problem by estimating the individual behavior of bacteria from population statistics obtained with flow cytometry. For this objective, a stochastic individual-based modeling framework is defined based on standard assumptions during division and exponential growth. The unknown single-cell parameters required for running the individual-based modeling simulations, such as cell size growth rate, are estimated from the flow cytometry data. Instead of using directly the individual-based model, we make use of a modified Fokker-Plank equation. This only equation simulates the population statistics in function of the unknown single-cell parameters. We test the validity of the approach by modeling the growth and division of Pediococcus acidilactici within the exponential phase. Estimations reveal the statistics of cell growth and division using only data from flow cytometry at a given time. From the relationship between the mother and daughter volumes, we also predict that P. acidilactici divide into two successive parallel planes.

Pollutant levels in discarded fish species by Spanish trawlers operating in the Great Sole Bank and the Atlantic coast of the Iberian Peninsula.

Organic and inorganic pollutant levels were determined for the most discarded species from trawlers operating in Great Sole and Spanish coastal fishing grounds. Results for heavy metals indicated that Cd can reach values higher than legal limits for some species and tissues, while Hg and Pb concentrations are below established values. No significant variation was noticed with fishing grounds, but both season influences in the case of Pb and interspecies variation for Hg and Cd have been detected. Valorization recommendations could be therefore established according to the levels found in the different species.

Dynamic optimization of bioprocesses: efficient and robust numerical strategies.

The dynamic optimization (open loop optimal control) of non-linear bioprocesses is considered in this contribution. These processes can be described by sets of non-linear differential and algebraic equations (DAEs), usually subject to constraints in the state and control variables. A review of the available solution techniques for this class of problems is presented, highlighting the numerical difficulties arising from the non-linear, constrained and often discontinuous nature of these systems. In order to surmount these difficulties, we present several alternative stochastic and hybrid techniques based on the control vector parameterization (CVP) approach. The CVP approach is a direct method which transforms the original problem into a non-linear programming (NLP) problem, which must be solved by a suitable (efficient and robust) solver. In particular, a hybrid technique uses a first global optimization phase followed by a fast second phase based on a local deterministic method, so it can handle the nonconvexity of many of these NLPs. The efficiency and robustness of these techniques is illustrated by solving several challenging case studies regarding the optimal control of fed-batch bioreactors and other bioprocesses. In order to fairly evaluate their advantages, a careful and critical comparison with several other direct approaches is provided. The results indicate that the two-phase hybrid approach presents the best compromise between robustness and efficiency.

Shaping protein distributions in stochastic self-regulated gene expression networks.

In this work, we study connections between dynamic behavior and network parameters, for self-regulatory networks. To that aim, a method to compute the regions in the space of parameters that sustain bimodal or binary protein distributions has been developed. Such regions are indicative of stochastic dynamics manifested either as transitions between absence and presence of protein or between two positive protein levels. The method is based on the continuous approximation of the chemical master equation, unlike other approaches that make use of a deterministic description, which as will be shown can be misleading. We find that bimodal behavior is a ubiquitous phenomenon in cooperative gene expression networks under positive feedback. It appears for any range of transcription and translation rate constants whenever leakage remains below a critical threshold. Above such a threshold, the region in the parameters space which sustains bimodality persists, although restricted to low transcription and high translation rate constants. Remarkably, such a threshold is independent of the transcription or translation rates or the proportion of an active or inactive promoter and depends only on the level of cooperativity. The proposed method can be employed to identify bimodal or binary distributions leading to stochastic dynamics with specific switching properties, by searching inside the parameter regions that sustain such behavior.

Valorisation of fish by-products against waste management treatments--Comparison of environmental impacts.

Reuse and valorisation of fish by-products is a key process for marine resources conservation. Usually, fishmeal and oil processing factories collect the by-products generated by fishing port and industry processing activities, producing an economical benefit to both parts. In the same way, different added-value products can be recovered by the valorisation industries whereas fishing companies save the costs associated with the management of those wastes. However, it is important to estimate the advantages of valorisation processes not only in terms of economic income, but also considering the environmental impacts. This would help to know if the valorisation of a residue provokes higher impact than other waste management options, which means that its advantages are probably not enough for guarantying a sustainable waste reuse. To that purpose, there are several methodologies to evaluate the environmental impacts of processes, including those of waste management, providing different indicators which give information on relevant environmental aspects. In the current study, a comparative environmental assessment between a valorisation process (fishmeal and oil production) and different waste management scenarios (composting, incineration and landfilling) was developed. This comparison is a necessary step for the development and industrial implementation of these processes as the best alternative treatment for fish by-products. The obtained results showed that both valorisation process and waste management treatments presented similar impacts. However, a significant benefit can be achieved through valorisation of fish by-products. Additionally, the implications of the possible presence of pollutants were discussed.

Quality and shelf-life prediction for retail fresh hake (Merluccius merluccius).

Fish quality has a direct impact on market price and its accurate assessment and prediction are of main importance to set prices, increase competitiveness, resolve conflicts of interest and prevent food wastage due to conservative product shelf-life estimations. In this work we present a general methodology to derive predictive models of fish freshness under different storage conditions. The approach makes use of the theory of optimal experimental design, to maximize data information and in this way reduce the number of experiments. The resulting growth model for specific spoilage microorganisms in hake (Merluccius merluccius) is sufficiently informative to estimate quality sensory indexes under time-varying temperature profiles. In addition it incorporates quantitative information of the uncertainty induced by fish variability. The model has been employed to test the effect of factors such as fishing gear or evisceration, on fish spoilage and therefore fish quality. Results show no significant differences in terms of microbial growth between hake fished by long-line or bottom-set nets, within the implicit uncertainty of the model. Similar conclusions can be drawn for gutted and un-gutted hake along the experiment horizon. In addition, whenever there is the possibility to carry out the necessary experiments, this approach is sufficiently general to be used in other fish species and under different stress variables.

Characterizing multistationarity regimes in biochemical reaction networks.

Switch like responses appear as common strategies in the regulation of cellular systems. Here we present a method to characterize bistable regimes in biochemical reaction networks that can be of use to both direct and reverse engineering of biological switches. In the design of a synthetic biological switch, it is important to study the capability for bistability of the underlying biochemical network structure. Chemical Reaction Network Theory (CRNT) may help at this level to decide whether a given network has the capacity for multiple positive equilibria, based on their structural properties. However, in order to build a working switch, we also need to ensure that the bistability property is robust, by studying the conditions leading to the existence of two different steady states. In the reverse engineering of biological switches, knowledge collected about the bistable regimes of the underlying potential model structures can contribute at the model identification stage to a drastic reduction of the feasible region in the parameter space of search. In this work, we make use and extend previous results of the CRNT, aiming not only to discriminate whether a biochemical reaction network can exhibit multiple steady states, but also to determine the regions within the whole space of parameters capable of producing multistationarity. To that purpose we present and justify a condition on the parameters of biochemical networks for the appearance of multistationarity, and propose an efficient and reliable computational method to check its satisfaction through the parameter space.

Fish discards management: pollution levels and best available removal techniques.

Fish discards and by-catch issues are highly topical subjects that are permanently under a social focus. Two main approaches are being considered to address this discard problem: reducing the by-catch and increasing by-catch utilization. Interest in increased by-catch valorization may arise from a greater demand for fish products, such as the development of new markets for previously discarded species, the use of low-value specimens for aquaculture or the creation of value-added fish products for the food, pharmaceutical or cosmetic industries. However, contaminants present in fish discards may be transferred to their valorized products, leading to possible long-term bioaccumulation and subsequent adverse health effects. In this valorization framework, the aim is to promote responsible and sustainable management of marine resources. The pollutant levels in catches from European fisheries and the best available decontamination techniques for marine valorized discards/by-products are compiled and analyzed in this work.

Dynamic optimization of distributed biological systems using robust and efficient numerical techniques.

BACKGROUND: Systems biology allows the analysis of biological systems behavior under different conditions through in silico experimentation. The possibility of perturbing biological systems in different manners calls for the design of perturbations to achieve particular goals. Examples would include, the design of a chemical stimulation to maximize the amplitude of a given cellular signal or to achieve a desired pattern in pattern formation systems, etc. Such design problems can be mathematically formulated as dynamic optimization problems which are particularly challenging when the system is described by partial differential equations.This work addresses the numerical solution of such dynamic optimization problems for spatially distributed biological systems. The usual nonlinear and large scale nature of the mathematical models related to this class of systems and the presence of constraints on the optimization problems, impose a number of difficulties, such as the presence of suboptimal solutions, which call for robust and efficient numerical techniques. RESULTS: Here, the use of a control vector parameterization approach combined with efficient and robust hybrid global optimization methods and a reduced order model methodology is proposed. The capabilities of this strategy are illustrated considering the solution of a two challenging problems: bacterial chemotaxis and the FitzHugh-Nagumo model. CONCLUSIONS: In the process of chemotaxis the objective was to efficiently compute the time-varying optimal concentration of chemotractant in one of the spatial boundaries in order to achieve predefined cell distribution profiles. Results are in agreement with those previously published in the literature. The FitzHugh-Nagumo problem is also efficiently solved and it illustrates very well how dynamic optimization may be used to force a system to evolve from an undesired to a desired pattern with a reduced number of actuators. The presented methodology can be used for the efficient dynamic optimization of generic distributed biological systems.

Inference of complex biological networks: distinguishability issues and optimization-based solutions.

BACKGROUND: The inference of biological networks from high-throughput data has received huge attention during the last decade and can be considered an important problem class in systems biology. However, it has been recognized that reliable network inference remains an unsolved problem. Most authors have identified lack of data and deficiencies in the inference algorithms as the main reasons for this situation. RESULTS: We claim that another major difficulty for solving these inference problems is the frequent lack of uniqueness of many of these networks, especially when prior assumptions have not been taken properly into account. Our contributions aid the distinguishability analysis of chemical reaction network (CRN) models with mass action dynamics. The novel methods are based on linear programming (LP), therefore they allow the efficient analysis of CRNs containing several hundred complexes and reactions. Using these new tools and also previously published ones to obtain the network structure of biological systems from the literature, we find that, often, a unique topology cannot be determined, even if the structure of the corresponding mathematical model is assumed to be known and all dynamical variables are measurable. In other words, certain mechanisms may remain undetected (or they are falsely detected) while the inferred model is fully consistent with the measured data. It is also shown that sparsity enforcing approaches for determining 'true' reaction structures are generally not enough without additional prior information. CONCLUSIONS: The inference of biological networks can be an extremely challenging problem even in the utopian case of perfect experimental information. Unfortunately, the practical situation is often more complex than that, since the measurements are typically incomplete, noisy and sometimes dynamically not rich enough, introducing further obstacles to the structure/parameter estimation process. In this paper, we show how the structural uniqueness and identifiability of the models can be guaranteed by carefully adding extra constraints, and that these important properties can be checked through appropriate computation methods.

An iterative identification procedure for dynamic modeling of biochemical networks.

BACKGROUND: Mathematical models provide abstract representations of the information gained from experimental observations on the structure and function of a particular biological system. Conferring a predictive character on a given mathematical formulation often relies on determining a number of non-measurable parameters that largely condition the model's response. These parameters can be identified by fitting the model to experimental data. However, this fit can only be accomplished when identifiability can be guaranteed. RESULTS: We propose a novel iterative identification procedure for detecting and dealing with the lack of identifiability. The procedure involves the following steps: 1) performing a structural identifiability analysis to detect identifiable parameters; 2) globally ranking the parameters to assist in the selection of the most relevant parameters; 3) calibrating the model using global optimization methods; 4) conducting a practical identifiability analysis consisting of two (a priori and a posteriori) phases aimed at evaluating the quality of given experimental designs and of the parameter estimates, respectively and 5) optimal experimental design so as to compute the scheme of experiments that maximizes the quality and quantity of information for fitting the model. CONCLUSIONS: The presented procedure was used to iteratively identify a mathematical model that describes the NF-kappaB regulatory module involving several unknown parameters. We demonstrated the lack of identifiability of the model under typical experimental conditions and computed optimal dynamic experiments that largely improved identifiability properties.