Search

This paper compares six different parameter estimation methods for shared frailty models via a series of simulation studies. A shared frailty model is a survival model that incorporates a random effect term, where the frailties are common or shared among individuals within specific groups. Several parameter estimation methods are available for fitting shared frailty models, such as penalized partial likelihood (PPL), expectation-maximization (EM), pseudo full likelihood (PFL), hierarchical likelihood (HL), maximum marginal likelihood (MML), and maximization penalized likelihood (MPL) algorithms. These estimation methods are implemented in various R packages, providing researchers with various options for analyzing clustered survival data using shared frailty models. However, there is a limited amount of research comparing the performance of these parameter estimation methods for fitting shared frailty models. Consequently, it can be challenging for users to determine the most appropriate method for analyzing clustered survival data. To address this gap, this paper aims to conduct a series of simulation studies to compare the performance of different parameter estimation methods implemented in R packages. We will evaluate several key aspects, including parameter estimation, bias and variance of the parameter estimates, rate of convergence, and computational time required by each package. Through this systematic evaluation, our goal is to provide a comprehensive understanding of the advantages and limitations associated with each estimation method.

Residual diagnostic methods play a critical role in assessing model assumptions and detecting outliers in statistical modelling. In the context of survival models with censored observations, Li et al. (2021) introduced the Z-residual, which follows an approximately normal distribution under the true model. This property makes it possible to use Z-residuals for diagnosing survival models in a way similar to how Pearson residuals are used in normal regression. However, computing residuals based on the full dataset can result in a conservative bias that reduces the power of detecting model mis-specification, as the same dataset is used for both model fitting and validation. Although cross-validation is a potential solution to this problem, it has not been commonly used in residual diagnostics due to computational challenges. In this paper, we propose a cross-validation approach for computing Z-residuals in the context of shared frailty models. Specifically, we develop a general function that calculates cross-validatory Z-residuals using the output from the \texttt{coxph} function in the \texttt{survival} package in R.Our simulation studies demonstrate that, for goodness-of-fit tests and outlier detection, cross-validatory Z-residuals are significantly more powerful and more discriminative than Z-residuals without cross-validation. We also compare the performance of Z-residuals with and without cross-validation in identifying outliers in a real application that models the recurrence time of kidney infection patients. Our findings suggest that cross-validatory Z-residuals can identify outliers that are missed by Z-residuals without cross-validation., Comment: 32 pages, 14 figures

In survival analysis, the hazard function often depends on a set of covariates. Martingale and deviance residual are most widely used for examining the validity of the function form of covariates by checking whether there is a discernible trend in their scatterplot against continuous covariates. However, visual inspection of martingale and deviance residuals is often subjective. In addition, these residuals lack a reference distribution due to censoring. It is therefore challenging to derive numerical statistical tests based on martingale or deviance residuals. In this paper, we extend the idea of randomized survival probability (Li et al. 2021) and develop a residual diagnostic tool that can provide both graphical and numerical tests for checking the covariate functional form in semi-parametric shared frailty models. We develop a general function that calculates Z-residuals for semi-parametric shared frailty models based on the output from the \texttt{coxph} function in the \texttt{survival} package in R. Our extensive simulation studies indicate that the derived numerical test based on Z-residuals has great power for checking the functional form of covariates. In a real data application on modelling the survival time of acute myeloid leukemia patients, the Z-residual diagnosis results show that a model with log-transformation is inappropriate for modelling the survival time, which could not be detected by other diagnostic methods., Comment: 21 pages, 7 figures

Residuals in normal regression are used to assess a model's goodness-of-fit (GOF) and discover directions for improving the model. However, there is a lack of residuals with a characterized reference distribution for censored regression. In this paper, we propose to diagnose censored regression with normalized randomized survival probabilities (RSP). The key idea of RSP is to replace the survival probability of a censored failure time with a uniform random number between 0 and the survival probability of the censored time. We prove that RSPs always have the uniform distribution on $(0,1)$ under the true model with the true generating parameters. Therefore, we can transform RSPs into normally-distributed residuals with the normal quantile function. We call such residuals by normalized RSP (NRSP residuals). We conduct simulation studies to investigate the sizes and powers of statistical tests based on NRSP residuals in detecting the incorrect choice of distribution family and non-linear effect in covariates. Our simulation studies show that, although the GOF tests with NRSP residuals are not as powerful as a traditional GOF test method, a non-linear test based on NRSP residuals has significantly higher power in detecting non-linearity. We also compared these model diagnostics methods with a breast-cancer recurrent-free time dataset. The results show that the NRSP residual diagnostics successfully captures a subtle non-linear relationship in the dataset, which is not detected by the graphical diagnostics with CS residuals and existing GOF tests., Comment: Accepted version. 12 pages

Examining residuals such as Pearson and deviance residuals, is a standard tool for assessing normal regression. However, for discrete response, these residuals cluster on lines corresponding to distinct response values. Their distributions are far from normality; graphical and quantitative inspection of these residuals provides little information for model diagnosis. Marshall and Spiegelhalter (2003) defined a cross-validatory predictive p-value for identifying outliers. Predictive p-values are uniformly distributed for continuous response but not for discrete response. We propose to use randomized predictive p-values (RPP) for diagnosing models with discrete responses. RPPs can be transformed to "residuals" with normal distribution, called NRPPs by us. NRPPs can be used to diagnose all regression models with scalar response using the same way for diagnosing normal regression. The NRPPs are nearly the same as the randomized quantile residuals (RQR), which are previously proposed by Dunn and Smyth (1996) but remain little known by statisticians. This paper provides an exposition of RQR using the RPP perspective. The contributions of this exposition include: (1) we give a rigorous proof of uniformity of RPP and illustrative examples to explain the uniformity under the true model; (2) we conduct extensive simulation studies to demonstrate the normality of NRPPs under the true model; (3) our simulation studies also show that the NRPP method is a versatile diagnostic tool for detecting many kinds of model inadequacies due to lack of complexity. The effectiveness of NRPP is further demonstrated with a health utilization dataset., Comment: 26 pages; we've updated some figures for better visualization, and fixed a few errors in text; R code for producing the results of this paper is available upon request

An important statistical task in disease mapping problems is to identify divergent regions with unusually high or low risk of disease. Leave-one-out cross-validatory (LOOCV) model assessment is the gold standard for estimating predictive p-values that can flag such divergent regions. However, actual LOOCV is time-consuming because one needs to rerun a Markov chain Monte Carlo analysis for each posterior distribution in which an observation is held out as a test case. This paper introduces a new method, called integrated importance sampling (iIS), for estimating LOOCV predictive p-values with only Markov chain samples drawn from the posterior based on a full data set. The key step in iIS is that we integrate away the latent variables associated the test observation with respect to their conditional distribution \textit{without} reference to the actual observation. By following the general theory for importance sampling, the formula used by iIS can be proved to be equivalent to the LOOCV predictive p-value. We compare iIS and other three existing methods in the literature with two disease mapping datasets. Our empirical results show that the predictive p-values estimated with iIS are almost identical to the predictive p-values estimated with actual LOOCV, and outperform those given by the existing three methods, namely, the posterior predictive checking, the ordinary importance sampling, and the ghosting method by Marshall and Spiegelhalter (2003)., Comment: 18 pages. Accepted version

A natural method for approximating out-of-sample predictive evaluation is leave-one-out cross-validation (LOOCV) --- we alternately hold out each case from a full data set and then train a Bayesian model using Markov chain Monte Carlo (MCMC) without the held-out; at last we evaluate the posterior predictive distribution of all cases with their actual observations. However, actual LOOCV is time-consuming. This paper introduces two methods, namely iIS and iWAIC, for approximating LOOCV with only Markov chain samples simulated from a posterior based on a \textit{full} data set. iIS and iWAIC aim at improving the approximations given by importance sampling (IS) and WAIC in Bayesian models with possibly correlated latent variables. In iIS and iWAIC, we first integrate the predictive density over the distribution of the latent variables associated with the held-out without reference to its observation, then apply IS and WAIC approximations to the integrated predictive density. We compare iIS and iWAIC with other approximation methods in three real data examples that respectively use mixture models, models with correlated spatial effects, and a random effect logistic model. Our empirical results show that iIS and iWAIC give substantially better approximates than non-integrated IS and WAIC and other methods., Comment: 38 pages