High temporal resolution RNA-seq time course data reveals widespread synchronous activation between mammalian lncRNAs and neighboring protein-coding genes.

The advent of massively parallel sequencing revealed extensive transcription beyond protein-coding genes, identifying tens of thousands of long noncoding RNAs (lncRNAs). Selected functional examples raised the possibility that lncRNAs, as a class, may maintain broad regulatory roles. Expression of lncRNAs is strongly linked with adjacent protein-coding gene expression, suggesting potential cis-regulatory functions. A more detailed understanding of these regulatory roles may be obtained through careful examination of the precise timing of lncRNA expression relative to adjacent protein-coding genes. Despite the diversity of reported lncRNA regulatory mechanisms, where causal cis-regulatory relationships exist, lncRNA transcription is expected to precede changes in target gene expression. Using a high temporal resolution RNA-seq time course, we profiled the expression dynamics of several thousand lncRNAs and protein-coding genes in synchronized, transitioning human cells. Our findings reveal that lncRNAs are expressed synchronously with adjacent protein-coding genes. Analysis of lipopolysaccharide-activated mouse dendritic cells revealed the same temporal relationship observed in transitioning human cells. Our findings suggest broad-scale cis-regulatory roles for lncRNAs are not common. The strong association between lncRNAs and adjacent genes may instead indicate an origin as transcriptional by-products from active protein-coding gene promoters and enhancers.

Multiscale modelling of desquamation in the interfollicular epidermis.

Maintenance of epidermal thickness is critical to the barrier function of the skin. Decreased tissue thickness, specifically in the stratum corneum (the outermost layer of the tissue), causes discomfort and inflammation, and is related to several severe diseases of the tissue. In order to maintain both stratum corneum thickness and overall tissue thickness it is necessary for the system to balance cell proliferation and cell loss. Cell proliferation in the epidermis occurs in the basal layer and causes constant upwards movement in the tissue. Cell loss occurs when dead cells at the top of the tissue are lost to the environment through a process called desquamation. Desquamation is thought to occur through a gradual reduction in adhesion between cells, due to the cleaving of adhesion proteins by enzymes, in the stratum corneum. In this paper we will investigate combining a (mass action) subcellular model of desquamation with a three dimensional (cell centre based) multicellular model of the interfollicular epidermis to better understand maintenance of epidermal thickness. Specifically, our aim is to determine if a hypothesised biological model for the degradation of cell-cell adhesion, from the literature, is sufficient to maintain a steady state tissue thickness. These investigations show the model is able to provide a consistent rate of cell loss in the multicellular model. This loss balances proliferation, and hence maintains a homeostatic tissue thickness. Moreover, we find that multiple proliferative cell populations in the basal layer can be represented by a single proliferative cell population, simplifying investigations with this model. The model is used to investigate a disorder (Netherton Syndrome) which disrupts desquamation. The model shows how biochemical changes can cause disruptions to the tissue, resulting in a reduced tissue thickness and consequently diminishing the protective role of the tissue. A hypothetical treatment result is also investigated: we compare the cases of a partially effective homogeneous treatment (where all cells partially recover) and a totally effective heterogeneous treatment (in which a proportion of the cells totally recover) with the aim to determine the difference in the response of the tissue to these different scenarios. Results show an increased benefit to corneum thickness from the heterogeneous treatment over the homogeneous treatment.

Modular assembly of dynamic models in systems biology.

It is widely acknowledged that the construction of large-scale dynamic models in systems biology requires complex modelling problems to be broken up into more manageable pieces. To this end, both modelling and software frameworks are required to enable modular modelling. While there has been consistent progress in the development of software tools to enhance model reusability, there has been a relative lack of consideration for how underlying biophysical principles can be applied to this space. Bond graphs combine the aspects of both modularity and physics-based modelling. In this paper, we argue that bond graphs are compatible with recent developments in modularity and abstraction in systems biology, and are thus a desirable framework for constructing large-scale models. We use two examples to illustrate the utility of bond graphs in this context: a model of a mitogen-activated protein kinase (MAPK) cascade to illustrate the reusability of modules and a model of glycolysis to illustrate the ability to modify the model granularity.

Hierarchical semantic composition of biosimulation models using bond graphs.

Simulating complex biological and physiological systems and predicting their behaviours under different conditions remains challenging. Breaking systems into smaller and more manageable modules can address this challenge, assisting both model development and simulation. Nevertheless, existing computational models in biology and physiology are often not modular and therefore difficult to assemble into larger models. Even when this is possible, the resulting model may not be useful due to inconsistencies either with the laws of physics or the physiological behaviour of the system. Here, we propose a general methodology for composing models, combining the energy-based bond graph approach with semantics-based annotations. This approach improves model composition and ensures that a composite model is physically plausible. As an example, we demonstrate this approach to automated model composition using a model of human arterial circulation. The major benefit is that modellers can spend more time on understanding the behaviour of complex biological and physiological systems and less time wrangling with model composition.

Maintaining the proliferative cell niche in multicellular models of epithelia.

The maintenance of the proliferative cell niche is critical to epithelial tissue morphology and function. In this paper we investigate how current modelling methods can result in the erroneous loss of proliferative cells from the proliferative cell niche. Using an established model of the inter-follicular epidermis we find there is a limit to the proliferative cell densities that can be maintained in the basal layer (the niche) if we do not include additional mechanisms to stop the loss of proliferative cells from the niche. We suggest a new methodology that enables maintenance of a desired homeostatic population of proliferative cells in the niche: a rotational force is applied to the two daughter cells during the mitotic phase of division to enforce a particular division direction. We demonstrate that this new methodology achieves this goal. This methodology reflects the regulation of the orientation of cell division.

Analysing and simulating energy-based models in biology using BondGraphTools.

Like all physical systems, biological systems are constrained by the laws of physics. However, mathematical models of biochemistry frequently neglect the conservation of energy, leading to unrealistic behaviour. Energy-based models that are consistent with conservation of mass, charge and energy have the potential to aid the understanding of complex interactions between biological components, and are becoming easier to develop with recent advances in experimental measurements and databases. In this paper, we motivate the use of bond graphs (a modelling tool from engineering) for energy-based modelling and introduce, BondGraphTools, a Python library for constructing and analysing bond graph models. We use examples from biochemistry to illustrate how BondGraphTools can be used to automate model construction in systems biology while maintaining consistency with the laws of physics.

Ca2+ Release via IP3 Receptors Shapes the Cardiac Ca2+ Transient for Hypertrophic Signaling.

Calcium (Ca2+) plays a central role in mediating both contractile function and hypertrophic signaling in ventricular cardiomyocytes. L-type Ca2+ channels trigger release of Ca2+ from ryanodine receptors for cellular contraction, whereas signaling downstream of G-protein-coupled receptors stimulates Ca2+ release via inositol 1,4,5-trisphosphate receptors (IP3Rs), engaging hypertrophic signaling pathways. Modulation of the amplitude, duration, and duty cycle of the cytosolic Ca2+ contraction signal and spatial localization have all been proposed to encode this hypertrophic signal. Given current knowledge of IP3Rs, we develop a model describing the effect of functional interaction (cross talk) between ryanodine receptor and IP3R channels on the Ca2+ transient and examine the sensitivity of the Ca2+ transient shape to properties of IP3R activation. A key result of our study is that IP3R activation increases Ca2+ transient duration for a broad range of IP3R properties, but the effect of IP3R activation on Ca2+ transient amplitude is dependent on IP3 concentration. Furthermore we demonstrate that IP3-mediated Ca2+ release in the cytosol increases the duty cycle of the Ca2+ transient, the fraction of the cycle for which [Ca2+] is elevated, across a broad range of parameter values and IP3 concentrations. When coupled to a model of downstream transcription factor (NFAT) activation, we demonstrate that there is a high correspondence between the Ca2+ transient duty cycle and the proportion of activated NFAT in the nucleus. These findings suggest increased cytosolic Ca2+ duty cycle as a plausible mechanism for IP3-dependent hypertrophic signaling via Ca2+-sensitive transcription factors such as NFAT in ventricular cardiomyocytes.

Physically-plausible modelling of biomolecular systems: A simplified, energy-based model of the mitochondrial electron transport chain.

Advances in systems biology and whole-cell modelling demand increasingly comprehensive mathematical models of cellular biochemistry. Such models require the development of simplified representations of specific processes which capture essential biophysical features but without unnecessarily complexity. Recently there has been renewed interest in thermodynamically-based modelling of cellular processes. Here we present an approach to developing of simplified yet thermodynamically consistent (hence physically plausible) models which can readily be incorporated into large scale biochemical descriptions but which do not require full mechanistic detail of the underlying processes. We illustrate the approach through development of a simplified, physically plausible model of the mitochondrial electron transport chain and show that the simplified model behaves like the full system.

A thermodynamic framework for modelling membrane transporters.

Membrane transporters contribute to the regulation of the internal environment of cells by translocating substrates across cell membranes. Like all physical systems, the behaviour of membrane transporters is constrained by the laws of thermodynamics. However, many mathematical models of transporters, especially those incorporated into whole-cell models, are not thermodynamically consistent, leading to unrealistic behaviour. In this paper we use a physics-based modelling framework, in which the transfer of energy is explicitly accounted for, to develop thermodynamically consistent models of transporters. We then apply this methodology to model two specific transporters: the cardiac sarcoplasmic/endoplasmic Ca2+ ATPase (SERCA) and the cardiac Na+/K+ ATPase.

Insights on the impact of mitochondrial organisation on bioenergetics in high-resolution computational models of cardiac cell architecture.

Recent electron microscopy data have revealed that cardiac mitochondria are not arranged in crystalline columns but are organised with several mitochondria aggregated into columns of varying sizes spanning the cell cross-section. This raises the question-how does the mitochondrial arrangement affect the metabolite distributions within cardiomyocytes and what is its impact on force dynamics? Here, we address this question by employing finite element modeling of cardiac bioenergetics on computational meshes derived from electron microscope images. Our results indicate that heterogeneous mitochondrial distributions can lead to significant spatial variation across the cell in concentrations of inorganic phosphate, creatine (Cr) and creatine phosphate (PCr). However, our model predicts that sufficient activity of the creatine kinase (CK) system, coupled with rapid diffusion of Cr and PCr, maintains near uniform ATP and ADP ratios across the cell cross sections. This homogenous distribution of ATP and ADP should also evenly distribute force production and twitch duration with contraction. These results suggest that the PCr shuttle and associated enzymatic reactions act to maintain uniform force dynamics in the cell despite the heterogeneous mitochondrial organization. However, our model also predicts that under hypoxia activity of mitochondrial CK enzymes and diffusion of high-energy phosphate compounds may be insufficient to sustain uniform ATP/ADP distribution and hence force generation.

Changes in mitochondrial morphology and organization can enhance energy supply from mitochondrial oxidative phosphorylation in diabetic cardiomyopathy.

Diabetic cardiomyopathy is accompanied by metabolic and ultrastructural alterations, but the impact of the structural changes on metabolism itself is yet to be determined. Morphometric analysis of mitochondrial shape and spatial organization within transverse sections of cardiomyocytes from control and streptozotocin-induced type I diabetic Sprague-Dawley rats revealed that mitochondria are 20% smaller in size while their spatial density increases by 53% in diabetic cells relative to control myocytes. Diabetic cells formed larger clusters of mitochondria (60% more mitochondria per cluster) and the effective surface-to-volume ratio of these clusters increased by 22.5%. Using a biophysical computational model we found that this increase can have a moderate compensatory effect by increasing the availability of ATP in the cytosol when ATP synthesis within the mitochondrial matrix is compromised.

Experimental and modelling evidence of shortening heat in cardiac muscle.

KEY POINTS: Heat associated with muscle shortening has been repeatedly demonstrated in skeletal muscle, but its existence in cardiac muscle remains contentious after five decades of study. By iterating between experiments and computational modelling, we show compelling evidence for the existence of shortening heat in cardiac muscle and reveal, mechanistically, the source of this excess heat. Our results clarify a long-standing uncertainty in the field of cardiac muscle energetics. We provide a revised partitioning of cardiac muscle energy expenditure to include this newly revealed thermal component. ABSTRACT: When a muscle shortens against an afterload, the heat that it liberates is greater than that produced by the same muscle contracting isometrically at the same level of force. This excess heat is defined as 'shortening heat', and has been repeatedly demonstrated in skeletal muscle but not in cardiac muscle. Given the micro-structural similarities between these two muscle types, and since we imagine that shortening heat is the thermal accompaniment of cross-bridge cycling, we have re-examined this issue. Using our flow-through microcalorimeter, we measured force and heat generated by isolated rat trabeculae undergoing isometric contractions at different muscle lengths and work-loop (shortening) contractions at different afterloads. We simulated these experimental protocols using a thermodynamically constrained model of cross-bridge cycling and probed the mechanisms underpinning shortening heat. Predictions generated by the model were subsequently validated by a further set of experiments. Both our experimental and modelling results show convincing evidence for the existence of shortening heat in cardiac muscle. Its magnitude is inversely related to the afterload or, equivalently, directly related to the extent of shortening. Computational simulations reveal that the heat of shortening arises from the cycling of cross-bridges, and that the rate of ATP hydrolysis is more sensitive to change of muscle length than to change of afterload. Our results clarify a long-standing uncertainty in the field of cardiac muscle energetics.

Predictive modelling of gene expression from transcriptional regulatory elements.

Predictive modelling of gene expression provides a powerful framework for exploring the regulatory logic underpinning transcriptional regulation. Recent studies have demonstrated the utility of such models in identifying dysregulation of gene and miRNA expression associated with abnormal patterns of transcription factor (TF) binding or nucleosomal histone modifications (HMs). Despite the growing popularity of such approaches, a comparative review of the various modelling algorithms and feature extraction methods is lacking. We define and compare three methods of quantifying pairwise gene-TF/HM interactions and discuss their suitability for integrating the heterogeneous chromatin immunoprecipitation (ChIP)-seq binding patterns exhibited by TFs and HMs. We then construct log-linear and ϵ-support vector regression models from various mouse embryonic stem cell (mESC) and human lymphoblastoid (GM12878) data sets, considering both ChIP-seq- and position weight matrix- (PWM)-derived in silico TF-binding. The two algorithms are evaluated both in terms of their modelling prediction accuracy and ability to identify the established regulatory roles of individual TFs and HMs. Our results demonstrate that TF-binding and HMs are highly predictive of gene expression as measured by mRNA transcript abundance, irrespective of algorithm or cell type selection and considering both ChIP-seq and PWM-derived TF-binding. As we encourage other researchers to explore and develop these results, our framework is implemented using open-source software and made available as a preconfigured bootable virtual environment.

Virtual Reference Environments: a simple way to make research reproducible.

'Reproducible research' has received increasing attention over the past few years as bioinformatics and computational biology methodologies become more complex. Although reproducible research is progressing in several valuable ways, we suggest that recent increases in internet bandwidth and disk space, along with the availability of open-source and free-software licences for tools, enable another simple step to make research reproducible. In this article, we urge the creation of minimal virtual reference environments implementing all the tools necessary to reproduce a result, as a standard part of publication. We address potential problems with this approach, and show an example environment from our own work.

Distributed gene expression modelling for exploring variability in epigenetic function.

BACKGROUND: Predictive gene expression modelling is an important tool in computational biology due to the volume of high-throughput sequencing data generated by recent consortia. However, the scope of previous studies has been restricted to a small set of cell-lines or experimental conditions due an inability to leverage distributed processing architectures for large, sharded data-sets. RESULTS: We present a distributed implementation of gene expression modelling using the MapReduce paradigm and prove that performance improves as a linear function of available processor cores. We then leverage the computational efficiency of this framework to explore the variability of epigenetic function across fifty histone modification data-sets from variety of cancerous and non-cancerous cell-lines. CONCLUSIONS: We demonstrate that the genome-wide relationships between histone modifications and mRNA transcription are lineage, tissue and karyotype-invariant, and that models trained on matched -omics data from non-cancerous cell-lines are able to predict cancerous expression with equivalent genome-wide fidelity.

Systems analysis identifies miR-29b regulation of invasiveness in melanoma.

BACKGROUND: In many cancers, microRNAs (miRs) contribute to metastatic progression by modulating phenotypic reprogramming processes such as epithelial-mesenchymal plasticity. This can be driven by miRs targeting multiple mRNA transcripts, inducing regulated changes across large sets of genes. The miR-target databases TargetScan and DIANA-microT predict putative relationships by examining sequence complementarity between miRs and mRNAs. However, it remains a challenge to identify which miR-mRNA interactions are active at endogenous expression levels, and of biological consequence. METHODS: We developed a workflow to integrate TargetScan and DIANA-microT predictions into the analysis of data-driven associations calculated from transcript abundance (RNASeq) data, specifically the mutual information and Pearson's correlation metrics. We use this workflow to identify putative relationships of miR-mediated mRNA repression with strong support from both lines of evidence. Applying this approach systematically to a large, published collection of unique melanoma cell lines - the Ludwig Melbourne melanoma (LM-MEL) cell line panel - we identified putative miR-mRNA interactions that may contribute to invasiveness. This guided the selection of interactions of interest for further in vitro validation studies. RESULTS: Several miR-mRNA regulatory relationships supported by TargetScan and DIANA-microT demonstrated differential activity across cell lines of varying matrigel invasiveness. Strong negative statistical associations for these putative regulatory relationships were consistent with target mRNA inhibition by the miR, and suggest that differential activity of such miR-mRNA relationships contribute to differences in melanoma invasiveness. Many of these relationships were reflected across the skin cutaneous melanoma TCGA dataset, indicating that these observations also show graded activity across clinical samples. Several of these miRs are implicated in cancer progression (miR-211, -340, -125b, -221, and -29b). The specific role for miR-29b-3p in melanoma has not been well studied. We experimentally validated the predicted miR-29b-3p regulation of LAMC1 and PPIC and LASP1, and show that dysregulation of miR-29b-3p or these mRNA targets can influence cellular invasiveness in vitro. CONCLUSIONS: This analytic strategy provides a comprehensive, systems-level approach to identify miR-mRNA regulation in high-throughput cancer data, identifies novel putative interactions with functional phenotypic relevance, and can be used to direct experimental resources for subsequent experimental validation. Computational scripts are available: http://github.com/uomsystemsbiology/LMMEL-miR-miner.

Myocardial energetics is not compromised during compensated hypertrophy in the Dahl salt-sensitive rat model of hypertension.

Salt-induced hypertension leads to development of left ventricular hypertrophy in the Dahl salt-sensitive (Dahl/SS) rat. Before progression to left ventricular failure, the heart initially undergoes a compensated hypertrophic response. We hypothesized that changes in myocardial energetics may be an early indicator of transition to failure. Dahl/SS rats and their salt-resistant consomic controls (SS-13(BN)) were placed on either a low- or high-salt diet to generate four cohorts: Dahl-SS rats on a low- (Dahl-LS) or high-salt diet (Dahl-HS), and SS-13(BN) rats on a low- (SSBN-LS) or high-salt diet (SSBN-HS). We isolated left ventricular trabeculae and characterized their mechanoenergetic performance. Our results show, at most, modest effects of salt-induced compensated hypertrophy on myocardial energetics. We found that the Dahl-HS cohort had a higher work-loop heat of activation (estimated from the intercept of the heat vs. relative afterload relationship generated from work-loop contractions) relative to the SSBN-HS cohort and a higher economy of contraction (inverse of the slope of the heat vs. active stress relation) relative to the Dahl-LS cohort. The maximum extent of shortening and maximum shortening velocity of the Dahl/SS groups were higher than those of the SS-13(BN) groups. Despite these differences, no significant effect of salt-induced hypertension was observed for either peak work output or peak mechanical efficiency during compensated hypertrophy.

NAIL, a software toolset for inferring, analyzing and visualizing regulatory networks.

UNLABELLED: The wide variety of published approaches for the problem of regulatory network inference makes using multiple inference algorithms complex and time-consuming. Network Analysis and Inference Library (NAIL) is a set of software tools to simplify the range of computational activities involved in regulatory network inference. It uses a modular approach to connect different network inference algorithms to the same visualization and network-based analyses. NAIL is technology-independent and includes an interface layer to allow easy integration of components into other applications. AVAILABILITY AND IMPLEMENTATION: NAIL is implemented in MATLAB, runs on Windows, Linux and OSX, and is available from SourceForge at https://sourceforge.net/projects/nailsystemsbiology/ for all researchers to use. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

A Framework to Account for Sedimentation and Diffusion in Particle-Cell Interactions.

In vitro experiments provide a solid basis for understanding the interactions between particles and biological systems. An important confounding variable for these studies is the difference between the amount of particles administered and that which reaches the surface of cells. Here, we engineer a hydrogel-based nanoparticle system and combine in situ characterization techniques, 3D-printed cell cultures, and computational modeling to evaluate and study particle-cell interactions of advanced particle systems. The framework presented demonstrates how sedimentation and diffusion can explain differences in particle-cell association, and provides a means to account for these effects. Finally, using in silico modeling, we predict the proportion of particles that reaches the cell surface using common experimental conditions for a wide range of inorganic and organic micro- and nanoparticles. This work can assist in the understanding and control of sedimentation and diffusion when investigating cellular interactions of engineered particles.

Examination of the Effects of Heterogeneous Organization of RyR Clusters, Myofibrils and Mitochondria on Ca2+ Release Patterns in Cardiomyocytes.

Spatio-temporal dynamics of intracellular calcium, [Ca2+]i, regulate the contractile function of cardiac muscle cells. Measuring [Ca2+]i flux is central to the study of mechanisms that underlie both normal cardiac function and calcium-dependent etiologies in heart disease. However, current imaging techniques are limited in the spatial resolution to which changes in [Ca2+]i can be detected. Using spatial point process statistics techniques we developed a novel method to simulate the spatial distribution of RyR clusters, which act as the major mediators of contractile Ca2+ release, upon a physiologically-realistic cellular landscape composed of tightly-packed mitochondria and myofibrils. We applied this method to computationally combine confocal-scale (~ 200 nm) data of RyR clusters with 3D electron microscopy data (~ 30 nm) of myofibrils and mitochondria, both collected from adult rat left ventricular myocytes. Using this hybrid-scale spatial model, we simulated reaction-diffusion of [Ca2+]i during the rising phase of the transient (first 30 ms after initiation). At 30 ms, the average peak of the simulated [Ca2+]i transient and of the simulated fluorescence intensity signal, F/F0, reached values similar to that found in the literature ([Ca2+]i ≈1 µM; F/F0≈5.5). However, our model predicted the variation in [Ca2+]i to be between 0.3 and 12.7 µM (~3 to 100 fold from resting value of 0.1 µM) and the corresponding F/F0 signal ranging from 3 to 9.5. We demonstrate in this study that: (i) heterogeneities in the [Ca2+]i transient are due not only to heterogeneous distribution and clustering of mitochondria; (ii) but also to heterogeneous local densities of RyR clusters. Further, we show that: (iii) these structure-induced heterogeneities in [Ca2+]i can appear in line scan data. Finally, using our unique method for generating RyR cluster distributions, we demonstrate the robustness in the [Ca2+]i transient to differences in RyR cluster distributions measured between rat and human cardiomyocytes.