Statistical inference for a quasi birth-death model of RNA transcription.

BACKGROUND: A birth-death process of which the births follow a hypoexponential distribution with L phases and are controlled by an on/off mechanism, is a population process which we call the on/off-seq-L process. It is a suitable model for the dynamics of a population of RNA molecules in a single living cell. Motivated by this biological application, our aim is to develop a statistical method to estimate the model parameters of the on/off-seq-L process, based on observations of the population size at discrete time points, and to apply this method to real RNA data. METHODS: It is shown that the on/off-seq-L process can be seen as a quasi birth-death process, and an Erlangization technique can be used to approximate the corresponding likelihood function. An extensive simulation-based numerical study is carried out to investigate the performance of the resulting estimation method. RESULTS AND CONCLUSION: A statistical method is presented to find maximum likelihood estimates of the model parameters for the on/off-seq-L process. Numerical complications related to the likelihood maximization are identified and analyzed, and solutions are presented. The proposed estimation method is a highly accurate method to find the parameter estimates. Based on real RNA data, the on/off-seq-3 process emerges as the best model to describe RNA transcription.

Analysis of Twitter data with the Bayesian fused graphical lasso.

We propose a method to simplify textual Twitter data into understandable networks of terms that can signify important events and their possible changes over time. The method allows for common characteristics of the networks across time periods and each period can comprise multiple unknown sub-networks. The networks are described by Gaussian graphical models and their parameter values are estimated through a Bayesian approach with a fused lasso-type prior on the precision matrices of the underlying mixtures of the sub-models. A flexible data allocation scheme is at the heart of an MCMC algorithm to recover mean and covariance parameters of the mixture components. Several implementations of the outlined estimation procedure are studied and compared based on simulated data. The procedure with the highest predictive power is used for mining tweets regarding the 2009 Iranian presidential election.

Bayesian mixture regression analysis for regulation of Pluripotency in ES cells.

BACKGROUND: Observed levels of gene expression strongly depend on both activity of DNA binding transcription factors (TFs) and chromatin state through different histone modifications (HMs). In order to recover the functional relationship between local chromatin state, TF binding and observed levels of gene expression, regression methods have proven to be useful tools. They have been successfully applied to predict mRNA levels from genome-wide experimental data and they provide insight into context-dependent gene regulatory mechanisms. However, heterogeneity arising from gene-set specific regulatory interactions is often overlooked. RESULTS: We show that regression models that predict gene expression by using experimentally derived ChIP-seq profiles of TFs can be significantly improved by mixture modelling. In order to find biologically relevant gene clusters, we employ a Bayesian allocation procedure which allows us to integrate additional biological information such as three-dimensional nuclear organization of chromosomes and gene function. The data integration procedure involves transforming the additional data into gene similarity values. We propose a generic similarity measure that is especially suitable for situations where the additional data are of both continuous and discrete type, and compare its performance with similar measures in the context of mixture modelling. CONCLUSIONS: We applied the proposed method on a data from mouse embryonic stem cells (ESC). We find that including additional data results in mixture components that exhibit biologically meaningful gene clusters, and provides valuable insight into the heterogeneity of the regulatory interactions.

Reconstruction of molecular network evolution from cross-sectional omics data.

Cross-sectional studies may shed light on the evolution of a disease like cancer through the comparison of patient traits among disease stages. This problem is especially challenging when a gene-gene interaction network needs to be reconstructed from omics data, and, in addition, the patients of each stage need not form a homogeneous group. Here, the problem is operationalized as the estimation of stage-wise mixtures of Gaussian graphical models (GGMs) from high-dimensional data. These mixtures are fitted by a (fused) ridge penalized EM algorithm. The fused ridge penalty shrinks GGMs of contiguous stages. The (fused) ridge penalty parameters are chosen through cross-validation. The proposed estimation procedures are shown to be consistent and their performance in other respects is studied in simulation. The down-stream exploitation of the fitted GGMs is outlined. In a data illustration the methodology is employed to identify gene-gene interaction network changes in the transition from normal to cancer prostate tissue.

Gene Network Reconstruction using Global-Local Shrinkage Priors.

Reconstructing a gene network from high-throughput molecular data is an important but challenging task, as the number of parameters to estimate easily is much larger than the sample size. A conventional remedy is to regularize or penalize the model likelihood. In network models, this is often done locally in the neighbourhood of each node or gene. However, estimation of the many regularization parameters is often difficult and can result in large statistical uncertainties. In this paper we propose to combine local regularization with global shrinkage of the regularization parameters to borrow strength between genes and improve inference. We employ a simple Bayesian model with non-sparse, conjugate priors to facilitate the use of fast variational approximations to posteriors. We discuss empirical Bayes estimation of hyper-parameters of the priors, and propose a novel approach to rank-based posterior thresholding. Using extensive model- and data-based simulations, we demonstrate that the proposed inference strategy outperforms popular (sparse) methods, yields more stable edges, and is more reproducible. The proposed method, termed ShrinkNet, is then applied to Glioblastoma to investigate the interactions between genes associated with patient survival.

A three domain covariance framework for EEG/MEG data.

In this paper we introduce a covariance framework for the analysis of single subject EEG and MEG data that takes into account observed temporal stationarity on small time scales and trial-to-trial variations. We formulate a model for the covariance matrix, which is a Kronecker product of three components that correspond to space, time and epochs/trials, and consider maximum likelihood estimation of the unknown parameter values. An iterative algorithm that finds approximations of the maximum likelihood estimates is proposed. Our covariance model is applicable in a variety of cases where spontaneous EEG or MEG acts as source of noise and realistic noise covariance estimates are needed, such as in evoked activity studies, or where the properties of spontaneous EEG or MEG are themselves the topic of interest, like in combined EEG-fMRI experiments in which the correlation between EEG and fMRI signals is investigated. We use a simulation study to assess the performance of the estimator and investigate the influence of different assumptions about the covariance factors on the estimated covariance matrix and on its components. We apply our method to real EEG and MEG data sets.

A morpho-density approach to estimating neural connectivity.

Neuronal signal integration and information processing in cortical neuronal networks critically depend on the organization of synaptic connectivity. Because of the challenges involved in measuring a large number of neurons, synaptic connectivity is difficult to determine experimentally. Current computational methods for estimating connectivity typically rely on the juxtaposition of experimentally available neurons and applying mathematical techniques to compute estimates of neural connectivity. However, since the number of available neurons is very limited, these connectivity estimates may be subject to large uncertainties. We use a morpho-density field approach applied to a vast ensemble of model-generated neurons. A morpho-density field (MDF) describes the distribution of neural mass in the space around the neural soma. The estimated axonal and dendritic MDFs are derived from 100,000 model neurons that are generated by a stochastic phenomenological model of neurite outgrowth. These MDFs are then used to estimate the connectivity between pairs of neurons as a function of their inter-soma displacement. Compared with other density-field methods, our approach to estimating synaptic connectivity uses fewer restricting assumptions and produces connectivity estimates with a lower standard deviation. An important requirement is that the model-generated neurons reflect accurately the morphology and variation in morphology of the experimental neurons used for optimizing the model parameters. As such, the method remains subject to the uncertainties caused by the limited number of neurons in the experimental data set and by the quality of the model and the assumptions used in creating the MDFs and in calculating estimating connectivity. In summary, MDFs are a powerful tool for visualizing the spatial distribution of axonal and dendritic densities, for estimating the number of potential synapses between neurons with low standard deviation, and for obtaining a greater understanding of the relationship between neural morphology and network connectivity.

Independently outgrowing neurons and geometry-based synapse formation produce networks with realistic synaptic connectivity.

Neuronal signal integration and information processing in cortical networks critically depend on the organization of synaptic connectivity. During development, neurons can form synaptic connections when their axonal and dendritic arborizations come within close proximity of each other. Although many signaling cues are thought to be involved in guiding neuronal extensions, the extent to which accidental appositions between axons and dendrites can already account for synaptic connectivity remains unclear. To investigate this, we generated a local network of cortical L2/3 neurons that grew out independently of each other and that were not guided by any extracellular cues. Synapses were formed when axonal and dendritic branches came by chance within a threshold distance of each other. Despite the absence of guidance cues, we found that the emerging synaptic connectivity showed a good agreement with available experimental data on spatial locations of synapses on dendrites and axons, number of synapses by which neurons are connected, connection probability between neurons, distance between connected neurons, and pattern of synaptic connectivity. The connectivity pattern had a small-world topology but was not scale free. Together, our results suggest that baseline synaptic connectivity in local cortical circuits may largely result from accidentally overlapping axonal and dendritic branches of independently outgrowing neurons.

Comparison of targeted maximum likelihood and shrinkage estimators of parameters in gene networks.

Gene regulatory networks, in which edges between nodes describe interactions between transcription factors (TFs) and their target genes, model regulatory interactions that determine the cell-type and condition-specific expression of genes. Regression methods can be used to identify TF-target gene interactions from gene expression and DNA sequence data. The response variable, i.e. observed gene expression, is modeled as a function of many predictor variables simultaneously. In practice, it is generally not possible to select a single model that clearly achieves the best fit to the observed experimental data and the selected models typically contain overlapping sets of predictor variables. Moreover, parameters that represent the marginal effect of the individual predictors are not always present. In this paper, we use the statistical framework of estimation of variable importance to define variable importance as a parameter of interest and study two different estimators of this parameter in the context of gene regulatory networks. On yeast data we show that the resulting parameter has a biologically appealing interpretation. We apply the proposed methodology on mammalian gene expression data to gain insight into the temporal activity of TFs that underly gene expression changes in F11 cells in response to Forskolin stimulation.

Identification of context-specific gene regulatory networks with GEMULA--gene expression modeling using LAsso.

MOTIVATION: Gene regulatory networks, in which edges between nodes describe interactions between transcriptional regulators and their target genes, determine the coordinated spatiotemporal expression of genes. Especially in higher organisms, context-specific combinatorial regulation by transcription factors (TFs) is believed to determine cellular states and fates. TF-target gene interactions can be studied using high-throughput techniques such as ChIP-chip or ChIP-Seq. These experiments are time and cost intensive, and further limited by, for instance, availability of high affinity TF antibodies. Hence, there is a practical need for methods that can predict TF-TF and TF-target gene interactions in silico, i.e. from gene expression and DNA sequence data alone. We propose GEMULA, a novel approach based on linear models to predict TF-gene expression associations and TF-TF interactions from experimental data. GEMULA is based on linear models, fast and considers a wide range of biologically plausible models that describe gene expression data as a function of predicted TF binding to gene promoters. RESULTS: We show that models inferred with GEMULA are able to explain roughly 70% of the observed variation in gene expression in the yeast heat shock response. The functional relevance of the inferred TF-TF interactions in these models are validated by different sources of independent experimental evidence. We also have applied GEMULA to an in vitro model of neuronal outgrowth. Our findings confirm existing knowledge on gene regulatory interactions underlying neuronal outgrowth, but importantly also generate new insights into the temporal dynamics of this gene regulatory network that can now be addressed experimentally. AVAILABILITY: The GEMULA R-package is available from http://www.few.vu.nl/~degunst/gemula_1.0.tar.gz.

Novel candidate genes associated with hippocampal oscillations.

The hippocampus is critical for a wide range of emotional and cognitive behaviors. Here, we performed the first genome-wide search for genes influencing hippocampal oscillations. We measured local field potentials (LFPs) using 64-channel multi-electrode arrays in acute hippocampal slices of 29 BXD recombinant inbred mouse strains. Spontaneous activity and carbachol-induced fast network oscillations were analyzed with spectral and cross-correlation methods and the resulting traits were used for mapping quantitative trait loci (QTLs), i.e., regions on the genome that may influence hippocampal function. Using genome-wide hippocampal gene expression data, we narrowed the QTLs to eight candidate genes, including Plcb1, a phospholipase that is known to influence hippocampal oscillations. We also identified two genes coding for calcium channels, Cacna1b and Cacna1e, which mediate presynaptic transmitter release and have not been shown to regulate hippocampal network activity previously. Furthermore, we showed that the amplitude of the hippocampal oscillations is genetically correlated with hippocampal volume and several measures of novel environment exploration.

LLM3D: a log-linear modeling-based method to predict functional gene regulatory interactions from genome-wide expression data.

All cellular processes are regulated by condition-specific and time-dependent interactions between transcription factors and their target genes. While in simple organisms, e.g. bacteria and yeast, a large amount of experimental data is available to support functional transcription regulatory interactions, in mammalian systems reconstruction of gene regulatory networks still heavily depends on the accurate prediction of transcription factor binding sites. Here, we present a new method, log-linear modeling of 3D contingency tables (LLM3D), to predict functional transcription factor binding sites. LLM3D combines gene expression data, gene ontology annotation and computationally predicted transcription factor binding sites in a single statistical analysis, and offers a methodological improvement over existing enrichment-based methods. We show that LLM3D successfully identifies novel transcriptional regulators of the yeast metabolic cycle, and correctly predicts key regulators of mouse embryonic stem cell self-renewal more accurately than existing enrichment-based methods. Moreover, in a clinically relevant in vivo injury model of mammalian neurons, LLM3D identified peroxisome proliferator-activated receptor γ (PPARγ) as a neuron-intrinsic transcriptional regulator of regenerative axon growth. In conclusion, LLM3D provides a significant improvement over existing methods in predicting functional transcription regulatory interactions in the absence of experimental transcription factor binding data.

Growth references for height, weight and body mass index of twins aged 0-2.5 years.

AIM: To determine the size of the growth deficit in Dutch monozygotic and dizygotic twins aged 0-2.5 years as compared to singletons and to construct reference growth charts for twins. METHODS: Growth of twins was studied using longitudinal data on over 4000 twins aged 0-2.5 years of the Netherlands Twin Register. The LMS method was used to obtain growth references for length/height, weight, and body mass index (BMI) for twins. RESULTS: During the first 2.5 years of age, differences in length/height and weight between twins and singletons decrease but do not disappear. BMI of twins deviates less than that of singletons. Approximately half of the growth retardation from birth until 1.5 years of age was attributable to gestational age. Between 1.5 years and 2.5 years of age, this difference was reduced to one-third. Thus, a substantial part of the growth difference could not be explained by gestational age. CONCLUSIONS: During the first 2.5 years of life, there is a difference in growth between twins and singletons. Twins catch up in their body size, i.e. they grow faster after birth, but do not yet achieve the same height and weight till they reach 2.5 years of age. We recommend the use of the growth references for twins.

Identification of candidate transcriptional modulators involved in successful regeneration after nerve injury.

Successful regeneration of injured neurons requires a complex molecular response that involves the expression, modification and transport of a large number of proteins. The identity of neuronal proteins responsible for the initiation of regenerative neurite outgrowth is largely unknown. Dorsal root ganglion (DRG) neurons display robust and successful regeneration following lesion of their peripheral neurite, whereas outgrowth of central neurites is weak and does not lead to functional recovery. We have utilized this differential response to gain insight in the early transcriptional events associated with successful regeneration. Surprisingly, our study shows that peripheral and central nerve crushes elicit very distinct transcriptional activation, revealing a large set of novel genes that are differentially regulated within the first 24 h after the lesion. Here we show that Ankrd1, a gene known to act as a transcriptional modulator, is involved in neurite outgrowth of a DRG neuron-derived cell line as well as in cultured adult DRG neurons. This gene, and others identified in this study, may be part of the transcriptional regulatory module that orchestrates the onset of successful regeneration.

On the power for linkage detection using a test based on scan statistics.

We analyze some aspects of scan statistics, which have been proposed to help for the detection of weak signals in genetic linkage analysis. We derive approximate expressions for the power of a test based on moving averages of the identity by descent allele sharing proportions for pairs of relatives at several contiguous markers. We confirm these approximate formulae by simulation. The results show that when there is a single trait-locus on a chromosome, the test based on the scan statistic is slightly less powerful than that based on the customary allele sharing statistic. On the other hand, if two genes having a moderate effect on a trait lie close to each other on the same chromosome, scan statistics improve power to detect linkage.

Genetic study of the height and weight process during infancy.

Longitudinal height and weight data from 4649 Dutch twin pairs between birth and 2.5 years of age were analyzed. The data were first summarized into parameters of a polynomial of degree 4 by a mixed-effects procedure. Next, the variation and covariation in the parameters of the growth curve (size at one year of age, growth velocity, deceleration of growth, rate of change in deceleration [i.e., jerk] and rate of change in jerk [i.e., snap]) were decomposed into genetic and nongenetic sources. Additionally, the variation in the estimated size at birth and at 2 years of age interpolated from the polynomial was decomposed into genetic and nongenetic components. Variation in growth was best characterized by a genetic model which included additive genetic, common environmental and specific environmental influences, plus effects of gestational age. The effect of gestational age was largest for size at birth, explaining 39% of the variance. The differences between monozygotic and dizygotic twin correlations were largest for size at 1 and 2 years of age and growth velocity of weight, which suggests that these parameters are more influenced by heritability than size at birth, deceleration and jerk. The percentage of variance explained by additive genetic influences for height at 2 years of age was 52% for females and 58% for males. For weight at 2 years of age, heritability was approximately 58% for both sexes. Variation in snap height for males was also mainly influenced by additive genetic factors, while snap for females was influenced by both additive genetic and common environmental factors. The correlations for the additive genetic and common environmental factors for deceleration and snap are large, indicating that these parameters are almost entirely under control of the same additive genetic and common environmental factors. Female jerk and snap, and also female height at birth and height at 2 years of age, are mostly under control of the same additive genetic factor.

Morphine exposure and abstinence define specific stages of gene expression in the rat nucleus accumbens.

Intermittent exposure to addictive drugs causes long-lasting changes in responsiveness to these substances due to persistent molecular and cellular alterations within the meso-corticolimbic system. In this report, we studied the expression profiles of 159 genes in the rat nucleus accumbens during morphine exposure (14 days, 10 mg/kg s.c.) and drug-abstinence (3 weeks). We used real-time quantitative PCR to monitor gene expression after establishing its sensitivity and resolution to resolve small changes in expression for genes in various abundance classes. Morphine-exposure (5 time points) and subsequent abstinence (6 time points) induced phase-specific temporal gene expression of distinct functional groups of genes, for example, short-term homeostatic responses. Opiate withdrawal appeared to be a new stimulus in terms of gene expression and mediates a marked wave of gene repression. Prolonged abstinence resulted in persistently changed expression levels of genes involved in neuronal outgrowth and re-wiring. Our findings substantiate the hypothesis that this new gene program, initiated upon morphine-withdrawal, may subserve long-term neuronal plasticity involved in the persistent behavioral consequences of repeated drug-exposure.

Search for haplotype interactions that influence susceptibility to type 1 diabetes, through use of unphased genotype data.

Type 1 diabetes is a T-cell-mediated chronic disease characterized by the autoimmune destruction of pancreatic insulin-producing beta cells and complete insulin deficiency. It is the result of a complex interrelation of genetic and environmental factors, most of which have yet to be identified. Simultaneous identification of these genetic factors, through use of unphased genotype data, has received increasing attention in the past few years. Several approaches have been described, such as the modified transmission/disequilibrium test procedure, the conditional extended transmission/disequilibrium test, and the stepwise logistic-regression procedure. These approaches are limited either by being restricted to family data or by ignoring so-called "haplotype interactions" between alleles. To overcome this limit, the present study provides a general method to identify, on the basis of unphased genotype data, the haplotype blocks that interact to define the risk for a complex disease. The principle underpinning the proposal is minimal entropy. The performance of our procedure is illustrated for both simulated and real data. In particular, for a set of Dutch type 1 diabetes data, our procedure suggests some novel evidence of the interactions between and within haplotype blocks that are across chromosomes 1, 2, 3, 4, 5, 6, 7, 8, 11, 12, 15, 16, 17, 19, and 21. The results demonstrate that, by considering interactions between potential disease haplotype blocks, we may succeed in identifying disease-predisposing genetic variants that might otherwise have remained undetected.

Exploring heterogeneity in tumour data using Markov chain Monte Carlo.

We describe a Bayesian approach to incorporate between-individual heterogeneity associated with parameters of complicated biological models. We emphasize the use of the Markov chain Monte Carlo (MCMC) method in this context and demonstrate the implementation and use of MCMC by analysis of simulated overdispersed Poisson counts and by analysis of an experimental data set on preneoplastic liver lesions (their number and sizes) in the presence of heterogeneity. These examples show that MCMC-based estimates, derived from the posterior distribution with uniform priors, may agree well with maximum likelihood estimates (if available). However, with heterogeneous parameters, maximum likelihood estimates can be difficult to obtain, involving many integrations. In this case, the MCMC method offers substantial computational advantages.