Monte Carlo EM Algorithm Research Articles

A Laplace-based model with flexible tail behavior

The proposed multiple scaled contaminated asymmetric Laplace (MSCAL) distribution is an extension of the multivariate asymmetric Laplace distribution to allow for a different excess kurtosis on each dimension and for more flexible shapes of the hyper-contours. These peculiarities are obtained by working on the principal component (PC) space. The structure of the MSCAL distribution has the further advantage of allowing for automatic PC-wise outlier detection – i.e., detection of outliers separately on each PC – when convenient constraints on the parameters are imposed. The MSCAL is fitted using a Monte Carlo expectation-maximization (MCEM) algorithm that uses a Monte Carlo method to estimate the orthogonal matrix of eigenvectors. A simulation study is used to assess the proposed MCEM in terms of computational efficiency and parameter recovery. In a real data application, the MSCAL is fitted to a real data set containing the anthropometric measurements of monozygotic/dizygotic twins. Both a skewed bivariate subset of the full data, perturbed by some outlying points, and the full data are considered.

A Proportional Intensity Model with Frailty for Missing Recurrent Failure Data

In some practical circumstances, data are recorded after the systems have begun operations, and data collection is stopped at a predetermined time or after a predetermined number of failures. In such circumstances, incompleteness of various types exists in the aspect of the missing number of failures and their occurrence times beyond the duration of the pilot study. Additionally, multiple repairable systems may present system-to-system variability caused by differences in the operating environments or working loads of individual systems. With respect to left-truncated and right-censored recurrent failure data from multiple repairable systems, we propose a reliability model based on a proportional intensity model with frailty. The frailty model explicitly models unobserved heterogeneity among systems. Covariates incorporated into the proportional intensity model additionally account for the heterogeneity between different operating conditions. To estimate the model parameters for the left-truncated and right-censored recurrent failure data, a Monte Carlo expectation maximization algorithm is proposed. Details of the estimation of the model parameters and the construction of their confidence intervals are examined. A real-world example and simulation studies under various scenarios show prominent applications of the proportional intensity model with frailty to left-truncated and right-censored multiple repairable systems for reliability prediction.

Statistical inference of multi-state transition model for longitudinal data with measurement error and heterogeneity

. Multi-state transition model is typically used to analyze longitudinal data in medicine and sociology. Moreover, variables in longitudinal studies usually are error-prone, and random effects are heterogeneous, which will result in biased estimates of the interest parameters. This article is intended to estimate the parameters of the multi-state transition model for longitudinal data with measurement error and heterogeneous random effects and further consider the covariate related to the covariance matrix of random effects is also error-prone when the covariate in the transition model is error-prone. We model the covariance matrix of random effects through the modified Cholesky decomposition and propose a pseudo-likelihood method based on the Monte Carlo expectation-maximization algorithm and the Bayesian method based on Markov Chain Monte Carlo to infer and calculate the whole estimates. Meanwhile, we obtain the asymptotic properties and evaluate the finite sample performance of the proposed method by simulation, which is well in terms of Bias, RMSE, and coverage rate of confidence intervals. In addition, we apply the proposed method to the MFUS data.

Misspecification analysis of gamma‐ and inverse Gaussian‐based perturbed degradation processes

AbstractAlbeit not equivalent, in many applications the gamma and the inverse Gaussian processes are treated as if they were. This circumstance makes the misspecification problem of these models interesting and important, especially when data are affected by measurement errors, since noisy/perturbed data do not allow to verify whether the selected model is actually able to adequately fit the real (hidden) degradation process. Motivated by the above considerations, in this paper we conduct a large Monte Carlo study to evaluate whether and how the presence of measurement errors affects this misspecification issue. The study is performed considering as reference models a perturbed gamma process recently proposed in the literature and a new perturbed inverse Gaussian process that share the same non‐Gaussian distributed error term. As an alternative option, we also analyze the more classical case where the error term is Gaussian distributed. We consider both the situation where the true model is the perturbed gamma and the one where it is the perturbed inverse Gaussian. Model parameters are estimated from perturbed data using the maximum likelihood method. Estimates are retrieved by using a new sequential Monte Carlo EM algorithm, which use allows to hugely mitigate the severe numerical issues posed by the direct maximization of the likelihood. The risk of incurring in a misspecification is evaluated as percentage of times the Akaike information criterion leads to select the wrong model. The severity of a misspecification is evaluated in terms of its impact on maximum likelihood estimate of the mean remaining useful life.

A composite Bayesian approach for quantile curve fitting with non-crossing constraints

To fit a set of quantile curves, Bayesian simultaneous quantile curve fitting methods face some challenges in properly specifying a feasible formulation and efficiently accommodating the non-crossing constraints. In this article, we propose a new minimization problem and develop its corresponding Bayesian analysis. The new minimization problem imposes two penalties to control not only the smoothness of fitted quantile curves but also the differences between quantile curves. This enables a direct inference on differences of quantile curves and facilitates improved information sharing among quantiles. After adopting B-spline approximation for the positive smoothing functions in the minimization problem, we specify the pseudo composite asymmetric Laplace likelihood and derive its priors. The computation algorithm, including partially collapsed Gibbs sampling for model parameters and Monte Carlo Expectation-Maximization algorithm for penalty parameters, are provided to carry out the proposed approach. The extensive simulation studies show that, compared with other candidate methods, the proposed approach yields more robust estimation. More advantages of the proposed approach are observed for the extreme quantiles, heavy-tailed random errors, and inference on the differences of quantiles. We also demonstrate the relative performances of the proposed approach and other competing methods through two real data analyses.

Spatio-Temporal Analysis and Prediction of Mass Telecommunication Base Station Failure Events

Large-scale telecommunication systems are lifeline infrastructure in modern society. A telecommunication system typically consists of a huge number of base stations with diverse geographical locations across a country, which highly complicates maintenance operations. To allocate maintenance resource properly, it is important to have a good understanding on the failure pattern of these base stations. Statistical inference of recurrent failures of these base stations is challenging because of the large number of base stations and the spatial correlation of their failure processes. Based on eight-month failure data of telecommunication base stations in Harbin, China, we propose a customized nonhomogeneous Poisson process (NHPP) model for recurrent failure data from telecommunication systems. The model consists of two layers, where the temporal layer applies an NHPP with station-specific frailty for failures of each base station, and the spatial layer uses a multivariate lognormal distribution to characterize the correlation among the frailties. The Monte Carlo EM (MCEM) algorithm is applied to estimate parameters included in the proposed model. We demonstrate the proposed model using the Harbin telecommunication system example with 7725 base stations and 4615 failure records.

Variable selection for joint models of multivariate skew-normal longitudinal and survival data.

Many joint models of multivariate skew-normal longitudinal and survival data have been presented to accommodate for the non-normality of longitudinal outcomes in recent years. But existing work did not consider variable selection. This article investigates simultaneous parameter estimation and variable selection in joint modeling of longitudinal and survival data. The penalized splines technique is used to estimate unknown log baseline hazard function, the rectangle integral method is adopted to approximate conditional survival function. Monte Carlo expectation-maximization algorithm is developed to estimate model parameters. Based on local linear approximations to conditional expectation of likelihood function and penalty function, a one-step sparse estimation procedure is proposed to circumvent the computationally challenge in optimizing the penalized conditional expectation of likelihood function, which is utilized to select significant covariates and trajectory functions, and identify the departure from normality of longitudinal data. The conditional expectation of likelihood function-based Bayesian information criterion is developed to select the optimal tuning parameter. Simulation studies and a real example from the clinical trial are used to illustrate the proposed methodologies.

Modeling and maximum likelihood based inference ofinterval-censored data with unknown upper limits andtime-dependent covariates.

Due to the nature of study design or other reasons, the upper limits of the interval-censored data with multiple visits are unknown. A naïve approach is to treat the last observed time as the exact event time, which may induce biased estimators of the model parameters. In this paper, we first develop a Cox model with time-dependent covariates for the event time and a proportional hazards model with frailty for the gap time. We then construct the upper limits using the latent gap times to resolve the issue of interval-censored event time data with unknown upper limits. A data-augmentation technique and a Monte Carlo EM (MCEM) algorithm are developed to facilitate computation. Theoretical properties of the computational algorithm are also investigated. Additionally, new model comparison criteria are developed to assess the fit of the gap time data as well as the fit of the event time data conditional on the gap time data. Our proposed method compares favorably with competing methods in both simulation study and real data analysis.

A general algorithm for error-in-variables regression modelling using Monte Carlo expectation maximization.

In regression modelling, measurement error models are often needed to correct for uncertainty arising from measurements of covariates/predictor variables. The literature on measurement error (or errors-in-variables) modelling is plentiful, however, general algorithms and software for maximum likelihood estimation of models with measurement error are not as readily available, in a form that they can be used by applied researchers without relatively advanced statistical expertise. In this study, we develop a novel algorithm for measurement error modelling, which could in principle take any regression model fitted by maximum likelihood, or penalised likelihood, and extend it to account for uncertainty in covariates. This is achieved by exploiting an interesting property of the Monte Carlo Expectation-Maximization (MCEM) algorithm, namely that it can be expressed as an iteratively reweighted maximisation of complete data likelihoods (formed by imputing the missing values). Thus we can take any regression model for which we have an algorithm for (penalised) likelihood estimation when covariates are error-free, nest it within our proposed iteratively reweighted MCEM algorithm, and thus account for uncertainty in covariates. The approach is demonstrated on examples involving generalized linear models, point process models, generalized additive models and capture-recapture models. Because the proposed method uses maximum (penalised) likelihood, it inherits advantageous optimality and inferential properties, as illustrated by simulation. We also study the model robustness of some violations in predictor distributional assumptions. Software is provided as the refitME package on R, whose key function behaves like a refit() function, taking a fitted regression model object and re-fitting with a pre-specified amount of measurement error.

A Hamiltonian Monte Carlo EM algorithm for generalized linear mixed models with spatial skew latent variables

A Hamiltonian Monte Carlo EM algorithm for generalized linear mixed models with spatial skew latent variables

Modeling left-truncated degradation data using random drift-diffusion Wiener processes

ABSTRACT For products whose performance characteristic (PC) gradually degrades with time, one usually observes its degradation levels repeatedly to predict its remaining useful life (RUL). Due to the limited storage space of the server and the low resolution of a measurement instrument, we seldom record the low-magnitude degradation values at the early degradation stage in applications. Such observation setting introduces left-truncated degradation data, in which the data collection starts later than the unit’s installation. This brings sampling biases and complicates the degradation data analysis. Moreover, due to the uncontrollable factors in applications, the degradation drift and the degradation diffusion may differ among various units. Motivated by an application of high-speed train bearings, we propose a Wiener process model for the left-truncated degradation data and jointly consider the drift-diffusion random effects. Closed-form formulas are available in the expectation-maximization (EM) algorithm for estimating the model parameters. We derive the RUL distribution in closed form. We also extend the proposed model to the multivariate degradation process. The parameters are estimated with the help of the Monte Carlo EM (MCEM) algorithm. An additional laser application illustrates the performance of the proposed model in RUL prediction, which may help to design a predictive maintenance strategy

Small area estimation of general finite-population parameters based on grouped data

This paper proposes a new model-based approach to small area estimation of general finite-population parameters based on grouped data or frequency data, often available from sample surveys. Grouped data contains information on frequencies of some pre-specified groups in each area, for example, the numbers of households in the income classes. Thus, grouped data provide more detailed insight into small areas than area-level aggregated data. A direct application of the widely used small area methods, such as the Fay–Herriot model for area-level data and nested error regression model for unit-level data, is not appropriate since they are not designed for grouped data. Our novel method adopts the multinomial likelihood function for the grouped data. In order to connect the group probabilities of the multinomial likelihood and the auxiliary variables within the framework of small area estimation, we introduce the unobserved unit-level quantities of interest. They follow a linear mixed model with random intercepts and dispersions after some transformation. Then the probabilities that a unit belongs to the groups can be derived and are used to construct the likelihood function for the grouped data given the random effects. The unknown model parameters (hyperparameters) are estimated by a newly developed Monte Carlo EM algorithm which uses an efficient importance sampling. The empirical best predicts (empirical Bayes estimates) of small area parameters are calculated by a simple Gibbs sampling algorithm. The numerical performance of the proposed method is illustrated based on the model-based and design-based simulations. In the application to the city-level grouped income data of Japan, we complete the patchy maps of the Gini coefficient as well as mean income across the country.

Multivariate claim count regression model with varying dispersion and dependence parameters

AbstractThe aim of this paper is to present a regression model for multivariate claim frequency data with dependence structures across the claim count responses, which may be of different sign and range, and overdispersion from the unobserved heterogeneity due to systematic effects in the data. For illustrative purposes, we consider the bivariate Poisson-lognormal regression model with varying dispersion. Maximum likelihood estimation of the model parameters is achieved through a novel Monte Carlo expectation–maximization algorithm, which is shown to have a satisfactory performance when we exemplify our approach to Local Government Property Insurance Fund data from the state of Wisconsin.

On robust probabilistic principal component analysis using multivariate t-distributions

Probabilistic principal component analysis (PPCA) is a probabilistic reformulation of principal component analysis (PCA), under the framework of a Gaussian latent variable model. To improve the robustness of PPCA, it has been proposed to change the underlying Gaussian distributions to multivariate t-distributions. Based on the representation of t-distribution as a scale mixture of Gaussian distributions, a hierarchical model is used for implementation. However, in the existing literature, the hierarchical model implemented does not yield the equivalent interpretation. In this paper, we present two sets of equivalent relationships between the high-level multivariate t-PPCA framework and the hierarchical model used for implementation. In doing so, we clarify a current misrepresentation in the literature, by specifying the correct correspondence. In addition, we discuss the performance of different multivariate t robust PPCA methods both in theory and simulation studies, and propose a new Monte Carlo expectation-maximization (MCEM) algorithm to implement one general type of such models.

Fast inference for robust nonlinear mixed-effects models

The interest for nonlinear mixed-effects models comes from application areas as pharmacokinetics, growth curves and HIV viral dynamics. However, the modeling procedure usually leads to many difficulties, as the inclusion of random effects, the estimation process and the model sensitivity to atypical or nonnormal data. The scale mixture of normal distributions include heavy-tailed models, as the Student-t, slash and contaminated normal distributions, and provide competitive alternatives to the usual models, enabling the obtention of robust estimates against outlying observations. Our proposal is to compare two estimation methods in nonlinear mixed-effects models where the random components follow a multivariate scale mixture of normal distributions. For this purpose, a Monte Carlo expectation-maximization algorithm (MCEM) and an efficient likelihood-based approximate method are developed. Results show that the approximate method is much faster and enables a fairly efficient likelihood maximization, although a slightly larger bias may be produced, especially for the fixed-effects parameters. A discussion on the robustness aspects of the proposed models are also provided. Two real nonlinear applications are discussed and a brief simulation study is presented.

Gaussian copula joint models for mixed longitudinal zero-inflated count and continuous responses

Gaussian copula joint models for mixed correlated longitudinal continuous and count responses with random effects are presented where the count responses have zero-inflated power series distribution. To account for associations between zero-inflated count and continuous responses, we use the Gaussian copula to indirectly specify their joint distributions. A full likelihood-based approach is applied to obtain IFM method to estimate marginal parameters marginally and share parameters jointly. In this method, we used the Monte Carlo EM algorithm to obtain the parameter estimates of Gaussian copula joint models. To illustrate the utility of the models, some simulation studies are performed. Finally, the proposed models are motivated by applying a medical data set. The data set is extracted from an observational study where the correlated responses are the continuous response of body mass index and the power series response of the number of joint damages.

Penalized Logistic Regression Analysis for Genetic Association Studies of Binary Phenotypes

Introduction: Increasingly, logistic regression methods for genetic association studies of binary phenotypes must be able to accommodate data sparsity, which arises from unbalanced case-control ratios and/or rare genetic variants. Sparseness leads to maximum likelihood estimators (MLEs) of log-OR parameters that are biased away from their null value of zero and tests with inflated type I errors. Different penalized likelihood methods have been developed to mitigate sparse data bias. We study penalized logistic regression using a class of log-F priors indexed by a shrinkage parameter m to shrink the biased MLE toward zero. Methods: We proposed a two-step approach to the analysis of a genetic association study: first, a set of variants that show evidence of association with the trait is used to estimate m; second, the estimated m is used for log-F-penalized logistic regression analyses of all variants using data augmentation with standard software. Our estimate of m is the maximizer of a marginal likelihood obtained by integrating the latent log-ORs out of the joint distribution of the parameters and observed data. We consider two approximate approaches to maximizing the marginal likelihood: (i) a Monte Carlo EM algorithm and (ii) a Laplace approximation to each integral, followed by derivative-free optimization of the approximation. Results: We evaluated the statistical properties of our proposed two-step method and compared its performance to other shrinkage methods by a simulation study. Our simulation studies suggest that the proposed log-F-penalized approach has lower bias and mean squared error than other methods considered. We also illustrated the approach on data from a study of genetic associations with “Super Senior” cases and middle-aged controls. Discussion/Conclusion: We have proposed a method for single rare variant analysis with binary phenotypes by logistic regression penalized by log-F priors. Our method has the advantage of being easily extended to correct for confounding due to population structure and genetic relatedness through a data augmentation approach.

Reducing subspace models for large-scale covariance regression.

We develop an envelope model for joint mean and covariance regression in the large p, small n setting. In contrast to existing envelope methods, which improve mean estimates by incorporating estimates of the covariance structure, we focus on identifying covariance heterogeneity by incorporating information about mean-level differences. We use a Monte Carlo EM algorithm to identify a low-dimensional subspace that explains differences in both means and covariances as a function of covariates, and then use MCMC to estimate the posterior uncertainty conditional on the inferred low-dimensional subspace. We demonstrate the utility of our model on a motivating application on the metabolomics of aging. We also provide R code that can be used to develop and test other generalizations of the response envelopemodel.

Censored Regression for Modelling Small Arms Trade Volumes and Its ‘Forensic’ Use for Exploring Unreported Trades

Abstract In this paper, we use a censored regression model to investigate data on the international trade of small arms and ammunition provided by the Norwegian Initiative on Small Arms Transfers. Taking a network-based view on the transfers, we do not only rely on exogenous covariates but also estimate endogenous network effects. We apply a spatial autocorrelation gravity model with multiple weight matrices. The likelihood is maximized employing the Monte Carlo expectation maximization algorithm. Our approach reveals strong and stable endogenous network effects. Furthermore, we find evidence for a substantial path dependence as well as a close connection between exports of civilian and military small arms. The model is then used in a ‘forensic’ manner to analyse latent network structures and thereby to identify countries with higher or lower tendency to export or import than reflected in the data. The approach is also validated using a simulation study.

Estimating Size and Number Density of Three-Dimensional Particles Using Truncated Cross-Sectional Data

Abstract The need for estimating three-dimensional (3D) information based on two-dimensional (2D) images has been increasing in numerous fields. It is essential in quality assessment, quality control, and process optimization. However, all the existing methods have not considered the data truncation issue, which is commonly faced in metrology. This paper proposes a new statistical approach to infer size distribution and volume number density (VND) of 3D particles based on 2D cross-sectional images with data truncation considered. In order to estimate the size distribution, a linkage is established between 3D particles and 2D observations with the existence of data truncation. Subsequently, this paper derives the likelihood function of 2D observations and an efficient Monte Carlo expectation-maximization algorithm is developed to estimate the parameters of size distribution. In addition, an explicit relationship between the 3D and 2D particle number densities is established and leveraged to estimate the VND and volume fraction. The effectiveness of the proposed method is demonstrated through both simulation study and real case studies in metal additive manufacturing and metal-matrix nanocomposites manufacturing.