Monte Carlo EM Algorithm Research Articles

Efficient inference for nonlinear state space models: An automatic sample size selection rule

This paper studies the maximum likelihood estimation of nonlinear state space models. Particle Markov chain Monte Carlo method is introduced to implement the Monte Carlo expectation maximization algorithm for more accurate and robust estimation. Under this framework, an automated sample size selection criterion is constructed via renewal theory. This criterion would increase the sample size when the relative likelihood indicates that the parameters are close to each other. The proposed methodology is applied to the stochastic volatility model and another nonlinear state space model for illustration, where the results show better estimation performance.

Joint modelling of longitudinal binary data and survival data

The medical costs in an ageing society substantially increase when the incidences of chronic diseases, disabilities and inability to live independently are high. Healthy lifestyles not only affect elderly individuals but also influence the entire community. When assessing treatment efficacy, survival and quality of life should be considered simultaneously. This paper proposes the joint likelihood approach for modelling survival and longitudinal binary covariates simultaneously. Because some unobservable information is present in the model, the Monte Carlo EM algorithm and Metropolis-Hastings algorithm are used to find the estimators. Monte Carlo simulations are performed to evaluate the performance of the proposed model based on the accuracy and precision of the estimates. Real data are used to demonstrate the feasibility of the proposed model.

A pathway for multivariate analysis of ecological communities using copulas.

We describe a new pathway for multivariate analysis of data consisting of counts of species abundances that includes two key components: copulas, to provide a flexible joint model of individual species, and dissimilarity‐based methods, to integrate information across species and provide a holistic view of the community. Individual species are characterized using suitable (marginal) statistical distributions, with the mean, the degree of over‐dispersion, and/or zero‐inflation being allowed to vary among a priori groups of sampling units. Associations among species are then modeled using copulas, which allow any pair of disparate types of variables to be coupled through their cumulative distribution function, while maintaining entirely the separate individual marginal distributions appropriate for each species. A Gaussian copula smoothly captures changes in an index of association that excludes joint absences in the space of the original species variables. A permutation‐based filter with exact family‐wise error can optionally be used a priori to reduce the dimensionality of the copula estimation problem. We describe in detail a Monte Carlo expectation maximization algorithm for efficient estimation of the copula correlation matrix with discrete marginal distributions (counts). The resulting fully parameterized copula models can be used to simulate realistic ecological community data under fully specified null or alternative hypotheses. Distributions of community centroids derived from simulated data can then be visualized in ordinations of ecologically meaningful dissimilarity spaces. Multinomial mixtures of data drawn from copula models also yield smooth power curves in dissimilarity‐based settings. Our proposed analysis pathway provides new opportunities to combine model‐based approaches with dissimilarity‐based methods to enhance understanding of ecological systems. We demonstrate implementation of the pathway through an ecological example, where associations among fish species were found to increase after the establishment of a marine reserve.

Parameter-driven state-space model for integer-valued time series with application

ABSTRACTTime series of counts occur in many different contexts, the counts being usually of certain events or objects in specified time intervals. In this paper we introduce a model called parameter-driven state-space model to analyse integer-valued time series data. A key property of such model is that the distribution of the observed count data is independent, conditional on the latent process, although the observations are correlated marginally. Our simulation shows that the Monte Carlo Expectation Maximization (MCEM) algorithm and the particle method are useful for the parameter estimation of the proposed model. In the application to Malaysia dengue data, our model fits better when compared with several other models including that of Yang et al. (2015)

Likelihood-based quantile mixed effects models for longitudinal data with multiple features via MCEM algorithm

In longitudinal data analysis, many data features are frequently encountered in practice. For example, variates are usually measured with censoring, substantial errors and non-normal feature, etc. In such case, it is impossible to give a precise result by conducting statistical inference of separate analysis. The joint modeling may be a considerable alternative. The superiority of joint modeling is that it can implement a simultaneous analysis for longitudinal mixed models with random effects, censoring, nor-normal errors, errors in covariates and heavy-tailed feature by pooling all data information together. However, most of traditional modeling methods aim at depicting the average variation of outcome variable conditionally on covariates, which may result in non-robust estimation results when suffering outliers or non-normal errors. Quantile regression provides an attractive alternative to model longitudinal data with multiple data features, which can draw a complete picture of the conditional distributions of outcome variable and bring out robust estimation results. In this paper, we conduct the likelihood-based joint quantile regression for longitudinal mixed models by accounting for the above multiple data features simultaneously. Based on the asymmetric Laplace distribution (ALD), the Monte Carlo Expectation-Maximization (MCEM) algorithm is employed to address the estimation problem. Finally, some Monte Carlo simulations are conducted and an AIDS data analysis is provided to illustrate the developed procedures.

Letter to the Editor.

Galarza, Lachos and Bandyopadhyay (2017) have recently proposed a method of estimating linear quantile mixed models (Geraci and Bottai, 2014) based on a Monte Carlo EM algorithm. They assert that their procedure represents an improvement over the numerical quadrature and non-smooth optimization approach implemented by Geraci (2014). The objective of this note is to demonstrate that this claim is incorrect. We also point out several inaccuracies and shortcomings in their paper which affect other results and conclusions that can be drawn.

A semi-supervised clustering-based approach for stratification identification using borehole and cone penetration test data

Borehole drilling and cone penetration test (CPT) are frequently employed site investigation methods for identifying subsurface stratification. However, these two methods have their respective pros and cons, and their corresponding soil type classification protocols are different. Therefore, an approach that can jointly interpret raw data from both investigation methods and provide unified soil classification results is in great demand. Motivated by the aforementioned point, this paper presents a novel semi-supervised clustering based stratification identification approach using information from both boreholes and CPT logs. The proposed approach is established on a hidden Markov random field (HMRF) framework so that the supervision constraints could be introduced by using borehole data during the clustering of CPT sounding samples. Further, the presented approach employs a Monte Carlo Expectation Maximization (MCEM) algorithm to perform the clustering process, which enables estimating the subsurface stratification in a probabilistic manner. The performances of the proposed approach are evaluated using real-world site investigation data. The test results indicate that the proposed approach is effective and robust for identifying subsurface stratification.

Spatial models for non‐Gaussian data with covariate measurement error

AbstractSpatial models have been widely used in the public health setup. In the case of continuous outcomes, the traditional approaches to model spatial data are based on the Gaussian distribution. This assumption might be overly restrictive to represent the data. The real data could be highly non‐Gaussian and may show features like heavy tails and/or skewness. In spatial data modeling, it is also commonly assumed that the covariates are observed without errors, but for various reasons, such as measurement techniques or instruments used, uncertainty is inherent in spatial (especially geostatistics) data, and so, these data are susceptible to measurement errors in the covariates of interest. In this paper, we introduce a general class of spatial models with covariate measurement error that can account for heavy tails, skewness, and uncertainty of the covariates. A likelihood method, which leads to the maximum likelihood estimation approach, is used for inference through the Monte Carlo expectation–maximization algorithm. The predictive distribution at nonsampled sites is approximated based on the Markov chain Monte Carlo algorithm. The proposed approach is evaluated through a simulation study and by a real application (particulate matter data set).

Maximum Likelihood Estimation in a Semicontinuous Survival Model with Covariates Subject to Detection Limits.

Semicontinuous data are common in biological studies, occurring when a variable is continuous over a region but has a point mass at one or more points. In the motivating Genetic and Inflammatory Markers of Sepsis (GenIMS) study, it was of interest to determine how several biomarkers subject to detection limits were related to survival for patients entering the hospital with community acquired pneumonia. While survival times were recorded for all individuals in the study, the primary endpoint of interest was the binary event of 90-day survival, and no patients were lost to follow-up prior to 90 days. In order to use all of the available survival information, we propose a two-part regression model where the probability of surviving to 90 days is modeled using logistic regression and the survival distribution for those experiencing the event prior to this time is modeled with a truncated accelerated failure time model. We assume a series of mixture of normal regression models to model the joint distribution of the censored biomarkers. To estimate the parameters in this model, we suggest a Monte Carlo EM algorithm where multiple imputations are generated for the censored covariates in order to estimate the expectation in the E-step and then weighted maximization is applied to the observed and imputed data in the M-step. We conduct simulations to assess the proposed model and maximization method, and we analyze the GenIMS data set.

Extending Tests of Random Effects to Assess for Measurement Invariance in Factor Models.

Factor analysis models are widely used in health research to summarize hard to measure predictor or outcome variable constructs. For example, in the ELEMENT study, factor models are used to summarize lead exposure biomarkers which are thought to indirectly measure prenatal exposure to lead. Classic latent factor models are fitted assuming that factor loadings are constant across all covariate levels (e.g., maternal age in ELEMENT); that is, measurement invariance (MI) is assumed. When the MI is not met, measurement bias is introduced. Traditionally, MI is examined by defining subgroups of the data based on covariates, fitting multi-group factor analysis, and testing differences in factor loadings across covariate groups. In this paper, we develop novel tests of measurement invariance by modeling the factor loadings as varying coeffcients, i.e., letting the factor loading vary across continuous covariate values instead of groups. These varying coeffcients are estimated using penalized splines, where spline coeffcients are penalized by treating them as random coeffcients. The test of MI is then carried out by conducting a likelihood ratio test for the null hypothesis that the variance of the random spline coeffcients equals zero. We use a Monte-Carlo EM algorithm for estimation, and obtain the likelihood using Monte-Carlo in tegration. Using simulations, we compare the Type I error and power of our testing approach and the multi-group testing method. We apply the proposed methods to to summarize data on prenatal biomarkers of lead exposure from the ELEMENT study and find violations of MI due to maternal age.

INVESTIGATING BARGAINING POWER OF FARMERS AND PROCESSORS IN IRAN'S DAIRY MARKET

AbstractThe farm-gate price of raw milk in Iran is determined annually in negotiations among representatives of dairy processors, milk producers, and government officials. This study estimates the average bargaining power of dairy farmers and processors, through applying the generalized axiomatic Nash approach in a bilateral bargaining model. We employ annual data from 1990 to 2013 to estimate econometric representation of a bilateral bargaining model using a Monte Carlo expectation maximization algorithm. Results imply a higher bargaining power of 0.69 for processors, compared with 0.31 for farmers. This asymmetry of bargaining power causes unequal allocation of gains in the milk market.

A coupled ETAS-I2GMM point process with applications to seismic fault detection

Epidemic-type aftershock sequence (ETAS) point process is a common model for the occurrence of earthquake events. The ETAS model consists of a stationary background Poisson process modeling spontaneous earthquakes and a triggering kernel representing the space–time-magnitude distribution of aftershocks. Popular nonparametric methods for estimation of the background intensity include histograms and kernel density estimators. While these methods are able to capture local spatial heterogeneity in the intensity of spontaneous events, they do not capture well patterns resulting from fault line structure over larger spatial scales. Here we propose a two-layer infinite Gaussian mixture model for clustering of earthquake events into fault-like groups over intermediate spatial scales. We introduce a Monte Carlo expectation-maximization (EM) algorithm for joint inference of the ETAS-I2GMM model and then apply the model to the Southern California Earthquake Catalog. We illustrate the advantages of the ETAS-I2GMM model in terms of both goodness of fit of the intensity and recovery of fault line clusters in the Community Fault Model 3.0 from earthquake occurrence data.

A joint quantile regression model for multiple longitudinal outcomes

Complexity of longitudinal data lies in the inherent dependence among measurements from same subject over different time points. For multiple longitudinal responses, the problem is challenging due to inter-trait and intra-trait dependence. While linear mixed models are popularly used for analysing such data, appropriate inference on the shape of the population cannot be drawn for non-normal data sets. We propose a linear mixed model for joint quantile regression of multiple longitudinal responses. We consider an asymmetric Laplace distribution for quantile regression and estimate model parameters by Monte Carlo EM algorithm. Nonparametric bootstrap resampling method is used for estimating confidence intervals of parameter estimates. Through extensive simulation studies, we investigate the operating characteristics of our proposed model and compare the performance to other traditional quantile regression models. We apply proposed model for analysing data from nutrition education programme on hypercholesterolemic children of the USA.

Joint modelling for organ transplantation outcomes for patients with diabetes and the end-stage renal disease.

This article is motivated by jointly modelling longitudinal and time-to-event clinical data of patients with diabetes and end-stage renal disease. All patients are on the waiting list for the pancreas transplant after kidney transplant, and some of them have a pancreas transplant before kidney transplant failure or death. Scant literature has studied the dynamical joint relationship of the estimated glomerular filtration rates trajectory, the effect of pancreas transplant, and time-to-event outcomes, although it remains an important clinical question. In an attempt to describe the association in the multiple outcomes, we propose a new joint model with a longitudinal submodel and an accelerated failure time submodel, which are linked by some latent variables. The accelerated failure time submodel is used to determine the relationship of the time-to-event outcome with all predictors. In addition, the piecewise linear function in the survival submodel is used to calculate the dynamic hazard ratio curve of a time-dependent side event, because the effect of the side event on the time-to-event outcome is non-proportional. The model parameters are estimated with a Monte Carlo EM algorithm. The finite sample performance of the proposed method is investigated in simulation studies. Our method is demonstrated by fitting the joint model for the clinical data of 13,635 patients with diabetes and the end-stage renal disease.

Latent mixture modeling for clustered data

This article proposes a mixture modeling approach to estimating cluster-wise conditional distributions in clustered (grouped) data. We adapt the mixture-of-experts model to the latent distributions, and propose a model in which each cluster-wise density is represented as a mixture of latent experts with cluster-wise mixing proportions distributed as Dirichlet distribution. The model parameters are estimated by maximizing the marginal likelihood function using a newly developed Monte Carlo Expectation-Maximization algorithm. We also extend the model such that the distribution of cluster-wise mixing proportions depends on some cluster-level covariates. The finite sample performance of the proposed model is compared with some existing mixture modeling approaches as well as linear mixed model through the simulation studies. The proposed model is also illustrated with the posted land price data in Japan.

A quantile variant of the expectation–maximization algorithm and its application to parameter estimation with interval data

The expectation–maximization algorithm is a powerful computational technique for finding the maximum likelihood estimates for parametric models when the data are not fully observed. The expectation–maximization is best suited for situations where the expectation in each E-step and the maximization in each M-step are straightforward. A difficulty with the implementation of the expectation–maximization algorithm is that each E-step requires the integration of the log-likelihood function in closed form. The explicit integration can be avoided by using what is known as the Monte Carlo expectation–maximization algorithm. The Monte Carlo expectation–maximization uses a random sample to estimate the integral at each E-step. But the problem with the Monte Carlo expectation–maximization is that it often converges to the integral quite slowly and the convergence behavior can also be unstable, which causes computational burden. In this paper, we propose what we refer to as the quantile variant of the expectation–maximization algorithm. We prove that the proposed method has an accuracy of [Formula: see text], while the Monte Carlo expectation–maximization method has an accuracy of [Formula: see text]. Thus, the proposed method possesses faster and more stable convergence properties when compared with the Monte Carlo expectation–maximization algorithm. The improved performance is illustrated through the numerical studies. Several practical examples illustrating its use in interval-censored data problems are also provided.

Survival analysis of DNA mutation motifs with penalized proportional hazards.

Antibodies, an essential part of our immune system, develop through an intricate process to bind a wide array of pathogens. This process involves randomly mutating DNA sequences encoding these antibodies to find variants with improved binding, though mutations are not distributed uniformly across sequence sites. Immunologists observe this nonuniformity to be consistent with "mutation motifs", which are short DNA subsequences that affect how likely a given site is to experience a mutation. Quantifying the effect of motifs on mutation rates is challenging: a large number of possible motifs makes this statistical problem high dimensional, while the unobserved history of the mutation process leads to a nontrivial missing data problem. We introduce an ℓ 1-penalized proportional hazards model to infer mutation motifs and their effects. In order to estimate model parameters, our method uses a Monte Carlo EM algorithm to marginalize over the unknown ordering of mutations. We show that our method performs better on simulated data compared to current methods and leads to more parsimonious models. The application of proportional hazards to mutation processes is, to our knowledge, novel and formalizes the current methods in a statistical framework that can be easily extended to analyze the effect of other biological features on mutation rates.

High-dimensional variable selection with the plaid mixture model for clustering

With high-dimensional data, the number of covariates is considerably larger than the sample size. We propose a sound method for analyzing these data. It performs simultaneously clustering and variable selection. The method is inspired by the plaid model. It may be seen as a multiplicative mixture model that allows for overlapping clustering. Unlike conventional clustering, within this model an observation may be explained by several clusters. This characteristic makes it specially suitable for gene expression data. Parameter estimation is performed with the Monte Carlo expectation maximization algorithm and importance sampling. Using extensive simulations and comparisons with competing methods, we show the advantages of our methodology, in terms of both variable selection and clustering. An application of our approach to the gene expression data of kidney renal cell carcinoma taken from The Cancer Genome Atlas validates some previously identified cancer biomarkers.

A Non-Linear Model for Censored and Mismeasured Time Varying Covariates in Survival Models, with Applications in Human Immunodeficiency Virus and Acquired Immune Deficiency Syndrome Studies

Summary In survival regression analysis, when the time-dependent covariates are censored and measured with errors, a joint model is often considered for the longitudinal covariate data and the survival data. Typically, an empirical linear (mixed) model is assumed for the time-dependent covariates. However, such an empirical linear covariate model may be inappropriate for the (unobserved) censored covariate values that may behave quite differently from the observed covariate process. In applications such as human immunodeficiency virus–acquired immune deficiency syndrome studies, a mechanistic non-linear model can be derived for the covariate process on the basis of the underlying data generation mechanisms and such a non-linear covariate model may provide better ‘predictions’ for the censored and mismeasured covariate values. We propose a joint Cox and non-linear mixed effect model to model survival data with censored and mismeasured time varying covariates. We use likelihood methods for inference, implemented by the Monte Carlo EM algorithm. The models and methods are evaluated by simulations. An acquired immune deficiency syndrome data set is analysed in detail, where the time-dependent covariate is a viral load which may be censored because of a lower detection limit and may also be measured with errors. The results based on linear and non-linear covariate models are compared and new insights are gained.

Joint analysis of left-censored longitudinal biomarker and binary outcome via latent class modeling.

Joint latentclass modeling is an appealing approach for evaluating the association between a longitudinal biomarker and clinical outcome when the study population is heterogeneous. The link between the biomarker trajectory and the risk of event is reflected by the latent classes, which accommodate the underlying population heterogeneity. The estimation of joint latent class models may be complicated by the censored data in the biomarker measurements due to detection limits. We propose a modified likelihood function under the parametric assumption of biomarker distribution and develop a Monte Carlo expectation-maximization algorithm for joint analysis of a biomarker and a binary outcome. We conduct simulation studies to demonstrate the satisfactory performance of our Monte Carlo expectation-maximization algorithm and the superiority of our method to the naive imputation method for handling censored biomarker data. In addition, we apply our method to the Genetic and Inflammatory Markers of Sepsis study to investigate the role of inflammatory biomarker profile in predicting 90-day mortality for patients hospitalized with community-acquired pneumonia.