Using fractional polynomials and restricted cubic splines to model non-proportional hazards or time-varying covariate effects in the Cox regression model.

The Cox proportional hazards model is used extensively in clinical and epidemiological research. A key assumption of this model is that of proportional hazards. A variable satisfies the proportional hazards assumption if the effect of that variable on the hazard function is constant over time. When the proportional hazards assumption is violated for a given variable, a common approach is to modify the model so that the regression coefficient associated with the given variable is assumed to be a linear function of time (or of log-time), rather than being constant or fixed. However, this is an unnecessarily restrictive assumption. We describe two different methods to allow a regression coefficient, and thus the hazard ratio, in a Cox model to vary as a flexible function of time. These methods use either fractional polynomials or restricted cubic splines to model the log-hazard ratio as a function of time. We illustrate the utility of these methods using data on 12 705 patients who presented to a hospital emergency department with a primary diagnosis of heart failure. We used a Cox model to assess the association between elevated cardiac troponin at presentation and the hazard of death after adjustment for an extensive set of covariates. SAS code for implementing the restricted cubic spline approach is provided, while an existing Stata function allows for the use of fractional polynomials.

Investigating treatment-effect modification by a continuous covariate in IPD meta-analysis: an approach using fractional polynomials.

BACKGROUND: In clinical trials, there is considerable interest in investigating whether a treatment effect is similar in all patients, or that one or more prognostic variables indicate a differential response to treatment. To examine this, a continuous predictor is usually categorised into groups according to one or more cutpoints. Several weaknesses of categorization are well known. To avoid the disadvantages of cutpoints and to retain full information, it is preferable to keep continuous variables continuous in the analysis. To handle this issue, the Subpopulation Treatment Effect Pattern Plot (STEPP) was proposed about two decades ago, followed by the multivariable fractional polynomial interaction (MFPI) approach. Provided individual patient data (IPD) from several studies are available, it is possible to investigate for treatment heterogeneity with meta-analysis techniques. Meta-STEPP was recently proposed and in patients with primary breast cancer an interaction of estrogen receptors with chemotherapy was investigated in eight randomized controlled trials (RCTs). METHODS: We use data from eight randomized controlled trials in breast cancer to illustrate issues from two main tasks. The first task is to derive a treatment effect function (TEF), that is, a measure of the treatment effect on the continuous scale of the covariate in the individual studies. The second is to conduct a meta-analysis of the continuous TEFs from the eight studies by applying pointwise averaging to obtain a mean function. We denote the method metaTEF. To improve reporting of available data and all steps of the analysis we introduce a three-part profile called MethProf-MA. RESULTS: Although there are considerable differences between the studies (populations with large differences in prognosis, sample size, effective sample size, length of follow up, proportion of patients with very low estrogen receptor values) our results provide clear evidence of an interaction, irrespective of the choice of the FP function and random or fixed effect models. CONCLUSIONS: In contrast to cutpoint-based analyses, metaTEF retains the full information from continuous covariates and avoids several critical issues when performing IPD meta-analyses of continuous effect modifiers in randomised trials. Early experience suggests it is a promising approach. TRIAL REGISTRATION: Not applicable.

Exposure assessment for Cox proportional hazards cure models with interval-censored survival data.

Mixture cure models have been developed as an effective tool to analyze failure time data with a cure fraction. Used in conjunction with the logistic regression model, this model allows covariate-adjusted inference of an exposure effect on the cured probability and the hazard of failure for the uncured subjects. However, the covariate-adjusted inference for the overall exposure effect is not directly provided. In this paper, we describe a Cox proportional hazards cure model to analyze interval-censored survival data in the presence of a cured fraction and then apply a post-estimation approach by using model-predicted estimates difference to assess the overall exposure effect on the restricted mean survival time scale. For baseline hazard/survival function estimation, simple parametric models as fractional polynomials or restricted cubic splines are utilized to approximate the baseline logarithm cumulative hazard function, or, alternatively, the full likelihood is specified through a piecewise linear approximation for the cumulative baseline hazard function. Simulation studies were conducted to demonstrate the unbiasedness of both estimation methods for the overall exposure effect estimates over various baseline hazard distribution shapes. The methods are applied to analyze the interval-censored relapse time data from a smoking cessation study.

Inference about time-dependent prognostic accuracy measures in the presence of competing risks.

BACKGROUND: Evaluating a candidate marker or developing a model for predicting risk of future conditions is one of the major goals in medicine. However, model development and assessment for a time-to-event outcome may be complicated in the presence of competing risks. In this manuscript, we propose a local and a global estimators of cause-specific AUC for right-censored survival times in the presence of competing risks. METHODS: The local estimator - cause-specific weighted mean rank (cWMR) - is a local average of time-specific observed cause-specific AUCs within a neighborhood of given time t. The global estimator - cause-specific fractional polynomials (cFPL) - is based on modelling the cause-specific AUC as a function of t through fractional polynomials. RESULTS: We investigated the performance of the proposed cWMR and cFPL estimators through simulation studies and real-life data analysis. The estimators perform well in small samples, have minimal bias and appropriate coverage. CONCLUSIONS: The local estimator cWMR and the global estimator cFPL will provide computationally efficient options for assessing the prognostic accuracy of markers for time-to-event outcome in the presence of competing risks in many practical settings.

Evaluating Flexible Modeling of Continuous Covariates in Inverse-Weighted Estimators.

Correct specification of the exposure model is essential for unbiased estimation in marginal structural models with inverse-probability-of-treatment weights. However, although flexible modeling is commonplace when estimating effects of continuous covariates in outcome models, its use is less frequent in estimation of inverse probability weights. Using simulations, we assess the accuracy of the treatment effect estimates and covariate balance obtained with different exposure model specifications when the true relationship between a continuous, possibly time-varying covariate Lt and the logit of the probability of exposure is nonlinear. Specifically, we compare 4 approaches to modeling the effect of Lt when estimating inverse probability weights: a linear function, the covariate-balancing propensity score, and 2 easy-to-implement flexible methods that relax the assumption of linearity: cubic regression splines and fractional polynomials. Using data from 2 empirical studies, we compare linear exposure models with flexible exposure models to estimate the effect of sustained virological response to hepatitis C virus treatment on the progression of liver fibrosis. Our simulation results demonstrate that ignoring important nonlinear relationships when fitting the exposure model may provide poorer covariate balance and induce substantial bias in the estimated exposure-outcome associations. Analysts should routinely consider flexible modeling of continuous covariates when estimating inverse-probability-of-treatment weights.

A comparative review of network meta-analysis models in longitudinal randomized controlled trial.

Network meta-analysis (NMA) technique extends the standard meta-analysis methods, allowing pairwise comparison of all treatments in a network in the absence of head-to-head comparisons. Traditional NMA models consider a single endpoint for each trial. However, in many cases, trials in the network have different durations and/or report data at multiple time points. Moreover, these time points are often not the same for all trials. In this work, we review the most relevant methods that incorporate multiple time points and allow indirect comparisons of treatment effects across different longitudinal studies. In particular, we focus on the mixed treatment comparison developed by Dakin et al,[10] on the Bayesian evidence synthesis techniques-integrated two-component prediction developed by Ding et al,[11] and on the more recent method based on fractional polynomials by Jansen et al.[12] We highlight the main features of each model and illustrate them in simulations and in a real data application. Our study shows that methods based on fractional polynomials offer a flexible modeling strategy in most applications.

Meta-analysis of non-linear exposure-outcome relationships using individual participant data: A comparison of two methods.

Non-linear exposure-outcome relationships such as between body mass index (BMI) and mortality are common. They are best explored as continuous functions using individual participant data from multiple studies. We explore two two-stage methods for meta-analysis of such relationships, where the confounder-adjusted relationship is first estimated in a non-linear regression model in each study, then combined across studies. The "metacurve" approach combines the estimated curves using multiple meta-analyses of the relative effect between a given exposure level and a reference level. The "mvmeta" approach combines the estimated model parameters in a single multivariate meta-analysis. Both methods allow the exposure-outcome relationship to differ across studies. Using theoretical arguments, we show that the methods differ most when covariate distributions differ across studies; using simulated data, we show that mvmeta gains precision but metacurve is more robust to model mis-specification. We then compare the two methods using data from the Emerging Risk Factors Collaboration on BMI, coronary heart disease events, and all-cause mortality (>80 cohorts, >18 000 events). For each outcome, we model BMI using fractional polynomials of degree 2 in each study, with adjustment for confounders. For metacurve, the powers defining the fractional polynomials may be study-specific or common across studies. For coronary heart disease, metacurve with common powers and mvmeta correctly identify a small increase in risk in the lowest levels of BMI, but metacurve with study-specific powers does not. For all-cause mortality, all methods identify a steep U-shape. The metacurve and mvmeta methods perform well in combining complex exposure-disease relationships across studies.

Correcting for measurement error in fractional polynomial models using Bayesian modelling and regression calibration, with an application to alcohol and mortality.

Exposure measurement error can result in a biased estimate of the association between an exposure and outcome. When the exposure-outcome relationship is linear on the appropriate scale (e.g. linear, logistic) and the measurement error is classical, that is the result of random noise, the result is attenuation of the effect. When the relationship is non-linear, measurement error distorts the true shape of the association. Regression calibration is a commonly used method for correcting for measurement error, in which each individual's unknown true exposure in the outcome regression model is replaced by its expectation conditional on the error-prone measure and any fully measured covariates. Regression calibration is simple to execute when the exposure is untransformed in the linear predictor of the outcome regression model, but less straightforward when non-linear transformations of the exposure are used. We describe a method for applying regression calibration in models in which a non-linear association is modelled by transforming the exposure using a fractional polynomial model. It is shown that taking a Bayesian estimation approach is advantageous. By use of Markov chain Monte Carlo algorithms, one can sample from the distribution of the true exposure for each individual. Transformations of the sampled values can then be performed directly and used to find the expectation of the transformed exposure required for regression calibration. A simulation study shows that the proposed approach performs well. We apply the method to investigate the relationship between usual alcohol intake and subsequent all-cause mortality using an error model that adjusts for the episodic nature of alcohol consumption.

Semiparametric methods for estimation of a nonlinear exposure-outcome relationship using instrumental variables with application to Mendelian randomization.

Mendelian randomization, the use of genetic variants as instrumental variables (IV), can test for and estimate the causal effect of an exposure on an outcome. Most IV methods assume that the function relating the exposure to the expected value of the outcome (the exposure-outcome relationship) is linear. However, in practice, this assumption may not hold. Indeed, often the primary question of interest is to assess the shape of this relationship. We present two novel IV methods for investigating the shape of the exposure-outcome relationship: a fractional polynomial method and a piecewise linear method. We divide the population into strata using the exposure distribution, and estimate a causal effect, referred to as a localized average causal effect (LACE), in each stratum of population. The fractional polynomial method performs metaregression on these LACE estimates. The piecewise linear method estimates a continuous piecewise linear function, the gradient of which is the LACE estimate in each stratum. Both methods were demonstrated in a simulation study to estimate the true exposure-outcome relationship well, particularly when the relationship was a fractional polynomial (for the fractional polynomial method) or was piecewise linear (for the piecewise linear method). The methods were used to investigate the shape of relationship of body mass index with systolic blood pressure and diastolic blood pressure.

Non-linear transformations of age at diagnosis, tumor size, and number of positive lymph nodes in prediction of clinical outcome in breast cancer.

BACKGROUND: Prognostic factors in breast cancer are often measured on a continuous scale, but categorized for clinical decision-making. The primary aim of this study was to evaluate if accounting for continuous non-linear effects of the three factors age at diagnosis, tumor size, and number of positive lymph nodes improves prognostication. These factors will most likely be included in the management of breast cancer patients also in the future, after an expected implementation of gene expression profiling for adjuvant treatment decision-making. METHODS: Four thousand four hundred forty seven and 1132 women with primary breast cancer constituted the derivation and validation set, respectively. Potential non-linear effects on the log hazard of distant recurrences of the three factors were evaluated during 10 years of follow-up. Cox-models of successively increasing complexity: dichotomized predictors, predictors categorized into three or four groups, and predictors transformed using fractional polynomials (FPs) or restricted cubic splines (RCS), were used. Predictive performance was evaluated by Harrell's C-index. RESULTS: Using FP-transformations, non-linear effects were detected for tumor size and number of positive lymph nodes in univariable analyses. For age, non-linear transformations did, however, not improve the model fit significantly compared to the linear identity transformation. As expected, the C-index increased with increasing model complexity for multivariable models including the three factors. By allowing more than one cut-point per factor, the C-index increased from 0.628 to 0.674. The additional gain, as measured by the C-index, when using FP- or RCS-transformations was modest (0.695 and 0.696, respectively). The corresponding C-indices for these four models in the validation set, based on the same transformations and parameter estimates from the derivation set, were 0.675, 0.700, 0.706, and 0.701. CONCLUSIONS: Categorization of each factor into three to four groups was found to improve prognostication compared to dichotomization. The additional gain by allowing continuous non-linear effects modeled by FPs or RCS was modest. However, the continuous nature of these transformations has the advantage of making it possible to form risk groups of any size.

Prediction of age at menopause in women with polycystic ovary syndrome.

OBJECTIVE: Considering the role of anti-Müllerian hormone (AMH) in female fertility and its high levels in women with polycystic ovary syndrome (PCOS), the longer reproductive span of these women is in doubt. In the present study, we aimed to improve earlier predictions using a non-linear model to substantiate the question as to whether PCOS women reach menopause later. METHODS: In total, 1162 women aged 20-50 years, comprising 378 PCOS cases and 784 eumenorrheic non-hirsute women, met the eligibility criteria. A scatterplot matrix was drawn to detect the association between age and AMH; this association was explored using a fractional polynomial regression model. Model assumptions were checked by examining the distribution of the residuals and plotting the standardized residuals against the functional form of AMH. RESULTS: The serum concentration of AMH among PCOS participants was significantly higher than in the controls (5.4 ng/ml (IQR 2.8-9.1 ng/ml) vs. 1.4 ng/ml (IQR 0.6-2.7 ng/ml), p < 0.001). The estimated mean age at menopause was 51.4 (95% CI 45-59) years and 49.7 (95% CI 45-55) years in PCOS cases and healthy controls, respectively. CONCLUSIONS: These findings provide the insight that, as reflected through significantly higher average levels of AMH in PCOS women, their predicted reproductive lifespan could be 2 years longer than their normo-ovulatory counterparts.

Modeling Variables With a Spike at Zero: Examples and Practical Recommendations.

In most epidemiologic studies and in clinical research generally, there are variables with a spike at zero, namely variables for which a proportion of individuals have zero exposure (e.g., never smokers) and among those exposed the variable has a continuous distribution. Different options exist for modeling such variables, such as categorization where the nonexposed form the reference group, or ignoring the spike by including the variable in the regression model with or without some transformation or modeling procedures. It has been shown that such situations can be analyzed by adding a binary indicator (exposed/nonexposed) to the regression model, and a method based on fractional polynomials with which to estimate a suitable functional form for the positive portion of the spike-at-zero variable distribution has been developed. In this paper, we compare different approaches using data from 3 case-control studies carried out in Germany: the Mammary Carcinoma Risk Factor Investigation (MARIE), a breast cancer study conducted in 2002-2005 (Flesch-Janys et al., Int J Cancer. 2008;123(4):933-941); the Rhein-Neckar Larynx Study, a study of laryngeal cancer conducted in 1998-2000 (Dietz et al., Int J Cancer. 2004;108(6):907-911); and a lung cancer study conducted in 1988-1993 (Jöckel et al., Int J Epidemiol. 1998;27(4):549-560). Strengths and limitations of different procedures are demonstrated, and some recommendations for practical use are given.

Incorporating nonlinearity into mediation analyses.

BACKGROUND: Mediation is an important issue considered in the behavioral, medical, and social sciences. It addresses situations where the effect of a predictor variable X on an outcome variable Y is explained to some extent by an intervening, mediator variable M. Methods for addressing mediation have been available for some time. While these methods continue to undergo refinement, the relationships underlying mediation are commonly treated as linear in the outcome Y, the predictor X, and the mediator M. These relationships, however, can be nonlinear. Methods are needed for assessing when mediation relationships can be treated as linear and for estimating them when they are nonlinear. METHODS: Existing adaptive regression methods based on fractional polynomials are extended here to address nonlinearity in mediation relationships, but assuming those relationships are monotonic as would be consistent with theories about directionality of such relationships. RESULTS: Example monotonic mediation analyses are provided assessing linear and monotonic mediation of the effect of family functioning (X) on a child's adaptation (Y) to a chronic condition by the difficulty (M) for the family in managing the child's condition. Example moderated monotonic mediation and simulation analyses are also presented. CONCLUSIONS: Adaptive methods provide an effective way to incorporate possibly nonlinear monotonicity into mediation relationships.

Adaptive Modeling: An Approach for Incorporating Nonlinearity in Regression Analyses.

Although regression relationships commonly are treated as linear, this often is not the case. An adaptive approach is described for identifying nonlinear relationships based on power transforms of predictor (or independent) variables and for assessing whether or not relationships are distinctly nonlinear. It is also possible to model adaptively both means and variances of continuous outcome (or dependent) variables and to adaptively power transform positive-valued continuous outcomes, along with their predictors. Example analyses are provided of data from parents in a nursing study on emotional-health-related quality of life for childhood brain tumor survivors as a function of the effort to manage the survivors' condition. These analyses demonstrate that relationships, including moderation relationships, can be distinctly nonlinear, that conclusions about means can be affected by accounting for non-constant variances, and that outcome transformation along with predictor transformation can provide distinct improvements and can resolve skewness problems.© 2017 Wiley Periodicals, Inc.

Model selection for univariable fractional polynomials.

Since Royston and Altman's 1994 publication (Journal of the Royal Statistical Society, Series C 43: 429-467), fractional polynomials have steadily gained popularity as a tool for flexible parametric modeling of regression relationships. In this article, I present fp_select, a postestimation tool for fp that allows the user to select a parsimonious fractional polynomial model according to a closed test procedure called the fractional polynomial selection procedure or function selection procedure. I also give a brief introduction to fractional polynomial models and provide examples of using fp and fp_select to select such models with real data.

A comparison of methods to adjust for continuous covariates in the analysis of randomised trials.

BACKGROUND: Although covariate adjustment in the analysis of randomised trials can be beneficial, adjustment for continuous covariates is complicated by the fact that the association between covariate and outcome must be specified. Misspecification of this association can lead to reduced power, and potentially incorrect conclusions regarding treatment efficacy. METHODS: We compared several methods of adjustment to determine which is best when the association between covariate and outcome is unknown. We assessed (a) dichotomisation or categorisation; (b) assuming a linear association with outcome; (c) using fractional polynomials with one (FP1) or two (FP2) polynomial terms; and (d) using restricted cubic splines with 3 or 5 knots. We evaluated each method using simulation and through a re-analysis of trial datasets. RESULTS: Methods which kept covariates as continuous typically had higher power than methods which used categorisation. Dichotomisation, categorisation, and assuming a linear association all led to large reductions in power when the true association was non-linear. FP2 models and restricted cubic splines with 3 or 5 knots performed best overall. CONCLUSIONS: For the analysis of randomised trials we recommend (1) adjusting for continuous covariates even if their association with outcome is unknown; (2) keeping covariates as continuous; and (3) using fractional polynomials with two polynomial terms or restricted cubic splines with 3 to 5 knots when a linear association is in doubt.

mfpa: Extension of mfp using the ACD covariate transformation for enhanced parametric multivariable modeling.

In a recent article, Royston (2015, Stata Journal 15: 275-291) introduced the approximate cumulative distribution (acd) transformation of a continuous covariate x as a route toward modeling a sigmoid relationship between x and an outcome variable. In this article, we extend the approach to multivariable modeling by modifying the standard Stata program mfp. The result is a new program, mfpa, that has all the features of mfp plus the ability to fit a new model for user-selected covariates that we call fp1(p1, p2). The fp1(p1, p2) model comprises the best-fitting combination of a dimension-one fractional polynomial (fp1) function of x and an fp1 function of acd (x). We describe a new model-selection algorithm called function-selection procedure with acd transformation, which uses significance testing to attempt to simplify an fp1(p1, p2) model to a submodel, an fp1 or linear model in x or in acd (x). The function-selection procedure with acd transformation is related in concept to the fsp (fp function-selection procedure), which is an integral part of mfp and which is used to simplify a dimension-two (fp2) function. We describe the mfpa command and give univariable and multivariable examples with real data to demonstrate its use.

Modeling continuous covariates with a "spike" at zero: Bivariate approaches.

In epidemiology and clinical research, predictors often take value zero for a large amount of observations while the distribution of the remaining observations is continuous. These predictors are called variables with a spike at zero. Examples include smoking or alcohol consumption. Recently, an extension of the fractional polynomial (FP) procedure, a technique for modeling nonlinear relationships, was proposed to deal with such situations. To indicate whether or not a value is zero, a binary variable is added to the model. In a two stage procedure, called FP-spike, the necessity of the binary variable and/or the continuous FP function for the positive part are assessed for a suitable fit. In univariate analyses, the FP-spike procedure usually leads to functional relationships that are easy to interpret. This paper introduces four approaches for dealing with two variables with a spike at zero (SAZ). The methods depend on the bivariate distribution of zero and nonzero values. Bi-Sep is the simplest of the four bivariate approaches. It uses the univariate FP-spike procedure separately for the two SAZ variables. In Bi-D3, Bi-D1, and Bi-Sub, proportions of zeros in both variables are considered simultaneously in the binary indicators. Therefore, these strategies can account for correlated variables. The methods can be used for arbitrary distributions of the covariates. For illustration and comparison of results, data from a case-control study on laryngeal cancer, with smoking and alcohol intake as two SAZ variables, is considered. In addition, a possible extension to three or more SAZ variables is outlined. A combination of log-linear models for the analysis of the correlation in combination with the bivariate approaches is proposed.

Combining fractional polynomial model building with multiple imputation.

Multivariable fractional polynomial (MFP) models are commonly used in medical research. The datasets in which MFP models are applied often contain covariates with missing values. To handle the missing values, we describe methods for combining multiple imputation with MFP modelling, considering in turn three issues: first, how to impute so that the imputation model does not favour certain fractional polynomial (FP) models over others; second, how to estimate the FP exponents in multiply imputed data; and third, how to choose between models of differing complexity. Two imputation methods are outlined for different settings. For model selection, methods based on Wald-type statistics and weighted likelihood-ratio tests are proposed and evaluated in simulation studies. The Wald-based method is very slightly better at estimating FP exponents. Type I error rates are very similar for both methods, although slightly less well controlled than analysis of complete records; however, there is potential for substantial gains in power over the analysis of complete records. We illustrate the two methods in a dataset from five trauma registries for which a prognostic model has previously been published, contrasting the selected models with that obtained by analysing the complete records only.

Network meta-analysis of longitudinal data using fractional polynomials.

Network meta-analysis of randomized controlled trials (RCTs) are often based on one treatment effect measure per study. However, many studies report data at multiple time points. Furthermore, not all studies measure the outcomes at the same time points. As an alternative to a network meta-analysis based on a synthesis of the results at one time point, a network meta-analysis method is presented that allows for the simultaneous analysis of outcomes at multiple time points. The development of outcomes over time of interventions compared in an RCT is modeled with fractional polynomials, and the differences between the parameters of these polynomials within a trial are synthesized across studies with a Bayesian network meta-analysis. The proposed models are illustrated with an analysis of RCTs evaluating interventions for osteoarthritis of the knee. Fixed and random effects second order fractional polynomials were applied to the case study. Network meta-analysis with models that represent the treatment effects in terms of several parameters using fractional polynomials can be considered a useful addition to models for network meta-analysis of repeated measures previously proposed. When RCTs report treatment effects at multiple follow-up times, these models can be used to synthesize the results even if reporting times differ across the studies.