High-dimensional Ordinary Differential Equations Research Articles

Stability and bifurcation of a two competing prey-one predator system with anti-predator behavior

This article considers the impact of competitive response to interfering time and anti-predator behavior of a three species system in which one predator consumes both the competing prey species. Here one of the competing species shows anti-predator behavior. We have shown that its solutions are non-negative and bounded. Further, we analyze the existence and stability of all the feasible equilibria. Conditions for uniform persistence of the system are derived. Applying Bendixson’s criterion for high-dimensional ordinary differential equations, we prove that the coexistence equilibrium point is globally stable under specific conditions. The system admits Hopf bifurcation when anti-predator behavior rate crosses a critical value. Analytical results are verified numerically.

Identifiability analysis of linear ordinary differential equation systems with a single trajectory

Ordinary differential equations (ODEs) are widely used to model dynamical behavior of systems. It is important to perform identifiability analysis prior to estimating unknown parameters in ODEs (a.k.a. inverse problem), because if a system is unidentifiable, the estimation procedure may fail or produce erroneous and misleading results. Although several qualitative identifiability measures have been proposed, much less effort has been given to developing quantitative (continuous) scores that are robust to uncertainties in the data, especially for those cases in which the data are presented as a single trajectory beginning with one initial value. In this paper, we first derived a closed-form representation of linear ODE systems that are not identifiable based on a single trajectory. This representation helps researchers design practical systems and choose the right prior structural information in practice. Next, we proposed several quantitative scores for identifiability analysis in practice. In simulation studies, the proposed measures outperformed the main competing method significantly, especially when noise was presented in the data. We also discussed the asymptotic properties of practical identifiability for high-dimensional ODE systems and conclude that, without additional prior information, many random ODE systems are practically unidentifiable when the dimension approaches infinity.

Dynamical modeling for non-Gaussian data with high-dimensional sparse ordinary differential equations

Ordinary differential equations (ODE) have been widely used for modeling dynamical complex systems. For high-dimensional ODE models where the number of differential equations is large, it remains challenging to estimate the ODE parameters and to identify the sparse structure of the ODE models. Most existing methods exploit the least-square based approach and are only applicable to Gaussian observations. However, as discrete data are ubiquitous in applications, it is of practical importance to develop dynamic modeling for non-Gaussian observations. New methods and algorithms are developed for both parameter estimation and sparse structure identification in high-dimensional linear ODE systems. First, the high-dimensional generalized profiling method is proposed as a likelihood-based approach with ODE fidelity and sparsity-inducing regularization, along with efficient computation based on parameter cascading. Second, two versions of the two-step collocation methods are extended to the non-Gaussian set-up by incorporating the iteratively reweighted least squares technique. Simulations show that the profiling procedure has excellent performance in latent process and derivative fitting and ODE parameter estimation, while the two-step collocation approach excels in identifying the sparse structure of the ODE system. The usefulness of the proposed methods is also demonstrated by analyzing three real datasets from Google trends, stock market sectors, and yeast cell cycle studies.

Real-time Human Response Prediction Using a Non-intrusive Data-driven Model Reduction Scheme

Recent research in non-intrusive data-driven model order reduction (MOR) enabled accurate and efficient approximation of parameterized ordinary differential equations (ODEs). However, previous studies have focused on constant parameters, whereas time-dependent parameters have been neglected. The purpose of this paper is to introduce a novel two-step MOR scheme to tackle this issue. In a first step, classic MOR approaches are applied to calculate a low-dimensional representation of high-dimensional ODE solutions, i.e., to extract the most important features of simulation data. Based on this representation, a long short-term memory (LSTM) is trained to predict the reduced dynamics iteratively in a second step considering the parameters at the respective time step. The potential of this approach is demonstrated on an occupant model within a car driving scenario. The reduced model's response to time-varying accelerations matches the reference data with high accuracy for a limited amount of time. Furthermore, real-time capability is achieved. Accordingly, it is concluded that the presented method is well suited to approximate parameterized ODEs and can handle time-dependent parameters in contrast to common methods.

DeepCME: A deep learning framework for computing solution statistics of the chemical master equation.

Stochastic models of biomolecular reaction networks are commonly employed in systems and synthetic biology to study the effects of stochastic fluctuations emanating from reactions involving species with low copy-numbers. For such models, the Kolmogorov's forward equation is called the chemical master equation (CME), and it is a fundamental system of linear ordinary differential equations (ODEs) that describes the evolution of the probability distribution of the random state-vector representing the copy-numbers of all the reacting species. The size of this system is given by the number of states that are accessible by the chemical system, and for most examples of interest this number is either very large or infinite. Moreover, approximations that reduce the size of the system by retaining only a finite number of important chemical states (e.g. those with non-negligible probability) result in high-dimensional ODE systems, even when the number of reacting species is small. Consequently, accurate numerical solution of the CME is very challenging, despite the linear nature of the underlying ODEs. One often resorts to estimating the solutions via computationally intensive stochastic simulations. The goal of the present paper is to develop a novel deep-learning approach for computing solution statistics of high-dimensional CMEs by reformulating the stochastic dynamics using Kolmogorov's backward equation. The proposed method leverages superior approximation properties of Deep Neural Networks (DNNs) to reliably estimate expectations under the CME solution for several user-defined functions of the state-vector. This method is algorithmically based on reinforcement learning and it only requires a moderate number of stochastic simulations (in comparison to typical simulation-based approaches) to train the "policy function". This allows not just the numerical approximation of various expectations for the CME solution but also of its sensitivities with respect to all the reaction network parameters (e.g. rate constants). We provide four examples to illustrate our methodology and provide several directions for future research.

A New Set of Stability Criteria Extending Lyapunov's Direct Method

A dynamical system is a mathematical model described by a high dimensional ordinary differential equation for a wide variety of real world phenomena, which can be as simple as a clock pendulum or as complex as a chaotic Lorenz system. Stability is an important topic in the studies of the dynamical system. A major challenge is that the analytical solution of a time-varying nonlinear dynamical system is in general not known. Lyapunov's direct method is a classical approach used for many decades to study stability without explicitly solving the dynamical system, and has been successfully employed in numerous applications ranging from aerospace guidance systems, chaos theory, to traffic assignment. Roughly speaking, an equilibrium is stable if an energy function monotonically decreases along the trajectory of the dynamical system. This paper extends Lyapunov's direct method by allowing the energy function to follow a rich set of dynamics. More precisely, the paper proves two theorems, one on globally uniformly asymptotic stability and the other on stability in the sense of Lyapunov, where stability is guaranteed provided that the evolution of the energy function satisfies an inequality of a non-negative Hurwitz polynomial differential operator, which uses not only the first-order but also high-order time derivatives of the energy function. The classical Lyapunov theorems are special cases of the extended theorems. the paper provides an example in which the new theorem successfully determines stability while the classical Lyapunov's direct method fails.

Robust Calibration of High Dimension Nonlinear Dynamical Models for Omics Data: An Application in Cancer Systems Biology

In computational mathematical modeling of biological systems, most model parameters, such as initial conditions, kinetics, and scale factors, are usually unknown because they cannot be directly measured. Therefore, key issues in system identification of nonlinear systems are model calibration and identifiability analysis. Currently, existing methodologies for parameter estimation are divided in two classes: frequentist and Bayesian methods. The first optimize a cost function, while the second estimate the posterior distribution of parameters through different sampling techniques. However, when dealing with high-dimensional models, these methodologies suffer from an increasing computational cost due to the important volume of -omic data necessary to carry out reliable and robust solutions. Here, we present an innovative Bayesian method, called conditional robust calibration (CRC), for model calibration and identifiability analysis. The algorithm is an iterative procedure based on a uniform and joint perturbation of the parameter space. At each step the algorithm returns the probability density functions of all parameters that progressively shrink toward specific points in the parameter space. These distributions are estimated on parameter samples that guarantee a certain level of agreement between each observable and the corresponding in silico measure. We apply CRC to a nonlinear high-dimensional ordinary differential equations model representing the signaling pathway of p38MAPK in multiple myeloma. The available data set consists of time courses of proteomic cancerous data. We test CRC performances in comparison with profile-likelihood and approximate Bayesian computation sequential Monte Carlo. We obtain a more precise and robust solution with a reduced computational cost.

Bayesian inference of a directional brain network model for intracranial EEG data

The human brain is a network system in which brain regions, as network nodes, constantly interact with each other. The directional effect exerted by one brain component on another is referred to as directional connectivity. Since the brain is also a continuous time dynamic system, it is natural to use ordinary differential equations (ODEs) to model directional connections among brain regions. The authors propose a high-dimensional ODE model to explore directional connectivity among many small brain regions recorded by intracranial EEG (iEEG). The new ODE model, motivated by the physical mechanism of the damped harmonic oscillator, is effective for approximating neural oscillation, a rhythmic or repetitive neural activity involved in many important brain functions. To produce scientifically meaningful network results, a cluster structure is assumed for the ODE model parameters that quantify directional connectivity among regions. The cluster structure is in line with the functional specialization of the human brain; the brain areas specialized in the same function tend to be in the same cluster. Two Bayesian methods are developed to estimate the model parameters of the proposed ODE model and to identify clusters of strongly connected brain regions. The proposed ODE model and Bayesian method are applied to iEEG data collected from a patient with medically intractable epilepsy and used to examine the patient’s brain networks before the seizure onset.

Local and global Hopf bifurcation analysis on simplified bidirectional associative memory neural networks with multiple delays

In this paper, a class of simplified bidirectional associative memory (BAM) neural networks with multiple delays are considered. By analyzing the associated characteristic transcendental equation, their linear stability is investigated and Hopf bifurcation is demonstrated. By applying Nyquist criterion, the length of delay which preserves the stability of the zero equilibrium is estimated. Some explicit results are derived for stability and direction of the bifurcating periodic orbit by using the normal form theory and center manifold arguments. Global existence of periodic orbits is also established by using a global Hopf bifurcation theorem for functional differential equations (FDE) and a Bendixson’s criterion for high-dimensional ordinary differential equations (ODE) due to Li and Muldowney. Finally, numerical simulations supporting the theoretical analysis are carried out.

Identifying the dynamic gene regulatory network during latent HIV-1 reactivation using high-dimensional ordinary differential equations.

Reactivation of latently infected cells has emerged as an important strategy for eradication of HIV. However, genetic mechanisms of regulation after reactivation remain unclear. We describe a five-step pipeline to study the dynamics of the gene regulatory network following a viral reactivation using high-dimensional ordinary differential equations. Our pipeline implements a combination of five different methods, by detecting temporally differentially expressed genes (step 1), clustering genes with similar temporal expression patterns into a small number of response modules (step2), performing a functional enrichment analysis within each gene response module (step 3), identifying a network structure based on the gene response modules using ordinary differential equations (ODE) and a high-dimensional variable selection technique (step 4), and obtaining a gene regulatory model based on refined parameter estimates using nonlinear least squares (step 5). We applied our pipeline to a time course gene expression data of latently infected T-cells following a latency-reversion.

Identifying the dynamic gene regulatory network during latent HIV-1 reactivation using high-dimensional ordinary differential equations

Reactivation of latently infected cells has emerged as an important strategy for eradication of HIV. However, genetic mechanisms of regulation after reactivation remain unclear. We describe a five-step pipeline to study the dynamics of the gene regulatory network following a viral reactivation using high-dimensional ordinary differential equations. Our pipeline implements a combination of five different methods, by detecting temporally differentially expressed genes (step 1), clustering genes with similar temporal expression patterns into a small number of response modules (step 2), performing a functional enrichment analysis within each gene response module (step 3), identifying a network structure based on the gene response modules using ordinary differential equations (ODE) and a high-dimensional variable selection technique (step 4), and obtaining a gene regulatory model based on refined parameter estimates using nonlinear least squares (step 5). We applied our pipeline to a time course gene expression data of latently infected T-cells following a latency-reversion.

Support vector machine classifier for predicting extubation outcome in intensive care unit

Ordinary differential equations (ODEs) are widely used to model dynamical behavior of systems. It is important to perform identifiability analysis prior to estimating unknown parameters in ODEs (a.k.a. inverse problem), because if a system is unidentifiable, the estimation procedure may fail or produce erroneous and misleading results. Although several qualitative identifiability measures have been proposed, much less effort has been given to developing quantitative (continuous) scores that are robust to uncertainties in the data, especially for those cases in which the data are presented as a single trajectory beginning with one initial value. In this paper, we first derived a closed-form representation of linear ODE systems that are not identifiable based on a single trajectory. This representation helps researchers design practical systems and choose the right prior structural information in practice. Next, we proposed several quantitative scores for identifiability analysis in practice. In simulation studies, the proposed measures outperformed the main competing method significantly, especially when noise was presented in the data. We also discussed the asymptotic properties of practical identifiability for high-dimensional ODE systems and conclude that, without additional prior information, many random ODE systems are practically unidentifiable when the dimension approaches infinity.

Network Reconstruction From High-Dimensional Ordinary Differential Equations

ABSTRACTWe consider the task of learning a dynamical system from high-dimensional time-course data. For instance, we might wish to estimate a gene regulatory network from gene expression data measured at discrete time points. We model the dynamical system nonparametrically as a system of additive ordinary differential equations. Most existing methods for parameter estimation in ordinary differential equations estimate the derivatives from noisy observations. This is known to be challenging and inefficient. We propose a novel approach that does not involve derivative estimation. We show that the proposed method can consistently recover the true network structure even in high dimensions, and we demonstrate empirical improvement over competing approaches. Supplementary materials for this article are available online.

Bayesian inference of high-dimensional, cluster-structured ordinary differential equation models with applications to brain connectivity studies

We build a new ordinary differential equation (ODE) model for the directional interaction, also called effective connectivity, among brain regions whose activities are measured by intracranial electrocorticography (ECoG) data. In contrast to existing ODE models that focus on effective connectivity among only a few large anatomic brain regions and that rely on strong prior belief of the existence and strength of the connectivity, the proposed high-dimensional ODE model, motivated by statistical considerations, can be used to explore connectivity among multiple small brain regions. The new model, called the modular and indicator-based dynamic directional model (MIDDM), features a cluster structure, which consists of modules of densely connected brain regions, and uses indicators to differentiate significant and void directional interactions among brain regions. We develop a unified Bayesian framework to quantify uncertainty in the assumed ODE model, identify clusters, select strongly connected brain regions, and make statistical comparison between brain networks across different experimental trials. The prior distributions in the Bayesian model for MIDDM parameters are carefully designed such that the ensuing joint posterior distributions for ODE state functions and the MIDDM parameters have well-defined and easy-to-simulate posterior conditional distributions. To further speed up the posterior simulation, we employ parallel computing schemes in Markov chain Monte Carlo steps. We show that the proposed Bayesian approach outperforms an existing optimization-based ODE estimation method. We apply the proposed method to an auditory electrocorticography dataset and evaluate brain auditory network changes across trials and different auditory stimuli.

Local Hopf bifurcation analysis and global existence of periodic solutions in a gene expression model with delays

In this paper, local Hopf bifurcation of a gene expression model with three delays is investigated by applying the frequency domain approach. It is shown that Hopf bifurcation will occur as the bifurcation parameter, the sum of all delays, passes through a sequence of critical values. The direction and the stability of bifurcating periodic solutions are determined by the Nyquist criterion and the graphical Hopf bifurcation theorem. By using the Bendixson’s criterion for high-dimensional ordinary differential equations and global Hopf bifurcation theorem for functional differential equations, the global existence of periodic solutions is established. Numerical results are also included to give a verification test of the theoretical analysis.

Global Existence of Periodic Solutions in a Nonlinear Delay-Coupling Chaos System

The dynamics of a unidirectional nonlinear delayed-coupling chaos system is investigated. Based on the local Hopf bifurcation at the zero equilibrium, we prove the global existence of periodic solutions using a global Hopf bifurcation result due to Wu and a Bendixson’s criterion for higher dimensional ordinary differential equations due to Li & Muldowney.

Solving systems of differential equations of block structure with nonseparated boundary conditions

We study the solution of a system of higher-dimensional ordinary differential equations of block structure. Some separate subsystems are connected with each other by the nonseparated boundary conditions caused by an arbitrary relation between the boundary values of the solutions to the subsystems. For numerical solution, we propose some scheme of the method of transferring the boundary conditions that takes into account the specific characteristics of the systems under consideration. The results of numerical experiments are given. DOI: 10.1134/S1990478915010019 In this article, taking into account the block structure of the system of differential equations and the weak but arbitrary filling of the boundary condition matrix, we propose a version of the method of transfer of the boundary conditions. The advantage of this approach in comparison with the direct use of transfer methods in general form is obvious since here the transfer is carried out only with respect to those variables whose coefficients in the boundary conditions are nonzero; moreover, the transfer is carried out with the use of the only subsystem of differential equations involving the variable that is being transferred. We show the results of the numerical experiments obtained by solving a test problem based on the problem of calculating an unsteady fluid motion in a fragment of a pipeline network of complex structure.

Stability and global Hopf bifurcation for neutral BAM neural network

In this paper, a simplified neutral BAM neural network model is considered. Some results of stability and Hopf bifurcations occurring at the zero equilibrium as the delay increases are exhibited. Global existence of periodic solutions is established by using a global Hopf bifurcation result of Wu and a Bendixson criterion for higher dimensional ordinary differential equations due to Li and Muldowney. Finally, some numerical simulations are carried out to support the analytical results.

Modeling genome-wide dynamic regulatory network in mouse lungs with influenza infection using high-dimensional ordinary differential equations.

The immune response to viral infection is regulated by an intricate network of many genes and their products. The reverse engineering of gene regulatory networks (GRNs) using mathematical models from time course gene expression data collected after influenza infection is key to our understanding of the mechanisms involved in controlling influenza infection within a host. A five-step pipeline: detection of temporally differentially expressed genes, clustering genes into co-expressed modules, identification of network structure, parameter estimate refinement, and functional enrichment analysis, is developed for reconstructing high-dimensional dynamic GRNs from genome-wide time course gene expression data. Applying the pipeline to the time course gene expression data from influenza-infected mouse lungs, we have identified 20 distinct temporal expression patterns in the differentially expressed genes and constructed a module-based dynamic network using a linear ODE model. Both intra-module and inter-module annotations and regulatory relationships of our inferred network show some interesting findings and are highly consistent with existing knowledge about the immune response in mice after influenza infection. The proposed method is a computationally efficient, data-driven pipeline bridging experimental data, mathematical modeling, and statistical analysis. The application to the influenza infection data elucidates the potentials of our pipeline in providing valuable insights into systematic modeling of complicated biological processes.

Global Hopf Bifurcation Analysis of a Nicholson’s Blowflies Equation of Neutral Type

We investigate Hopf bifurcations in a delayed Nicholson’s blowflies equation of neutral type, derived from the Gurtin–MacCamy model. A key parameter that determines the direction of the Hopf bifurcation and the stability of the bifurcating periodic solutions is derived. Global extension of local Hopf branches is established by combining a global Hopf bifurcation theorem with a Bendixson criterion for higher dimensional ordinary differential equations. We show that a branch of slowly varying periodic solutions and a branch of fast oscillating periodic solutions coexist for all large delays.