The sensitivity analysis of stochastic computer codes poses several challenges. One of them is the computation of Sobol indices of the expected output of a code of interest, commonly done from a Monte Carlo sample with two nested loops. The number of iterations of each one of them need to go to infinity for the estimator to converge, which is a major hindrance in practice. Here it is shown that Sobol indices can be estimated from a Monte Carlo sample with a single loop. The estimator is asymptotically normal and unbiased as the number of iterations goes to infinity.
Reconstructing gene regulatory networks from large-scale heterogeneous data is a key challenge in biology. In multi-omics data analysis, networks based on pairwise statistical association measures remain popular, as they are easy to build and understand. In the presence of mixed-type (discrete and continuous) data, however, the choice of good association measures remains an important issue. It is proposed a novel approach based on the Gaussian copula, the parameters of which represent the links of the network. Novel properties of the model are obtained to guide the interpretation of the network. To estimate the copula parameters, a semiparametric pairwise likelihood for mixed data was calculated. An extensive simulation study showed that the proposed estimation procedure was able to accurately estimate the copula correlation matrix. The proposed methodology was also applied to a real ICGC dataset on breast cancer, and is implemented in a freely available R package heterocop.
BACKGROUND:Inferring partial correlation networks is essential in systems biology to uncover direct interactions between biological entities. Traditional Gaussian graphical models rely on the assumption of normally distributed data; this assumption is not satisfied when dealing with multi-omics datasets comprising heterogeneous data types such as continuous and discrete variables. RESULTS:We propose a novel likelihood-based approach for network inference using a Gaussian copula model with semiparametric pairwise-likelihood estimation of the latent correlation matrix. The inferred correlation structure is then inverted and regularized via the graphical lasso to recover latent partial correlations. Compared to a moment-based approach employing bridge functions, our method demonstrates significantly improved computational efficiency and estimation accuracy, particularly for discrete data with many categories and/or large values, such as count data. This result is important for biological applications, especially for the integration of RNA-seq count data. An application to a breast cancer data set from the International Cancer Genome Consortium (ICGC) successfully identified biologically relevant interactions. CONCLUSIONS:The proposed approach, based on the Gaussian copula and likelihood-based estimation, provides a novel, effective and computationally efficient mathematical framework for integrative multi-omics data analysis and network inference.
Although there is a plethora of methods to estimate sensitivity indices associated with individual inputs, there is much less work on interaction effects of every order, especially when it comes to make inferences about the true underlying values of the indices. To fill this gap, a method that allows one to make such inferences simultaneously from a Monte Carlo sample is given. One advantage of this method is its simplicity: it leverages the fact that Shapley effects and Sobol indices are only linear transformations of total indices, so that standard asymptotic theory suffices to get confidence intervals and to carry out statistical tests. To perform the numerical computations efficiently, M & ouml;bius inversion formulas are used, and linked to the fast M & ouml;bius transform algorithm. The method is illustrated on two dynamical systems, both with an application in life sciences: a Boolean network modeling a cellular decision-making process involving 12 inputs, and a system of ordinary differential equations modeling some population dynamics involving 10 inputs.
Sobol' sensitivity index estimators for stochastic models are functions of nested Monte Carlo estimators, which are estimators built from two nested Monte Carlo loops. The outer loop explores the input space and, for each of the explorations, the inner loop repeats model runs to estimate conditional expectations. Although the optimal allocation between explorations and repetitions of one's computational budget is well-known for nested Monte Carlo estimators, it is less clear how to deal with functions of nested Monte Carlo estimators, especially when those functions have unbounded Hessian matrices, as it is the case for Sobol' index estimators. To address this problem, a regularization method is introduced to bound the mean squared error of functions of nested Monte Carlo estimators. Based on a heuristic, an allocation strategy that seeks to minimize a bias-variance trade-off is proposed. The method is applied to Sobol' index estimators for stochastic models. A practical algorithm that adapts to the level of intrinsic randomness in the models is given and illustrated on numerical experiments.
Reconstructing gene regulatory networks from large-scale heterogeneous data is a key challenge in biology. In multi-omics data analysis, networks based on pairwise statistical association measures remain popular, as they are easy to build and understand. In the presence of mixed-type (discrete and continuous) data, however, the choice of good association measures remains an important issue. We propose here a novel approach based on the Gaussian copula, the parameters of which represent the links of the network. Novel properties of the model are obtained to guide the interpretation of the network. To estimate the copula parameters, we calculated a semiparametric pairwise likelihood for mixed data. In an extensive simulation study, we showed that the proposed estimation procedure was able to accurately estimate the copula correlation matrix. The proposed methodology was also applied to a real ICGC dataset on breast cancer, and is implemented in a freely available R package heterocop.
It is well-known that Sobol indices, which count among the most popular sensitivity indices, are based on the Sobol decomposition. Here we challenge this construction by redefining Sobol indices without the Sobol decomposition. In fact, we show that Sobol indices are a particular instance of a more general concept which we call sensitivity measures. A sensitivity measure of a system taking inputs and returning outputs is a set function that is null at a subset of inputs if and only if, with probability one, the output actually does not depend on those inputs. A sensitivity measure evaluated at the whole set of inputs represents the uncertainty about the output. We show that measuring sensitivity to a particular subset is akin to measuring the expected output's uncertainty conditionally on the fact that the inputs belonging to that subset have been fixed to random values. By considering all of the possible combinations of inputs, sensitivity measures induce an implicit symmetric factorial experiment with two levels, the factorial effects of which can be calculated. This new paradigm generalizes many known sensitivity indices, can create new ones, and defines interaction effects independently of the choice of the sensitivity measure. No assumption about the distribution of the inputs is required.
In this manuscript, we consider a finite multivariate nonparametric mixture model where the dependence between the marginal densities is modeled using the copula device. Pseudo expectation–maximization (EM) stochastic algorithms were recently proposed to estimate all of the components of this model under a location-scale constraint on the marginals. Here, we introduce a deterministic algorithm that seeks to maximize a smoothed semiparametric likelihood. No location-scale assumption is made about the marginals. The algorithm is monotonic in one special case, and, in another, leads to "approximate monotonicity"—whereby the difference between successive values of the objective function becomes non-negative up to an additive term that becomes negligible after a sufficiently large number of iterations. The behavior of this algorithm is illustrated on several simulated and real datasets. The results suggest that, under suitable conditions, the proposed algorithm may indeed be monotonic in general. A discussion of the results and some possible future research directions round out our presentation.
In this paper we apply a methodology introduced in Navarro Jimenez et al (2016) in the framework of chemical reaction networks to perform a global sensitivity analysis on simulations of a continuous-time Markov chain model motivated by epidemiology. Our goal is to quantify not only the effects of uncertain parameters such as epidemic parameters (transmission rate, mean sojourn duration in compartments), but also those of intrinsic randomness and interactions between epidemic parameters and intrinsic randomness. For that purpose, following what was proposed in Navarro Jimenez et al, we leverage three exact simulation algorithms for continuous-time Markov chains from the state of the art which we combine with common tools from variance-based sensitivity analysis as introduced in Sobol (1993). Also, we discuss the impact of the choice of the simulation algorithm used for the simulations on the results of sensitivity analysis. Such a discussion is new, at least to our knowledge. In a numerical section, we implement and compare three sensitivity analyses based on simulations obtained from different exact simulation algorithms of a SARS-CoV-2 epidemic model.
Pairwise likelihood methods are commonly used for inference in parametric statistical models in cases where the full likelihood is too complex to be used, such as multivariate count data. Although pairwise likelihood methods represent a useful solution to perform inference for intractable likelihoods, several computational challenges remain. The pairwise likelihood function still requires the computation of a sum over all pairs of variables and all observations, which may be prohibitive in high dimensions. Moreover, it may be difficult to calculate confidence intervals of the resulting estimators, as they involve summing all pairs of pairs and all of the four-dimensional marginals. To alleviate these issues, we consider a randomized pairwise likelihood approach, where only summands randomly sampled across observations and pairs are used for the estimation. In addition to the usual tradeoff between statistical and computational efficiency, it is shown that, under a condition on the sampling parameter, this two-way random sampling mechanism makes the individual bivariate likelihood scores become asymptotically independent, allowing more computationally efficient confidence intervals to be constructed. The proposed approach is illustrated in tandem with copula-based models for multivariate count data in simulations, and in real data from a transcriptome study. Supplementary materials for this article are available online.
The fitness coefficient, introduced in this paper, results from a competition between parametric and nonparametric density estimators within the likelihood of the data. As illustrated on several real datasets, the fitness coefficient generally agrees with p-values but is easier to compute and interpret. Namely, the fitness coefficient can be interpreted as the proportion of data coming from the parametric model. Moreover, the fitness coefficient can be used to build a semiparamteric compromise which improves inference over the parametric and nonparametric approaches. From a theoretical perspective, the fitness coefficient is shown to converge in probability to one if the model is true and to zero if the model is false. From a practical perspective, the utility of the fitness coefficient is illustrated on real and simulated datasets.
A Markov tree is a probabilistic graphical model for a random vector indexed by the nodes of an undirected tree encoding conditional independence relations between variables. One possible limit distribution of partial maxima of samples from such a Markov tree is a max-stable Hüsler–Reiss distribution whose parameter matrix inherits its structure from the tree, each edge contributing one free dependence parameter. Our central assumption is that, upon marginal standardization, the data-generating distribution is in the max-domain of attraction of the said Hüsler–Reiss distribution, an assumption much weaker than the one that data are generated according to a graphical model. Even if some of the variables are unobservable (latent), we show that the underlying model parameters are still identifiable if and only if every node corresponding to a latent variable has degree at least three. Three estimation procedures, based on the method of moments, maximum composite likelihood, and pairwise extremal coefficients, are proposed for usage on multivariate peaks over thresholds data when some variables are latent. A typical application is a river network in the form of a tree where, on some locations, no data are available. We illustrate the model and the identifiability criterion on a data set of high water levels on the Seine, France, with two latent variables. The structured Hüsler–Reiss distribution is found to fit the observed extremal dependence patterns well. The parameters being identifiable we are able to quantify tail dependence between locations for which there are no data.
Sobol sensitivity indices assess how the output of a given mathematical model is sensitive to its inputs. If the model is stochastic, then it cannot be represented as a function of the inputs, thus raising questions about how to do a sensitivity analysis in those models. Practitioners have been using an approach that exploits the availability of methods for deterministic models. For each input, the stochastic model is repeated and the outputs are averaged. These averages are seen as if they came from a deterministic model and hence Sobol's method can be used. We show that the estimator so obtained is asymptotically biased if the number of repetitions goes to infinity too slowly. With limited computational resources, the number of repetitions of the stochastic model and the number of explorations of the input space cannot be large together and hence some balance must be found. We find the pair of numbers that minimizes a bound on some rank-based error criterion, penalizing bad rankings of the inputs' sensitivities. Also, under minimal distributional assumptions, we derive a functional relationship between the output, the input, and some random noise; the Sobol--Hoeffding decomposition can be applied to it to define a new sensitivity index, which asymptotically is estimated without bias even though the number of repetitions remains fixed. The theory is illustrated on numerical experiments.
Sensitivity analysis often accompanies computer modeling to understand what are the important factors of a model of interest. In particular , Sobol indices, naturally estimated by Monte-Carlo sampling, permit to quantify the contribution of the inputs to the variability of the output. However, when the model is stochastic, the problem of carrying out a sensitivity analysis remains open. There is no unique definition of Sobol indices and their estimation is more difficult because a good balance between repetitions of the computer code and explorations of the input space must be found. The problem of performing a sensitivity analysis for stochastic computer models with the Sobol method is addressed. Two Sobol indices are considered, their estimators constructed and their asymptotic properties established. An optimal balance between repetitions and explorations is proposed under a limited computing budget. A two-stage procedure is built: the first stage permits to find the optimal balance and the second stage produces the Sobol estimates based on the balance obtained in the first stage. The procedure is asymptotically oracle and the optimal convergence rates are derived. The theoretical results are tested with numerical experiments.
We present a general construction principle for copulas that is inspired by the celebrated Marshall–Olkin exponential model. From this general construction method, we derive special subclasses of copulas that could be useful in different situations and recall their main properties. Moreover, we discuss possible estimation strategy for the proposed copulas. The presented results are expected to be useful in the construction of stochastic models for lifetimes (e.g., in reliability theory) or in credit risk models.
A novel algorithm for performing inference and/or clustering in semiparametric copula-based mixture models is presented. The standard kernel density estimator is replaced by a weighted version that permits to take into account the constraints put on the underlying marginal densities. Lower misclassification error rates and better estimates are obtained on simulations. The pointwise consistency of the weighted kernel density estimator is established under an assumption on the rate of convergence of the sample maximum.