The vast majority of studies have employed a non-parametric approach, by treating BMI as a categorical variable. The mortality risk of individuals in different BMI groups is computed relative to a reference BMI category. World Health Organization (WHO) BMI classifications [16] are typically used. For example, numerous studies have computed the excess mortality due to being overweight (BMI = 25-30) and obese (BMI ≥ 30) [17-22]. Categorizing continuous variables has been a popular approach, particularly because of a priori knowledge that the relationship between the two measures is nonlinear. However, Royston, et al. (2006) [23] pointed out a number of disadvantages with this approach. The most significant drawback is the loss of information and power through what is equivalent to rounding. For example, studies of excess mortality using WHO BMI classification implicitly assume all individuals that are in the normal category (BMI = 18.5-25) exhibit the same mortality risk. However, the normal category is heterogenous and includes individuals who are healthy and those who are chronically ill. Categorization is particularly problematic if groups are large. In addition, the number of cutoff points and where to place cutoff points is arbitrary. Finally, results are not necessarily robust across the choice of reference category.

A number of studies have used other approaches that maintain BMI as a continuous variable. Estimation of the BMI-mortality curve using a continuous measure is problematic because the relationship is nonlinear and the distribution of BMI is right skewed. Schauer et al. (2010) [24] included linear and squared terms to account for nonlinearities, but truncated their sample to respondents with a BMI of 25 or greater to address the skewness in BMI. Durazo et al. (1998) [6] transformed BMI into a normally distributed variable using Tukey's "ladder of powers" method. Gronniger et al. (2006) [25] treated BMI non-parametrically without categorization.

We investigated the relationship between mortality and obesity through BMI using the logistic regression model, stratified by gender and adjusting for age and smoking history. The logistic model was chosen over the Cox proportional hazards model because the proportional hazards condition did not hold for the BMI fractional polynomial (FP) terms. We used 5-year all-cause mortality as the dependent outcome because of the very low incidence of death annually. To check for robustness, we also estimated all models using 3-year mortality as the outcome, which produced similar results. Analyses were stratified by gender because the biological process by which men and women gain and maintain weight is different [29]. We adjusted for smoking status because it confounded the BMI-mortality relationship, which if ignored may result in overestimation of the BMI associated with minimum mortality [30]. Sample adult weights from the NHIS, which denoted the inverse probability of inclusion into the sample were used within the logistic regression model to correct for potential biases resulting from the NHIS sampling design. Because data were pooled, sampling weights were divided by the number of years to generate a sample that is representative of the U.S. population on average from 1997-2000.

We maintained BMI as a continuous variable in our analysis. To account for the nonlinear and asymmetric relationship between BMI and mortality, we first applied the fractional polynomials [31,32] method. To allow for flexibility in fitting a curve with a single turning point, we considered second degree polynomial transformations for BMI. We used the closed test procedure [33] which first determined the best fitting second degree polynomial by choosing power transformations from the set {-2, -1, -0.5, 0, 0.5, 1, 2, 3}, where 0 denotes the log transformation. The best fitting second-degree FP was then compared against the null model using a deviance difference test with four degrees of freedom to determine whether BMI should be included in the model. If the first test was statistically significant, a second deviance difference test with three degrees of freedom was applied to compare the best fitting second degree FP against the linear model. If the second test was significant, a final deviance difference test with two degrees of freedom compared the best fitting second degree FP with the best fitting first degree FP. If the final test was significant, the second degree FP was included, otherwise the first degree FP was chosen. To prevent collinearity and model overfit, the best fitting first degree polynomial was chosen for age. The selection of powers for BMI and age was computed simultaneously using the multivariable fractional polynomials (MFP) method [26,27], which combined backward elimination to select the best fitting model. The regression model we estimated was

where π_i was the 5-year death probability for individual i, p1 and p2 were the fractional powers for BMI, and q1 was the fractional power for age. The MFP method also scaled and centered variables in model selection process to improve numerical stability and to provide a model intercept that was easier to interpret. A nominal p-value of 0.05 was used to test all hypotheses. To evaluate the validity of the FP model for BMI, we graphically compared three models with the main FP model defined by (1). First, we estimated the model which categorized BMI into 30 narrow bins (1 bin for each BMI unit between 18 and 40, for every two BMI units between 40 and 54 and a single bin for BMI above 54) while also adjusting for age and smoking status. We then estimated separate FP models after omitting subjects with early death (< 1 year from baseline) and extreme BMI values (> 50).

Interactions between adjustment variables were tested to address the possibility of differences in the BMI-mortality curve across the age distribution, and by smoking history. The multivariable fractional polynomial interaction (MFPI) algorithm [26] was used to assess interactions, which first determined the best fitting polynomial functions for BMI and age using MFP and then tested for significant interactions between fractionally transformed variables and smoking history using a deviance difference test. We then verified interactions found by the MFPI algorithm graphically using Lowess smoothed curves. The use of FPs when fitting models using BMI as a continuous variable avoided inclusion of spurious interactions in a strictly linear model.

The BMI associated with minimum mortality was calculated by first estimating the final MFP model (including interaction terms) using logistic regression. To derive the optimal BMI, we set the first derivative of the estimated FP model equal to 0 and solved for BMI. As an example, the optimal BMI for the model with linear and quadratic BMI terms without interactions is , where and are the logistic regression coefficients for the linear and quadratic BMI terms, respectively. Confidence intervals were based on standard errors computed using the delta method.

The MFP approach provides a robust alternative to other commonly used methods for addressing the nonlinear and asymmetric relationship between BMI and mortality by allowing the data itself to determine the functional form for BMI and other adjustment variables. The closed test algorithm used in this study determined the best fitting model from a predefined set of candidate models based on power transformations of BMI. Improvements in model fit were found relative to other commonly used models including the linear-quadratic BMI model, which imposes symmetry on the BMI-mortality curve. This in conjunction with low variation in mortality among those in the 21-30 BMI range and high mortality among those with BMI over 40 resulted in an overly flat curve in the center. Allowing the estimated curve to have a flexible shape made the MFP model more sensitive to variation in death rates in the data. Improvements in model fit were generally accompanied by smaller estimates of the BMI associated with minimum mortality and narrower confidence intervals. While the differences in optimal BMI were small for the representative 50-year old female, there was a large discrepancy in the male sample with the nadir extending into the class I obese range (BMI 30.0-34.9).

The other common approach to modelling the nonlinear functional form for BMI is a nonparametric approach incorporating categorical variables defined by WHO BMI classifications. Assessing the risk of a BMI category relative to the normal category is a convenient method to account for the nonlinear form, but assumes mortality is uniform across a BMI category, which is problematic when a category is heterogenous. In particular, the normal category consists of a mix of healthy and sick lean. Studies employing the categorical approach also typically take the mortality risk of obese individuals beyond a given threshold as constant. We addressed this difficulty by allowing the data from across the entire BMI distribution to predict mortality risk at extreme obesity levels, where fewer observations exist. Because our results showed that mortality increased exponentially for extremely obese individuals, categorization can drastically underestimate mortality at the right tail of the BMI distribution. The use of wide BMI categories are also inadequate for the purposes of prediction. Categorizing BMI using finer intervals can alleviate some of these difficulties, but the decision of which categories to add is generally an arbitrary choice. Moreover, the additional categories increase the variance of parameter estimates, particularly in high BMI categories where the sample size is small.

Our study also employed methods that differentiate the effect of BMI across gender and age. The identification of distinct FP terms across gender samples point to important differences in the BMI-mortality relationship. Further differences were identified within the female sample through interactions. The finding of significant interactions between age and BMI is not common in literature, but has been identified in a select number of studies [39,40]. The inclusion of age-BMI interactions resulted in optimal BMI estimates that varied with age and predicted the nadir to exist in the overweight range for women at age 71 and above.

This paper was written for the Bariatric Obesity Outcome Modeling Collaborative. Members are listed at the end of the manuscript. This work was supported by the HQ Air Force Surgeon General under Award No. FA7014-08-2-0002. Opinions, interpretations, conclusions and recommendations are those of the author and are not necessarily endorsed by the U.S. Air Force. Dr. Wong was supported by the Department of Veterans Affairs, Health Services Research and Development Post-Doctoral Fellowship TPP 61-024.