Sparse Additive Ordinary Differential Equations for Dynamic Gene Regulatory Network Modeling

Exploring gene regulation and biological processes in insects: Insights from omics data using gene regulatory network models

Gene regulatory network (GRN) comprises complicated yet intertwined gene-regulator relationships. Understanding the GRN dynamics will unravel the complexity behind the observed gene expressions. Insect gene regulation is often complicated due to their complex life cycles and diverse ecological adaptations. The main interest of this review is to have an update on the current mathematical modelling methods of GRNs to explain insect science. Several popular GRN architecture models are discussed, together with examples of applications in insect science. In the last part of this review, each model is compared from different aspects, including network scalability, computation complexity, robustness to noise and biological relevancy.

Identifying temporal pathways using biomarkers in the presence of latent non-Gaussian components

Time-series data collected from a network of random variables are useful for identifying temporal pathways among the network nodes. Observed measurements may contain multiple sources of signals and noises, including Gaussian signals of interest and non-Gaussian noises, including artifacts, structured noise, and other unobserved factors (eg, genetic risk factors, disease susceptibility). Existing methods, including vector autoregression (VAR) and dynamic causal modeling do not account for unobserved non-Gaussian components. Furthermore, existing methods cannot effectively distinguish contemporaneous relationships from temporal relations. In this work, we propose a novel method to identify latent temporal pathways using time-series biomarker data collected from multiple subjects. The model adjusts for the non-Gaussian components and separates the temporal network from the contemporaneous network. Specifically, an independent component analysis (ICA) is used to extract the unobserved non-Gaussian components, and residuals are used to estimate the contemporaneous and temporal networks among the node variables based on method of moments. The algorithm is fast and can easily scale up. We derive the identifiability and the asymptotic properties of the temporal and contemporaneous networks. We demonstrate superior performance of our method by extensive simulations and an application to a study of attention-deficit/hyperactivity disorder (ADHD), where we analyze the temporal relationships between brain regional biomarkers. We find that temporal network edges were across different brain regions, while most contemporaneous network edges were bilateral between the same regions and belong to a subset of the functional connectivity network.

The WRKY46-MYC2 module plays a critical role in E-2-hexenal-induced anti-herbivore responses by promoting flavonoid accumulation

Volatile organic compounds (VOCs) play key roles in plant-plant communication, especially in response to pest attack. E-2-hexenal is an important component of VOCs, but it is unclear whether it can induce endogenous plant resistance to insects. Here, we show that E-2-hexenal activates early signaling events in Arabidopsis (Arabidopsis thaliana) mesophyll cells, including an H2O2 burst at the plasma membrane, the directed flow of calcium ions, and an increase in cytosolic calcium concentration. Treatment of wild-type Arabidopsis plants with E-2-hexenal increases their resistance when challenged with the diamondback moth Plutella xylostella L., and this phenomenon is largely lost in the wrky46 mutant. Mechanistically, E-2-hexenal induces the expression of WRKY46 and MYC2, and the physical interaction of their encoded proteins was verified by yeast two-hybrid, firefly luciferase complementation imaging, and in vitro pull-down assays. The WRKY46-MYC2 complex directly binds to the promoter of RBOHD to promote its expression, as demonstrated by luciferase reporter, yeast one-hybrid, chromatin immunoprecipitation, and electrophoretic mobility shift assays. This module also positively regulates the expression of E-2-hexenal-induced naringenin biosynthesis genes (TT4 and CHIL) and the accumulation of total flavonoids, thereby modulating plant tolerance to insects. Together, our results highlight an important role for the WRKY46-MYC2 module in the E-2-hexenal-induced defense response of Arabidopsis, providing new insights into the mechanisms by which VOCs trigger plant defense responses.

The identification of piecewise non-linear dynamical system without understanding the mechanism

This paper designs an algorithm to distill the piecewise non-linear dynamical system from the data without prior knowledge. The system to be identified does not have to be written as a known model term or be thoroughly understood. We exploit the fact that an unknown piecewise non-linear system can be decomposed into the Fourier series as long as its equations of motion are Riemann integrable. Based on this property, we reduce the challenge of finding the correct model to discovering the Fourier series approximation. However, the Fourier series approximation of the piecewise function is inaccurate. The new method takes advantage of this weakness to determine whether the model has piecewise features and to find a way to discover the discontinuity set. Then, the dynamical system on each segment is identified as a pure Fourier series. Identification of intricate models can be achieved in simple steps. The results show that the method can accurately discover the equation of motion and precisely capture the non-smooth characteristic. Next, the prediction and further detailed analysis can be carried out.

A Network-Based Approach for Improving Annotation of Transcription Factor Functions and Binding Sites in Arabidopsis thaliana

Transcription factors are an integral component of the cellular machinery responsible for regulating many biological processes, and they recognize distinct DNA sequence patterns as well as internal/external signals to mediate target gene expression. The functional roles of an individual transcription factor can be traced back to the functions of its target genes. While such functional associations can be inferred through the use of binding evidence from high-throughput sequencing technologies available today, including chromatin immunoprecipitation sequencing, such experiments can be resource-consuming. On the other hand, exploratory analysis driven by computational techniques can alleviate this burden by narrowing the search scope, but the results are often deemed low-quality or non-specific by biologists. In this paper, we introduce a data-driven, statistics-based strategy to predict novel functional associations for transcription factors in the model plant Arabidopsis thaliana. To achieve this, we leverage one of the largest available gene expression compendia to build a genome-wide transcriptional regulatory network and infer regulatory relationships among transcription factors and their targets. We then use this network to build a pool of likely downstream targets for each transcription factor and query each target pool for functionally enriched gene ontology terms. The results exhibited sufficient statistical significance to annotate most of the transcription factors in Arabidopsis with highly specific biological processes. We also perform DNA binding motif discovery for transcription factors based on their target pool. We show that the predicted functions and motifs strongly agree with curated databases constructed from experimental evidence. In addition, statistical analysis of the network revealed interesting patterns and connections between network topology and system-level transcriptional regulation properties. We believe that the methods demonstrated in this work can be extended to other species to improve the annotation of transcription factors and understand transcriptional regulation on a system level.

Vaginal microbiota networks as a mechanistic predictor of aerobic vaginitis

Aerobic vaginitis (AV) is a complex vaginal dysbiosis that is thought to be caused by the micro-ecological change of the vaginal microbiota. While most studies have focused on how changes in the abundance of individual microbes are associated with the emergence of AV, we still do not have a complete mechanistic atlas of the microbe-AV link. Network modeling is central to understanding the structure and function of any microbial community assembly. By encapsulating the abundance of microbes as nodes and ecological interactions among microbes as edges, microbial networks can reveal how each microbe functions and how one microbe cooperate or compete with other microbes to mediate the dynamics of microbial communities. However, existing approaches can only estimate either the strength of microbe-microbe link or the direction of this link, failing to capture full topological characteristics of a network, especially from high-dimensional microbial data. We combine allometry scaling law and evolutionary game theory to derive a functional graph theory that can characterize bidirectional, signed, and weighted interaction networks from any data domain. We apply our theory to characterize the causal interdependence between microbial interactions and AV. From functional networks arising from different functional modules, we find that, as the only favorable genus from Firmicutes among all identified genera, the role of Lactobacillus in maintaining vaginal microbial symbiosis is enabled by upregulation from other microbes, rather than through any intrinsic capacity. Among Lactobacillus species, the proportion of L. crispatus to L. iners is positively associated with more healthy acid vaginal ecosystems. In a less healthy alkaline ecosystem, L. crispatus establishes a contradictory relationship with other microbes, leading to population decrease relative to L. iners. We identify topological changes of vaginal microbiota networks when the menstrual cycle of women changes from the follicular to luteal phases. Our network tool provides a mechanistic approach to disentangle the internal workings of the microbiota assembly and predict its causal relationships with human diseases including AV.

A Joint estimation approach to sparse additive ordinary differential equations

Ordinary differential equations (ODEs) are widely used to characterize the dynamics of complex systems in real applications. In this article, we propose a novel joint estimation approach for generalized sparse additive ODEs where observations are allowed to be non-Gaussian. The new method is unified with existing collocation methods by considering the likelihood, ODE fidelity and sparse regularization simultaneously. We design a block coordinate descent algorithm for optimizing the non-convex and non-differentiable objective function. The global convergence of the algorithm is established. The simulation study and two applications demonstrate the superior performance of the proposed method in estimation and improved performance of identifying the sparse structure.

A statistical framework for recovering pseudo-dynamic networks from static data

MOTIVATION

The collection of temporal or perturbed data is often a prerequisite for reconstructing dynamic networks in most cases. However, these types of data are seldom available for genomic studies in medicine, thus significantly limiting the use of dynamic networks to characterize the biological principles underlying human health and diseases.

RESULTS

We proposed a statistical framework to recover disease risk-associated pseudo-dynamic networks (DRDNet) from steady-state data. We incorporated a varying coefficient model with multiple ordinary differential equations to learn a series of networks. We analyzed the publicly available Genotype-Tissue Expression data to construct networks associated with hypertension risk, and biological findings showed that key genes constituting these networks had pivotal and biologically relevant roles associated with the vascular system. We also provided the selection consistency of the proposed learning procedure and evaluated its utility through extensive simulations.

AVAILABILITY AND IMPLEMENTATION

DRDNet is implemented in the R language, and the source codes are available at https://github.com/chencxxy28/DRDnet/.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

A nonlinear sparse neural ordinary differential equation model for multiple functional processes

In this article, we propose a new sparse neural ordinary differential equation (ODE) model to characterize flexible relations among multiple functional processes. We characterize the latent states of the functions via a set of ordinary differential equations. We then model the dynamic changes of the latent states using a deep neural network (DNN) with a specially designed architecture and a sparsity-inducing regularization. The new model is able to capture both nonlinear and sparse dependent relations among multivariate functions. We develop an efficient optimization algorithm to estimate the unknown weights for the DNN under the sparsity constraint. We establish both the algorithmic convergence and selection consistency, which constitute the theoretical guarantees of the proposed method. We illustrate the efficacy of the method through simulations and a gene regulatory network example.

Granger Causality: A Review and Recent Advances

Introduced more than a half-century ago, Granger causality has become a popular tool for analyzing time series data in many application domains, from economics and finance to genomics and neuroscience. Despite this popularity, the validity of this framework for inferring causal relationships among time series has remained the topic of continuous debate. Moreover, while the original definition was general, limitations in computational tools have constrained the applications of Granger causality to primarily simple bivariate vector autoregressive processes. Starting with a review of early developments and debates, this article discusses recent advances that address various shortcomings of the earlier approaches, from models for high-dimensional time series to more recent developments that account for nonlinear and non-Gaussian observations and allow for subsampled and mixed-frequency time series.

Kernel Ordinary Differential Equations

Ordinary differential equation (ODE) is widely used in modeling biological and physical processes in science. In this article, we propose a new reproducing kernel-based approach for estimation and inference of ODE given noisy observations. We do not assume the functional forms in ODE to be known, or restrict them to be linear or additive, and we allow pairwise interactions. We perform sparse estimation to select individual functionals, and construct confidence intervals for the estimated signal trajectories. We establish the estimation optimality and selection consistency of kernel ODE under both the low-dimensional and high-dimensional settings, where the number of unknown functionals can be smaller or larger than the sample size. Our proposal builds upon the smoothing spline analysis of variance (SS-ANOVA) framework, but tackles several important problems that are not yet fully addressed, and thus extends the scope of existing SS-ANOVA as well. We demonstrate the efficacy of our method through numerous ODE examples.

Rejuvenation remodels transcriptional network to improve rhizogenesis in mature Juglans tree

Adventitious rooting of walnut species (Juglans L.) is known to be rather difficult, especially for mature trees. The adventitious root formation (ARF) capacities of mature trees can be significantly improved by rejuvenation. However, the underlying gene regulatory networks (GRNs) of rejuvenation remain largely unknown. To characterize such regulatory networks, we carried out the transcriptomic study using RNA samples of the cambia and peripheral tissues on the bottom of rejuvenated and mature walnut (Juglans hindsii × J. regia) cuttings during the ARF. The RNA sequencing data suggested that zeatin biosynthesis, energy metabolism and substance metabolism were activated by rejuvenation, whereas photosynthesis, fatty acid biosynthesis and the synthesis pathways for secondary metabolites were inhibited. The inter- and intra-module GRNs were constructed using differentially expressed genes. We identified 35 hub genes involved in five modules associated with ARF. Among these hub genes, particularly, beta-glucosidase-like (BGLs) family members involved in auxin metabolism were overexpressed at the early stage of the ARF. Furthermore, BGL12 from the cuttings of Juglans was overexpressed in Populus alba × P. glandulosa. Accelerated ARF and increased number of ARs were observed in the transgenic poplars. These results provide a high-resolution atlas of gene activity during ARF and help to uncover the regulatory modules associated with the ARF promoted by rejuvenation.

Reverse engineering gene regulatory network based on complex-valued ordinary differential equation model

BACKGROUND

The growing researches of molecular biology reveal that complex life phenomena have the ability to demonstrating various types of interactions in the level of genomics. To establish the interactions between genes or proteins and understand the intrinsic mechanisms of biological systems have become an urgent need and study hotspot.

RESULTS

In order to forecast gene expression data and identify more accurate gene regulatory network, complex-valued version of ordinary differential equation (CVODE) is proposed in this paper. In order to optimize CVODE model, a complex-valued hybrid evolutionary method based on Grammar-guided genetic programming and complex-valued firefly algorithm is presented.

CONCLUSIONS

When tested on three real gene expression datasets from E. coli and Human Cell, the experiment results suggest that CVODE model could improve 20-50% prediction accuracy of gene expression data, which could also infer more true-positive regulatory relationships and less false-positive regulations than ordinary differential equation.

Modeling genome-wide by environment interactions through omnigenic interactome networks

How genes interact with the environment to shape phenotypic variation and evolution is a fundamental question intriguing to biologists from various fields. Existing linear models built on single genes are inadequate to reveal the complexity of genotype-environment (G-E) interactions. Here, we develop a conceptual model for mechanistically dissecting G-E interplay by integrating previously disconnected theories and methods. Under this integration, evolutionary game theory, developmental modularity theory, and a variable selection method allow us to reconstruct environment-induced, maximally informative, sparse, and casual multilayer genetic networks. We design and conduct two mapping experiments by using a desert-adapted tree species to validate the biological application of the model proposed. The model identifies previously uncharacterized molecular mechanisms that mediate trees' response to saline stress. Our model provides a tool to comprehend the genetic architecture of trait variation and evolution and trace the information flow of each gene toward phenotypes within omnigenic networks.

Modeling Temporal Biomarkers With Semiparametric Nonlinear Dynamical Systems

Dynamical systems based on differential equations are useful for modeling the temporal evolution of biomarkers. These systems can characterize the temporal patterns of biomarkers and inform the detection of interactions among biomarkers. Existing statistical methods for dynamical systems mostly target single time-course data based on a linear model or generalized additive model. Hence, they cannot adequately capture the complex interactions among biomarkers; neither can they take into account the heterogeneity between systems or subjects. in this work, we propose a semiparametric dynamical system based on multi-index models for multiple subjects time-course data. Our model accounts for between-subject heterogeneity by introducing system-level or subject-level covariates to dynamic systems, and it allows for nonlinear relationship and interaction between the combined biomarkers and the temporal rate of each biomarker. For estimation and inference, we consider a two-step procedure based on integral equations from the proposed model. We propose an algorithm that iterates between the estimation of the link function through splines and the estimation of index parameters and that allows for regularization to achieve sparsity. We prove model identifiability and derive the asymptotic properties of the estimated model parameters. A benefit of our approach is to pool information from multiple subjects to identify the interaction among biomarkers. We apply the method to analyze electroencephalogram (EEG) data for patients affected by alcohol dependence. The results reveal new insight on patients' brain activities and demonstrate differential interaction patterns in patients compared to health control subjects.

Genome-wide association studies of callus differentiation for the desert tree, Populus euphratica

Callus differentiation is a key developmental process in plant regeneration from cells. A better understanding of the genetic architecture of callus differentiation timing can help improve tissue transformation and the efficiency of artificial propagation. In this study, we investigated genotypic variation in callus differentiation capacity among 297 diverse P. euphratica trees sampled from a natural population. We employed a genome-wide association study (GWAS) of binary and growth-based parameters to identify loci and characterize the genetic architecture and genetic network underlying regulation of callus differentiation in P. euphratica. The results of this GWAS experiment suggested potential associations controlling whether the callus could differentiate and the process of callus differentiation. We identified multiple significant quantitative trait loci (QTLs), including the genes LOG1 and LOG7 and a locus containing WOX1. We reconstructed a genetic network that visualizes how each QTL interacts uniquely with other variants, and several core QTLs were detected that are involved in the degree of callus differentiation, providing potential targets for selection. This study represents one of the first to identify genetic variants affecting callus differentiation in a forest tree. Our results suggest that callus differentiation may be a typical qualitative-quantitative trait controlled by a major gene as well as polygenes across the genome of P. euphratica. This GWAS will help to design more complex and specific molecular tools for systematically manipulating organ regeneration.

SEGN: Inferring real-time gene networks mediating phenotypic plasticity

The capacity of an organism to alter its phenotype in response to environmental perturbations changes over developmental time and is a process determined by multiple genes that are co-expressed in intricate but organized networks. Characterizing the spatiotemporal change of such gene networks can offer insight into the genomic signatures underlying organismic adaptation, but it represents a major methodological challenge. Here, we integrate the holistic view of systems biology and the interactive notion of evolutionary game theory to reconstruct so-called systems evolutionary game networks (SEGN) that can autonomously detect, track, and visualize environment-induced gene networks along the time axis. The SEGN overcomes the limitations of traditional approaches by inferring context-specific networks, encapsulating bidirectional, signed, and weighted gene-gene interactions into fully informative networks, and monitoring the process of how networks topologically alter across environmental and developmental cues. Based on the design principle of SEGN, we perform a transcriptional plasticity study by culturing Euphrates poplar, a tree that can grow in the saline desert, in saline-free and saline-stress conditions. SEGN characterize previously unknown gene co-regulation that modulates the time trajectories of the trees' response to salt stress. As a marriage of multiple disciplines, SEGN shows its potential to interpret gene interdependence, predict how transcriptional co-regulation responds to various regimes, and provides a hint for exploring the mass, energetic, or signal basis that drives various types of gene interactions.

Gene Regulatory Network Inference: Connecting Plant Biology and Mathematical Modeling

Plant responses to environmental and intrinsic signals are tightly controlled by multiple transcription factors (TFs). These TFs and their regulatory connections form gene regulatory networks (GRNs), which provide a blueprint of the transcriptional regulations underlying plant development and environmental responses. This review provides examples of experimental methodologies commonly used to identify regulatory interactions and generate GRNs. Additionally, this review describes network inference techniques that leverage gene expression data to predict regulatory interactions. These computational and experimental methodologies yield complex networks that can identify new regulatory interactions, driving novel hypotheses. Biological properties that contribute to the complexity of GRNs are also described in this review. These include network topology, network size, transient binding of TFs to DNA, and competition between multiple upstream regulators. Finally, this review highlights the potential of machine learning approaches to leverage gene expression data to predict phenotypic outputs.

Mapping epileptic directional brain networks using intracranial EEG data

The human brain is a directional network system, in which brain regions are network nodes and the influence exerted by one region on another is a network edge. We refer to this directional information flow from one region to another as directional connectivity. Seizures arise from an epileptic directional network; abnormal neuronal activities start from a seizure onset zone and propagate via a network to otherwise healthy brain regions. As such, effective epilepsy diagnosis and treatment require accurate identification of directional connections among regions, i.e., mapping of epileptic patients' brain networks. This article aims to understand the epileptic brain network using intracranial electroencephalographic data-recordings of epileptic patients' brain activities in many regions. The most popular models for directional connectivity use ordinary differential equations (ODE). However, ODE models are sensitive to data noise and computationally costly. To address these issues, we propose a high-dimensional state-space multivariate autoregression (SSMAR) model for the brain's directional connectivity. Different from standard multivariate autoregression and SSMAR models, the proposed SSMAR features a cluster structure, where the brain network consists of several clusters of densely connected brain regions. We develop an expectation-maximization algorithm to estimate the proposed model and use it to map the interregional networks of epileptic patients in different seizure stages. Our method reveals the evolution of brain networks during seizure development.

The geometry of masking in neural populations

The normalization model provides an elegant account of contextual modulation in individual neurons of primary visual cortex. Understanding the implications of normalization at the population level is hindered by the heterogeneity of cortical neurons, which differ in the composition of their normalization pools and semi-saturation constants. Here we introduce a geometric approach to investigate contextual modulation in neural populations and study how the representation of stimulus orientation is transformed by the presence of a mask. We find that population responses can be embedded in a low-dimensional space and that an affine transform can account for the effects of masking. The geometric analysis further reveals a link between changes in discriminability and bias induced by the mask. We propose the geometric approach can yield new insights into the image processing computations taking place in early visual cortex at the population level while coping with the heterogeneity of single cell behavior.

An omnidirectional visualization model of personalized gene regulatory networks

Gene regulatory networks (GRNs) have been widely used as a fundamental tool to reveal the genomic mechanisms that underlie the individual's response to environmental and developmental cues. Standard approaches infer GRNs as holistic graphs of gene co-expression, but such graphs cannot quantify how gene-gene interactions vary among individuals and how they alter structurally across spatiotemporal gradients. Here, we develop a general framework for inferring informative, dynamic, omnidirectional, and personalized networks (idopNetworks) from routine transcriptional experiments. This framework is constructed by a system of quasi-dynamic ordinary differential equations (qdODEs) derived from the combination of ecological and evolutionary theories. We reconstruct idopNetworks using genomic data from a surgical experiment and illustrate how network structure is associated with surgical response to infrainguinal vein bypass grafting and the outcome of grafting. idopNetworks may shed light on genotype-phenotype relationships and provide valuable information for personalized medicine.

Co-emergence and Collapse: The Mesoscopic Approach for Conceptualizing and Investigating the Functional Integration of Organisms

The fall of reductionist approaches to explanation leaves biology with an unescapable challenge: how to decipher complex systems. This entails a number of very critical questions, the most basic ones being: "What do we mean by 'complex'?" and "What is the system we should look for?" In complex systems, constraints belong to a higher level that the molecular one and their effect reduces and constrains the manifold of the accessible internal states of the system itself. Function is related but not deterministically imposed by the underlying structure. It is quite unlikely that such kind of complexity could be grasped by current approaches focusing on a single organization scale. The natural co-emergence of systems, parts and properties can be adopted as a hypothesis-free conceptual framework to understand functional integration of organisms, including their hierarchical or multilevel patterns, and including the way scientific practice proceeds in approaching such complexity. External, "driving" factors - order parameters and control parameters provided by the surrounding microenvironment - are always required to "push" the components' fate into well-defined developmental directions. In the negative, we see that in pathological processes such as cancer, organizational fluidity, collapse of levels and dynamic heterogeneity make it hard to even find a level of observation for a stable explanandum to persist in scientific practice. Parts and the system both lose their properties once the system is destabilized. The mesoscopic approach is our proposal to conceptualizing, investigating and explaining in biology. "Mesoscopic way of thinking" is increasingly popular in the epistemology of biology and corresponds to looking for an explanation (and possibly a prediction) where "non-trivial determinism is maximal": the "most microscopic" level of organization is not necessarily the place where "the most relevant facts do happen." A fundamental re-thinking of the concept of causality is also due for order parameters to be carefully and correctly identified. In the biological realm, entities have relational properties only, as they depend ontologically on the context they happen to be in. The basic idea of a relational ontology is that, in our inventory of the world, relations are somehow prior to the relata (i.e., entities).

Ensemble Method of Feature Selection and Reverse Construction of Gene Logical Network Based on Information Entropy

In this paper, we propose a novel ensemble gene selection method to obtain a gene subset. Then we provide a reverse construction method of gene network derived from expression profile data of the gene subset. The uncertainty coefficient based on information entropy are used to define the existence of logical relations among these genes. If the uncertainty coefficient between some genes exceeds predefined thresholds, the gene nodes will be connected by directed edges. Thus, a gene network is generated, which we define as gene logical network. This method is applied to the breast cancer data including control group and experimental group, with comparisons of the 2nd-order logic type distribution, average degree as well as average path length of the networks. It is found that these structures with different networks are quite distinct. By the comparison of the degree difference between control group and experimental group, the key genes are picked up. By defining the dynamics evolution rules of state transition based on the logical regulation among the key genes in the network, the dynamic behaviors for normal breast cells and cells with cancer of different stages are simulated numerically. Some of them are highly related to the development of breast cancer through literature inquiry. The study may provide a useful revelation to the biological mechanism in the formation and development of cancer.

Boolean gene regulatory network model of centromere function in Saccharomyces cerevisiae

Centromeres, a highly conserved locus of eukaryotic chromosomes, have critical function for genome stability and integrity. Because their centromeric DNA sequences are necessary and sufficient for kinetochore recruitment and DNA segregation, point centromeres of Saccharomyces cerevisiae chromosomes provide an attractive system for the study of the regulation of centromere function. Using the mathematical model of Boolean gene regulatory networks, the gene regulatory dynamics of centromere region of S. cerevisiae (budding yeast), which is actively involved in the cell-cycle, has been examined. A gene regulatory network containing the relevant centromere genes of the model organism from biological databases was established and all possible cellular phenotypes subjected to a synchronous gene regulation and attracted to several basins. Gene expression in the largest attractor was compared with the biological data by obtaining changes in the cell-cycle. We show that the model for centromere function recovers a single cyclic attractor. The trajectory flow diagram plotted over all initial conditions of the system also shows good correspondence with the cell-cycle phases. Although other upstream signals are possibly involved in the regulation of centromere genes, proposed interactions with selected cell-cycle genes were sufficient to recover whole cell-cycle process. To truly clarify these proposed regulatory interactions of candidate genes for centromere function, profiling and analyzing their expression levels over time with expanded nodes/edges are required. Moreover, a previously modeled gene knock-down mechanism applied to the network and robustness versus knock-down was interpreted based on the obtained consequences.

Parameter Estimation and Variable Selection for Big Systems of Linear Ordinary Differential Equations: A Matrix-Based Approach

Ordinary differential equations (ODEs) are widely used to model the dynamic behavior of a complex system. Parameter estimation and variable selection for a "Big System" with linear ODEs are very challenging due to the need of nonlinear optimization in an ultra-high dimensional parameter space. In this article, we develop a parameter estimation and variable selection method based on the ideas of similarity transformation and separable least squares (SLS). Simulation studies demonstrate that the proposed matrix-based SLS method could be used to estimate the coefficient matrix more accurately and perform variable selection for a linear ODE system with thousands of dimensions and millions of parameters much better than the direct least squares (LS) method and the vector-based two-stage method that are currently available. We applied this new method to two real data sets: a yeast cell cycle gene expression data set with 30 dimensions and 930 unknown parameters and the Standard & Poor 1500 index stock price data with 1250 dimensions and 1,563,750 unknown parameters, to illustrate the utility and numerical performance of the proposed parameter estimation and variable selection method for big systems in practice.

Restructured GEO: restructuring Gene Expression Omnibus metadata for genome dynamics analysis

Motivation

Gene Expression Omnibus (GEO) and other publicly available data store their metadata in the format of unstructured English text, which is very difficult for automated reuse.

Results

We employed text mining techniques to analyze the metadata of GEO and developed Restructured GEO database (ReGEO). ReGEO reorganizes and categorizes GEO series and makes them searchable by two new attributes extracted automatically from each series' metadata. These attributes are the number of time points tested in the experiment and the disease being investigated. ReGEO also makes series searchable by other attributes available in GEO, such as platform organism, experiment type, associated PubMed ID as well as general keywords in the study's description. Our approach greatly expands the usability of GEO data, demonstrating a credible approach to improve the utility of vast amount of publicly available data in the era of Big Data research.

Inference of Large-scale Time-delayed Gene Regulatory Network with Parallel MapReduce Cloud Platform

Inference of gene regulatory network (GRN) is crucial to understand intracellular physiological activity and function of biology. The identification of large-scale GRN has been a difficult and hot topic of system biology in recent years. In order to reduce the computation load for large-scale GRN identification, a parallel algorithm based on restricted gene expression programming (RGEP), namely MPRGEP, is proposed to infer instantaneous and time-delayed regulatory relationships between transcription factors and target genes. In MPRGEP, the structure and parameters of time-delayed S-system (TDSS) model are encoded into one chromosome. An original hybrid optimization approach based on genetic algorithm (GA) and gene expression programming (GEP) is proposed to optimize TDSS model with MapReduce framework. Time-delayed GRNs (TDGRN) with hundreds of genes are utilized to test the performance of MPRGEP. The experiment results reveal that MPRGEP could infer more accurately gene regulatory network than other state-of-art methods, and obtain the convincing speedup.

Independence screening for high dimensional nonlinear additive ODE models with applications to dynamic gene regulatory networks

Mechanism-driven low-dimensional ordinary differential equation (ODE) models are often used to model viral dynamics at cellular levels and epidemics of infectious diseases. However, low-dimensional mechanism-based ODE models are limited for modeling infectious diseases at molecular levels such as transcriptomic or proteomic levels, which is critical to understand pathogenesis of diseases. Although linear ODE models have been proposed for gene regulatory networks (GRNs), nonlinear regulations are common in GRNs. The reconstruction of large-scale nonlinear networks from time-course gene expression data remains an unresolved issue. Here, we use high-dimensional nonlinear additive ODEs to model GRNs and propose a 4-step procedure to efficiently perform variable selection for nonlinear ODEs. To tackle the challenge of high dimensionality, we couple the 2-stage smoothing-based estimation method for ODEs and a nonlinear independence screening method to perform variable selection for the nonlinear ODE models. We have shown that our method possesses the sure screening property and it can handle problems with non-polynomial dimensionality. Numerical performance of the proposed method is illustrated with simulated data and a real data example for identifying the dynamic GRN of Saccharomyces cerevisiae.

fastBMA: scalable network inference and transitive reduction

Inferring genetic networks from genome-wide expression data is extremely demanding computationally. We have developed fastBMA, a distributed, parallel, and scalable implementation of Bayesian model averaging (BMA) for this purpose. fastBMA also includes a computationally efficient module for eliminating redundant indirect edges in the network by mapping the transitive reduction to an easily solved shortest-path problem. We evaluated the performance of fastBMA on synthetic data and experimental genome-wide time series yeast and human datasets. When using a single CPU core, fastBMA is up to 100 times faster than the next fastest method, LASSO, with increased accuracy. It is a memory-efficient, parallel, and distributed application that scales to human genome-wide expression data. A 10 000-gene regulation network can be obtained in a matter of hours using a 32-core cloud cluster (2 nodes of 16 cores). fastBMA is a significant improvement over its predecessor ScanBMA. It is more accurate and orders of magnitude faster than other fast network inference methods such as the 1 based on LASSO. The improved scalability allows it to calculate networks from genome scale data in a reasonable time frame. The transitive reduction method can improve accuracy in denser networks. fastBMA is available as code (M.I.T. license) from GitHub (https://github.com/lhhunghimself/fastBMA), as part of the updated networkBMA Bioconductor package (https://www.bioconductor.org/packages/release/bioc/html/networkBMA.html) and as ready-to-deploy Docker images (https://hub.docker.com/r/biodepot/fastbma/).

Inferring a nonlinear biochemical network model from a heterogeneous single-cell time course data

Mathematical modeling and analysis of biochemical reaction networks are key routines in computational systems biology and biophysics; however, it remains difficult to choose the most valid model. Here, we propose a computational framework for data-driven and systematic inference of a nonlinear biochemical network model. The framework is based on the expectation-maximization algorithm combined with particle smoother and sparse regularization techniques. In this method, a “redundant” model consisting of an excessive number of nodes and regulatory paths is iteratively updated by eliminating unnecessary paths, resulting in an inference of the most likely model. Using artificial single-cell time-course data showing heterogeneous oscillatory behaviors, we demonstrated that this algorithm successfully inferred the true network without any prior knowledge of network topology or parameter values. Furthermore, we showed that both the regulatory paths among nodes and the optimal number of nodes in the network could be systematically determined. The method presented in this study provides a general framework for inferring a nonlinear biochemical network model from heterogeneous single-cell time-course data.

Physical constraints in cell fate specification. A case in point: Microgravity and phenotypes differentiation

Data obtained by studying mammalian cells in absence of gravity strongly support the notion that cell fate specification cannot be understood according to the current molecular model. A paradigmatic case in point is provided by studying cell populations growing in absence of gravity. When the physical constraint (gravity) is 'experimentally removed', cells spontaneously allocate into two morphologically different phenotypes. Such phenomenon is likely enacted by the intrinsic stochasticity, which, in turn, is successively 'canalized' by a specific gene regulatory network. Both phenotypes are thermodynamically and functionally 'compatibles' with the new, modified environment. However, when the two cell subsets are reseeded into the 1g gravity field the two phenotypes collapse into one. Gravity constraints the system in adopting only one phenotype, not by selecting a pre-existing configuration, but more precisely shaping it de-novo through the modification of the cytoskeleton three-dimensional structure. Overall, those findings highlight how macro-scale features are irreducible to lower-scale explanations. The identification of macroscale control parameters - as those depending on the field (gravity, electromagnetic fields) or emerging from the cooperativity among the field's components (tissue stiffness, cell-to-cell connectivity) - are mandatory for assessing boundary conditions for models at lower scales, thus providing a concrete instantiation of top-down effects.

Statistical inference in mechanistic models: time warping for improved gradient matching

Inference in mechanistic models of non-linear differential equations is a challenging problem in current computational statistics. Due to the high computational costs of numerically solving the differential equations in every step of an iterative parameter adaptation scheme, approximate methods based on gradient matching have become popular. However, these methods critically depend on the smoothing scheme for function interpolation. The present article adapts an idea from manifold learning and demonstrates that a time warping approach aiming to homogenize intrinsic length scales can lead to a significant improvement in parameter estimation accuracy. We demonstrate the effectiveness of this scheme on noisy data from two dynamical systems with periodic limit cycle, a biopathway, and an application from soft-tissue mechanics. Our study also provides a comparative evaluation on a wide range of signal-to-noise ratios.

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.

Estimating time-varying directed gene regulation networks

The problem of modeling the dynamical regulation process within a gene network has been of great interest for a long time. We propose to model this dynamical system with a large number of nonlinear ordinary differential equations (ODEs), in which the regulation function is estimated directly from data without any parametric assumption. Most current research assumes the gene regulation network is static, but in reality, the connection and regulation function of the network may change with time or environment. This change is reflected in our dynamical model by allowing the regulation function varying with the gene expression and forcing this regulation function to be zero if no regulation happens. We introduce a statistical method called functional SCAD to estimate a time-varying sparse and directed gene regulation network, and simultaneously, to provide a smooth estimation of the regulation function and identify the interval in which no regulation effect exists. The finite sample performance of the proposed method is investigated in a Monte Carlo simulation study. Our method is demonstrated by estimating a time-varying directed gene regulation network of 20 genes involved in muscle development during the embryonic stage of Drosophila melanogaster.

Integrative Network-based Analysis of Colonic Detoxification Gene Expression in Ulcerative Colitis According to Smoking Status

BACKGROUNDS AND AIMS

The effect of cigarette smoking [CS] is ambivalent since smoking improves ulcerative colitis [UC] while it worsens Crohn's disease [CD]. Although this clinical relationship between inflammatory bowel disease [IBD] and tobacco is well established, only a few experimental works have investigated the effect of smoking on the colonic barrier homeostasis focusing on xenobiotic detoxification genes.

METHODS

A comprehensive and integrated comparative analysis of the global xenobiotic detoxification capacity of the normal colonic mucosa of healthy smokers [n = 8] and non-smokers [n = 9] versus the non-affected colonic mucosa of UC patients [n = 19] was performed by quantitative real-time polymerase chain reaction [qRT PCR]. The detoxification gene expression profile was analysed in CD patients [n = 18], in smoking UC patients [n = 5], and in biopsies from non-smoking UC patients cultured or not with cigarette smoke extract [n = 8].

RESULTS

Of the 244 detoxification genes investigated, 65 were dysregulated in UC patients in comparison with healthy controls or CD patients. The expression of ≥ 45/65 genes was inversed by CS in biopsies of smoking UC patients in remission and in colonic explants of UC patients exposed to cigarette smoke extract. We devised a network-based data analysis approach for differentially assessing changes in genetic interactions, allowing identification of unexpected regulatory detoxification genes that may play a major role in the beneficial effect of smoking on UC.

CONCLUSIONS

Non-inflamed colonic mucosa in UC is characterised by a specifically altered detoxification gene network, which is partially restored by tobacco. These mucosal signatures could be useful for developing new therapeutic strategies and biomarkers of drug response in UC.

Using Triplet Ordering Preferences for Estimating Causal Effects in the Analysis of Gene Expression Data

Triplet ordering preferences are used to perform Monte Carlo sampling of the posterior causal orderings originating from the analysis of gene-expression experiments involving observation as well as, usually few, interventions, like knock-outs. The performance of this sampling approach is compared to a previously used sampling via pairwise ordering preference as well as to the sampling of the full posterior distribution. For a fair comparison, the latter approach is restricted to twice the numerical effort of the triplet-based approach. This is done for artificially generated causal, i.e., directed acyclic graphs (DAGs) and for actual experimental data taken from the ROSETTA challenge. The sampling using the triplets ordering turns out to be superior to both other approaches.

Estimating varying coefficients for partial differential equation models

Partial differential equations (PDEs) are used to model complex dynamical systems in multiple dimensions, and their parameters often have important scientific interpretations. In some applications, PDE parameters are not constant but can change depending on the values of covariates, a feature that we call varying coefficients. We propose a parameter cascading method to estimate varying coefficients in PDE models from noisy data. Our estimates of the varying coefficients are shown to be consistent and asymptotically normally distributed. The performance of our method is evaluated by a simulation study and by an empirical study estimating three varying coefficients in a PDE model arising from LIDAR data.

Network Reconstruction From High-Dimensional Ordinary Differential Equations

We 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.