# Multiple regression or anova or bestglm or forestplot or Boruta

I have data on a continuous health variable and following others: age, gender, height, weight, waist, city and season.

I applied multiple regression and got following output: (age, gender, height, waist and city were significant)

> summary(lm(y~., data=mydf))

Call:
lm(formula = y ~ ., data = mydf)

Residuals:
Min      1Q  Median      3Q     Max
-73.111  -9.528  -0.897   8.907  78.653

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept)  107.20300    2.83286  37.843  < 2e-16
age           -0.87090    0.12356  -7.048 1.97e-12  # SIGNIFICANT
genderM       -6.34184    0.33625 -18.861  < 2e-16  # SIGNIFICANT
htcm          -0.05992    0.02657  -2.255  0.02415  # SIGNIFICANT
wtkg           0.01247    0.04037   0.309  0.75745
waistcm        0.08095    0.03434   2.358  0.01842  # SIGNIFICANT
cityP          1.18070    0.38454   3.070  0.00214  # SIGNIFICANT
seasonsummer   0.28349    0.66278   0.428  0.66886
seasonwinter  -1.25711    0.67247  -1.869  0.06161

Residual standard error: 14.32 on 7767 degrees of freedom
(396 observations deleted due to missingness)
Multiple R-squared:  0.08514,   Adjusted R-squared:  0.08419
F-statistic: 90.35 on 8 and 7767 DF,  p-value: < 2.2e-16


With anova, I get following result: (age, gender, height, weight, city and season were significant)

> summary(aov(y~., data=mydf))
Df  Sum Sq Mean Sq F value  Pr(>F)
age            1   68902   68902 335.992 < 2e-16    # SIGNIFICANT
gender         1   72243   72243 352.280 < 2e-16    # SIGNIFICANT
htcm           1     149     149   0.726 0.39409
wtkg           1    1592    1592   7.762 0.00535    # SIGNIFICANT
waistcm        1     767     767   3.738 0.05323
city           1     829     829   4.043 0.04440    # SIGNIFICANT
season         2    3742    1871   9.124 0.00011    # SIGNIFICANT
Residuals   7767 1592791     205
396 observations deleted due to missingness


I applied bestglm and got following output: (only age and gender were significant)

> bestglm(mydf)
Morgan-Tatar search since factors present with more than 2 levels.
BIC
Best Model:
Df  Sum Sq Mean Sq F value Pr(>F)
age            1   68902   68902   334.8 <2e-16 # SIGNIFICANT
gender         1   72243   72243   351.0 <2e-16 # SIGNIFICANT
Residuals   7773 1599869     206
396 observations deleted due to missingness


Using randomforest, following is the importance: (in decreasing order: height, waist, weight, age, gender, season and city)

> library(randomForest)
> fit <- randomForest(y~., data=mydf, importance=TRUE)
> print(fit)

Call:
randomForest(formula = y ~ ., data = mydf)
Type of random forest: regression
Number of trees: 500
No. of variables tried at each split: 2

Mean of squared residuals: 207.2199
% Var explained: 7.45

# FOLLOWING IS FROM fit\$importance:

IncNodePurity
htcm       219809.13
waistcm    196753.10
wtkg       181179.19
age        119446.90
gender      83154.71
season      42938.42
city        27040.10

%IncMSE
htcm       72.663197
wtkg       68.040321
age        48.075415
waistcm    33.267517
gender     26.680004
season      5.932131
city        3.905936


It is especially surprising since this gives low importance to age and gender, while bestglm selected only these 2 variables for the model.

Using Boruta, importance table is:

var        importance

gender     55.4005861
waistcm    34.4082250
age        32.3720673
htcm       28.6817975
wtkg       26.7268140
season      8.0689392
city        7.9994742


Gender is highest here while it was much lower in randomForest.

Which one should I use?

Edit: The goal of this analysis is to find out which variables (out of age, gender, height, weight, waist, city and season) are independent predictors of y variable.

• What is the goal of this analysis? Dec 11 '14 at 1:38
• @AndrewM : I have added the goal in question above.
– rnso
Dec 11 '14 at 1:51
• In all cases it looks they are awful predictors, so maybe it doesn't matter which of a bunch of bad models you choose... .
– whuber
Dec 11 '14 at 3:11
• Why do you say they are awful predictors? Based on coefficients or R^2?
– rnso
Dec 11 '14 at 3:41
• Perhaps based on % Var explained: 7.45 in RF and Adjusted R-squared: 0.08419 in MR. But I don't think any of these methods is ideal for finding out which independent predictors are best. Dec 11 '14 at 18:01