We develop a new regression algorithm, cMIKANA, for inference of gene regulatory networks from combinations of steady-state and time-series gene expression data. Using simulated gene expression datasets to assess the accuracy of reconstructing gene regulatory networks, we show that steady-state and time-series data sets can successfully be combined to identify gene regulatory interactions using the new algorithm. Inferring gene networks from combined data sets was found to be advantageous when using noisy measurements collected with either lower sampling rates or a limited number of experimental replicates. We illustrate our method by applying it to a microarray gene expression dataset from human umbilical vein endothelial cells (HUVECs) which combines time series data from treatment with growth factor TNF and steady state data from siRNA knockdown treatments. Our results suggest that the combination of steady-state and time-series datasets may provide better prediction of RNA-to-RNA interactions, and may also reveal biological features that cannot be identified from dynamic or steady state information alone. Finally, we consider the experimental design of genomics experiments for gene regulatory network inference and show that network inference can be improved by incorporating steady-state measurements with time-series data.
An understanding of Ca2+ signalling in saliva-secreting acinar cells is important, as Ca2+ is the second messenger linking stimulation of cells to production of saliva. Ca2+ signals effect secretion via the ion channels located both apically and basolaterally in the cell. By approximating Ca2+ waves with periodic functions on the apical and basolateral membranes, we isolate individual wave properties and investigate them for their effect on fluid secretion in a mathematical model of the acinar cell. Mean Ca2+ concentration is found to be the most significant property in signalling secretion. Wave speed was found to encode a range of secretion rates. Ca2+ oscillation frequency and amplitude had little effect on fluid secretion.
mathematical model; parotid acinar cell; oscillation frequency; Ca2+ wave speed; calcium signalling
Many studies have revealed correlations between breast tumour phenotypes, variations in gene expression, and patient survival outcomes. The molecular heterogeneity between breast tumours revealed by these studies has allowed prediction of prognosis and has underpinned stratified therapy, where groups of patients with particular tumour types receive specific treatments. The molecular tests used to predict prognosis and stratify treatment usually utilise fixed sets of genomic biomarkers, with the same biomarker sets being used to test all patients. In this paper we suggest that instead of fixed sets of genomic biomarkers, it may be more effective to use a stratified biomarker approach, where optimal biomarker sets are automatically chosen for particular patient groups, analogous to the choice of optimal treatments for groups of similar patients in stratified therapy. We illustrate the effectiveness of a biclustering approach to select optimal gene sets for determining the prognosis of specific strata of patients, based on potentially overlapping, non-discrete molecular characteristics of tumours.
Biclustering identified tightly co-expressed gene sets in the tumours of restricted subgroups of breast cancer patients. The co-expressed genes in these biclusters were significantly enriched for particular biological annotations and gene regulatory modules associated with breast cancer biology. Tumours identified within the same bicluster were more likely to present with similar clinical features. Bicluster membership combined with clinical information could predict patient prognosis in conditional inference tree and ridge regression class prediction models.
The increasing clinical use of genomic profiling demands identification of more effective methods to segregate patients into prognostic and treatment groups. We have shown that biclustering can be used to select optimal gene sets for determining the prognosis of specific strata of patients.
Biclustering; Gene expression profiles; Tumour classification; Survival prediction; Breast cancer
Secretion from the salivary glands is driven by osmosis following the establishment of osmotic gradients between the lumen, the cell and the interstitium by active ion transport. We consider a dynamic model of osmotically-driven primary saliva secretion, and use singular perturbation approaches and scaling assumptions to reduce the model. Our analysis shows that isosmotic secretion is the most efficient secretion regime, and that this holds for single isolated cells and for multiple cells assembled into an acinus. For typical parameter variations, we rule out any significant synergistic effect on total water secretion of an acinar arrangement of cells about a single shared lumen. Conditions for the attainment of isosmotic secretion are considered, and we derive an expression for how the concentration gradient between the interstitium and the lumen scales with water and chloride transport parameters. Aquaporin knockout studies are interpreted in the context of our analysis and further investigated using simulations of transport efficiency with different membrane water permeabilities. We conclude that recent claims that aquaporin knockout studies can be interpreted as evidence against a simple osmotic mechanism are not supported by our work. Many of the results that we obtain are independent of specific transporter details, and our analysis can be easily extended to apply to models that use other proposed ionic mechanisms of saliva secretion.
The cell division cycle and the circadian clock represent two major cellular rhythms. These two periodic processes are coupled in multiple ways, given that several molecular components of the cell cycle network are controlled in a circadian manner. For example, in the network of cyclin-dependent kinases (Cdks) that governs progression along the successive phases of the cell cycle, the synthesis of the kinase Wee1, which inhibits the G2/M transition, is enhanced by the complex CLOCK-BMAL1 that plays a central role in the circadian clock network. Another component of the latter network, REV-ERBα, inhibits the synthesis of the Cdk inhibitor p21. Moreover, the synthesis of the oncogene c-Myc, which promotes G1 cyclin synthesis, is repressed by CLOCK-BMAL1. Using detailed computational models for the two networks we investigate the conditions in which the mammalian cell cycle can be entrained by the circadian clock. We show that the cell cycle can be brought to oscillate at a period of 24 h or 48 h when its autonomous period prior to coupling is in an appropriate range. The model indicates that the combination of multiple modes of coupling does not necessarily facilitate entrainment of the cell cycle by the circadian clock. Entrainment can also occur as a result of circadian variations in the level of a growth factor controlling entry into G1. Outside the range of entrainment, the coupling to the circadian clock may lead to disconnected oscillations in the cell cycle and the circadian system, or to complex oscillatory dynamics of the cell cycle in the form of endoreplication, complex periodic oscillations or chaos. The model predicts that the transition from entrainment to 24 h or 48 h might occur when the strength of coupling to the circadian clock or the level of growth factor decrease below critical values.
The cell cycle and the circadian clock are two major cellular rhythms. These two periodic processes are tightly coupled through multiple regulatory interactions; several components of the cell cycle machinery are indeed controlled by the circadian network. By using detailed computational models for the cell cycle and circadian networks we investigate the conditions in which the mammalian cell cycle can be entrained by the circadian clock. We show that entrainment to a circadian period can occur when the period of the cell cycle prior to coupling is either smaller or larger than 24 h. Entrainment to 48 h can also be observed. The presence of multiple modes of coupling does not enlarge the domain of entrainment. Coupling to the circadian clock may also lead to complex oscillatory dynamics of the cell cycle in the form of endoreplication, complex periodic oscillations, or chaotic oscillations. The model predicts that entrainment of the cell cycle could also result from the circadian variation of a growth factor gating entry into G1, and that the transition from an entrained period of 24 h to 48 h might result from a decrease in coupling strength or in the level of growth factor.
Many organs of higher organisms, such as the vascular system, lung, kidney, pancreas, liver and glands, are heavily branched structures. The branching process during lung development has been studied in great detail and is remarkably stereotyped. The branched tree is generated by the sequential, non-random use of three geometrically simple modes of branching (domain branching, planar and orthogonal bifurcation). While many regulatory components and local interactions have been defined an integrated understanding of the regulatory network that controls the branching process is lacking. We have developed a deterministic, spatio-temporal differential-equation based model of the core signaling network that governs lung branching morphogenesis. The model focuses on the two key signaling factors that have been identified in experiments, fibroblast growth factor (FGF10) and sonic hedgehog (SHH) as well as the SHH receptor patched (Ptc). We show that the reported biochemical interactions give rise to a Schnakenberg-type Turing patterning mechanisms that allows us to reproduce experimental observations in wildtype and mutant mice. The kinetic parameters as well as the domain shape are based on experimental data where available. The developed model is robust to small absolute and large relative changes in the parameter values. At the same time there is a strong regulatory potential in that the switching between branching modes can be achieved by targeted changes in the parameter values. We note that the sequence of different branching events may also be the result of different growth speeds: fast growth triggers lateral branching while slow growth favours bifurcations in our model. We conclude that the FGF10-SHH-Ptc1 module is sufficient to generate pattern that correspond to the observed branching modes.
Most organs of higher organisms, such as the vascular system, lung, kidney, pancreas, liver and glands, are heavily branched structures. The branching process during lung development has been studied in great detail and is remarkably stereotyped. The branched tree is generated by the sequential, non-random use of three geometrically simple modes of branching. While the branching sequence is identical in mice of identical genetic background it differs between mouse strains. This suggests that the positioning of branch points and the type of branching sensitively depends on information encoded in the genome. Encoding every branching point independently in the genome would require a large number of genes, and it is more likely that a recursive, self-organized process exists that determines the patterning. While many regulatory molecules have been identified an integrated understanding of the regulatory network (program) is missing. Based on available experimental data we have developed a model for lung branching. The model correctly predicts branching phenotypes in mutants and suggests that also the growth speed of the lung tip can affect the positioning and type of the next branching event.
One of the main challenges in the biomedical sciences is the determination of reaction mechanisms that constitute a biochemical pathway. During the last decades, advances have been made in building complex diagrams showing the static interactions of proteins. The challenge for systems biologists is to build realistic models of the dynamical behavior of reactants, intermediates and products. For this purpose, several methods have been recently proposed to deduce the reaction mechanisms or to estimate the kinetic parameters of the elementary reactions that constitute the pathway. One such method is MIKANA: Method to Infer Kinetics And Network Architecture. MIKANA is a computational method to infer both reaction mechanisms and estimate the kinetic parameters of biochemical pathways from time course data. To make it available to the scientific community, we developed a Graphical User Interface (GUI) for MIKANA. Among other features, the GUI validates and processes an input time course data, displays the inferred reactions, generates the differential equations for the chemical species in the pathway and plots the prediction curves on top of the input time course data. We also added a new feature to MIKANA that allows the user to exclude a priori known reactions from the inferred mechanism. This addition improves the performance of the method. In this article, we illustrate the GUI for MIKANA with three examples: an irreversible Michaelis–Menten reaction mechanism; the interaction map of chemical species of the muscle glycolytic pathway; and the glycolytic pathway of Lactococcus lactis. We also describe the code and methods in sufficient detail to allow researchers to further develop the code or reproduce the experiments described. The code for MIKANA is open source, free for academic and non-academic use and is available for download (Information S1).
We construct a mathematical model of the parotid acinar cell with the aim of investigating how the distribution of K+ and Cl− channels affects saliva production. Secretion of fluid is initiated by Ca2+ signals acting the Ca2+ dependent K+ and Cl− channels. The opening of these channels facilitates the movement of Cl− ions into the lumen which water follows by osmosis. We use recent results into both the release of Ca2+ from internal stores via the inositol (1,4,5)-trisphosphate receptor (IP3R) and IP3 dynamics to create a physiologically realistic Ca2+ model which is able to recreate important experimentally observed behaviours seen in parotid acinar cells. We formulate an equivalent electrical circuit diagram for the movement of ions responsible for water flow which enables us to calculate and include distinct apical and basal membrane potentials to the model. We show that maximum saliva production occurs when a small amount of K+ conductance is located at the apical membrane, with the majority in the basal membrane. The maximum fluid output is found to coincide with a minimum in the apical membrane potential. The traditional model whereby all Cl− channels are located in the apical membrane is shown to be the most efficient Cl− channel distribution.
mathematical model; parotid acinar cell; K+ channel; IP3 receptor; calcium oscillations
Identifying transcription factor (TF) binding sites (TFBSs) is an important step towards understanding transcriptional regulation. A common approach is to use gaplessly aligned, experimentally supported TFBSs for a particular TF, and algorithmically search for more occurrences of the same TFBSs. The largest publicly available databases of TF binding specificities contain models which are represented as position weight matrices (PWM). There are other methods using more sophisticated representations, but these have more limited databases, or aren't publicly available. Therefore, this paper focuses on methods that search using one PWM per TF. An algorithm, MATCHTM, for identifying TFBSs corresponding to a particular PWM is available, but is not based on a rigorous statistical model of TF binding, making it difficult to interpret or adjust the parameters and output of the algorithm. Furthermore, there is no public description of the algorithm sufficient to exactly reproduce it. Another algorithm, MAST, computes a p-value for the presence of a TFBS using true probabilities of finding each base at each offset from that position. We developed a statistical model, BaSeTraM, for the binding of TFs to TFBSs, taking into account random variation in the base present at each position within a TFBS. Treating the counts in the matrices and the sequences of sites as random variables, we combine this TFBS composition model with a background model to obtain a Bayesian classifier. We implemented our classifier in a package (SBaSeTraM). We tested SBaSeTraM against a MATCHTM implementation by searching all probes used in an experimental Saccharomyces cerevisiae TF binding dataset, and comparing our predictions to the data. We found no statistically significant differences in sensitivity between the algorithms (at fixed selectivity), indicating that SBaSeTraM's performance is at least comparable to the leading currently available algorithm. Our software is freely available at: http://wiki.github.com/A1kmm/sbasetram/building-the-tools.
Predicting information about human physiology and pathophysiology from genomic data is a compelling, but unfulfilled goal of post-genomic biology. This is the aim of the so-called Physiome Project and is, undeniably, an ambitious goal. Yet if we can exploit even a small proportion of the rich and varied experimental data currently available, significant insights into clinically important aspects of human physiology will follow. To achieve this requires the integration of data from disparate sources into a common framework. Extrapolation of available data across species, laboratory techniques and conditions requires a quantitative approach. Mathematical models allow us to integrate molecular information into cellular, tissue and organ-level, and ultimately clinically relevant scales. In this paper we argue that biophysically detailed computational modelling provides the essential tool for this process and, furthermore, that an appropriate framework for annotating, databasing and critiquing these models will be essential for the development of integrative computational biology.
physiome; mathematical modelling; cardiac; multi-scale
Salivary fluid secretion is crucial for preventing problems such as dryness of mouth, difficulty with mastication and swallowing, as well as oral pain and dental cavities. Fluid flow is driven primarily by the transepithelial movement of chloride and sodium ions into the parotid acinus lumen. The activation of Cl− channels is calcium dependent, with the average elevated calcium concentration during calcium oscillations increasing the conductance of the channels, leading to an outflow of Cl−. The accumulation of NaCl in the lumen drives water flow by osmosis. We construct a mathematical model of the calcium concentration oscillations and couple this to a model for Cl− efflux. We also construct a model governing fluid flow in an isolated parotid acinar cell, which includes a description of the rate of change of intracellular ion concentrations, cell volume, membrane potential and water flow rate. We find that [Ca2+] oscillations lead to oscillations in fluid flow, and that the rate of fluid flow is regulated by the average calcium concentration and not the frequency of the oscillations.
Mathematical model; Salivary fluid secretion; Parotid acinar cells; intracellular calcium; Cl− channels and fluxes
Complexes of chemoreceptors in the bacterial cytoplasmic membrane allow for the sensing of ligands with remarkable sensitivity. Despite the excellent characterization of the chemotaxis signaling network, very little is known about what controls receptor complex size. Here we use in vitro signaling data to model the distribution of complex sizes. In particular, we model Tar receptors in membranes as an ensemble of different sized oligomer complexes, i.e., receptor dimers, dimers of dimers, and trimers of dimers, where the relative free energies, including receptor modification, ligand binding, and interaction with the kinase CheA determine the size distribution. Our model compares favorably with a variety of signaling data, including dose-response curves of receptor activity and the dependence of activity on receptor density in the membrane. We propose that the kinetics of complex assembly can be measured in vitro from the temporal response to a perturbation of the complex free energies, e.g., by addition of ligand.
Chemotaxis allows bacteria to sense and swim toward nutrients and away from toxins. The remarkable sensing properties of the chemotaxis network, such as high sensitivity to small changes in the chemical environment, are thought to originate from receptor complexes in the membrane, which act as antennas to magnify weak signals. To adapt to persistent stimulation, receptors are covalently modified. While the individual protein components of the chemotaxis network are well characterized, making the system well suited for quantitative and computational analysis, direct experimental visualization of receptors and receptor complexes is difficult within the current limits of fluorescence and electron microscopy. To address questions such as how large are complexes and why do they assemble, we analyze in vitro signaling data using a previously developed model of signaling by receptor complexes. Based on the data, we propose a statistical physics model for the distribution of complex sizes in the membrane. Within this model, complex size depends on the receptor free energy with contributions from receptor modification level, ligand binding, receptor–receptor coupling, and binding to accessory proteins. Our model results compare favorably with a variety of different signaling data, and suggest new experiments to measure the kinetics of assembly of receptor complexes.
The growth of scleractinian corals is strongly influenced by the effect of water motion. Corals are known to have a high level of phenotypic variation and exhibit a diverse range of growth forms, which often contain a high level of geometric complexity. Due to their complex shape, simulation models represent an important option to complement experimental studies of growth and flow. In this work, we analyzed the impact of flow on coral's morphology by an accretive growth model coupled with advection-diffusion equations. We performed simulations under no-flow and uni-directional flow setup with the Reynolds number constant. The relevant importance of diffusion to advection was investigated by varying the diffusion coefficient, rather than the flow speed in Péclet number. The flow and transport equations were coupled and solved using COMSOL Multiphysics. We then compared the simulated morphologies with a series of Computed Tomography (CT) scans of scleractinian corals Pocillopora verrucosa exposed to various flow conditions in the in situ controlled flume setup. As a result, we found a similar trend associated with the increasing Péclet for both simulated forms and in situ corals; that is uni-directional current tends to facilitate asymmetrical growth response resulting in colonies with branches predominantly developed in the upstream direction. A closer look at the morphological traits yielded an interesting property about colony symmetry and plasticity induced by uni-directional flow. Both simulated and in situ corals exhibit a tendency where the degree of symmetry decreases and compactification increases in conjunction with the augmented Péclet thus indicates the significant importance of hydrodynamics.
A long-standing question in marine biology and coral biology is the morphological plasticity of corals, sponges and other marine sessile organisms and the influence of water movement. Usually branching species tend to develop symmetrical colonies where branches are being formed in all directions. There is a long standing discussion if this process in which colonies develop symmetrical colonies is controlled by genes or by the environment. In this work, we address this question for the scleractinian coral Pocillopora verrucosa. We first have acquired coral colonies from a controlled in-situ flow experiment where the coral was growing under uni-directional flow conditions. The corals colonies were scanned using a Computed Tomography (CT) technique used for medical imaging and industrial imaging. We have developed a simulation for the growth and form of corals and the influence of water movement. We have compared the simulated morphologies to the three dimensional images obtained with the CT scanner. We have found that coral's branches predominantly develop in the upstream part of the colony and an asymmetrical colony is being formed under uni-directional flow conditions. Our results confirm that growth of the coral is strongly influenced by the flow conditions.
To characterize intracellular energy transfer in the heart, two organ-level methods have frequently been employed: inversion and saturation transfer, and dynamic labeling. Creatine kinase (CK) fluxes obtained by following oxygen labeling have been considerably smaller than the fluxes determined by saturation transfer. It has been proposed that dynamic labeling determines net flux through CK shuttle, whereas saturation transfer measures total unidirectional flux. However, to our knowledge, no sensitivity analysis of flux determination by oxygen labeling has been performed, limiting our ability to compare flux distributions predicted by different methods. Here we analyze oxygen labeling in a physiological heart phosphotransfer network with active CK and adenylate kinase (AdK) shuttles and establish which fluxes determine the labeling state. A mathematical model consisting of a system of ordinary differential equations was composed describing enrichment in each phosphoryl group and inorganic phosphate. By varying flux distributions in the model and calculating the labeling, we analyzed labeling sensitivity to different fluxes in the heart. We observed that the labeling state is predominantly sensitive to total unidirectional CK and AdK fluxes and not to net fluxes. We conclude that measuring dynamic incorporation of into the high-energy phosphotransfer network in heart does not permit unambiguous determination of energetic fluxes with a higher magnitude than the ATP synthase rate when the bidirectionality of fluxes is taken into account. Our analysis suggests that the flux distributions obtained using dynamic labeling, after removing the net flux assumption, are comparable with those from inversion and saturation transfer.
In heart, the movement of energy metabolites between force-producing myosin, other ATPases, and mitochondria is vital for its function and closely related to heart pathologies. In addition to diffusion, transport of ATP, ADP, Pi, and phosphocreatine occurs along parallel pathways such as the adenylate kinase and creatine kinase shuttles. Two organ-level methods have been developed to study the relative flux through these pathways. However, their results differ. It was recently demonstrated that studies often suffer from the exclusion of compartmentation from their metabolic models. One study overcame this limitation by using compartmental models and statistical methods on multiple experiments. Here, we analyzed the sensitivity of the other method - dynamic labeling of phosphoryl groups and inorganic phosphate. For that, we composed a mathematical model tracking enrichment of the metabolites and evaluated sensitivity of labeling to different flux distribution scenarios. Our study shows that the dynamic method provides a measure of total flux, and not net flux as presumed previously, making the fluxes predicted from both methods consistent. Importantly, conclusions derived on the basis of labeling analysis, particularly those regarding the net flux through the shuttles in control and pathological cases, need to be reevaluated.
Mitochondrial dynamics and mitophagy play a key role in ensuring mitochondrial quality control. Impairment thereof was proposed to be causative to neurodegenerative diseases, diabetes, and cancer. Accumulation of mitochondrial dysfunction was further linked to aging. Here we applied a probabilistic modeling approach integrating our current knowledge on mitochondrial biology allowing us to simulate mitochondrial function and quality control during aging in silico. We demonstrate that cycles of fusion and fission and mitophagy indeed are essential for ensuring a high average quality of mitochondria, even under conditions in which random molecular damage is present. Prompted by earlier observations that mitochondrial fission itself can cause a partial drop in mitochondrial membrane potential, we tested the consequences of mitochondrial dynamics being harmful on its own. Next to directly impairing mitochondrial function, pre-existing molecular damage may be propagated and enhanced across the mitochondrial population by content mixing. In this situation, such an infection-like phenomenon impairs mitochondrial quality control progressively. However, when imposing an age-dependent deceleration of cycles of fusion and fission, we observe a delay in the loss of average quality of mitochondria. This provides a rational why fusion and fission rates are reduced during aging and why loss of a mitochondrial fission factor can extend life span in fungi. We propose the ‘mitochondrial infectious damage adaptation’ (MIDA) model according to which a deceleration of fusion–fission cycles reflects a systemic adaptation increasing life span.
Mitochondria are organelles that play a central role as ‘cellular power plants’. The cellular organization of these organelles involves a dynamic spatial network where mitochondria constantly undergo fusion and fission associated with the mixing of their molecular content. Together with the processes of mitophagy and biogenesis of mitochondrial mass, this results into a cellular surveillance system for maintaining their bioenergetic quality. The accumulation of molecular damage in mitochondria is associated with various human disorders and with aging. However, how these processes affect aging and how they can be reconciled with existing aging theories is just at the beginning to be considered. Mathematical modeling allows simulating the dynamics of mitochondrial quality control during aging in silico and leads to the ‘mitochondrial infectious damage adaptation’ (MIDA) model of aging. It reconciles a number of counterintuitive observations obtained during the last decade including infection-like processes of molecular damage spread, the reduction of fusion and fission rates during cellular aging, and observed life span extension for reduced mitochondrial fission. Interestingly, the MIDA model suggests that a reduction in mitochondrial dynamics rather than being merely a sign or even cause of aging, may actually reflect a systemic adaptation to prolong organismic life span.
Choroidal neovascularization (CNV) of the macular area of the retina is the major cause of severe vision loss in adults. In CNV, after choriocapillaries initially penetrate Bruch's membrane (BrM), invading vessels may regress or expand (CNV initiation). Next, during Early and Late CNV, the expanding vasculature usually spreads in one of three distinct patterns: in a layer between BrM and the retinal pigment epithelium (sub-RPE or Type 1 CNV), in a layer between the RPE and the photoreceptors (sub-retinal or Type 2 CNV) or in both loci simultaneously (combined pattern or Type 3 CNV). While most studies hypothesize that CNV primarily results from growth-factor effects or holes in BrM, our three-dimensional simulations of multi-cell model of the normal and pathological maculae recapitulate the three growth patterns, under the hypothesis that CNV results from combinations of impairment of: 1) RPE-RPE epithelial junctional adhesion, 2) Adhesion of the RPE basement membrane complex to BrM (RPE-BrM adhesion), and 3) Adhesion of the RPE to the photoreceptor outer segments (RPE-POS adhesion). Our key findings are that when an endothelial tip cell penetrates BrM: 1) RPE with normal epithelial junctions, basal attachment to BrM and apical attachment to POS resists CNV. 2) Small holes in BrM do not, by themselves, initiate CNV. 3) RPE with normal epithelial junctions and normal apical RPE-POS adhesion, but weak adhesion to BrM (e.g. due to lipid accumulation in BrM) results in Early sub-RPE CNV. 4) Normal adhesion of RBaM to BrM, but reduced apical RPE-POS or epithelial RPE-RPE adhesion (e.g. due to inflammation) results in Early sub-retinal CNV. 5) Simultaneous reduction in RPE-RPE epithelial binding and RPE-BrM adhesion results in either sub-RPE or sub-retinal CNV which often progresses to combined pattern CNV. These findings suggest that defects in adhesion dominate CNV initiation and progression.
This paper tests hypotheses for the mechanisms of choroidal neovascularization (CNV), the pathological growth of capillaries in response to physical defects in a structured tissue, the retina, showing that previously neglected cell-cell, cell-ECM and ECM-ECM adhesion failures suffice to determine the loci and progression of neovascularization. Surprisingly, a simple theory based on classes of adhesion failures, which involve variation of only five parameters, can coherently explain the heterogeneous range of CNV growth patterns and dynamics. Our results are generally applicable to other types of tissues where capillaries are close to an epithelium, e.g., lung and gut.
Our understanding of the molecular pathways that underlie melanoma remains incomplete. Although several published microarray studies of clinical melanomas have provided valuable information, we found only limited concordance between these studies. Therefore, we took an in vitro functional genomics approach to understand melanoma molecular pathways.
Affymetrix microarray data were generated from A375 melanoma cells treated in vitro with siRNAs against 45 transcription factors and signaling molecules. Analysis of this data using unsupervised hierarchical clustering and Bayesian gene networks identified proliferation-association RNA clusters, which were co-ordinately expressed across the A375 cells and also across melanomas from patients. The abundance in metastatic melanomas of these cellular proliferation clusters and their putative upstream regulators was significantly associated with patient prognosis. An 8-gene classifier derived from gene network hub genes correctly classified the prognosis of 23/26 metastatic melanoma patients in a cross-validation study. Unlike the RNA clusters associated with cellular proliferation described above, co-ordinately expressed RNA clusters associated with immune response were clearly identified across melanoma tumours from patients but not across the siRNA-treated A375 cells, in which immune responses are not active. Three uncharacterised genes, which the gene networks predicted to be upstream of apoptosis- or cellular proliferation-associated RNAs, were found to significantly alter apoptosis and cell number when over-expressed in vitro.
This analysis identified co-expression of RNAs that encode functionally-related proteins, in particular, proliferation-associated RNA clusters that are linked to melanoma patient prognosis. Our analysis suggests that A375 cells in vitro may be valid models in which to study the gene expression modules that underlie some melanoma biological processes (e.g., proliferation) but not others (e.g., immune response). The gene expression modules identified here, and the RNAs predicted by Bayesian network inference to be upstream of these modules, are potential prognostic biomarkers and drug targets.
Gene regulatory networks inferred from RNA abundance data have generated significant interest, but despite this, gene network approaches are used infrequently and often require input from bioinformaticians. We have assembled a suite of tools for analysing regulatory networks, and we illustrate their use with microarray datasets generated in human endothelial cells. We infer a range of regulatory networks, and based on this analysis discuss the strengths and limitations of network inference from RNA abundance data. We welcome contact from researchers interested in using our inference and visualization tools to answer biological questions.
The motion of ions, molecules or proteins in dendrites is restricted by cytoplasmic obstacles such as organelles, microtubules and actin network. To account for molecular crowding, we study the effect of diffusion barriers on local calcium spread in a dendrite. We first present a model based on a dimension reduction approach to approximate a three dimensional diffusion in a cylindrical dendrite by a one-dimensional effective diffusion process. By comparing uncaging experiments of an inert dye in a spiny dendrite and in a thin glass tube, we quantify the change in diffusion constants due to molecular crowding as Dcyto/Dwater = 1/20. We validate our approach by reconstructing the uncaging experiments using Brownian simulations in a realistic 3D model dendrite. Finally, we construct a reduced reaction-diffusion equation to model calcium spread in a dendrite under the presence of additional buffers, pumps and synaptic input. We find that for moderate crowding, calcium dynamics is mainly regulated by the buffer concentration, but not by the cytoplasmic crowding, dendritic spines or synaptic inputs. Following high frequency stimulations, we predict that calcium spread in dendrites is limited to small microdomains of the order of a few microns (<5 μm).
Diffusion is one of the main transport phenomena involved in signaling mechanisms of ions and molecules in living cells, such as neurons. As the cell cytoplasmic medium is highly heterogeneous and filled with many organelles, the motion of a diffusing particle is affected by many interactions with its environment. Interestingly, the functional consequences of these interactions cannot be directly quantified. Thus, in parallel with experimental methods, we have developed a computational approach to decipher the role of crowding from binding. We first study here the diffusion of a fluorescent marker in dendrites by a one-dimensional effective diffusion equation and obtained an effective diffusion constant that accounts for the presence heterogeneity in the medium. Furthermore, comparing our experimental data with simulations of diffusion in a crowded environment, we estimate the intracellular calcium spread in dendrites after injection of calcium transients. We confirm that calcium spread is mainly regulated by fixed buffer molecules, that bind temporarily to calcium, and less by the heterogeneous structure of the surrounding medium. Finally, we find that after synaptic inputs, calcium remains restricted to a domain of 2.5 µm to each side of the input location independent of the input frequency.
Somitogenesis, the formation of the body's primary segmental structure common to all vertebrate development, requires coordination between biological mechanisms at several scales. Explaining how these mechanisms interact across scales and how events are coordinated in space and time is necessary for a complete understanding of somitogenesis and its evolutionary flexibility. So far, mechanisms of somitogenesis have been studied independently. To test the consistency, integrability and combined explanatory power of current prevailing hypotheses, we built an integrated clock-and-wavefront model including submodels of the intracellular segmentation clock, intercellular segmentation-clock coupling via Delta/Notch signaling, an FGF8 determination front, delayed differentiation, clock-wavefront readout, and differential-cell-cell-adhesion-driven cell sorting. We identify inconsistencies between existing submodels and gaps in the current understanding of somitogenesis mechanisms, and propose novel submodels and extensions of existing submodels where necessary. For reasonable initial conditions, 2D simulations of our model robustly generate spatially and temporally regular somites, realistic dynamic morphologies and spontaneous emergence of anterior-traveling stripes of Lfng. We show that these traveling stripes are pseudo-waves rather than true propagating waves. Our model is flexible enough to generate interspecies-like variation in somite size in response to changes in the PSM growth rate and segmentation-clock period, and in the number and width of Lfng stripes in response to changes in the PSM growth rate, segmentation-clock period and PSM length.
Recent decades have seen a revolution in experimental techniques that has shifted the focus of experimental biology from behaviors at the micron (cell) scale to those at the nanometer (molecular) scale. An ever-increasing number of studies detail subcellular behaviors, genetic pathways and protein interactions that relate to specific cell functions. This progress, while welcome, sometimes leads us to forget that these components do not exist or function in isolation. To understand their biological importance, in addition to exploring individual components in more detail, we must integrate them into comprehensive models of cells, tissues, organs and organisms. This integration has been incomplete for somitogenesis, an early developmental process that establishes the first signs of segmentation in all vertebrates, patterning the precursors of the vertebrae, ribs, and skeletal muscles of the back and limbs. In this paper, we make significant progress towards a comprehensive model of somitogenesis by combining specialized hypotheses for specific subcomponent mechanisms of somitogenesis into a unified multi-scale model that successfully reproduces many characteristic events seen in the embryo.
Tumor necrosis factor α (TNF-α) is a key regulator of inflammation and rheumatoid arthritis (RA). TNF-α blocker therapies can be very effective for a substantial number of patients, but fail to work in one third of patients who show no or minimal response. It is therefore necessary to discover new molecular intervention points involved in TNF-α blocker treatment of rheumatoid arthritis patients. We describe a data analysis strategy for predicting gene expression measures that are critical for rheumatoid arthritis using a combination of comprehensive genotyping, whole blood gene expression profiles and the component clinical measures of the arthritis Disease Activity Score 28 (DAS28) score. Two separate network ensembles, each comprised of 1024 networks, were built from molecular measures from subjects before and 14 weeks after treatment with TNF-α blocker. The network ensemble built from pre-treated data captures TNF-α dependent mechanistic information, while the ensemble built from data collected under TNF-α blocker treatment captures TNF-α independent mechanisms. In silico simulations of targeted, personalized perturbations of gene expression measures from both network ensembles identify transcripts in three broad categories. Firstly, 22 transcripts are identified to have new roles in modulating the DAS28 score; secondly, there are 6 transcripts that could be alternative targets to TNF-α blocker therapies, including CD86 - a component of the signaling axis targeted by Abatacept (CTLA4-Ig), and finally, 59 transcripts that are predicted to modulate the count of tender or swollen joints but not sufficiently enough to have a significant impact on DAS28.
The collection and analysis of clinical data has played a key role in providing insights into the diagnosis, prognosis and treatment of disease. However, it is imperative that molecular and genetic data also be collected and integrated into the creation of network models, which capture underlying mechanisms of disease and can be interrogated to elucidate previously unknown biology. Bringing data from the clinic to the bench completes the cycle of translational research, which we demonstrate with this work. We built disease models from genetics, whole blood gene expression profiles and the component clinical measures of rheumatoid arthritis using a data-driven approach that leverages supercomputing. Genetic factors can be utilized as a source of perturbation to the system such that causal connections between genetics, molecular entities and clinical outcomes can be inferred. The existing TNF-α blocker treatments for rheumatoid arthritis are only effective for approximately 2/3 of the affected population. We identified novel therapeutic intervention points that may lead to the development of alternatives to TNF-α blocker treatments. We believe this approach will provide improved drug discovery programs, new insights into disease progression, increased drug efficacy and novel biomarkers for chronic and complex diseases.
The organization of muscle is the product of functional adaptation over several length scales spanning from the sarcomere to the muscle bundle. One possible strategy for solving this multiscale coupling problem is to physically constrain the muscle cells in microenvironments that potentiate the organization of their intracellular space. We hypothesized that boundary conditions in the extracellular space potentiate the organization of cytoskeletal scaffolds for directed sarcomeregenesis. We developed a quantitative model of how the cytoskeleton of neonatal rat ventricular myocytes organizes with respect to geometric cues in the extracellular matrix. Numerical results and in vitro assays to control myocyte shape indicated that distinct cytoskeletal architectures arise from two temporally-ordered, organizational processes: the interaction between actin fibers, premyofibrils and focal adhesions, as well as cooperative alignment and parallel bundling of nascent myofibrils. Our results suggest that a hierarchy of mechanisms regulate the self-organization of the contractile cytoskeleton and that a positive feedback loop is responsible for initiating the break in symmetry, potentiated by extracellular boundary conditions, is required to polarize the contractile cytoskeleton.
How muscle is organized impacts its function. However, understanding how muscle organizes is challenging, as the process occurs over several length scales. We approach this multiscale coupling problem by constraining the overall shapes of muscle cells to indirectly control the organization of their intracellular space. We hypothesized the cellular boundary conditions direct the organization of cytoskeletal scaffolds. We developed a model of how the cytoskeleton of cardiomyocytes organizes with respect to boundary cues. Our computational and experimental results to control myocyte shape indicated that distinct muscle architectures arise from two main organizational mechanisms: the interaction between actin fibers, premyofibrils and focal adhesions, as well as cooperative alignment and parallel bundling of more mature myofibrils. We show that a hierarchy of processes regulate the self-organization of cardiomyocytes. Our results suggest that a symmetry break, due to the boundary conditions imposed on the cell, is responsible for polarization of the contractile cytoskeletal organization.
Calmodulin (CaM) is a ubiquitous Ca2+ buffer and second messenger that affects cellular function as diverse as cardiac excitability, synaptic plasticity, and gene transcription. In CA1 pyramidal neurons, CaM regulates two opposing Ca2+-dependent processes that underlie memory formation: long-term potentiation (LTP) and long-term depression (LTD). Induction of LTP and LTD require activation of Ca2+-CaM-dependent enzymes: Ca2+/CaM-dependent kinase II (CaMKII) and calcineurin, respectively. Yet, it remains unclear as to how Ca2+ and CaM produce these two opposing effects, LTP and LTD. CaM binds 4 Ca2+ ions: two in its N-terminal lobe and two in its C-terminal lobe. Experimental studies have shown that the N- and C-terminal lobes of CaM have different binding kinetics toward Ca2+ and its downstream targets. This may suggest that each lobe of CaM differentially responds to Ca2+ signal patterns. Here, we use a novel event-driven particle-based Monte Carlo simulation and statistical point pattern analysis to explore the spatial and temporal dynamics of lobe-specific Ca2+-CaM interaction at the single molecule level. We show that the N-lobe of CaM, but not the C-lobe, exhibits a nano-scale domain of activation that is highly sensitive to the location of Ca2+ channels, and to the microscopic injection rate of Ca2+ ions. We also demonstrate that Ca2+ saturation takes place via two different pathways depending on the Ca2+ injection rate, one dominated by the N-terminal lobe, and the other one by the C-terminal lobe. Taken together, these results suggest that the two lobes of CaM function as distinct Ca2+ sensors that can differentially transduce Ca2+ influx to downstream targets. We discuss a possible role of the N-terminal lobe-specific Ca2+-CaM nano-domain in CaMKII activation required for the induction of synaptic plasticity.
Calmodulin is a versatile Ca2+ signal mediator and a buffer in a wide variety of body organs including the heart and brain. In the brain, calmodulin regulates intracellular molecular processes that change the strength of connectivity between neurons, thus contributing to various brain functions including memory formation. The exact molecular mechanism as to how calmodulin regulates these processes is not yet known. Interestingly, in other excitable tissues, including the heart, each of two lobes of calmodulin responds differentially toward Ca2+ influx and toward its target molecules (e.g., ion channels). This way, calmodulin precisely controls the Ca2+ dynamics of the cell. We wish to test if a similar mechanism may be operational in neurons so that two lobes of calmodulin interact differentially with Ca2+ ions to activate different downstream molecules that control the strength of connections between neurons. We constructed a detailed simulation of calmodulin that allows us to keep track of its interactions with Ca2+ ions and target proteins at the single molecule level. The simulation predicts that two lobes of calmodulin respond differentially to Ca2+ influx both in space and in time. This work opens a door to future experimental testing of the lobe-specific control of neural function by calmodulin.
The kinetics of label uptake and dilution in dividing stem cells, e.g., using Bromodeoxyuridine (BrdU) as a labeling substance, are a common way to assess the cellular turnover of all hematopoietic stem cells (HSCs) in vivo. The assumption that HSCs form a homogeneous population of cells which regularly undergo cell division has recently been challenged by new experimental results. For a consistent functional explanation of heterogeneity among HSCs, we propose a concept in which stem cells flexibly and reversibly adapt their cycling state according to systemic needs. Applying a mathematical model analysis, we demonstrate that different experimentally observed label dilution kinetics are consistently explained by the proposed model. The dynamically stabilized equilibrium between quiescent and activated cells leads to a biphasic label dilution kinetic in which an initial and pronounced decline of label retaining cells is attributed to faster turnover of activated cells, whereas a secondary, decelerated decline results from the slow turnover of quiescent cells. These results, which support our previous model prediction of a reversible activation/deactivation of HSCs, are also consistent with recent findings that use GFP-conjugated histones as a label instead of BrdU. Based on our findings we interpret HSC organization as an adaptive and regulated process in which the slow activation of quiescent cells and their possible return into quiescence after division are sufficient to explain the simultaneous occurrence of self-renewal and differentiation. Furthermore, we suggest an experimental strategy which is suited to demonstrate that the repopulation ability among the population of label retaining cells changes during the course of dilution.
Hematopoietic stem cells (HSCs) are among the best-studied populations of adult stem cells. Commonly HSCs are considered to be in a so called quiescent state of reduced cellular turnover. However, it appears that even quiescent HSCs are activated into the cell cycle from time to time to support the continuous production of peripheral blood. Previous reports on cell kinetic studies using specific chromosomal labeling techniques suggest that all HSCs regularly undergo cell divisions on average once in two weeks. However, these results are challenged by recent experimental findings supporting the idea that HSC populations are heterogeneous with respect to their cell cycle activity. Applying an established model of HSC organization, we demonstrate that different data sets can be consistently explained by the assumption that HSCs are reversibly switching between a quiescent and an activated state. Based on this assumption, we provide a functional explanation for the experimentally observed biphasic label dilution kinetics of HSCs. We conclude that the ability of these cells to reversibly change between proliferation and quiescence is a critical and inherent property of the HSC system necessary for the sustained balance between stem cell self-renewal and the maintenance of functional blood cells.