Search

Identification of mutations of the genes that give cancer a selective advantage is an important step towards research and clinical objectives. As such, there has been a growing interest in developing methods for identification of driver genes and their temporal order within a single patient (intra-tumor) as well as across a cohort of patients (inter-tumor). In this paper, we develop a probabilistic model for tumor progression, in which the driver genes are clustered into several ordered driver pathways. We develop an efficient inference algorithm that exhibits favorable scalability to the number of genes and samples compared to a previously introduced ILP-based method. Adopting a probabilistic approach also allows principled approaches to model selection and uncertainty quantification. Using a large set of experiments on synthetic datasets, we demonstrate our superior performance compared to the ILP-based method. We also analyze two biological datasets of colorectal and glioblastoma cancers. We emphasize that while the ILP-based method puts many seemingly passenger genes in the driver pathways, our algorithm keeps focused on truly driver genes and outputs more accurate models for cancer progression., Author summary Cancer is a disease caused by the accumulation of somatic mutations in the genome. This process is mainly driven by mutations in certain genes that give the harboring cells some selective advantage. The rather few driver genes are usually masked amongst an abundance of so-called passenger mutations. Identification of the driver genes and the temporal order in which the mutations occur is of great importance towards research and clinical objectives. In this paper, we introduce a probabilistic model for cancer progression and devise an efficient inference algorithm to train the model. We show that our method scales favorably to large datasets and provides superior performance compared to an ILP-based counterpart on a wide set of synthetic data simulations. Our Bayesian approach also allows for systematic model selection and confidence quantification procedures in contrast to the previous non-probabilistic progression models. We also study two large datasets on colorectal and glioblastoma cancers and validate our inferred model in comparison to the ILP-based method.

Phylogenetic networks generalize phylogenetic trees by allowing the modelization of events of reticulate evolution. Among the different kinds of phylogenetic networks that have been proposed in the literature, the subclass of binary tree-child networks is one of the most studied ones. However, very little is known about the combinatorial structure of these networks. In this paper we address the problem of generating all possible binary tree-child (BTC) networks with a given number of leaves in an efficient way via reduction/augmentation operations that extend and generalize analogous operations for phylogenetic trees, and are biologically relevant. Since our solution is recursive, this also provides us with a recurrence relation giving an upper bound on the number of such networks. We also show how the operations introduced in this paper can be employed to extend the evolutive history of a set of sequences, represented by a BTC network, to include a new sequence. An implementation in python of the algorithms described in this paper, along with some computational experiments, can be downloaded from https://github.com/bielcardona/TCGenerators., Research of GC and JCP has been partially supported by the Spanish Ministry of Science, Innovation and Universities (http://www.ciencia.gob.es/) and European Regional Development Fund (https://ec.europa.eu/regional_policy/es/funding/erdf/) projects DPI2015-67082-P and PGC2018-096956-B-C43. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

The presence of treatment-resistant cells is an important factor that limits the efficacy of cancer therapy, and the prospect of resistance is considered the major cause of the treatment strategy. Several recent studies have employed mathematical models to elucidate the dynamics of generating resistant cancer cells and attempted to predict the probability of emerging resistant cells. The purpose of this paper is to present numerical approach to compute the number of resistant cells and the emerging probability of resistance. Stochastic model was designed and developed a method to approximately but efficiently compute the number of resistant cells and the probability of resistance. To model the progression of cancer, a discrete-state, two-dimensional Markov process whose states are the total number of cells and the number of resistant cells was employed. Then exact analysis and approximate aggregation approaches were proposed to calculate the number of resistant cells and the probability of resistance when the cell population reaches detection size. To confirm the accuracy of computed results of approximation, relative errors between exact analysis and approximation were computed. The numerical values of our approximation method were very close to those of exact analysis calculated in the range of small detection size M = 500, 100, and 1500. Then computer simulation was performed to confirm the accuracy of computed results of approximation when the detection size was M = 10000,30000,50000,100000 and 1000000. All the numerical results of approximation fell between the upper level and the lower level of 95% confidential intervals and our method took less time to compute over a broad range of cell size. The effects of parameter change on emerging probabilities of resistance were also investigated by computed values using approximation method. The results showed that the number of divisions until the cell population reached the detection size is important for emerging the probability of resistance. The next step of numerical approach is to compute the emerging probabilities of resistance under drug administration and with multiple mutation. Another effective approximation would be necessary for the analysis of the latter case., Author summary Drug therapies for cancer have dramatically succeeded since molecular-targeted drugs have been introduced in medical practice; however, drug treatment often fails owing to the emergence of drug-resistant cells. A variety of approaches, including mathematical modeling, has been undertaken to clarify the mechanism of resistance and subsequently avoid resistance to therapy. This paper proposes one of the mathematical approaches that uses a stochastic model and provides the emerging probabilities of resistance at detection size.

Integrated Information Theory (IIT) has become nowadays the most sensible general theory of consciousness. In addition to very important statements, it opens the door for an abstract (mathematical) formulation of the theory. Given a mechanism in a particular state, IIT identifies a conscious experience with a conceptual structure, an informational object which exists, is composed of identified parts, is informative, integrated and maximally irreducible. This paper introduces a space-time continuous version of the concept of integrated information. To this aim, a graph and a dynamical systems treatment is used to define, for a given mechanism in a state for which a dynamics is settled, an Informational Structure, which is associated to the global attractor at each time of the system. By definition, the informational structure determines all the past and future behavior of the system, possesses an informational nature and, moreover, enriches all the points of the phase space with cause-effect power by means of its associated Informational Field. A detailed description of its inner structure by invariants and connections between them allows to associate a transition probability matrix to each informational structure and to develop a measure for the level of integrated information of the system., Author summary In this paper we introduce a space-time continuous version for the level of integrated information of a network on which a dynamics is defined. The concept of integrated information comes from the IIT of consciousness. By a strict mathematical formulation, we complement the existing IIT theoretical framework from a dynamical systems perspective. In other words, we develop the bases for a continuous mathematical approach to IIT introducing a dynamical system as the driving rule of a given mechanism. We also introduce and define the concepts of Informational Structure and Informational Field as the complex network with the power to ascertain the dynamics (past and future scenarios) of the studied phenomena. The detailed description of an informational structure is showing the cause-effect power of a mechanism in a state and thus, a characterization of the quantity and quality of information, and the way this is integrated. We firstly introduce how network patterns arise from dynamic phenomena on networks, leading to the concept of informational structure. Then, we formally introduce the mathematical objects supporting the theory, from graphs to informational structures, throughout the integration of dynamics on graphs with a global model of differential equations. After this, we formally present some of the IIT’s postulates associated to a given mechanism. Finally, we provide the quantitative and qualitative characterization of the integrated information, and how it depends on the geometry of the mechanism.

Phylogenetic networks are rooted, directed, acyclic graphs that model reticulate evolutionary histories. Recently, statistical methods were devised for inferring such networks from either gene tree estimates or the sequence alignments of multiple unlinked loci. Bi-allelic markers, most notably single nucleotide polymorphisms (SNPs) and amplified fragment length polymorphisms (AFLPs), provide a powerful source of genome-wide data. In a recent paper, a method called SNAPP was introduced for statistical inference of species trees from unlinked bi-allelic markers. The generative process assumed by the method combined both a model of evolution for the bi-allelic markers, as well as the multispecies coalescent. A novel component of the method was a polynomial-time algorithm for exact computation of the likelihood of a fixed species tree via integration over all possible gene trees for a given marker. Here we report on a method for Bayesian inference of phylogenetic networks from bi-allelic markers. Our method significantly extends the algorithm for exact computation of phylogenetic network likelihood via integration over all possible gene trees. Unlike the case of species trees, the algorithm is no longer polynomial-time on all instances of phylogenetic networks. Furthermore, the method utilizes a reversible-jump MCMC technique to sample the posterior of phylogenetic networks given bi-allelic marker data. Our method has a very good performance in terms of accuracy and robustness as we demonstrate on simulated data, as well as a data set of multiple New Zealand species of the plant genus Ourisia (Plantaginaceae). We implemented the method in the publicly available, open-source PhyloNet software package., Author summary The availability of genomic data has revolutionized the study of evolutionary histories and phylogeny inference. Inferring evolutionary histories from genomic data requires, in most cases, accounting for the fact that different genomic regions could have evolutionary histories that differ from each other as well as from that of the species from which the genomes were sampled. In this paper, we introduce a method for inferring evolutionary histories while accounting for two processes that could give rise to such differences across the genomes, namely incomplete lineage sorting and hybridization. We introduce a novel algorithm for computing the likelihood of phylogenetic networks from bi-allelic genetic markers and use it in a Bayesian inference method. Analyses of synthetic and empirical data sets show a very good performance of the method in terms of the estimates it obtains.

Type 2 diabetes leads to premature death and reduced quality of life for 8% of Americans. Nutrition management is critical to maintaining glycemic control, yet it is difficult to achieve due to the high individual differences in glycemic response to nutrition. Anticipating glycemic impact of different meals can be challenging not only for individuals with diabetes, but also for expert diabetes educators. Personalized computational models that can accurately forecast an impact of a given meal on an individual’s blood glucose levels can serve as the engine for a new generation of decision support tools for individuals with diabetes. However, to be useful in practice, these computational engines need to generate accurate forecasts based on limited datasets consistent with typical self-monitoring practices of individuals with type 2 diabetes. This paper uses three forecasting machines: (i) data assimilation, a technique borrowed from atmospheric physics and engineering that uses Bayesian modeling to infuse data with human knowledge represented in a mechanistic model, to generate real-time, personalized, adaptable glucose forecasts; (ii) model averaging of data assimilation output; and (iii) dynamical Gaussian process model regression. The proposed data assimilation machine, the primary focus of the paper, uses a modified dual unscented Kalman filter to estimate states and parameters, personalizing the mechanistic models. Model selection is used to make a personalized model selection for the individual and their measurement characteristics. The data assimilation forecasts are empirically evaluated against actual postprandial glucose measurements captured by individuals with type 2 diabetes, and against predictions generated by experienced diabetes educators after reviewing a set of historical nutritional records and glucose measurements for the same individual. The evaluation suggests that the data assimilation forecasts compare well with specific glucose measurements and match or exceed in accuracy expert forecasts. We conclude by examining ways to present predictions as forecast-derived range quantities and evaluate the comparative advantages of these ranges., Author summary Type 2 diabetes is a devastating disease that requires constant patient self-management of glucose, insulin, nutrition and exercise. Nevertheless, glucose and insulin dynamics are complicated, nonstationary, nonlinear, and individual-dependent, making self-management of diabetes a complex task. To help alleviate some of the difficulty for patients, we develop a method for personalized, real-time, glucose forecasting based on nutrition. Specifically, we create and evaluate the computational machinery based on both Gaussian process models and data assimilation that leverages the physiologic knowledge of two mechanistic models to produce a personalized, nutrition-based glucose forecast for individuals with type 2 diabetes in real time that is robust to sparse data and nonstationary patients. Our computational engine was conceived to be of potential use for diabetes self-management.

Neural Mass Models provide a compact description of the dynamical activity of cell populations in neocortical regions. Moreover, models of regional activity can be connected together into networks, and inferences made about the strength of connections, using M/EEG data and Bayesian inference. To date, however, Bayesian methods have been largely restricted to the Variational Laplace (VL) algorithm which assumes that the posterior distribution is Gaussian and finds model parameters that are only locally optimal. This paper explores the use of Annealed Importance Sampling (AIS) to address these restrictions. We implement AIS using proposals derived from Langevin Monte Carlo (LMC) which uses local gradient and curvature information for efficient exploration of parameter space. In terms of the estimation of Bayes factors, VL and AIS agree about which model is best but report different degrees of belief. Additionally, AIS finds better model parameters and we find evidence of non-Gaussianity in their posterior distribution., Author Summary The activity of populations of neurons in the human brain can be described using a set of differential equations known as a neural mass model. These models can then be connected to describe activity in multiple brain regions and, by fitting them to human brain imaging data, statistical inferences can be made about changes in macroscopic connectivity among brain regions. For example, the strength of a connection from one region to another may be more strongly engaged in a particular patient population or during a specific cognitive task. Current statistical inference approaches use a Bayesian algorithm based on principles of local optimization and the assumption that uncertainty about model parameters (e.g. connectivity), having seen the data, follows a Gaussian distribution. This paper evaluates current methods against a global Bayesian optimization algorithm and finds that the two approaches (local/global) agree about which model is best, but finds that the global approach produces better parameter estimates.

System identification techniques -- projection pursuit regression models (PPRs) and convolutional neural networks (CNNs) -- provide state-of-the-art performance in predicting visual cortical neurons' responses to arbitrary input stimuli. However, the constituent kernels recovered by these methods are often noisy and lack coherent structure, making it difficult to understand the underlying component features of a neuron's receptive field. In this paper, we show that using a dictionary of diverse kernels with complex shapes learned from natural scenes based on efficient coding theory, as the front-end for PPRs and CNNs can improve their performance in neuronal response prediction as well as algorithmic data efficiency and convergence speed. Extensive experimental results also indicate that these sparse-code kernels provide important information on the component features of a neuron's receptive field. In addition, we find that models with the complex-shaped sparse code front-end are significantly better than models with a standard orientation-selective Gabor filter front-end for modeling V1 neurons that have been found to exhibit complex pattern selectivity. We show that the relative performance difference due to these two front-ends can be used to produce a sensitive metric for detecting complex selectivity in V1 neurons., 22 pages

Working memory is a cognitive function involving the storage and manipulation of latent information over brief intervals of time, thus making it crucial for context-dependent computation. Here, we use a top-down modeling approach to examine network-level mechanisms of working memory, an enigmatic issue and central topic of study in neuroscience. We optimize thousands of recurrent rate-based neural networks on a working memory task and then perform dynamical systems analysis on the ensuing optimized networks, wherein we find that four distinct dynamical mechanisms can emerge. In particular, we show the prevalence of a mechanism in which memories are encoded along slow stable manifolds in the network state space, leading to a phasic neuronal activation profile during memory periods. In contrast to mechanisms in which memories are directly encoded at stable attractors, these networks naturally forget stimuli over time. Despite this seeming functional disadvantage, they are more efficient in terms of how they leverage their attractor landscape and paradoxically, are considerably more robust to noise. Our results provide new hypotheses regarding how working memory function may be encoded within the dynamics of neural circuits., Author summary The ability to remember information for brief periods of time before using it is a key human ability. For example, retaining a phone number for a few moments prior to entering it into a keypad. Such ability, known as working memory, enables many more complex functions such as planning and reasoning. In this paper, we use theory and computational modeling approaches to try and better understand how circuits and networks in the brain might be achieving working memory. Specifically, we construct hypothetical network models to perform a task that embodies essential aspects of working memory. We then dissect our model to reveal how its components—simulated neural units—interact with each other to represent and maintain information. It turns out that some models can achieve memory by maintaining a fixed representation over time, i.e., the units remain ‘still’ during memory. However, other models fluctuate their activity during memory in a particular way that is seemingly quite efficient and less susceptible to distraction. In total, our computational study provides new theory for how this form of memory might be implemented in brain networks.

Gene function annotation is important for a variety of downstream analyses of genetic data. But experimental characterization of function remains costly and slow, making computational prediction an important endeavor. Phylogenetic approaches to prediction have been developed, but implementation of a practical Bayesian framework for parameter estimation remains an outstanding challenge. We have developed a computationally efficient model of evolution of gene annotations using phylogenies based on a Bayesian framework using Markov Chain Monte Carlo for parameter estimation. Unlike previous approaches, our method is able to estimate parameters over many different phylogenetic trees and functions. The resulting parameters agree with biological intuition, such as the increased probability of function change following gene duplication. The method performs well on leave-one-out cross-validation, and we further validated some of the predictions in the experimental scientific literature., Author summary Understanding the individual roles that genes play in life is a key issue in biomedical science. While information regarding gene functions is continuously growing, the number of genes with uncharacterized biological functions is still greater. Because of this, scientists have dedicated much of their time to build and design tools that automatically infer gene functions. One of the most promising approaches (sometimes called “phylogenomics”) attempts to construct a model of inheritance and divergence of function along branches of the phylogenetic tree that relates different members of a gene family. If the functions of one or more of the family members has been characterized experimentally, the presence or absence of these functions for other family members can be predicted, in a probabilistic framework, based on the evolutionary relationships. Previously proposed Bayesian approaches to parameter estimation have proved to be computationally intractable, preventing development of such a probabilistic model. In this paper, we present a simple model of gene-function evolution that is highly scalable, which means that it is possible to perform parameter estimation not only on one family, but simultaneously for hundreds of gene families, comprising thousands of genes. The parameter estimates we obtain coherently agree with what theory dictates regarding how gene-functions evolved. Finally, notwithstanding its simplicity, the model’s prediction quality is comparable to other more complex alternatives. Although we believe further improvements can be made to our model, even this simple model makes verifiable predictions, and suggests areas in which existing annotations show inconsistencies that may indicate errors or controversies.

Tumors contain multiple subpopulations of genetically distinct cancer cells. Reconstructing their evolutionary history can improve our understanding of how cancers develop and respond to treatment. Subclonal reconstruction methods cluster mutations into groups that co-occur within the same subpopulations, estimate the frequency of cells belonging to each subpopulation, and infer the ancestral relationships among the subpopulations by constructing a clone tree. However, often multiple clone trees are consistent with the data and current methods do not efficiently capture this uncertainty; nor can these methods scale to clone trees with a large number of subclonal populations. Here, we formalize the notion of a partially-defined clone tree (partial clone tree for short) that defines a subset of the pairwise ancestral relationships in a clone tree, thereby implicitly representing the set of all clone trees that have these defined pairwise relationships. Also, we introduce a special partial clone tree, the Maximally-Constrained Ancestral Reconstruction (MAR), which summarizes all clone trees fitting the input data equally well. Finally, we extend commonly used clone tree validity conditions to apply to partial clone trees and describe SubMARine, a polynomial-time algorithm producing the subMAR, which approximates the MAR and guarantees that its defined relationships are a subset of those present in the MAR. We also extend SubMARine to work with subclonal copy number aberrations and define equivalence constraints for this purpose. Further, we extend SubMARine to permit noise in the estimates of the subclonal frequencies while retaining its validity conditions and guarantees. In contrast to other clone tree reconstruction methods, SubMARine runs in time and space that scale polynomially in the number of subclones. We show through extensive noise-free simulation, a large lung cancer dataset and a prostate cancer dataset that the subMAR equals the MAR in all cases where only a single clone tree exists and that it is a perfect match to the MAR in most of the other cases. Notably, SubMARine runs in less than 70 seconds on a single thread with less than one Gb of memory on all datasets presented in this paper, including ones with 50 nodes in a clone tree. On the real-world data, SubMARine almost perfectly recovers the previously reported trees and identifies minor errors made in the expert-driven reconstructions of those trees. The freely-available open-source code implementing SubMARine can be downloaded at https://github.com/morrislab/submarine. Copyright: © 2021 Sundermann et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited., PLoS Computational Biology, 17 (1), ISSN:1553-734X, ISSN:1553-7358

When we face ambiguous images, the brain cannot commit to a single percept; instead, it switches between mutually exclusive interpretations every few seconds, a phenomenon known as bistable perception. While neuromechanistic models, e.g., adapting neural populations with lateral inhibition, may account for the dynamics of bistability, a larger question remains unresolved: how this phenomenon informs us on generic perceptual processes in less artificial contexts. Here, we propose that bistable perception is due to our prior beliefs being reverberated in the cortical hierarchy and corrupting the sensory evidence, a phenomenon known as “circular inference”. Such circularity could occur in a hierarchical brain where sensory responses trigger activity in higher-level areas but are also modulated by feedback projections from these same areas. We show that in the face of ambiguous sensory stimuli, circular inference can change the dynamics of the perceptual system and turn what should be an integrator of inputs into a bistable attractor switching between two highly trusted interpretations. The model captures various aspects of bistability, including Levelt’s laws and the stabilizing effects of intermittent presentation of the stimulus. Since it is related to the generic perceptual inference and belief updating mechanisms, this approach can be used to predict the tendency of individuals to form aberrant beliefs from their bistable perception behavior. Overall, we suggest that feedforward/feedback information loops in hierarchical neural networks, a phenomenon that could lead to psychotic symptoms when overly strong, could also underlie perception in nonclinical populations., Author summary In cases of high ambiguity, our perceptual system cannot commit to a single percept and switches between different interpretations, giving rise to bistable perception. In this paper we outline a computational model of bistability based on the notion of circular inference, i.e. a form of suboptimal hierarchical inference in which priors and / or sensory inputs are reverberated and over-counted. We suggest that descending loops (i.e. reverberated priors) transform our perceptual system from a simple accumulator of sensory inputs into a bistable attractor, that switches between two highly-trusted interpretations. Using analytical methods we derive the necessary conditions for bistable perception to occur. We show that our dynamical circular inference model is able to capture many features of bistability, such as Levelt’s laws and the stabilizing effects of intermittent presentation of the stimulus. Finally we make novel predictions about the behavior of psychotic patients.

There has been considerable debate and concern as to whether there is a replication crisis in the scientific literature. A likely cause of poor replication is the multiple comparisons problem. An important way in which this problem can manifest in the M/EEG context is through post hoc tailoring of analysis windows (a.k.a. regions-of-interest, ROIs) to landmarks in the collected data. Post hoc tailoring of ROIs is used because it allows researchers to adapt to inter-experiment variability and discover novel differences that fall outside of windows defined by prior precedent, thereby reducing Type II errors. However, this approach can dramatically inflate Type I error rates. One way to avoid this problem is to tailor windows according to a contrast that is orthogonal (strictly parametrically orthogonal) to the contrast being tested. A key approach of this kind is to identify windows on a fully flattened average. On the basis of simulations, this approach has been argued to be safe for post hoc tailoring of analysis windows under many conditions. Here, we present further simulations and mathematical proofs to show exactly why the Fully Flattened Average approach is unbiased, providing a formal grounding to the approach, clarifying the limits of its applicability and resolving published misconceptions about the method. We also provide a statistical power analysis, which shows that, in specific contexts, the fully flattened average approach provides higher statistical power than Fieldtrip cluster inference. This suggests that the Fully Flattened Average approach will enable researchers to identify more effects from their data without incurring an inflation of the false positive rate., Author summary It is clear from recent replicability studies that the replication rate in psychology and cognitive neuroscience is not high. One reason for this is that the noise in high dimensional neuroimaging data sets can “look-like” signal. A classic manifestation would be selecting a region in the data volume where an effect is biggest and then specifically reporting results on that region. There is a key trade-off in the selection of such regions of interest: liberal selection will inflate false positive rates, but conservative selection (e.g. strictly on the basis of prior precedent in the literature) can reduce statistical power, causing real effects to be missed. We propose a means to reconcile these two possibilities, by which regions of interest can be tailored to the pattern in the collected data, while not inflating false-positive rates. This is based upon generating what we call the Flattened Average. Critically, we validate the correctness of this method both in (ground-truth) simulations and with formal mathematical proofs. Given the replication “crisis”, there may be no more important issue in psychology and cognitive neuroscience than improving the application of methods. This paper makes a valuable contribution to this improvement.

Gene regulatory networks (GRNs) are complex systems in which many genes regulate mutually to adapt the cell state to environmental conditions. In addition to function, the GRNs possess several kinds of robustness. This robustness means that systems do not lose their functionality when exposed to disturbances such as mutations or noise, and is widely observed at many levels in living systems. Both function and robustness have been acquired through evolution. In this respect, GRNs utilized in living systems are rare among all possible GRNs. In this study, we explored the fitness landscape of GRNs and investigated how robustness emerged in highly-fit GRNs. We considered a toy model of GRNs with one input gene and one output gene. The difference in the expression level of the output gene between two input states, “on” and “off”, was considered as fitness. Thus, the determination of the fitness of a GRN was based on how sensitively it responded to the input. We employed the multicanonical Monte Carlo method, which can sample GRNs randomly in a wide range of fitness levels, and classified the GRNs according to their fitness. As a result, the following properties were found: (1) Highly-fit GRNs exhibited bistability for intermediate input between “on” and “off”. This means that such GRNs responded to two input states by using different fixed points of dynamics. This bistability emerges necessarily as fitness increases. (2) These highly-fit GRNs were robust against noise because of their bistability. In other words, noise robustness is a byproduct of high fitness. (3) GRNs that were robust against mutations were not extremely rare among the highly-fit GRNs. This implies that mutational robustness is readily acquired through the evolutionary process. These properties are universal irrespective of the evolutionary pathway, because the results do not rely on evolutionary simulation., Author summary Living systems have developed through a long history of Darwinian evolution. They acquired characteristic properties distinct from other physical systems; one is biological function. Another important property, which is overlooked by non-experts, is robustness to noise and mutation. Here, robustness means that a system does not lose its functionality when exposed to disturbances. Then, how do they relate to each other? In this paper, we explored this question using a toy model of gene regulatory networks (GRNs). While evolutionary simulations are usually used for such purposes, we instead generated GRNs randomly and classified them according to functionality. By requiring sensitive responses to environmental change as a function, we found that bistability emerges as a common property of highly-functional GRNs. Since this property does not depend on a particular evolutionary pathway, if the evolution was rewound and repeated over and over again, phenotypes with the same property would always evolve. At the same time, such bistable GRNs were robust to noise. We also found that GRNs robust to mutation were not extremely rare among the highly-functional GRNs. This implies that mutational robustness would be readily acquired through evolution.

The molecular mechanisms and functions in complex biological systems currently remain elusive. Recent high-throughput techniques, such as next-generation sequencing, have generated a wide variety of multiomics datasets that enable the identification of biological functions and mechanisms via multiple facets. However, integrating these large-scale multiomics data and discovering functional insights are, nevertheless, challenging tasks. To address these challenges, machine learning has been broadly applied to analyze multiomics. This review introduces multiview learning-an emerging machine learning field-and envisions its potentially powerful applications to multiomics. In particular, multiview learning is more effective than previous integrative methods for learning data's heterogeneity and revealing cross-talk patterns. Although it has been applied to various contexts, such as computer vision and speech recognition, multiview learning has not yet been widely applied to biological data-specifically, multiomics data. Therefore, this paper firstly reviews recent multiview learning methods and unifies them in a framework called multiview empirical risk minimization (MV-ERM). We further discuss the potential applications of each method to multiomics, including genomics, transcriptomics, and epigenomics, in an aim to discover the functional and mechanistic interpretations across omics. Secondly, we explore possible applications to different biological systems, including human diseases (e.g., brain disorders and cancers), plants, and single-cell analysis, and discuss both the benefits and caveats of using multiview learning to discover the molecular mechanisms and functions of these systems.

Predominantly heterosexual HIV-1 epidemics like those in sub-Saharan Africa continue to have high HIV incidence in young people. We used a stochastic, agent-based model for age-disparate networks to test the hypothesis that focusing uptake and retention of ART among youth could enhance the efficiency of treatment as prevention (TasP) campaigns. We used the model to identify strategies that reduce incidence to negligible levels (i.e., < 0.1 cases/100 person-years) 20–25 years after initiation of a targeted TasP campaign. The model was parameterized using behavioral, demographic, and clinical data from published papers and national reports. To keep a focus on the underlying age effects we model a generalized heterosexual population with average risks (i.e., no MSM, no PWIDs, no sex workers) and no entry of HIV+ people from other regions. The model assumes that most people (default 95%, range in variant simulations 60–95%) are “linkable”; i.e., could get linked to effective care given sufficient resources. To simplify the accounting, we assume a rapid jump in the number of people receiving treatment at the start of the TasP campaign, followed by a 2% annual increase that continues until all linkable HIV+ people have been treated. Under historical scenarios of CD4-based targeted ART allocation and current policies of untargeted (random) ART allocation, our model predicts that viral replication would need to be suppressed in 60–85% of infected people at the start of the TasP campaign to drive incidence to negligible levels. Under age-based strategies, by contrast, this percentage dropped by 18–54%, depending on the strength of the epidemic and the age target. For our baseline model, targeting those under age 30 halved the number of people who need to be treated. Age-based targeting also minimized total and time-discounted AIDS deaths over 25 years. Age-based targeting yielded benefits without being highly exclusive; in a model in which 60% of infected people were treated, ~87% and ~58% of those initiating therapy during a campaign targeting those, Author summary We used an agent-based model for heterosexual epidemics to evaluate the effects of focusing treatment-as-prevention (TasP) efforts on young people. To keep the focus on age effects, we modeled populations with average risks (i.e., no MSM, no PWIDs, no sex workers). For our baseline model, we found that a long-term TasP program that focuses on HIV+ people under the age of 30 can halve the number of people who need to be treated in order to drive incidence to zero. Within the context of an overall expansion of care, we were able to identify multiple scenarios in which a successful youth-focused TasP campaign also expanded care for older people. Interestingly, we continued to see strong advantages to age-based TasP in sensitivity tests in which we removed key age-related risks (e.g., the tendency for young people to have shorter relationships). Youth-focused TasP succeeds, in large part, because it protects HIV- adolescents entering the sexually active population; i.e., a population that would otherwise be at high risk not only for getting infected but also for spreading the virus to others. This protection of adolescents leads to a steadily expanding AIDS-free generation that drives the virus to extinction.

Phenotypic plasticity is a vital strategy for plants to deal with changing conditions by inducing phenotypes favourable in different environments. Understanding how natural selection acts on variation in phenotypic plasticity in plants is therefore a central question in ecology, but is often ignored in modelling studies. Here we present a new modelling approach that allows for the analysis of selection for variation in phenotypic plasticity as a response strategy. We assess selection for shade avoidance strategies of Arabidopsis thaliana in response to future neighbour shading signalled through a decrease in red:far-red (R:FR) ratio. For this, we used a spatially explicit 3D virtual plant model that simulates individual Arabidopsis plants competing for light in different planting densities. Plant structure and growth were determined by the organ-specific interactions with the light environment created by the vegetation structure itself. Shade avoidance plastic responses were defined by a plastic response curve relating petiole elongation and lamina growth to R:FR perceived locally. Different plasticity strategies were represented by different shapes of the response curve that expressed different levels of R:FR sensitivity. Our analyses show that the shape of the selected shade avoidance strategy varies with planting density. At higher planting densities, more sensitive response curves are selected for than at lower densities. In addition, the balance between lamina and petiole responses influences the sensitivity of the response curves selected for. Combining computational virtual plant modelling with a game theoretical analysis represents a new step towards analysing how natural selection could have acted upon variation in shade avoidance as a response strategy, which can be linked to genetic variation and underlying physiological processes., Author summary Plants are able to respond to changes in the environment. Particularly, plants can show different structural traits (e.g. stem height and leaf size) in different planting densities. These trait changes are the result of so-called plastic responses that can be induced by changes in the light spectrum. Although a great part of the physiological processes underlying these plastic responses have been identified, it remains unclear how these plastic responses, and variations therein, can be the result of selection. In this paper we analyse selection on different plastic responses within different dynamic competitive environments. We use a 3D virtual plant model that simulates realistic plant growth based on light absorption, photosynthesis and specific light signals that induce changes in leaf growth. The 3D model simulated various competitive vegetation stands consisting of plants with different plastic responses at different planting densities. We conclude that selection in different densities results in different plastic responses. We advocate that our modelling approach allows for analyses related to selection on plastic responses itself, instead of on specific trait values in different environments, and is therefore an important new step forward in understanding the role of plastic responses for plant performance.

Rare variants are of increasing interest to genetic association studies because of their etiological contributions to human complex diseases. Due to the rarity of the mutant events, rare variants are routinely analyzed on an aggregate level. While aggregation analyses improve the detection of global-level signal, they are not able to pinpoint causal variants within a variant set. To perform inference on a localized level, additional information, e.g., biological annotation, is often needed to boost the information content of a rare variant. Following the observation that important variants are likely to cluster together on functional domains, we propose a protein structure guided local test (POINT) to provide variant-specific association information using structure-guided aggregation of signal. Constructed under a kernel machine framework, POINT performs local association testing by borrowing information from neighboring variants in the 3-dimensional protein space in a data-adaptive fashion. Besides merely providing a list of promising variants, POINT assigns each variant a p-value to permit variant ranking and prioritization. We assess the selection performance of POINT using simulations and illustrate how it can be used to prioritize individual rare variants in PCSK9, ANGPTL4 and CETP in the Action to Control Cardiovascular Risk in Diabetes (ACCORD) clinical trial data., Author summary While it is known that rare variants play an important role in understanding associations between genotype and complex diseases, pinpointing individual rare variants likely to be responsible for association is still a daunting task. Due to their low frequency in the population and reduced signal, localizing causal rare variants often requires additional information, such as type of DNA change or location of variant along the sequence, to be incorporated in a biologically meaningful fashion that does not overpower the genotype data. In this paper, we use the observation that important variants tend to cluster together on functional domains to propose a new approach for prioritizing rare variants: the protein structure guided local test (POINT). POINT uses a gene’s 3-dimensional protein folding structure to guide aggregation of information from neighboring variants in the protein in a robust manner. We show how POINT improves selection performance over existing methods. We further illustrate how it can be used to prioritize individual rare variants using the Action to Control Cardiovascular Risk in Diabetes (ACCORD) clinical trial data, finding promising variants within genes in association with lipoprotein-related outcomes.

Cellular morphology and associated morphodynamics are widely used for qualitative and quantitative assessments of cell state. Here we implement a framework to profile cellular morphodynamics based on an adaptive decomposition of local cell boundary motion into instantaneous frequency spectra defined by the Hilbert-Huang transform (HHT). Our approach revealed that spontaneously migrating cells with approximately homogeneous molecular makeup show remarkably consistent instantaneous frequency distributions, though they have markedly heterogeneous mobility. Distinctions in cell edge motion between these cells are captured predominantly by differences in the magnitude of the frequencies. We found that acute photo-inhibition of Vav2 guanine exchange factor, an activator of the Rho family of signaling proteins coordinating cell motility, produces significant shifts in the frequency distribution, but does not affect frequency magnitude. We therefore concluded that the frequency spectrum encodes the wiring of the molecular circuitry that regulates cell boundary movements, whereas the magnitude captures the activation level of the circuitry. We also used HHT spectra as multi-scale spatiotemporal features in statistical region merging to identify subcellular regions of distinct motion behavior. In line with our conclusion that different HHT spectra relate to different signaling regimes, we found that subcellular regions with different morphodynamics indeed exhibit distinct Rac1 activities. This algorithm thus can serve as an accurate and sensitive classifier of cellular morphodynamics to pinpoint spatial and temporal boundaries between signaling regimes., Author summary Many studies in cell biology employ global shape descriptors to probe mechanisms of cell morphogenesis. Here, we implement a framework in this paper to profile cellular morphodynamics very locally. We employ the Hilbert-Huang transform (HHT) to extract along the entire cell edge spectra of instantaneous edge motion frequency and magnitude and use them to classify overall cell behavior as well as subcellular edge sectors of distinct dynamics. We find in fibroblast-like COS7 cells that the marked heterogeneity in mobility of an unstimulated population is fully captured by differences in the magnitude spectra, while the frequency spectra are conserved between cells. Using optogenetics to acutely inhibit morphogenetic signaling pathways we find that these molecular shifts are reflected by changes in the frequency spectra but not in the magnitude spectra. After clustering cell edge sectors with distinct morphodynamics we observe in cells expressing a Rac1 activity biosensor that the sectors with different frequency spectra associate with different signaling intensity and dynamics. Together, these observations let us conclude that the frequency spectrum encodes the wiring of the molecular circuitry that regulates edge movements, whereas the magnitude captures the activation level of the circuitry.

Many diverse biological systems are described by randomly moving particles that can be captured by traps in their environment. Examples include neurotransmitters diffusing in the synaptic cleft before binding to receptors and prey roaming an environment before capture by predators. In most cases, the traps cannot capture particles continuously. Rather, each trap must wait a transitory “recharge” time after capturing a particle before additional captures. This recharge time is often overlooked. In the case of instant recharge, the average number of particles captured before they escape grows linearly in the total number of particles. In stark contrast, we prove that for any nonzero recharge time, the average number of captured particles grows at most logarithmically in the total particle number. This is a fundamental effect of recharge, as it holds under very general assumptions on particle motion and spatial domain. Furthermore, we characterize the parameter regime in which a given recharge time will dramatically affect a system, allowing researchers to easily verify if they need to account for recharge in their specific system. Finally, we consider a few examples, including a neural system in which recharge reduces neurotransmitter bindings by several orders of magnitude., Author summary Consider particles that are released into an environment (think diffusing molecules or plankton), and suppose that there are traps in the environment. How many particles will be captured by the traps before they escape? In a standard model, the number of captured particles is proportional to the initial number released. In this paper, we show that for a more realistic model of a trap (one in which traps must recharge after every capture), the number of captures is proportional to the logarithm of the initial number released. That means that if 106 particles are released, only about 6 will be captured. We prove this result mathematically, and then consider a number of applications, including neuronal synapses and ambush predators.

The study of processes evolving on networks has recently become a very popular research field, not only because of the rich mathematical theory that underpins it, but also because of its many possible applications, a number of them in the field of biology. Indeed, molecular signaling pathways, gene regulation, predator-prey interactions and the communication between neurons in the brain can be seen as examples of networks with complex dynamics. The properties of such dynamics depend largely on the topology of the underlying network graph. In this work, we want to answer the following question: Knowing network connectivity, what can be said about the level of third-order correlations that will characterize the network dynamics? We consider a linear point process as a model for pulse-coded, or spiking activity in a neuronal network. Using recent results from theory of such processes, we study third-order correlations between spike trains in such a system and explain which features of the network graph (i.e. which topological motifs) are responsible for their emergence. Comparing two different models of network topology—random networks of Erdős-Rényi type and networks with highly interconnected hubs—we find that, in random networks, the average measure of third-order correlations does not depend on the local connectivity properties, but rather on global parameters, such as the connection probability. This, however, ceases to be the case in networks with a geometric out-degree distribution, where topological specificities have a strong impact on average correlations., Author Summary Many biological phenomena can be viewed as dynamical processes on a graph. Understanding coordinated activity of nodes in such a network is of some importance, as it helps to characterize the behavior of the complex system. Of course, the topology of a network plays a pivotal role in determining the level of coordination among its different vertices. In particular, correlations between triplets of events (here: action potentials generated by neurons) have recently garnered some interest in the theoretical neuroscience community. In this paper, we present a decomposition of an average measure of third-order coordinated activity of neurons in a spiking neuronal network in terms of the relevant topological motifs present in the underlying graph. We study different network topologies and show, in particular, that the presence of certain tree motifs in the synaptic connectivity graph greatly affects the strength of third-order correlations between spike trains of different neurons.