Volume 18, Issue 2

Effect modification occurs when the effect of a treatment on an outcome differsaccording to the level of some pre-treatment variable (the effect modifier). Assessing an effect modifier is not a straight-forward task even for a subject matter expert. In this paper, we propose a two-stageprocedure to automatically selecteffect modifying variables in a Marginal Structural Model (MSM) with a single time point exposure based on the two nuisance quantities (the conditionaloutcome expectation and propensity score). We highlight the performance of our proposal in a simulation study. Finally, to illustrate tractability of our proposed methods, we apply them to analyze a set of pregnancy data. We estimate the conditional expected difference in the counterfactual birth weight if all women were exposed to inhaled corticosteroids during pregnancy versus the counterfactual birthweight if all women were not, using data from asthma medications during pregnancy.

Estimating causal effects from randomized experiments is central to clinical research. Reducing the statistical uncertainty in these analyses is an important objective for statisticians. Registries, prior trials, and health records constitute a growing compendium of historical data on patients under standard-of-care that may be exploitable to this end. However, most methods for historical borrowing achieve reductions in variance by sacrificing strict type-I error rate control. Here, we propose a use of historical data that exploits linear covariate adjustment to improve the efficiency of trial analyses without incurring bias. Specifically, we train a prognostic model on the historical data, then estimate the treatment effect using a linear regression while adjusting for the trial subjects’ predicted outcomes (their prognostic scores ). We prove that, under certain conditions, this prognostic covariate adjustment procedure attains the minimum variance possible among a large class of estimators. When those conditions are not met, prognostic covariate adjustment is still more efficient than raw covariate adjustment and the gain in efficiency is proportional to a measure of the predictive accuracy of the prognostic model above and beyond the linear relationship with the raw covariates. We demonstrate the approach using simulations and a reanalysis of an Alzheimer’s disease clinical trial and observe meaningful reductions in mean-squared error and the estimated variance. Lastly, we provide a simplified formula for asymptotic variance that enables power calculations that account for these gains. Sample size reductions between 10% and 30% are attainable when using prognostic models that explain a clinically realistic percentage of the outcome variance.

In finding effects of a binary treatment, practitioners use mostly either propensity score matching (PSM) or inverse probability weighting (IPW). However, many new treatment effect estimators are available now using propensity score and “prognostic score”, and some of these estimators are much better than PSM and IPW in several aspects. In this paper, we review those recent treatment effect estimators to show how they are related to one another, and why they are better than PSM and IPW. We compare 26 estimators in total through extensive simulation and empirical studies. Based on these, we recommend recent treatment effect estimators using “overlap weight”, and “targeted MLE” using statistical/machine learning, as well as a simple regression imputation/adjustment estimator using linear prognostic score models.

In the context of finite mixture models one considers the problem of classifying as many observations as possible in the classes of interest while controlling the classification error rate in these same classes. Similar to what is done in the framework of statistical test theory, different type I and type II-like classification error rates can be defined, along with their associated optimal rules, where optimality is defined as minimizing type II error rate while controlling type I error rate at some nominal level. It is first shown that finding an optimal classification rule boils down to searching an optimal region in the observation space where to apply the classical Maximum A Posteriori (MAP) rule. Depending on the misclassification rate to be controlled, the shape of the optimal region is provided, along with a heuristic to compute the optimal classification rule in practice. In particular, a multiclass FDR-like optimal rule is defined and compared to the thresholded MAP rules that is used in most applications. It is shown on both simulated and real datasets that the FDR-like optimal rule may be significantly less conservative than the thresholded MAP rule.

The use of cumulative incidence functions for characterizing the risk of one type of event in the presence of others has become increasingly popular over the past two decades. The problems of modeling, estimation and inference have been treated using parametric, nonparametric and semi-parametric methods. Efforts to develop suitable extensions of machine learning methods, such as regression trees and ensemble methods, have begun comparatively recently. In this paper, we propose a novel approach to estimating cumulative incidence curves in a competing risks setting using regression trees and associated ensemble estimators. The proposed methods use augmented estimators of the Brier score risk as the primary basis for building and pruning trees, and lead to methods that are easily implemented using existing R packages. Data from the Radiation Therapy Oncology Group (trial 9410) is used to illustrate these new methods.

Many causal models of interest in epidemiology involve longitudinal exposures, confounders and mediators. However, repeated measurements are not always available or used in practice, leading analysts to overlook the time-varying nature of exposures and work under over-simplified causal models. Our objective is to assess whether – and how – causal effects identified under such misspecified causal models relates to true causal effects of interest. We derive sufficient conditions ensuring that the quantities estimated in practice under over-simplified causal models can be expressed as weighted averages of longitudinal causal effects of interest. Unsurprisingly, these sufficient conditions are very restrictive, and our results state that the quantities estimated in practice should be interpreted with caution in general, as they usually do not relate to any longitudinal causal effect of interest. Our simulations further illustrate that the bias between the quantities estimated in practice and the weighted averages of longitudinal causal effects of interest can be substantial. Overall, our results confirm the need for repeated measurements to conduct proper analyses and/or the development of sensitivity analyses when they are not available.

Group therapy is a common treatment modality for behavioral health conditions. Patients often enter and exit groups on an ongoing basis, leading to dynamic therapy groups. Examining the effect of high versus low session attendance on patient outcomes is a research question of interest. However, there are several challenges to identifying causal effects in this setting, including the lack of randomization, interference among patients, and the interrelatedness of patient participation. Dynamic therapy groups motivate a unique causal inference scenario, as the treatment statuses are completely defined by the patient attendance record for the therapy session, which is also the structure inducing interference. We adopt the Rubin causal model framework to define the causal effect of high versus low session attendance of group therapy at both the individual patient and peer levels. We propose a strategy to identify individual, peer, and total effects of high attendance versus low attendance on patient outcomes by the prognostic score stratification. We examine performance of our approach via simulation and apply it to data from a group cognitive behavioral therapy trial for treating depression among patients in a substance use disorders treatment setting.

Since the discovery of the human immunodeficiency virus (HIV) 35 years ago, the epidemic is still ongoing in France. To monitor the dynamics of HIV transmission and assess the impact of prevention campaigns, the main indicator is the incidence. One method to estimate the HIV incidence is based on biomarker values at diagnosis and their dynamics over time. Estimating the HIV incidence from biomarkers first requires modeling their dynamics since infection using external longitudinal data. The objective of the work presented here is to estimate the joint dynamics of two biomarkers from the PRIMO cohort. We thus jointly modeled the dynamics of two biomarkers (TM and V3) using a multi-response nonlinear mixed-effect model. The parameters were estimated using Bayesian Hamiltonian Monte Carlo inference. This procedure was first applied to the real data of the PRIMO cohort. In a simulation study, we then evaluated the performance of the Bayesian procedure for estimating the parameters of multi-response nonlinear mixed-effect models.

The accelerated failure time mixture cure (AFTMC) model is widely used for survival data when a portion of patients can be cured. In this paper, a Bayesian semiparametric method is proposed to obtain the estimation of parameters and density distribution for both the cure probability and the survival distribution of the uncured patients in the AFTMC model. Specifically, the baseline error distribution of the uncured patients is nonparametrically modeled by a mixture of Dirichlet process. Based on the stick-breaking formulation of the Dirichlet process, the techniques of retrospective and slice sampling, an efficient and easy-to-implement Gibbs sampler is developed for the posterior calculation. The proposed approach can be easily implemented in commonly used statistical softwares, and its performance is comparable to fully parametric method via comprehensive simulation studies. Besides, the proposed approach is adopted to the analysis of a colorectal cancer clinical trial data.

The timings of visits in observational longitudinal data may depend on the study outcome, and this can result in bias if ignored. Assessing the extent of visit irregularity is important because it can help determine whether visits can be treated as repeated measures or as irregular data. We propose plotting the mean proportions of individuals with 0 visits per bin against the mean proportions of individuals with >1 visit per bin as bin width is varied and using the area under the curve (AUC) to assess the extent of irregularity. The AUC is a single score which can be used to quantify the extent of irregularity and assess how closely visits resemble repeated measures. Simulation results confirm that the AUC increases with increasing irregularity while being invariant to sample size and the number of scheduled measurement occasions. A demonstration of the AUC was performed on the TARGet Kids! study which enrolls healthy children aged 0–5 years with the aim of investigating the relationship between early life exposures and later health problems. The quality of statistical analyses can be improved by using the AUC as a guide to select the appropriate analytic outcome approach and minimize the potential for biased results.

Interrupted time series (ITS) design is commonly used to evaluate the impact of interventions in healthcare settings. Segmented regression (SR) is the most commonly used statistical method and has been shown to be useful in practical applications involving ITS designs. Nevertheless, SR is prone to aggregation bias, which leads to imprecision and loss of power to detect clinically meaningful differences. The objective of this article is to present a weighted SR method, where variability across patients within the healthcare facility and across time points is incorporated through weights. We present the methodological framework, provide optimal weights associated with data at each time point and discuss relevant statistical inference. We conduct extensive simulations to evaluate performance of our method and provide comparative analysis with the traditional SR using established performance criteria such as bias, mean square error and statistical power. Illustrations using real data is also provided. In most simulation scenarios considered, the weighted SR method produced estimators that are uniformly more precise and relatively less biased compared to the traditional SR. The weighted approach also associated with higher statistical power in the scenarios considered. The performance difference is much larger for data with high variability across patients within healthcare facilities. The weighted method proposed here allows us to account for the heterogeneity in the patient population, leading to increased accuracy and power across all scenarios. We recommend researchers to carefully design their studies and determine their sample size by incorporating heterogeneity in the patient population.

Hazard ratios (HR) associated with the well-known proportional hazard Cox regression models are routinely used for measuring the impact of one factor of interest on a time-to-event outcome. However, if the underlying real model does not fit with the theoretical requirements, the interpretation of those HRs is not clear. We propose a new index, gHR, which generalizes the HR beyond the underlying survival model. We consider the case in which the study factor is a binary variable and we are interested in both the unadjusted and adjusted effect of this factor on a time-to-event variable, potentially, observed in a right-censored scenario. We propose non-parametric estimations for unadjusted gHR and semi-parametric regression-induced techniques for the adjusted case. The behavior of those estimators is studied in both large and finite sample situations. Monte Carlo simulations reveal that both estimators provide good approximations of their respective inferential targets. Data from the Health and Lifestyle Study are used for studying the relationship of the tobacco use and the age of death and illustrate the practical application of the proposed technique. gHR is a promising index which can help facilitate better understanding of the association of one study factor on a time-dependent outcome.

This paper considers the problem of semi-parametric proportional hazards model fitting where observed survival times contain event times and also interval, left and right censoring times. Although this is not a new topic, many existing methods suffer from poor computational performance. In this paper, we adopt a more versatile penalized likelihood method to estimate the baseline hazard and the regression coefficients simultaneously. The baseline hazard is approximated using basis functions such as M-splines. A penalty is introduced to regularize the baseline hazard estimate and also to ease dependence of the estimates on the knots of the basis functions. We propose a Newton–MI (multiplicative iterative) algorithm to fit this model. We also present novel asymptotic properties of our estimates, allowing for the possibility that some parameters of the approximate baseline hazard may lie on the parameter space boundary. Comparisons of our method against other similar approaches are made through an intensive simulation study. Results demonstrate that our method is very stable and encounters virtually no numerical issues. A real data application involving melanoma recurrence is presented and an R package ‘survivalMPL’ implementing the method is available on R CRAN.

With known cause of death (CoD), competing risk survival methods are applicable in estimating disease-specific survival. Relative survival analysis may be used to estimate disease-specific survival when cause of death is either unknown or subject to misspecification and not reliable for practical usage. This method is popular for population-based cancer survival studies using registry data and does not require CoD information. The standard estimator is the ratio of all-cause survival in the cancer cohort group to the known expected survival from a general reference population. Disease-specific death competes with other causes of mortality, potentially creating dependence among the CoD. The standard ratio estimate is only valid when death from disease and death from other causes are independent. To relax the independence assumption, we formulate dependence using a copula-based model. Likelihood-based parametric method is used to fit the distribution of disease-specific death without CoD information, where the copula is assumed known and the distribution of other cause of mortality is derived from the reference population. We propose a sensitivity analysis, where the analysis is conducted across a range of assumed dependence structures. We demonstrate the utility of our method through simulation studies and an application to French breast cancer data.

The choice of the number m of response categories is a crucial issue in categorization of a continuous response. The paper exploits the Proportional Odds Models’ property which allows to generate ordinal responses with a different number of categories from the same underlying variable. It investigates the asymptotic efficiency of the estimators of the regression coefficients and the accuracy of the derived inferential procedures when m varies. The analysis is based on models with closed-form information matrices so that the asymptotic efficiency can be analytically evaluated without need of simulations. The paper proves that a finer categorization augments the information content of the data and consequently shows that the asymptotic efficiency and the power of the tests on the regression coefficients increase with m . The impact of the loss of information produced by merging categories on the efficiency of the estimators is also considered, highlighting its risks especially when performed in its extreme form of dichotomization. Furthermore, the appropriate value of m for various sample sizes is explored, pointing out that a large number of categories can offset the limited amount of information of a small sample by a better quality of the data. Finally, two case studies on the quality of life of chemotherapy patients and on the perception of pain, based on discretized continuous scales, illustrate the main findings of the paper.

Dichotomization is often used on clinical and diagnostic settings to simplify interpretation. For example, a person with systolic and diastolic blood pressure above 140 over 90 may be prescribed medication. Blood pressure as well as other factors such as age and cholesterol and their interactions may lead to increased risk of certain diseases. When using a dichotomized variable to determine a diagnosis, if the interactions with other variables are not considered, then an incorrect threshold for the continuous variable may be selected. In this paper, we compare single dichotomization with joint dichotomization; the process of simultaneously optimizing cutpoints for multiple variables. A simulation study shows that simultaneous dichotomization of continuous variables is more accurate in recovering both ‘true’ thresholds given they exist.

We present in this paper a global methodology for the spike detection in a biological context of fluorescence recording of GnRH-neurons calcium activity. For this purpose we first propose a simple stochastic model that could mimic experimental time series by considering an autoregressive AR(1) process with a linear trend and specific innovations involving spiking times. Estimators of parameters with asymptotic normality are established and used to set up a statistical test on estimated innovations in order to detect spikes. We compare several procedures and illustrate on biological data the performance of our procedure.