8  Model selection

NoteIn this chapter

In this chapter, you will learn

  • why a model can be too simple or too complex, and how adding variables affects bias and variance (Section 8.1);
  • how to select models with likelihood-ratio tests, AIC and shrinkage, and what each of them answers (Section 8.2.1, Section 8.2.2, Section 8.2.3);
  • why model selection is a problem when you want to interpret effects and p-values, and how to evaluate predictions with cross-validation (Section 8.3, Section 8.4).

Workflow steps (Chapter 12): step 1 (question and goal), step 3 (choice of predictors), step 6 (inference after model selection).

Prerequisites: causal inference (Chapter 7), the linear model (Chapter 2).

8.1 The bias-variance trade-off

Apart from causality (see Chapter 7), the arguably most fundamental idea about modeling choice is the bias-variance trade-off, which applies regardless of whether we are interested in causal effects or predictive models. The idea is the following:

  • The more variables / complexity we include in the model, the better it can (in principle) adjust to the true relationship, thus reducing model error from bias.
  • The more variables / complexity we include in the model, the larger our error (variance) on the fitted coefficients, thus increasing model error from variance. This means, the model adopts to the given data but no longer to the underlying relationship.

If we sum both terms up, we see that the total error of a model that is too simple will be dominated by bias (underfitting), and the total error of a model that is too complex will be dominated by variance (overfitting):

Let’s confirm this with a small simulation - let’s assume I have a simple relationship between x and y:

set.seed(123)

x = runif(100)
y = 0.25 * x + rnorm(100, sd = 0.3)

coef(summary(lm(y~x)))["x", ]
  Estimate Std. Error    t value   Pr(>|t|) 
0.22305059 0.10254636 2.17511950 0.03202667

As you see, the effect is significant. Now, I add 80 new variables to the model which are just noise

xNoise = matrix(runif(8000), ncol = 80)
datNoise = data.frame(y=y,x=x, xNoise)

noiseModel = lm(y~., data = datNoise)
coef(summary(noiseModel))["x", ]
  Estimate Std. Error    t value   Pr(>|t|) 
 0.1872455  0.2756263  0.6793454  0.5055618

The effect estimates are relatively unchanged, but the CI has increased, and the p-values are n.s.

Tip 8.1: Precision variables

We have learned that adding variables increases the complexity of the model and thus the variance/uncertainty of the model. However, there are so-called precision variables that, when added to the model, actually reduce the variance/uncertainty of the model! Technically, they explain the response variable so well that they help to improve the overall fit of the model and thus help to reduce the overall variance of the model.

Here’s a small example based on toy data:

set.seed(42)
X = runif(100)
P = runif(100)
Y = 0.8*X + 10*P + rnorm(100, sd = 0.5)

Without the precision variable (P):

summary(lm(Y~X))

Call:
lm(formula = Y ~ X)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.3935 -2.7579  0.3701  2.2876  4.8379 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   5.3465     0.5826   9.177 7.39e-15 ***
X             0.4364     0.9639   0.453    0.652    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.896 on 98 degrees of freedom
Multiple R-squared:  0.002088,  Adjusted R-squared:  -0.008095 
F-statistic: 0.205 on 1 and 98 DF,  p-value: 0.6517

Effect of X is n.s. (standard error = 0.96)

With the precision variable (P):

summary(lm(Y~X+P))

Call:
lm(formula = Y ~ X + P)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.97608 -0.28268  0.00618  0.27548  1.38684 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.04947    0.12646  -0.391    0.697    
X            0.75999    0.15197   5.001 2.54e-06 ***
P           10.05137    0.16200  62.045  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.4563 on 97 degrees of freedom
Multiple R-squared:  0.9755,    Adjusted R-squared:  0.975 
F-statistic:  1929 on 2 and 97 DF,  p-value: < 2.2e-16

The effect of X on Y is now significant (standard error = 0.15)!

8.2 Model selection methods

Because of the bias-variance trade-off (Section 8.1), we cannot just fit the most complex model that we can imagine. Of course, such a model would have the lowest possible bias, but we would lose all power to see effects. Therefore, we must put bounds on model complexity. There are many methods to do so, but I would argue that three of them stand out (a fourth, cross-validation, follows at the end of the chapter, Section 8.4):

  1. Test whether there is evidence for going to a more complex model -> Likelihood Ratio Tests (Section 8.2.1)
  2. Optimize a penalized fit to the data -> AIC selection (Section 8.2.2)
  3. Specify a preference for parameter estimates to be small -> shrinkage estimation (Section 8.2.3)

Let’s go through them one by one

8.2.1 Likelihood-ratio tests

A likelihood-ratio test (LRT) compares 2 nested models, i.e. the simpler of the 2 models is a special case of the more complex one. The more complex model will always fit the data better, i.e. have a higher likelihood. Hence, fit alone cannot be used for model selection. The LRT gives a p-value for whether this improvement in likelihood is larger than one would expect if the simpler model were the correct model. The test statistic and its null distribution are explained in the section on testing in the ANOVA chapter (Section 3.2.2); here, we use the test to decide whether the more complex model is warranted.

In R, the anova function performs the test between 2 fitted models. For linear models with normal residuals, it reports an F-test, which is the exact version of the LRT for this case (Tip 3.1):

# Model 1
m1 = lm(Ozone ~ Wind , data = airquality)

# Model 2
m2 = lm(Ozone ~ Wind + Temp, data = airquality)

# LRT
anova(m1, m2)
Analysis of Variance Table

Model 1: Ozone ~ Wind
Model 2: Ozone ~ Wind + Temp
  Res.Df   RSS Df Sum of Sq      F    Pr(>F)    
1    114 79859                                  
2    113 53973  1     25886 54.196 3.149e-11 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

8.2.2 AIC model selection

Another method for model selection, and probably the most widely used, also because it does not require that models are nested, is the AIC = Akaike Information Criterion.

The AIC is defined as \(-2 \ln(\text{likelihood}) + 2k\), where \(k\) = number of parameters.

Essentially, this means AIC = misfit + penalty for complexity: the first term (\(-2 \ln(\text{likelihood})\)) gets smaller the better the model fits, while the second term (\(2k\)) grows with the number of parameters.

Lower AIC is better!

For the two ozone models m1 and m2 from the LRT example:

AIC(m1)
[1] 1093.187
AIC(m2)
[1] 1049.741

Note 1: It can be shown that AIC is asymptotically equivalent to leave-one-out cross-validation (under regularity conditions, for likelihood-based loss; see Section 8.4), so what AIC is optimizing is essentially the predictive error of the model on new data.

Note 2: AIC values can only be compared between models fitted to exactly the same rows of the data (Section 6.1). As a rule of thumb, models with a difference in AIC of less than about 2 are essentially equivalent. For a small n relative to the number of parameters, use the corrected AICc (MuMIn::AICc).

Note 8.1: AIC or BIC?

The most common alternative to the AIC is the BIC (Bayesian information criterion), \(-2 \ln(\text{likelihood}) + \ln(n) \, k\). Instead of 2 per parameter, BIC adds \(\ln(n)\) per parameter, which is more than 2 for any data set with at least 8 observations, so BIC selects smaller models, and the more so the larger the data. The two criteria answer different questions: AIC aims at the model that predicts best, BIC at the true model, assuming it is among the candidates (Aho, Derryberry and Peterson, 2014). If they disagree, the extra parameters improve the fit, but not enough to pass BIC’s stronger penalty: the larger model probably predicts better, but whether the extra structure really exists remains uncertain. State which of the two questions you ask. Note also that in mixed models, \(n\) is ambiguous (number of observations or number of groups?); R uses the number of observations (for an example where this matters, see the case study on species cover, Section 11.6.1). Other criteria (WAIC, DIC) exist, but are beyond the scope of this course; Bayesian model selection is covered in the chapter on Bayesian model selection of Hartig (2024).

Hartig, F. (2024) “Introduction to Bayesian statistics with R.” Available at: https://florianhartig.github.io/LearningBayes/.
Aho, K., Derryberry, D. and Peterson, T. (2014) “Model selection for ecologists: The worldviews of AIC and BIC,” Ecology, 95(3), pp. 631–636.
CautionExercise

Compare results of AIC with likelihood-ratio tests. Discuss: When to use one or the other?

The two answer different questions. A likelihood-ratio test is a hypothesis test: it requires nested models, returns a p-value, and controls the type-I error rate at your chosen \(\alpha\) - it answers “is there significant evidence to prefer the more complex model?”. AIC estimates the model’s out-of-sample predictive error (it is asymptotically equivalent to leave-one-out cross-validation), does not require nesting, and returns no p-value or error guarantee - it answers “which model is expected to predict best?”.

For the two ozone models above, both agree: the F-test (an LRT) is highly significant, and the AIC of the model with Temp is lower by more than 40. Often, however, they disagree. For adding a single parameter, AIC prefers the larger model whenever the deviance drops by more than 2, which corresponds to an LRT p-value of about 0.157 - so AIC is considerably more liberal about including predictors than an LRT at \(\alpha = 0.05\).

Rule of thumb: use the LRT (or F-test) for a confirmatory, pre-specified nested comparison where you want error control; use AIC when the goal is prediction, or to rank many (possibly non-nested) candidate models. For neither should you read the resulting coefficients as unbiased causal effects (see the section on problems with model selection for inference, Section 8.3, and the chapter on causal inference, Chapter 7).

8.2.3 Shrinkage estimation

A third option for model selection are shrinkage estimators. These include the LASSO and ridge, but also random-effects (see Section 4.3.2) could be seen as shrinkage estimators.

The basic idea behind these estimators is not to reduce the number of parameters, but to reduce the flexibility of the model by introducing a penalty on the regression coefficients that code a preference for smaller or zero coefficient values. Effectively, this can either amount to model selection (because some coefficients are shrunk directly to zero), or it can mean that we can fit very large models while still being able to do good predictions, or avoid overfitting.

To put a ridge penalty on the standard lm, we can use

library(MASS)
lm.ridge(Ozone ~ Wind + Temp + Solar.R, data = airquality, lambda = 2)
                     Wind         Temp      Solar.R 
-62.73376169  -3.30622990   1.62842247   0.05961015

Here, lambda is the strength of the penalty (how to choose it is shown at the end of this section). We can see how the regression estimates vary for different penalties via

plot( lm.ridge( Ozone ~ Wind + Temp + Solar.R, data = airquality,
              lambda = seq(0, 200, 0.1) ) )

We can use glmmTMB for shrinkage estimation by changing the slopes into random slopes with only one level in a dummy grouping variable:

  1. Constant shrinkage
library(glmmTMB)
airquality$group = as.factor(rep(1, nrow(airquality)))
ridgeGlmmTMB01 = glmmTMB(Ozone~(0+ Wind + Temp + Solar.R || group ),
                       data = airquality,
                       start = list(theta = rep(log(0.1), 3)),
                       map = list(theta = factor(c(NA, NA, NA))))

ridgeGlmmTMB10 = glmmTMB(Ozone~(0+ Wind + Temp + Solar.R || group ),
                       data = airquality,
                       start = list(theta = rep(log(10), 3)),
                       map = list(theta = factor(c(NA, NA, NA))))
summary(ridgeGlmmTMB01)
 Family: gaussian  ( identity )
Formula:          Ozone ~ (0 + Wind + Temp + Solar.R || group)
Data: airquality

      AIC       BIC    logLik -2*log(L)  df.resid 
   1082.8    1088.2    -539.4    1078.8       109 

Random effects:

Conditional model:
 Groups   Name    Variance Std.Dev.
 group    Wind      0.01    0.10   
          Temp      0.01    0.10   
          Solar.R   0.01    0.10   
 Residual         891.61   29.86   
Number of obs: 111, groups:  group, 1

Dispersion estimate for gaussian family (sigma^2):  892 

Conditional model:
            Estimate Std. Error z value Pr(>|z|)
(Intercept)    6.045      9.589    0.63    0.528
summary(ridgeGlmmTMB10)
 Family: gaussian  ( identity )
Formula:          Ozone ~ (0 + Wind + Temp + Solar.R || group)
Data: airquality

      AIC       BIC    logLik -2*log(L)  df.resid 
   1023.1    1028.5    -509.6    1019.1       109 

Random effects:

Conditional model:
 Groups   Name    Variance Std.Dev.
 group    Wind    100.0    10.00   
          Temp    100.0    10.00   
          Solar.R 100.0    10.00   
 Residual         444.5    21.08   
Number of obs: 111, groups:  group, 1

Dispersion estimate for gaussian family (sigma^2):  444 

Conditional model:
            Estimate Std. Error z value Pr(>|z|)   
(Intercept)   -64.60      22.92  -2.819  0.00482 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Note: The (…||…) syntax will model the individual random slopes individually (the default is to fit a covariance matrix between the random effects).

We can use the map and the start arguments to a) mark the variances of the random slopes as constant and b) to set them so specific values (smaller values mean higher regularization). Note that glmmTMB parameterizes theta as the log of the random effect standard deviation, so log(0.1) fixes the SD of the slopes to 0.1 (strong shrinkage) and log(10) to 10 (weak shrinkage).

Random effects are now our slopes:

ranef(ridgeGlmmTMB01)[[1]]$group
         Wind      Temp   Solar.R
1 -0.07822852 0.2132513 0.1095346
ranef(ridgeGlmmTMB10)[[1]]$group
      Wind     Temp    Solar.R
1 -3.32082 1.653688 0.05983465

For lower variances they are more strongly biased towards 0 (stronger shrinkage); for higher variances, shrinkage is reduced and estimates are less biased towards 0. With SD = 10, the estimates are practically identical to those of the unpenalized lm(Ozone ~ Wind + Temp + Solar.R, data = airquality) (the glmmTMB model still includes a fixed intercept).

  1. Adaptive shrinkage

The next obvious step would be to let the model and the data decide on the variance estimate (the amount of shrinkage), which is called adaptive shrinkage.

Small simulation example:

set.seed(123)
n = 30
p = 50
X = matrix(runif(p*n), n, p)
head(X)
          [,1]       [,2]       [,3]      [,4]      [,5]      [,6]      [,7]
[1,] 0.2875775 0.96302423 0.66511519 0.1306957 0.6478935 0.8474532 0.8397678
[2,] 0.7883051 0.90229905 0.09484066 0.6531019 0.3198206 0.4975273 0.3124482
[3,] 0.4089769 0.69070528 0.38396964 0.3435165 0.3077200 0.3879090 0.7082903
[4,] 0.8830174 0.79546742 0.27438364 0.6567581 0.2197676 0.2464490 0.2650178
[5,] 0.9404673 0.02461368 0.81464004 0.3203732 0.3694889 0.1110965 0.5943432
[6,] 0.0455565 0.47779597 0.44851634 0.1876911 0.9842192 0.3899944 0.4812898
           [,8]      [,9]       [,10]       [,11]     [,12]     [,13]     [,14]
[1,] 0.48204261 0.5468262 0.871043414 0.784575267 0.6172353 0.1438170 0.9770990
[2,] 0.25296493 0.6623176 0.006300784 0.009429905 0.2862853 0.1928159 0.2963022
[3,] 0.21625479 0.1716985 0.072057124 0.779065883 0.7377974 0.8967387 0.7259830
[4,] 0.67437639 0.6330554 0.164211225 0.729390652 0.8340543 0.3081196 0.7856878
[5,] 0.04766363 0.3118697 0.770334074 0.630131853 0.3142708 0.3633005 0.1054177
[6,] 0.70085309 0.7245543 0.735184306 0.480910830 0.4925665 0.7839465 0.2395946
         [,15]     [,16]     [,17]      [,18]      [,19]     [,20]     [,21]
[1,] 0.5474595 0.5238226 0.2151665 0.46471377 0.80094821 0.5533140 0.2372297
[2,] 0.6442402 0.3498018 0.8139337 0.08253118 0.09164469 0.9060481 0.6864904
[3,] 0.5962635 0.2405307 0.3077639 0.86010685 0.83210777 0.5874614 0.2258184
[4,] 0.3219374 0.0581918 0.6877427 0.39566064 0.27685507 0.4234637 0.3184946
[5,] 0.8911143 0.2366197 0.9326811 0.73589935 0.75311010 0.9495853 0.1739838
[6,] 0.6262569 0.8900779 0.1157797 0.17174341 0.96415257 0.7090379 0.8014296
         [,22]      [,23]      [,24]     [,25]      [,26]      [,27]     [,28]
[1,] 0.7591353 0.84555612 0.03978064 0.4186692 0.27776932 0.25880120 0.0494156
[2,] 0.8448346 0.80643516 0.63480051 0.7178845 0.24940659 0.81946714 0.6951052
[3,] 0.4579216 0.11733123 0.53956091 0.7425391 0.15943548 0.01547767 0.3632616
[4,] 0.7296317 0.71268655 0.14010355 0.8719988 0.07163006 0.65440077 0.8841336
[5,] 0.1040786 0.23526886 0.28371178 0.6078680 0.55743758 0.81171958 0.7752972
[6,] 0.2199832 0.07495675 0.58303060 0.7562034 0.48983705 0.46136893 0.1392036
         [,29]     [,30]     [,31]     [,32]     [,33]     [,34]     [,35]
[1,] 0.9624213 0.1406462 0.9236992 0.8868626 0.5081581 0.0872193 0.4612166
[2,] 0.7363340 0.1697691 0.5425984 0.2040959 0.4495999 0.2237800 0.7927383
[3,] 0.6127235 0.7619853 0.8523646 0.7706230 0.6232614 0.5721487 0.5855656
[4,] 0.1199289 0.5273949 0.5835629 0.5963630 0.1399779 0.4001692 0.3947590
[5,] 0.5502590 0.8609894 0.6683236 0.9576697 0.9079464 0.5654654 0.8112373
[6,] 0.2627563 0.6735550 0.5113146 0.1586884 0.5694433 0.8296239 0.3185631
         [,36]     [,37]      [,38]      [,39]     [,40]     [,41]       [,42]
[1,] 0.7088554 0.3923313 0.07228658 0.83179032 0.4105144 0.8585914 0.004638151
[2,] 0.8267155 0.1719368 0.82633158 0.40536388 0.4018139 0.8873848 0.277560080
[3,] 0.8918963 0.3196421 0.47954408 0.07263193 0.6485563 0.4890915 0.325203143
[4,] 0.5394107 0.4135551 0.39122878 0.95375862 0.9163057 0.7180918 0.588706277
[5,] 0.5553136 0.7593357 0.29912278 0.21718378 0.2166621 0.4867056 0.249684701
[6,] 0.2764399 0.2718902 0.24430825 0.40562527 0.5479989 0.9887089 0.043117281
          [,43]     [,44]     [,45]      [,46]      [,47]     [,48]      [,49]
[1,] 0.24872636 0.8157679 0.8061095 0.02489827 0.46357376 0.7223989 0.14878158
[2,] 0.61495395 0.5262726 0.3377528 0.00484174 0.35472962 0.9220046 0.61276590
[3,] 0.03175732 0.2030400 0.9162859 0.13137825 0.88242303 0.6138409 0.04196682
[4,] 0.14642354 0.8481077 0.3571064 0.02481710 0.07580286 0.1609665 0.49835899
[5,] 0.70307297 0.3704942 0.9753650 0.09212945 0.88329119 0.1829407 0.73950822
[6,] 0.06556071 0.3033252 0.3843785 0.47555822 0.18822184 0.7301836 0.75389978
         [,50]
[1,] 0.7001786
[2,] 0.0655749
[3,] 0.7330400
[4,] 0.3959316
[5,] 0.4614878
[6,] 0.7252994
effects = runif(p, -1,1)
Y = X%*%effects + rnorm(n)

30 observations, 50 predictors, each with its own true effect drawn from a Uniform(-1, 1) distribution, let’s start with a fixed-effect model:

m1 = glmmTMB(Y~X)
Warning in glmmTMB(Y ~ X): use of the 'data' argument is recommended
dropping columns from rank-deficient conditional model: X30, X31, X32, X33, X34, X35, X36, X37, X38, X39, X40, X41, X42, X43, X44, X45, X46, X47, X48, X49, X50
Warning in finalizeTMB(TMBStruc, obj, fit, h, data.tmb.old): Model convergence
problem; non-positive-definite Hessian matrix. See vignette('troubleshooting')
Warning in finalizeTMB(TMBStruc, obj, fit, h, data.tmb.old): Model convergence
problem; function evaluation limit reached without convergence (9). See
vignette('troubleshooting'), help('diagnose')
summary(m1)
 Family: gaussian  ( identity )
Formula:          Y ~ X

      AIC       BIC    logLik -2*log(L)  df.resid 
       NA        NA        NA        NA        -1 


Dispersion estimate for gaussian family (sigma^2): 1.03e-06 

Conditional model:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)  6.292147   0.008442   745.4   <2e-16 ***
X1           2.318114   0.002618   885.5   <2e-16 ***
X2          -2.781089   0.003650  -762.0   <2e-16 ***
X3          -2.164066   0.002563  -844.2   <2e-16 ***
X4          -2.296509   0.004205  -546.2   <2e-16 ***
X5          -3.596781   0.003102 -1159.4   <2e-16 ***
X6          -2.243081   0.002638  -850.4   <2e-16 ***
X7          -3.018425   0.003666  -823.3   <2e-16 ***
X8           1.626393   0.002524   644.3   <2e-16 ***
X9          -1.283200   0.002433  -527.4   <2e-16 ***
X10          0.599021   0.002334   256.6   <2e-16 ***
X11         -5.081301   0.004303 -1181.0   <2e-16 ***
X12         -0.160191   0.002116   -75.7   <2e-16 ***
X13          2.671564   0.004333   616.6   <2e-16 ***
X14          6.783907   0.004730  1434.2   <2e-16 ***
X15         -3.782583   0.002070 -1827.7   <2e-16 ***
X16          1.125239   0.004558   246.9   <2e-16 ***
X17         -0.438207   0.003366  -130.2   <2e-16 ***
X18         -2.433692   0.003201  -760.2   <2e-16 ***
X19          4.510137   0.002949  1529.2   <2e-16 ***
X20          1.182427   0.002699   438.0   <2e-16 ***
X21         -2.414001   0.002371 -1018.2   <2e-16 ***
X22         -5.986762   0.005489 -1090.7   <2e-16 ***
X23          6.962840   0.004361  1596.7   <2e-16 ***
X24         -1.733277   0.002813  -616.1   <2e-16 ***
X25          2.946605   0.003565   826.6   <2e-16 ***
X26         -0.373676   0.001769  -211.2   <2e-16 ***
X27          0.648173   0.001296   500.1   <2e-16 ***
X28         -6.191007   0.002481 -2495.1   <2e-16 ***
X29          0.796629   0.002039   390.8   <2e-16 ***
X30                NA         NA      NA       NA    
X31                NA         NA      NA       NA    
X32                NA         NA      NA       NA    
X33                NA         NA      NA       NA    
X34                NA         NA      NA       NA    
X35                NA         NA      NA       NA    
X36                NA         NA      NA       NA    
X37                NA         NA      NA       NA    
X38                NA         NA      NA       NA    
X39                NA         NA      NA       NA    
X40                NA         NA      NA       NA    
X41                NA         NA      NA       NA    
X42                NA         NA      NA       NA    
X43                NA         NA      NA       NA    
X44                NA         NA      NA       NA    
X45                NA         NA      NA       NA    
X46                NA         NA      NA       NA    
X47                NA         NA      NA       NA    
X48                NA         NA      NA       NA    
X49                NA         NA      NA       NA    
X50                NA         NA      NA       NA    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Not surprising, we get many warnings about convergence issues and NAs, because we have more predictors (p = 50) than observations (n = 30) - the fixed-effect model is not even identifiable. Also, the effects that can be estimated are not even close to their true values.

Let’s try adaptive shrinkage. As before, we turn each predictor’s slope into a random slope within a dummy group that has only one level:

data = data.frame(Y = Y, X, group = as.factor(rep(1, n)))
form = paste0("Y~ 0+ (0+",paste0(paste0(  "X", 1:p ),  collapse = "+"),"||group)" )
print(form)
[1] "Y~ 0+ (0+X1+X2+X3+X4+X5+X6+X7+X8+X9+X10+X11+X12+X13+X14+X15+X16+X17+X18+X19+X20+X21+X22+X23+X24+X25+X26+X27+X28+X29+X30+X31+X32+X33+X34+X35+X36+X37+X38+X39+X40+X41+X42+X43+X44+X45+X46+X47+X48+X49+X50||group)"

What we want for a proper (ridge-like) adaptive shrinkage is a single, shared variance for all p random slopes, estimated jointly by pooling information across all p coefficients at once - this is the classic “ridge regression as a mixed model” trick. We achieve this with the map argument. By default, the || syntax would give every predictor its own, independently estimated variance (p = 50 separate variance parameters) - but map lets us tell glmmTMB which of its internal parameters should be tied together (estimated as one shared value) instead of independently. Giving all p entries of theta the same factor level ties them into a single, freely estimated parameter. This is needed here because, with only one level in group, each individually estimated variance would be informed by essentially a single realized slope - nowhere near enough information to estimate 50 separate variances reliably. Tying them together pools all p coefficients into the estimation of that one shared variance, which is far better identified. Note that, unlike the constant-shrinkage example above, where NA fixed the variances at chosen constants, a shared (non-NA) level still lets the model estimate their common value from the data:

m2 = glmmTMB(as.formula(form), data = data,
             start = list(theta = rep(0, p)),
             map   = list(theta = factor(rep(1, p))))
summary(m2)
 Family: gaussian  ( identity )
Formula:          
Y ~ 0 + (0 + X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10 +  
    X11 + X12 + X13 + X14 + X15 + X16 + X17 + X18 + X19 + X20 +  
    X21 + X22 + X23 + X24 + X25 + X26 + X27 + X28 + X29 + X30 +  
    X31 + X32 + X33 + X34 + X35 + X36 + X37 + X38 + X39 + X40 +  
    X41 + X42 + X43 + X44 + X45 + X46 + X47 + X48 + X49 + X50 ||      group)
Data: data

      AIC       BIC    logLik -2*log(L)  df.resid 
    107.5     110.3     -51.8     103.5        28 

Random effects:

Conditional model:
 Groups   Name Variance Std.Dev.
 group    X1   0.2910   0.5395  
          X2   0.2910   0.5395  
          X3   0.2910   0.5395  
          X4   0.2910   0.5395  
          X5   0.2910   0.5395  
          X6   0.2910   0.5395  
          X7   0.2910   0.5395  
          X8   0.2910   0.5395  
          X9   0.2910   0.5395  
          X10  0.2910   0.5395  
          ...  ...      ...     
 Residual      0.6375   0.7984  
Number of obs: 30, groups:  group, 1

Dispersion estimate for gaussian family (sigma^2): 0.638

Now, extract the estimates and compare them to the true values

shrinkageEstimates = as.numeric(ranef(m2)$cond$group)

plot(effects, shrinkageEstimates, xlab = "True effects", ylab = "Estimated shrinkage effects")
abline(0,1)

Note that there is quite a bit of variance, but this is to be expected when we estimate 50 parameters with 30 observations. However, there is a clear signal, in the sense that the correlation between true and estimated effects is significantly positive:

cor.test(effects, shrinkageEstimates)

    Pearson's product-moment correlation

data:  effects and shrinkageEstimates
t = 4.3793, df = 48, p-value = 6.425e-05
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.3006812 0.7074441
sample estimates:
      cor 
0.5343104

Choosing the amount of shrinkage. In lm.ridge() above, we simply set the penalty lambda by hand. In practice, the penalty is usually chosen by cross-validation (see Section 8.4), i.e. we pick the value that gives the best out-of-sample predictions. The glmnet package does this automatically with cv.glmnet(), for ridge (alpha = 0), the LASSO (alpha = 1), and the elastic net (values in between). Applied to the simulation with 30 observations and 50 predictors from the box above:

library(glmnet)
set.seed(1) # the cross-validation folds are random
cvRidge = cv.glmnet(X, Y, alpha = 0) # ridge
cvLasso = cv.glmnet(X, Y, alpha = 1) # LASSO
plot(cvLasso) # prediction error against log(lambda)

ridgeEstimates = as.numeric(coef(cvRidge, s = "lambda.min"))[-1] # [-1] removes the intercept
lassoEstimates = as.numeric(coef(cvLasso, s = "lambda.min"))[-1]
cor(effects, ridgeEstimates)
[1] 0.4982708
cor(effects, lassoEstimates)
[1] 0.4310109
sum(lassoEstimates != 0) # number of predictors kept by the LASSO
[1] 8

Which method to use when?

  • Ridge shrinks all coefficients towards zero, but keeps all predictors in the model. It works well if there are many predictors with small effects, and with collinear predictors (it spreads the effect over them).
  • LASSO shrinks some coefficients exactly to zero, i.e. it also selects variables. It works well if only a few predictors have relevant effects. With collinear predictors, it tends to pick one of them more or less arbitrarily.
  • Elastic net (0 < alpha < 1) is a compromise between the two.
  • Random effects (the glmmTMB trick in the box above) give ridge-type shrinkage, where the amount of shrinkage is estimated from the data instead of by cross-validation. This is convenient if you need other features of mixed models (other distributions, further random effects).

Note that all shrinkage estimators are biased towards zero, and the standard errors and p-values of a normal regression do not apply to them. Use them for prediction, not to test hypotheses (see Section 8.3).

8.2.4 P-hacking

The most dubious model selection strategy, a questionable research practice that can amount to scientific misconduct, is p-hacking. The purpose of this exercise is to show you how not to do model selection, i.e, that by playing around with the variables, you can make any outcome significant. That is why your hypothesis needs to be fixed before looking at the data, ideally through pre-registration, based on an experimental plan or a causal analysis. Here is the example (why this works is explained in Section 8.3.2):

CautionExercise — p-hacking

We have (simulated) measurements of plant performance. The goal of the analysis was to find out if Gen1 has an effect on Performance. Various other variables are measured. As you can see, the way I simulated the data, none of the variables has an effect on the response, so this is pure noise

set.seed(1)
dat = data.frame(matrix(rnorm(300), ncol = 10))
colnames(dat) = c("Performance", "Gen1", "Gen2", "soilC", "soilP", "Temp",
                  "Humidity", "xPos", "yPos", "Water")
fullModel <- lm(Performance ~ ., data = dat)

When you run summary(fullModel), you will see that there is no significant effect of Gen1. Task for you: P-hack the analysis, i.e. make an effect appear, by trying around (systematically, e.g. with selecting with data, model selection, or by hand to find a model combination that has an effect). Popular strategies for p-hacking include changing the model, but also sub-selecting particular data. The group who finds the model with the highest significance for Gen1 wins!

summary(lm(Performance ~ Gen1 * Humidity, data = dat[20:30,]))

Call:
lm(formula = Performance ~ Gen1 * Humidity, data = dat[20:30, 
    ])

Residuals:
     Min       1Q   Median       3Q      Max 
-0.71665 -0.39627 -0.05915  0.28044  0.91257 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)   
(Intercept)    -0.5248     0.2277  -2.304  0.05465 . 
Gen1            0.8657     0.2276   3.804  0.00668 **
Humidity        0.6738     0.2544   2.649  0.03298 * 
Gen1:Humidity  -0.5480     0.1756  -3.122  0.01680 * 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.6102 on 7 degrees of freedom
Multiple R-squared:  0.7004,    Adjusted R-squared:  0.572 
F-statistic: 5.454 on 3 and 7 DF,  p-value: 0.03

We tried several subsets and models before finding this one; a single significant result among many tests is expected by chance (see Section 2.6.3).

Here some more inspiration on p-hacking:

  1. Hack Your Way To Scientific Glory: https://projects.fivethirtyeight.com/p-hacking/
  2. False-Positive Psychology: Undisclosed Flexibility in Data Collection and Analysis Allows Presenting Anything as Significant: Simmons, Nelson and Simonsohn (2011)
  3. Sixty seconds on … P-hacking: Hawkes (2018)
Simmons, J.P., Nelson, L.D. and Simonsohn, U. (2011) “False-positive psychology: Undisclosed flexibility in data collection and analysis allows presenting anything as significant,” Psychological Science, 22(11), pp. 1359–1366.
Hawkes, N. (2018) “Sixty seconds on … P-hacking,” BMJ, 362, p. k4039.

John Oliver about p-hacking:

8.3 Problems with model selection for inference

Model selection works great to generate predictive models. The big problem with model selection is when we want interpret effect estimates. Here, we have two problems

  1. Model selection doesn’t respect causal relationships
  2. Model selection methods often lead to wrong p-values, CIs etc.

Let’s look at them 1 by one:

8.3.1 Causality

In general, model selection does NOT solve the problem of estimating causal effects. Quite the contrary: most model selection methods act against estimating causal effects. Consider the following example, where we create a response that causally depends on two collinear predictors, both with a direct effect of 1. Moreover, x2 is partly caused by x1, so x1 also affects the response indirectly via x2 (x1 -> x2 -> y)

set.seed(123)
x1 = runif(100)
x2 = 0.8 * x1 + 0.2 *runif(100)
y = x1 + x2 + rnorm(100)

Given the structure of the data, we should run a multiple regression, and the multiple regression will get the two direct effects roughly right: both effects are n.s. (the CIs are very wide because of the collinearity), but the estimates are roughly right and the true values are in the 95% CI.

m1 = lm(y ~ x1 + x2)
summary(m1)

Call:
lm(formula = y ~ x1 + x2)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.8994 -0.6821 -0.1086  0.5749  3.3663 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)
(Intercept)  -0.1408     0.2862  -0.492    0.624
x1            1.2158     1.5037   0.809    0.421
x2            0.8518     1.8674   0.456    0.649

Residual standard error: 0.9765 on 97 degrees of freedom
Multiple R-squared:  0.237, Adjusted R-squared:  0.2212 
F-statistic: 15.06 on 2 and 97 DF,  p-value: 2.009e-06

A model selection will remove one of the variables. The remaining x1 then absorbs the effect that runs through x2, i.e. it no longer estimates the direct effect of x1, but its total effect (direct effect 1 plus indirect effect 0.8 × 1 = 1.8). The estimate of about 1.9 is consistent with this total effect, but the true direct effect of 1 is far outside the 95% CI.

m2 = MASS::stepAIC(m1)
summary(m2)

Call:
lm(formula = y ~ x1)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.9047 -0.6292 -0.1019  0.6077  3.3394 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.04633    0.19670  -0.236    0.814    
x1           1.88350    0.34295   5.492 3.13e-07 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.9725 on 98 degrees of freedom
Multiple R-squared:  0.2353,    Adjusted R-squared:  0.2275 
F-statistic: 30.16 on 1 and 98 DF,  p-value: 3.134e-07

Here, x2 is a mediator, and whether to adjust for it depends on whether we want the direct or the total effect of x1 (Section 7.3). The problem is that the model selection takes this decision for us, based on how much x2 improves the fit and not on causal reasoning, and thereby silently changes the meaning of the remaining estimate.

If you understand what AIC model selection (Section 8.2.2) is doing, you wouldn’t be concerned - this is actually not a bug, but rather, the method is doing exactly what it is designed for. However, it was not designed to estimate causal effects. Let’s do another simulation: here, x1,x2 and x3 all have the same effect on the response, but x1 and x2 are highly collinear (x1 and x2 are generated as in the collinearity simulation of the linear regression chapter, Section 2.8), while x3 is independent

set.seed(123)
x1 = runif(100)
x2 = 0.95 * x1 + 0.05 *runif(100)
x3 = runif(100)
y = x1 + x2 + x3 + rnorm(100)
m1 = lm(y ~ x1)
AIC(m1)
[1] 286.2746

Comparing a base model with only x1 to x1 + x2

m2 = lm(y ~ x1 + x2)
anova(m1, m2)
Analysis of Variance Table

Model 1: y ~ x1
Model 2: y ~ x1 + x2
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1     98 96.548                           
2     97 94.066  1    2.4817 2.5591 0.1129
AIC(m2)
[1] 285.6705

The same for the comparison of x1 to x1 + x3

m3 = lm(y ~ x1 + x3)
anova(m1, m3)
Analysis of Variance Table

Model 1: y ~ x1
Model 2: y ~ x1 + x3
  Res.Df    RSS Df Sum of Sq      F   Pr(>F)   
1     98 96.548                                
2     97 86.848  1    9.6998 10.834 0.001391 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AIC(m3)
[1] 277.6868

What we see is that the model selection methods are far more enthusiastic about including x3. And with good reason: x3 is an independent predictor, so adding x3 improves the model a lot. The effect of x2, however, is mostly already included in x1, so adding x2 will only slightly improve the model. Thus, from the point of view of how good the data is explained, it doesn’t pay off to add x2. From the causal viewpoint, however, adding x3 is irrelevant for the effect of x1, because it is not a confounder, while adding x2 changes the meaning of the effect of x1 (total instead of direct effect, as in the first example). If x2 were instead a common cause of x1 and the response (a confounder), leaving it out would bias the estimate, and the model selection would make exactly the wrong choice.

8.3.2 P-values

The second problem with model selection is the calibration of p-values. Let’s revisit our simulation example from the bias-variance trade-off (Section 8.1), where we added 80 noisy predictors with no effect to a situation where we had one variable with an effect (noiseModel). We saw before that the y~x effect is n.s. after adding all 80 predictors. Let’s use AIC in a stepwise model selection procedure to reduce model complexity

Tip 8.2: Stepwise vs. global model selection

Systematic LRT or AIC model selections are often used stepwise or global. Stepwise means that we start with the most simple (forward) or the most complex (backward) model, and then run a chain of model selection steps (AIC or LRT) adding (forward) or removing (backward) complexity until we arrive at an optimum. Global means that we run immediately all possible models and compare their AIC. The main function for stepwise selection is MASS::stepAIC, for global selection MuMIn::dredge. There was a lot of discussion about step-wise vs. global selection in the literature, that mostly revolved around the fact that a stepwise selection is faster, but will not always find the global optimum. However, compared to the other problems discussed here (causality, p-values), I do not consider this a serious problem. Thus, if you can, run a global selection by all means, but if this is computationally prohibitive, stepwise selections are also fine. Note that dredge() only runs with options(na.action = "na.fail"), which makes sure that all models use the same rows.

library(MASS)
reduced = stepAIC(noiseModel)
Note

When you inspect the output of stepAIC, you can see that the function calculates at each step the AIC improvement for each predictor that could be removed, and then chooses to remove the predictor that leads to the strongest AIC improvement first.

Here is the selected model

summary(reduced)

Call:
lm(formula = y ~ x + X2 + X3 + X4 + X7 + X8 + X9 + X10 + X11 + 
    X12 + X13 + X15 + X16 + X17 + X19 + X23 + X24 + X25 + X26 + 
    X28 + X32 + X34 + X36 + X37 + X38 + X40 + X41 + X42 + X44 + 
    X45 + X46 + X47 + X48 + X49 + X51 + X54 + X55 + X56 + X57 + 
    X58 + X60 + X61 + X65 + X66 + X67 + X69 + X71 + X72 + X74 + 
    X75 + X76 + X79 + X80, data = datNoise)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.29437 -0.11119 -0.00176  0.10598  0.33330 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.66297    0.37833   1.752 0.086372 .  
x            0.23719    0.10004   2.371 0.021990 *  
X2          -0.16068    0.10439  -1.539 0.130590    
X3          -0.11884    0.10389  -1.144 0.258568    
X4          -0.13278    0.10011  -1.326 0.191278    
X7           0.17556    0.11021   1.593 0.118026    
X8           0.17163    0.12746   1.347 0.184710    
X9          -0.21560    0.10822  -1.992 0.052301 .  
X10          0.24336    0.09903   2.457 0.017823 *  
X11         -0.14428    0.11623  -1.241 0.220757    
X12         -0.16108    0.10125  -1.591 0.118472    
X13         -0.19011    0.14470  -1.314 0.195407    
X15         -0.13111    0.10703  -1.225 0.226796    
X16         -0.47202    0.10938  -4.315 8.37e-05 ***
X17          0.25903    0.11279   2.297 0.026243 *  
X19         -0.15402    0.10664  -1.444 0.155427    
X23          0.17998    0.11368   1.583 0.120224    
X24          0.23410    0.10586   2.212 0.032009 *  
X25         -0.23611    0.11903  -1.984 0.053285 .  
X26          0.29341    0.10141   2.893 0.005811 ** 
X28          0.22012    0.11068   1.989 0.052687 .  
X32          0.10311    0.10640   0.969 0.337543    
X34         -0.21144    0.12567  -1.682 0.099252 .  
X36          0.15986    0.10697   1.494 0.141896    
X37          0.24533    0.10131   2.422 0.019452 *  
X38          0.29514    0.10844   2.722 0.009144 ** 
X40          0.15977    0.11203   1.426 0.160575    
X41          0.30199    0.11911   2.535 0.014695 *  
X42         -0.18966    0.11599  -1.635 0.108864    
X44          0.18595    0.10397   1.788 0.080285 .  
X45          0.20115    0.09857   2.041 0.047051 *  
X46         -0.32878    0.12243  -2.685 0.010040 *  
X47          0.17192    0.10575   1.626 0.110837    
X48         -0.36013    0.10129  -3.556 0.000886 ***
X49         -0.13339    0.11330  -1.177 0.245137    
X51         -0.15802    0.10369  -1.524 0.134378    
X54         -0.17854    0.11057  -1.615 0.113206    
X55         -0.23024    0.11682  -1.971 0.054766 .  
X56         -0.25930    0.10889  -2.381 0.021447 *  
X57         -0.24371    0.10668  -2.284 0.027009 *  
X58          0.31375    0.10379   3.023 0.004084 ** 
X60          0.28972    0.10872   2.665 0.010586 *  
X61          0.14802    0.10430   1.419 0.162616    
X65         -0.47949    0.11591  -4.137 0.000148 ***
X66          0.11926    0.10296   1.158 0.252723    
X67         -0.10825    0.09985  -1.084 0.283944    
X69          0.29208    0.10506   2.780 0.007850 ** 
X71         -0.30513    0.11891  -2.566 0.013607 *  
X72          0.17224    0.11238   1.533 0.132204    
X74         -0.23912    0.11049  -2.164 0.035683 *  
X75          0.25758    0.11616   2.217 0.031570 *  
X76         -0.39372    0.11185  -3.520 0.000985 ***
X79         -0.14147    0.11483  -1.232 0.224221    
X80         -0.37070    0.11865  -3.124 0.003082 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.2132 on 46 degrees of freedom
Multiple R-squared:  0.7594,    Adjusted R-squared:  0.4821 
F-statistic: 2.739 on 53 and 46 DF,  p-value: 0.0003334

The result is good and bad. Good is that we now get the effect of y ~ x significant. Bad is that a lot of the other noisy variables are also significant, and the rate at which this occurs is much higher than we should expect (22 of the 52 retained noise predictors are significant, while about 4 would be expected for 80 noise predictors at alpha = 0.05).

The phenomenon is well-known in the stats literature, and the reason is that performing a stepwise / global selection + calculating a regression table for the selected model is hidden multiple testing and has inflated Type I error rates! Remember, you are implicitly trying out hundreds or thousand of models, and are taking the one that is showing the strongest effects.

There are methods to correct for the problem (keyword: post-selection inference), e.g. the selectiveInference package for lasso and forward stepwise selection, but not for stepAIC or dredge, and mostly those corrected p-values have lower power than the full model. The simplest valid strategy is to split the data: select the model on one half, test it on the other.

8.4 Cross-validation

Cross-validation estimates how well a model predicts new data. The data are split into \(K\) parts (folds); the model is fitted to \(K-1\) folds and predicts the left-out fold, and this is repeated until each fold was left out once. Comparing these out-of-sample predictions with the observed values gives the prediction error on new data.

Compare this to the R² (Section 2.2): it measures how well the model fits the data it was fitted to, and it always increases when we add predictors. If we calculate an R² from the cross-validated predictions instead, we get an estimate of how well the model would fit new data (it can even become negative if the model predicts worse than the mean). The gap between the two shows how much the model overfits. Like AIC (Section 8.2.2), which approximates leave-one-out cross-validation, cross-validation can therefore be used to compare models, but it does not rely on a likelihood or on asymptotic approximations.

For a glm, the cv.glm() function from the boot package does the work:

library(boot)
data(mammals, package = "MASS")
mammals.glm = glm(log(brain) ~ log(body), data = mammals)

set.seed(123)
cv.err = cv.glm(mammals, mammals.glm, K = 5)$delta
cv.err
[1] 0.5079227 0.5034182
y = log(mammals$brain)
1 - mean(residuals(mammals.glm)^2) / mean((y - mean(y))^2)  # R2 on the fitted data
[1] 0.9207837
1 - cv.err[1] / mean((y - mean(y))^2)                       # cross-validated R2
[1] 0.9137489

delta contains the cross-validated mean squared error (the second value is adjusted for the fact that each fit uses only \((K-1)/K\) of the data; the two values are close here). The R² on the fitted data is 0.92, the cross-validated R² 0.91: with one predictor, the model hardly overfits. With many predictors, as in the simulation in Section 8.1, the gap becomes large.

Choosing a performance metric. The squared error (or the R² derived from it) is a natural metric for continuous responses, but not for binary ones. A simple metric for binary data is the misclassification rate: each predicted probability is turned into 0 or 1 with a threshold of 0.5, and we count how often this is wrong. This depends on the arbitrary threshold, and with rare events, a model that always predicts 0 already looks good. The AUC, which we used in the elephant example (Section 9.7.4), avoids both problems because it summarizes all thresholds. (cv.glm() accepts such metrics via its cost argument, but calculates them separately for each fold and then averages.)

For metrics like the AUC, it is better to collect the out-of-sample predictions of all folds first and to calculate the metric once on all of them: with few observations per fold, the AUC of a single fold is very noisy, and with leave-one-out it is not even defined. For most models, you will have to program this by hand, which is only a few lines. As an example, a logistic regression for the nodal data from the boot package (lymph node involvement of 53 prostate cancer patients):

library(pROC)
nodal.glm = glm(r ~ stage + xray + acid, family = binomial, data = nodal)

cvAUC = function(K = 5){
  fold = sample(rep(1:K, length.out = nrow(nodal)))  # random assignment to folds
  predCV = numeric(nrow(nodal))
  for (k in 1:K) {
    fit = glm(r ~ stage + xray + acid, family = binomial, data = nodal[fold != k, ])
    predCV[fold == k] = predict(fit, newdata = nodal[fold == k, ], type = "response")
  }
  as.numeric(auc(nodal$r, predCV, quiet = TRUE))
}

auc(nodal$r, predict(nodal.glm, type = "response"), quiet = TRUE) # on the fitted data
Area under the curve: 0.8341
set.seed(123)
cvAUC()                                                           # cross-validated
[1] 0.7643939
summary(replicate(20, cvAUC()))                                   # 20 different splits
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.7417  0.7625  0.7742  0.7725  0.7850  0.7909

On the fitted data, the AUC is 0.83; on new data, we expect about 0.77 (0.74 to 0.79 across splits). The last line shows that the result depends on the random split into folds. It is therefore good practice to repeat the cross-validation with several splits and to report the mean.

Note that neither the spread over repetitions nor the spread of the metric over the single folds is a confidence interval for the performance: the training sets of different folds share most of their data, so these estimates are not independent, and their spread underestimates the true uncertainty (Bengio and Grandvalet, 2004; Bates, Hastie and Tibshirani, 2024). Valid intervals require more elaborate methods such as nested cross-validation (Bates, Hastie and Tibshirani, 2024), which are beyond the scope of this book.

Bengio, Y. and Grandvalet, Y. (2004) “No unbiased estimator of the variance of k-fold cross-validation,” Journal of Machine Learning Research, 5, pp. 1089–1105. Available at: https://www.jmlr.org/papers/v5/grandvalet04a.html.
Bates, S., Hastie, T. and Tibshirani, R. (2024) “Cross-validation: What does it estimate and how well does it do it?” Journal of the American Statistical Association, 119(546), pp. 1434–1445.
Roberts, D.R. et al. (2017) “Cross-validation strategies for data with temporal, spatial, hierarchical, or phylogenetic structure,” Ecography, 40(8), pp. 913–929.

Note that cross-validation requires independence of data points. For non-independent data, it is possible to block the cross-validation, see Roberts et al. (2017), methods implemented in package blockCV, see https://cran.r-project.org/web/packages/blockCV/vignettes/BlockCV_for_SDM.html.

8.5 Take-home messages

Tip
  1. Do not fit the most complex model you can imagine: added variables usually increase the variance of the estimates (bias-variance trade-off). The exception are precision variables that explain much of the response and thereby reduce the uncertainty (Tip 8.1).
  2. Choose the method by the question. An LRT tests a pre-specified comparison of nested models with error control; AIC ranks models (also non-nested) by expected predictive performance and is more liberal than an LRT at alpha = 0.05; shrinkage and cross-validation serve prediction.
  3. If you want to estimate an effect, choose the variables from causal reasoning (Chapter 7), not by model selection: it can silently change the meaning of the remaining estimates.
  4. Do not interpret p-values and CIs of a model that was selected on the same data: stepwise or global selection is hidden multiple testing and inflates the type I error. Fix your hypothesis before looking at the data to avoid p-hacking.
  5. Evaluate predictive models on data that were not used for fitting (cross-validation); repeat it with several random splits, and do not treat the spread over folds as a confidence interval. For non-independent data (e.g. repeated measures, space), block the cross-validation, otherwise the error is too optimistic.

8.6 Case studies

8.6.1 Exercise — Global plant trait analysis #3

CautionExercise

Revisit our previous analyses of the dataset plantHeight (linear regression: Caution 2.2, Caution 2.3; ANOVA: Section 3.5.1; random effects: Caution 4.1), and discuss / analyze:

Which would be the appropriate model, if we want to get a predictive model for plant height, based on all the variables in the data set? Note: some text-based variables may need to be included, so probably it’s the easiest if you start with a large model that you specify by hand. You can also include interactions. The syntax:

fit <- lm(y ~ (x1 + x2 + x3)^2, data = dat)

includes all possible 2nd-order interactions between the variables in your model. You can extend this to x3, x4 but I would not recommend it, your model will get too large.

Note that we use all data here. In the earlier exercises (Caution 2.2), we removed the two growth forms with only one observation (Fern and Herb/Shrub). For model selection this does not matter, because we are interested in prediction and not in the contrasts between growth forms.

library(EcoData)
library(MASS)
head(plantHeight)
##   sort_number site         Genus_species        Family growthform height
## 1        1402  193    Acer_macryophyllum   Sapindaceae       Tree   28.0
## 2       25246  103    Quararibea_cordata     Malvaceae       Tree   26.6
## 3       11648   54    Eragrostis_dielsii       Poaceae       Herb    0.3
## 4        8168  144    Cistus_salvifolius     Cistaceae      Shrub    1.6
## 5       22422  178          Phlox_bifida Polemoniaceae       Herb    0.2
## 6       15925   59 Homalium_betulifolium    Salicaceae      Shrub    1.7
##        loght       Country              Site    lat     long entered.by alt
## 1  1.4471580           USA    Oregon - McDun 44.600 -123.334     Angela 179
## 2  1.4248816          Peru              Manu 12.183  -70.550     Angela 386
## 3 -0.5228787     Australia Central Australia 23.800  133.833   Michelle 553
## 4  0.2041200        Israel           Hanadiv 32.555   34.938     Angela 115
## 5 -0.6989700           USA     Indiana Dunes 41.617  -86.950   Michelle 200
## 6  0.2304489 New Caledonia              <NA> 21.500  165.500      Laura  95
##   temp diurn.temp isotherm temp.seas temp.max.warm temp.min.cold temp.ann.range
## 1 10.8       11.8      4.4       5.2          27.0           0.3           26.7
## 2 24.5       10.8      7.4       0.9          31.2          16.7           14.5
## 3 20.9       16.3      4.8       6.0          37.0           3.6           33.4
## 4 19.9        9.7      4.4       4.9          30.7           8.7           22.0
## 5  9.7       10.7      2.8       9.7          28.6          -9.5           38.1
## 6 22.6        7.4      5.4       2.2          29.0          15.5           13.5
##   temp.mean.wetqr temp.mean.dryqr temp.mean.warmqr temp.mean.coldqr rain
## 1             4.9            17.4             17.6              4.5 1208
## 2            25.1            23.2             25.3             23.1 3015
## 3            28.1            14.8             28.1             12.8  278
## 4            13.6            25.3             25.7             13.6  598
## 5            21.6            -3.3             21.6             -3.3  976
## 6            25.4            20.4             25.4             19.7 1387
##   rain.wetm rain.drym rain.seas rain.wetqr rain.dryqr rain.warmqr rain.coldqr
## 1       217        13        69        601         68          75         560
## 2       416        99        45       1177        340         928         359
## 3        37         9        42        109         35         109          42
## 4       159         0       115        408          0           2         408
## 5       104        44        23        299        165         299         165
## 6       216        59        46        600        186         600         212
##    LAI  NPP hemisphere
## 1 2.51  572          1
## 2 4.26 1405         -1
## 3 1.32  756         -1
## 4 1.01  359          1
## 5 3.26 1131          1
## 6 6.99 1552         -1
fullModel <- lm(loght ~ (growthform + Family + lat + long + alt + temp + NPP )^2, data = plantHeight)
selectedModel = stepAIC(fullModel, trace = 0)
summary(selectedModel)
## 
## Call:
## lm(formula = loght ~ growthform + Family + lat + long + alt + 
##     temp + NPP + growthform:Family + growthform:lat + growthform:long + 
##     growthform:NPP + Family:lat + Family:long + Family:alt + 
##     Family:temp + lat:long + lat:temp + lat:NPP + long:alt + 
##     long:temp + long:NPP + alt:temp + temp:NPP, data = plantHeight)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.1897  0.0000  0.0000  0.0000  0.1998 
## 
## Coefficients: (521 not defined because of singularities)
##                                                 Estimate Std. Error t value
## (Intercept)                                   -1.612e+01  2.152e+01  -0.749
## growthformHerb                                 1.598e+02  1.822e+02   0.877
## growthformHerb/Shrub                          -2.013e-01  3.854e+00  -0.052
## growthformShrub                               -1.992e+02  2.590e+02  -0.769
## growthformShrub/Tree                          -1.126e+02  1.952e+02  -0.577
## growthformTree                                 2.053e+01  1.530e+01   1.342
## FamilyAsteraceae                              -3.562e+02  1.352e+02  -2.636
## FamilyAtherospermataceae                      -2.866e+01  2.295e+01  -1.248
## FamilyBalsaminaceae                            1.530e+02  1.927e+02   0.794
## FamilyBetulaceae                              -3.010e+01  1.056e+02  -0.285
## FamilyBrassicaceae                            -1.271e+02  1.910e+02  -0.665
## FamilyCactaceae                                2.620e+02  2.820e+02   0.929
## FamilyCasuarinaceae                            1.995e+02  2.582e+02   0.773
## FamilyChloranthaceae                           7.419e-01  9.460e-01   0.784
## FamilyChrysobalanaceae                         5.392e+00  8.314e+00   0.649
## FamilyCistaceae                                1.999e+02  2.586e+02   0.773
## FamilyCornaceae                                3.361e+02  4.220e+02   0.796
## FamilyCrassulaceae                             1.954e+02  2.551e+02   0.766
## FamilyCunoniaceae                              1.716e-01  3.008e-01   0.571
## FamilyCupressaceae                            -3.435e+01  3.185e+01  -1.079
## FamilyCyperaceae                              -1.814e+01  2.054e+01  -0.883
## FamilyDennstaedtiaceae                                NA         NA      NA
## FamilyDicksoniaceae                            1.878e+02  2.490e+02   0.754
## FamilyDipterocarpaceae                         4.334e+00  4.919e+00   0.881
## FamilyEbenaceae                                5.281e+00  5.975e+00   0.884
## FamilyElaeocarpaceae                           5.440e-01  5.242e-01   1.038
## FamilyEricaceae                                2.186e+02  2.792e+02   0.783
## FamilyEuphorbiaceae                            3.472e+00  4.301e+00   0.807
## FamilyFabaceae - C                             2.644e+00  3.387e+00   0.781
## FamilyFabaceae - M                            -2.768e+00  2.238e+00  -1.237
## FamilyFabaceae - P                             5.128e+01  5.128e+01   1.000
## FamilyFagaceae                                -2.036e+00  7.802e+00  -0.261
## FamilyGentianaceae                            -1.963e+02  2.264e+02  -0.867
## FamilyHeliconiaceae                           -1.193e+02  1.413e+02  -0.845
## FamilyJuglandaceae                            -2.223e+01  2.069e+01  -1.074
## FamilyJuncaginaceae                            2.855e+02  3.596e+02   0.794
## FamilyLamiaceae                                1.876e+02  2.522e+02   0.744
## FamilyLauraceae                               -2.595e+00  4.184e+00  -0.620
## FamilyMaesaceae                                2.097e+02  2.652e+02   0.791
## FamilyMalvaceae                                3.153e+00  4.933e+00   0.639
## FamilyMelastomataceae                          2.135e+02  2.701e+02   0.790
## FamilyMoraceae                                -5.805e+00  6.329e+00  -0.917
## FamilyMyristicaceae                            8.819e-01  8.987e-01   0.981
## FamilyMyrsinaceae                             -9.893e+00  1.109e+01  -0.892
## FamilyMyrtaceae                               -3.193e+01  1.040e+01  -3.070
## FamilyOchnaceae                                2.125e+02  2.692e+02   0.789
## FamilyOnagraceae                                      NA         NA      NA
## FamilyOrobanchaceae                           -1.075e+02  1.752e+02  -0.613
## FamilyPhyllanthaceae                          -1.211e+01  1.357e+01  -0.892
## FamilyPicrodendraceae                         -1.256e+01  1.307e+01  -0.961
## FamilyPinaceae                                 1.522e+01  8.029e+00   1.895
## FamilyPoaceae                                  1.059e+02  1.625e+02   0.652
## FamilyPolemoniaceae                           -1.369e+02  1.613e+02  -0.849
## FamilyPolygonaceae                            -1.292e+02  1.622e+02  -0.796
## FamilyProteaceae                               1.685e+01  9.511e+00   1.772
## FamilyRanunculaceae                            1.136e+02  1.438e+02   0.790
## FamilyRhamnaceae                               1.944e+02  2.570e+02   0.756
## FamilyRosaceae                                 1.148e+01  3.983e+00   2.883
## FamilyRubiaceae                                2.165e+02  2.702e+02   0.801
## FamilyRutaceae                                 1.979e+02  2.555e+02   0.774
## FamilySalicaceae                               2.326e+02  2.849e+02   0.816
## FamilySapindaceae                              4.753e+00  8.974e+00   0.530
## FamilySapotaceae                              -1.511e+01  1.313e+01  -1.151
## FamilyScrophulariaceae                         2.099e+02  2.660e+02   0.789
## FamilyThymelaeaceae                           -9.324e+00  9.227e+00  -1.010
## FamilyUlmaceae                                 1.783e+01  1.401e+01   1.273
## FamilyUrticaceae                               8.939e+00  8.901e+00   1.004
## FamilyViolaceae                                2.145e+02  2.712e+02   0.791
## FamilyXanthorrhoeaceae                         2.016e+02  2.596e+02   0.777
## FamilyZygophyllaceae                           1.871e+02  2.474e+02   0.756
## lat                                            6.678e-01  6.451e-01   1.035
## long                                           6.846e-03  3.808e-02   0.180
## alt                                            5.440e-03  8.816e-03   0.617
## temp                                          -1.851e-01  2.343e-01  -0.790
## NPP                                            2.978e-03  2.716e-03   1.097
## growthformHerb:FamilyAsteraceae               -1.535e+01  4.924e+01  -0.312
## growthformHerb/Shrub:FamilyAsteraceae                 NA         NA      NA
## growthformShrub:FamilyAsteraceae                      NA         NA      NA
## growthformShrub/Tree:FamilyAsteraceae                 NA         NA      NA
## growthformTree:FamilyAsteraceae                       NA         NA      NA
## growthformHerb:FamilyAtherospermataceae               NA         NA      NA
## growthformHerb/Shrub:FamilyAtherospermataceae         NA         NA      NA
## growthformShrub:FamilyAtherospermataceae              NA         NA      NA
## growthformShrub/Tree:FamilyAtherospermataceae         NA         NA      NA
## growthformTree:FamilyAtherospermataceae               NA         NA      NA
## growthformHerb:FamilyBalsaminaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyBalsaminaceae              NA         NA      NA
## growthformShrub:FamilyBalsaminaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyBalsaminaceae              NA         NA      NA
## growthformTree:FamilyBalsaminaceae                    NA         NA      NA
## growthformHerb:FamilyBetulaceae                       NA         NA      NA
## growthformHerb/Shrub:FamilyBetulaceae                 NA         NA      NA
## growthformShrub:FamilyBetulaceae               2.035e+02  2.795e+02   0.728
## growthformShrub/Tree:FamilyBetulaceae                 NA         NA      NA
## growthformTree:FamilyBetulaceae                       NA         NA      NA
## growthformHerb:FamilyBrassicaceae                     NA         NA      NA
## growthformHerb/Shrub:FamilyBrassicaceae               NA         NA      NA
## growthformShrub:FamilyBrassicaceae                    NA         NA      NA
## growthformShrub/Tree:FamilyBrassicaceae               NA         NA      NA
## growthformTree:FamilyBrassicaceae                     NA         NA      NA
## growthformHerb:FamilyCactaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyCactaceae                  NA         NA      NA
## growthformShrub:FamilyCactaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyCactaceae                  NA         NA      NA
## growthformTree:FamilyCactaceae                        NA         NA      NA
## growthformHerb:FamilyCasuarinaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyCasuarinaceae              NA         NA      NA
## growthformShrub:FamilyCasuarinaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyCasuarinaceae              NA         NA      NA
## growthformTree:FamilyCasuarinaceae                    NA         NA      NA
## growthformHerb:FamilyChloranthaceae                   NA         NA      NA
## growthformHerb/Shrub:FamilyChloranthaceae             NA         NA      NA
## growthformShrub:FamilyChloranthaceae                  NA         NA      NA
## growthformShrub/Tree:FamilyChloranthaceae             NA         NA      NA
## growthformTree:FamilyChloranthaceae                   NA         NA      NA
## growthformHerb:FamilyChrysobalanaceae                 NA         NA      NA
## growthformHerb/Shrub:FamilyChrysobalanaceae           NA         NA      NA
## growthformShrub:FamilyChrysobalanaceae                NA         NA      NA
## growthformShrub/Tree:FamilyChrysobalanaceae           NA         NA      NA
## growthformTree:FamilyChrysobalanaceae                 NA         NA      NA
## growthformHerb:FamilyCistaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyCistaceae                  NA         NA      NA
## growthformShrub:FamilyCistaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyCistaceae                  NA         NA      NA
## growthformTree:FamilyCistaceae                        NA         NA      NA
## growthformHerb:FamilyCornaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyCornaceae                  NA         NA      NA
## growthformShrub:FamilyCornaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyCornaceae                  NA         NA      NA
## growthformTree:FamilyCornaceae                        NA         NA      NA
## growthformHerb:FamilyCrassulaceae                     NA         NA      NA
## growthformHerb/Shrub:FamilyCrassulaceae               NA         NA      NA
## growthformShrub:FamilyCrassulaceae                    NA         NA      NA
## growthformShrub/Tree:FamilyCrassulaceae               NA         NA      NA
## growthformTree:FamilyCrassulaceae                     NA         NA      NA
## growthformHerb:FamilyCunoniaceae                      NA         NA      NA
## growthformHerb/Shrub:FamilyCunoniaceae                NA         NA      NA
## growthformShrub:FamilyCunoniaceae                     NA         NA      NA
## growthformShrub/Tree:FamilyCunoniaceae         2.124e+02  2.694e+02   0.789
## growthformTree:FamilyCunoniaceae                      NA         NA      NA
## growthformHerb:FamilyCupressaceae                     NA         NA      NA
## growthformHerb/Shrub:FamilyCupressaceae               NA         NA      NA
## growthformShrub:FamilyCupressaceae             2.214e+02  2.802e+02   0.790
## growthformShrub/Tree:FamilyCupressaceae               NA         NA      NA
## growthformTree:FamilyCupressaceae                     NA         NA      NA
## growthformHerb:FamilyCyperaceae                       NA         NA      NA
## growthformHerb/Shrub:FamilyCyperaceae                 NA         NA      NA
## growthformShrub:FamilyCyperaceae                      NA         NA      NA
## growthformShrub/Tree:FamilyCyperaceae                 NA         NA      NA
## growthformTree:FamilyCyperaceae                       NA         NA      NA
## growthformHerb:FamilyDennstaedtiaceae                 NA         NA      NA
## growthformHerb/Shrub:FamilyDennstaedtiaceae           NA         NA      NA
## growthformShrub:FamilyDennstaedtiaceae                NA         NA      NA
## growthformShrub/Tree:FamilyDennstaedtiaceae           NA         NA      NA
## growthformTree:FamilyDennstaedtiaceae                 NA         NA      NA
## growthformHerb:FamilyDicksoniaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyDicksoniaceae              NA         NA      NA
## growthformShrub:FamilyDicksoniaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyDicksoniaceae              NA         NA      NA
## growthformTree:FamilyDicksoniaceae                    NA         NA      NA
## growthformHerb:FamilyDipterocarpaceae                 NA         NA      NA
## growthformHerb/Shrub:FamilyDipterocarpaceae           NA         NA      NA
## growthformShrub:FamilyDipterocarpaceae                NA         NA      NA
## growthformShrub/Tree:FamilyDipterocarpaceae           NA         NA      NA
## growthformTree:FamilyDipterocarpaceae                 NA         NA      NA
## growthformHerb:FamilyEbenaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyEbenaceae                  NA         NA      NA
## growthformShrub:FamilyEbenaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyEbenaceae                  NA         NA      NA
## growthformTree:FamilyEbenaceae                        NA         NA      NA
## growthformHerb:FamilyElaeocarpaceae                   NA         NA      NA
## growthformHerb/Shrub:FamilyElaeocarpaceae             NA         NA      NA
## growthformShrub:FamilyElaeocarpaceae                  NA         NA      NA
## growthformShrub/Tree:FamilyElaeocarpaceae             NA         NA      NA
## growthformTree:FamilyElaeocarpaceae                   NA         NA      NA
## growthformHerb:FamilyEricaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyEricaceae                  NA         NA      NA
## growthformShrub:FamilyEricaceae               -1.611e+00  4.441e+00  -0.363
## growthformShrub/Tree:FamilyEricaceae                  NA         NA      NA
## growthformTree:FamilyEricaceae                        NA         NA      NA
## growthformHerb:FamilyEuphorbiaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyEuphorbiaceae              NA         NA      NA
## growthformShrub:FamilyEuphorbiaceae            2.201e+02  2.766e+02   0.796
## growthformShrub/Tree:FamilyEuphorbiaceae              NA         NA      NA
## growthformTree:FamilyEuphorbiaceae                    NA         NA      NA
## growthformHerb:FamilyFabaceae - C                     NA         NA      NA
## growthformHerb/Shrub:FamilyFabaceae - C               NA         NA      NA
## growthformShrub:FamilyFabaceae - C                    NA         NA      NA
## growthformShrub/Tree:FamilyFabaceae - C               NA         NA      NA
## growthformTree:FamilyFabaceae - C                     NA         NA      NA
## growthformHerb:FamilyFabaceae - M                     NA         NA      NA
## growthformHerb/Shrub:FamilyFabaceae - M               NA         NA      NA
## growthformShrub:FamilyFabaceae - M                    NA         NA      NA
## growthformShrub/Tree:FamilyFabaceae - M        2.338e+02  2.857e+02   0.818
## growthformTree:FamilyFabaceae - M                     NA         NA      NA
## growthformHerb:FamilyFabaceae - P              1.161e+02  1.357e+02   0.855
## growthformHerb/Shrub:FamilyFabaceae - P               NA         NA      NA
## growthformShrub:FamilyFabaceae - P             2.457e+02  3.020e+02   0.814
## growthformShrub/Tree:FamilyFabaceae - P               NA         NA      NA
## growthformTree:FamilyFabaceae - P                     NA         NA      NA
## growthformHerb:FamilyFagaceae                         NA         NA      NA
## growthformHerb/Shrub:FamilyFagaceae                   NA         NA      NA
## growthformShrub:FamilyFagaceae                 2.144e+02  2.735e+02   0.784
## growthformShrub/Tree:FamilyFagaceae                   NA         NA      NA
## growthformTree:FamilyFagaceae                         NA         NA      NA
## growthformHerb:FamilyGentianaceae                     NA         NA      NA
## growthformHerb/Shrub:FamilyGentianaceae               NA         NA      NA
## growthformShrub:FamilyGentianaceae                    NA         NA      NA
## growthformShrub/Tree:FamilyGentianaceae               NA         NA      NA
## growthformTree:FamilyGentianaceae                     NA         NA      NA
## growthformHerb:FamilyHeliconiaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyHeliconiaceae              NA         NA      NA
## growthformShrub:FamilyHeliconiaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyHeliconiaceae              NA         NA      NA
## growthformTree:FamilyHeliconiaceae                    NA         NA      NA
## growthformHerb:FamilyJuglandaceae                     NA         NA      NA
## growthformHerb/Shrub:FamilyJuglandaceae               NA         NA      NA
## growthformShrub:FamilyJuglandaceae                    NA         NA      NA
## growthformShrub/Tree:FamilyJuglandaceae               NA         NA      NA
## growthformTree:FamilyJuglandaceae                     NA         NA      NA
## growthformHerb:FamilyJuncaginaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyJuncaginaceae              NA         NA      NA
## growthformShrub:FamilyJuncaginaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyJuncaginaceae              NA         NA      NA
## growthformTree:FamilyJuncaginaceae                    NA         NA      NA
## growthformHerb:FamilyLamiaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyLamiaceae                  NA         NA      NA
## growthformShrub:FamilyLamiaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyLamiaceae                  NA         NA      NA
## growthformTree:FamilyLamiaceae                        NA         NA      NA
## growthformHerb:FamilyLauraceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyLauraceae                  NA         NA      NA
## growthformShrub:FamilyLauraceae                       NA         NA      NA
## growthformShrub/Tree:FamilyLauraceae                  NA         NA      NA
## growthformTree:FamilyLauraceae                        NA         NA      NA
## growthformHerb:FamilyMaesaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyMaesaceae                  NA         NA      NA
## growthformShrub:FamilyMaesaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyMaesaceae                  NA         NA      NA
## growthformTree:FamilyMaesaceae                        NA         NA      NA
## growthformHerb:FamilyMalvaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyMalvaceae                  NA         NA      NA
## growthformShrub:FamilyMalvaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyMalvaceae                  NA         NA      NA
## growthformTree:FamilyMalvaceae                        NA         NA      NA
## growthformHerb:FamilyMelastomataceae                  NA         NA      NA
## growthformHerb/Shrub:FamilyMelastomataceae            NA         NA      NA
## growthformShrub:FamilyMelastomataceae                 NA         NA      NA
## growthformShrub/Tree:FamilyMelastomataceae            NA         NA      NA
## growthformTree:FamilyMelastomataceae                  NA         NA      NA
## growthformHerb:FamilyMoraceae                         NA         NA      NA
## growthformHerb/Shrub:FamilyMoraceae                   NA         NA      NA
## growthformShrub:FamilyMoraceae                        NA         NA      NA
## growthformShrub/Tree:FamilyMoraceae                   NA         NA      NA
## growthformTree:FamilyMoraceae                         NA         NA      NA
## growthformHerb:FamilyMyristicaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyMyristicaceae              NA         NA      NA
## growthformShrub:FamilyMyristicaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyMyristicaceae              NA         NA      NA
## growthformTree:FamilyMyristicaceae                    NA         NA      NA
## growthformHerb:FamilyMyrsinaceae                      NA         NA      NA
## growthformHerb/Shrub:FamilyMyrsinaceae                NA         NA      NA
## growthformShrub:FamilyMyrsinaceae                     NA         NA      NA
## growthformShrub/Tree:FamilyMyrsinaceae         2.252e+02  2.829e+02   0.796
## growthformTree:FamilyMyrsinaceae                      NA         NA      NA
## growthformHerb:FamilyMyrtaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyMyrtaceae                  NA         NA      NA
## growthformShrub:FamilyMyrtaceae                2.171e+02  2.710e+02   0.801
## growthformShrub/Tree:FamilyMyrtaceae                  NA         NA      NA
## growthformTree:FamilyMyrtaceae                        NA         NA      NA
## growthformHerb:FamilyOchnaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyOchnaceae                  NA         NA      NA
## growthformShrub:FamilyOchnaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyOchnaceae                  NA         NA      NA
## growthformTree:FamilyOchnaceae                        NA         NA      NA
## growthformHerb:FamilyOnagraceae                       NA         NA      NA
## growthformHerb/Shrub:FamilyOnagraceae                 NA         NA      NA
## growthformShrub:FamilyOnagraceae                      NA         NA      NA
## growthformShrub/Tree:FamilyOnagraceae                 NA         NA      NA
## growthformTree:FamilyOnagraceae                       NA         NA      NA
## growthformHerb:FamilyOrobanchaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyOrobanchaceae              NA         NA      NA
## growthformShrub:FamilyOrobanchaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyOrobanchaceae              NA         NA      NA
## growthformTree:FamilyOrobanchaceae                    NA         NA      NA
## growthformHerb:FamilyPhyllanthaceae                   NA         NA      NA
## growthformHerb/Shrub:FamilyPhyllanthaceae             NA         NA      NA
## growthformShrub:FamilyPhyllanthaceae                  NA         NA      NA
## growthformShrub/Tree:FamilyPhyllanthaceae             NA         NA      NA
## growthformTree:FamilyPhyllanthaceae                   NA         NA      NA
## growthformHerb:FamilyPicrodendraceae                  NA         NA      NA
## growthformHerb/Shrub:FamilyPicrodendraceae            NA         NA      NA
## growthformShrub:FamilyPicrodendraceae                 NA         NA      NA
## growthformShrub/Tree:FamilyPicrodendraceae            NA         NA      NA
## growthformTree:FamilyPicrodendraceae                  NA         NA      NA
## growthformHerb:FamilyPinaceae                         NA         NA      NA
## growthformHerb/Shrub:FamilyPinaceae                   NA         NA      NA
## growthformShrub:FamilyPinaceae                        NA         NA      NA
## growthformShrub/Tree:FamilyPinaceae                   NA         NA      NA
## growthformTree:FamilyPinaceae                         NA         NA      NA
## growthformHerb:FamilyPoaceae                  -2.601e+02  3.246e+02  -0.801
## growthformHerb/Shrub:FamilyPoaceae                    NA         NA      NA
## growthformShrub:FamilyPoaceae                         NA         NA      NA
## growthformShrub/Tree:FamilyPoaceae                    NA         NA      NA
## growthformTree:FamilyPoaceae                          NA         NA      NA
## growthformHerb:FamilyPolemoniaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyPolemoniaceae              NA         NA      NA
## growthformShrub:FamilyPolemoniaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyPolemoniaceae              NA         NA      NA
## growthformTree:FamilyPolemoniaceae                    NA         NA      NA
## growthformHerb:FamilyPolygonaceae                     NA         NA      NA
## growthformHerb/Shrub:FamilyPolygonaceae               NA         NA      NA
## growthformShrub:FamilyPolygonaceae                    NA         NA      NA
## growthformShrub/Tree:FamilyPolygonaceae               NA         NA      NA
## growthformTree:FamilyPolygonaceae                     NA         NA      NA
## growthformHerb:FamilyProteaceae                       NA         NA      NA
## growthformHerb/Shrub:FamilyProteaceae                 NA         NA      NA
## growthformShrub:FamilyProteaceae               2.174e+02  2.714e+02   0.801
## growthformShrub/Tree:FamilyProteaceae                 NA         NA      NA
## growthformTree:FamilyProteaceae                       NA         NA      NA
## growthformHerb:FamilyRanunculaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyRanunculaceae              NA         NA      NA
## growthformShrub:FamilyRanunculaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyRanunculaceae              NA         NA      NA
## growthformTree:FamilyRanunculaceae                    NA         NA      NA
## growthformHerb:FamilyRhamnaceae                       NA         NA      NA
## growthformHerb/Shrub:FamilyRhamnaceae                 NA         NA      NA
## growthformShrub:FamilyRhamnaceae                      NA         NA      NA
## growthformShrub/Tree:FamilyRhamnaceae                 NA         NA      NA
## growthformTree:FamilyRhamnaceae                       NA         NA      NA
## growthformHerb:FamilyRosaceae                  3.723e+02  4.458e+02   0.835
## growthformHerb/Shrub:FamilyRosaceae                   NA         NA      NA
## growthformShrub:FamilyRosaceae                 2.126e+02  2.715e+02   0.783
## growthformShrub/Tree:FamilyRosaceae                   NA         NA      NA
## growthformTree:FamilyRosaceae                         NA         NA      NA
## growthformHerb:FamilyRubiaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyRubiaceae                  NA         NA      NA
## growthformShrub:FamilyRubiaceae                6.243e+00  3.379e+00   1.848
## growthformShrub/Tree:FamilyRubiaceae                  NA         NA      NA
## growthformTree:FamilyRubiaceae                        NA         NA      NA
## growthformHerb:FamilyRutaceae                         NA         NA      NA
## growthformHerb/Shrub:FamilyRutaceae                   NA         NA      NA
## growthformShrub:FamilyRutaceae                        NA         NA      NA
## growthformShrub/Tree:FamilyRutaceae                   NA         NA      NA
## growthformTree:FamilyRutaceae                         NA         NA      NA
## growthformHerb:FamilySalicaceae                       NA         NA      NA
## growthformHerb/Shrub:FamilySalicaceae                 NA         NA      NA
## growthformShrub:FamilySalicaceae               9.252e+00  6.665e+00   1.388
## growthformShrub/Tree:FamilySalicaceae                 NA         NA      NA
## growthformTree:FamilySalicaceae                       NA         NA      NA
## growthformHerb:FamilySapindaceae                      NA         NA      NA
## growthformHerb/Shrub:FamilySapindaceae                NA         NA      NA
## growthformShrub:FamilySapindaceae                     NA         NA      NA
## growthformShrub/Tree:FamilySapindaceae                NA         NA      NA
## growthformTree:FamilySapindaceae                      NA         NA      NA
## growthformHerb:FamilySapotaceae                       NA         NA      NA
## growthformHerb/Shrub:FamilySapotaceae                 NA         NA      NA
## growthformShrub:FamilySapotaceae                      NA         NA      NA
## growthformShrub/Tree:FamilySapotaceae                 NA         NA      NA
## growthformTree:FamilySapotaceae                       NA         NA      NA
## growthformHerb:FamilyScrophulariaceae                 NA         NA      NA
## growthformHerb/Shrub:FamilyScrophulariaceae           NA         NA      NA
## growthformShrub:FamilyScrophulariaceae                NA         NA      NA
## growthformShrub/Tree:FamilyScrophulariaceae           NA         NA      NA
## growthformTree:FamilyScrophulariaceae                 NA         NA      NA
## growthformHerb:FamilyThymelaeaceae                    NA         NA      NA
## growthformHerb/Shrub:FamilyThymelaeaceae              NA         NA      NA
## growthformShrub:FamilyThymelaeaceae                   NA         NA      NA
## growthformShrub/Tree:FamilyThymelaeaceae              NA         NA      NA
## growthformTree:FamilyThymelaeaceae                    NA         NA      NA
## growthformHerb:FamilyUlmaceae                         NA         NA      NA
## growthformHerb/Shrub:FamilyUlmaceae                   NA         NA      NA
## growthformShrub:FamilyUlmaceae                        NA         NA      NA
## growthformShrub/Tree:FamilyUlmaceae                   NA         NA      NA
## growthformTree:FamilyUlmaceae                         NA         NA      NA
## growthformHerb:FamilyUrticaceae                       NA         NA      NA
## growthformHerb/Shrub:FamilyUrticaceae                 NA         NA      NA
## growthformShrub:FamilyUrticaceae                      NA         NA      NA
## growthformShrub/Tree:FamilyUrticaceae                 NA         NA      NA
## growthformTree:FamilyUrticaceae                       NA         NA      NA
## growthformHerb:FamilyViolaceae                        NA         NA      NA
## growthformHerb/Shrub:FamilyViolaceae                  NA         NA      NA
## growthformShrub:FamilyViolaceae                       NA         NA      NA
## growthformShrub/Tree:FamilyViolaceae                  NA         NA      NA
## growthformTree:FamilyViolaceae                        NA         NA      NA
## growthformHerb:FamilyXanthorrhoeaceae                 NA         NA      NA
## growthformHerb/Shrub:FamilyXanthorrhoeaceae           NA         NA      NA
## growthformShrub:FamilyXanthorrhoeaceae                NA         NA      NA
## growthformShrub/Tree:FamilyXanthorrhoeaceae           NA         NA      NA
## growthformTree:FamilyXanthorrhoeaceae                 NA         NA      NA
## growthformHerb:FamilyZygophyllaceae                   NA         NA      NA
## growthformHerb/Shrub:FamilyZygophyllaceae             NA         NA      NA
## growthformShrub:FamilyZygophyllaceae                  NA         NA      NA
## growthformShrub/Tree:FamilyZygophyllaceae             NA         NA      NA
## growthformTree:FamilyZygophyllaceae                   NA         NA      NA
## growthformHerb:lat                            -7.667e+00  9.125e+00  -0.840
## growthformHerb/Shrub:lat                              NA         NA      NA
## growthformShrub:lat                            7.514e-02  3.299e-02   2.278
## growthformShrub/Tree:lat                      -4.602e+00  4.806e+00  -0.957
## growthformTree:lat                                    NA         NA      NA
## growthformHerb:long                            3.209e-02  3.417e-02   0.939
## growthformHerb/Shrub:long                             NA         NA      NA
## growthformShrub:long                          -6.040e-03  1.358e-02  -0.445
## growthformShrub/Tree:long                             NA         NA      NA
## growthformTree:long                                   NA         NA      NA
## growthformHerb:NPP                             4.933e-04  4.811e-04   1.025
## growthformHerb/Shrub:NPP                              NA         NA      NA
## growthformShrub:NPP                            6.473e-04  3.846e-04   1.683
## growthformShrub/Tree:NPP                              NA         NA      NA
## growthformTree:NPP                                    NA         NA      NA
## FamilyAsteraceae:lat                           9.729e+00  8.577e+00   1.134
## FamilyAtherospermataceae:lat                          NA         NA      NA
## FamilyBalsaminaceae:lat                               NA         NA      NA
## FamilyBetulaceae:lat                          -5.249e-01  1.875e+00  -0.280
## FamilyBrassicaceae:lat                         6.681e+00  9.160e+00   0.729
## FamilyCactaceae:lat                           -1.966e+00  8.672e-01  -2.267
## FamilyCasuarinaceae:lat                               NA         NA      NA
## FamilyChloranthaceae:lat                              NA         NA      NA
## FamilyChrysobalanaceae:lat                    -7.988e-01  1.363e+00  -0.586
## FamilyCistaceae:lat                                   NA         NA      NA
## FamilyCornaceae:lat                                   NA         NA      NA
## FamilyCrassulaceae:lat                                NA         NA      NA
## FamilyCunoniaceae:lat                                 NA         NA      NA
## FamilyCupressaceae:lat                                NA         NA      NA
## FamilyCyperaceae:lat                                  NA         NA      NA
## FamilyDennstaedtiaceae:lat                            NA         NA      NA
## FamilyDicksoniaceae:lat                               NA         NA      NA
## FamilyDipterocarpaceae:lat                            NA         NA      NA
## FamilyEbenaceae:lat                           -2.703e-01  3.089e-01  -0.875
## FamilyElaeocarpaceae:lat                              NA         NA      NA
## FamilyEricaceae:lat                           -8.632e-01  6.646e-01  -1.299
## FamilyEuphorbiaceae:lat                       -7.697e-01  7.758e-01  -0.992
## FamilyFabaceae - C:lat                        -6.439e-01  7.319e-01  -0.880
## FamilyFabaceae - M:lat                                NA         NA      NA
## FamilyFabaceae - P:lat                        -2.849e+00  2.838e+00  -1.004
## FamilyFagaceae:lat                            -6.391e-01  6.767e-01  -0.944
## FamilyGentianaceae:lat                         7.953e+00  9.600e+00   0.828
## FamilyHeliconiaceae:lat                               NA         NA      NA
## FamilyJuglandaceae:lat                                NA         NA      NA
## FamilyJuncaginaceae:lat                               NA         NA      NA
## FamilyLamiaceae:lat                                   NA         NA      NA
## FamilyLauraceae:lat                                   NA         NA      NA
## FamilyMaesaceae:lat                                   NA         NA      NA
## FamilyMalvaceae:lat                           -5.759e-01  7.352e-01  -0.783
## FamilyMelastomataceae:lat                             NA         NA      NA
## FamilyMoraceae:lat                                    NA         NA      NA
## FamilyMyristicaceae:lat                               NA         NA      NA
## FamilyMyrsinaceae:lat                                 NA         NA      NA
## FamilyMyrtaceae:lat                           -2.630e-01  6.091e-01  -0.432
## FamilyOchnaceae:lat                                   NA         NA      NA
## FamilyOnagraceae:lat                                  NA         NA      NA
## FamilyOrobanchaceae:lat                        6.419e+00  8.683e+00   0.739
## FamilyPhyllanthaceae:lat                              NA         NA      NA
## FamilyPicrodendraceae:lat                             NA         NA      NA
## FamilyPinaceae:lat                            -1.017e+00  6.698e-01  -1.518
## FamilyPoaceae:lat                              7.161e+00  8.523e+00   0.840
## FamilyPolemoniaceae:lat                        6.911e+00  8.461e+00   0.817
## FamilyPolygonaceae:lat                         6.776e+00  8.505e+00   0.797
## FamilyProteaceae:lat                          -1.078e+00  6.906e-01  -1.560
## FamilyRanunculaceae:lat                               NA         NA      NA
## FamilyRhamnaceae:lat                                  NA         NA      NA
## FamilyRosaceae:lat                            -8.875e-01  5.851e-01  -1.517
## FamilyRubiaceae:lat                           -1.008e+00  7.723e-01  -1.305
## FamilyRutaceae:lat                                    NA         NA      NA
## FamilySalicaceae:lat                          -1.327e+00  1.062e+00  -1.249
## FamilySapindaceae:lat                         -7.032e-01  7.657e-01  -0.918
## FamilySapotaceae:lat                                  NA         NA      NA
## FamilyScrophulariaceae:lat                            NA         NA      NA
## FamilyThymelaeaceae:lat                               NA         NA      NA
## FamilyUlmaceae:lat                            -1.490e+00  1.311e+00  -1.137
## FamilyUrticaceae:lat                          -1.170e+00  1.217e+00  -0.962
## FamilyViolaceae:lat                                   NA         NA      NA
## FamilyXanthorrhoeaceae:lat                            NA         NA      NA
## FamilyZygophyllaceae:lat                              NA         NA      NA
## FamilyAsteraceae:long                          3.953e-02  8.917e-03   4.433
## FamilyAtherospermataceae:long                         NA         NA      NA
## FamilyBalsaminaceae:long                              NA         NA      NA
## FamilyBetulaceae:long                          1.284e+00  1.382e+00   0.929
## FamilyBrassicaceae:long                               NA         NA      NA
## FamilyCactaceae:long                           3.726e-02  2.427e-02   1.535
## FamilyCasuarinaceae:long                              NA         NA      NA
## FamilyChloranthaceae:long                             NA         NA      NA
## FamilyChrysobalanaceae:long                           NA         NA      NA
## FamilyCistaceae:long                                  NA         NA      NA
## FamilyCornaceae:long                                  NA         NA      NA
## FamilyCrassulaceae:long                               NA         NA      NA
## FamilyCunoniaceae:long                                NA         NA      NA
## FamilyCupressaceae:long                               NA         NA      NA
## FamilyCyperaceae:long                                 NA         NA      NA
## FamilyDennstaedtiaceae:long                           NA         NA      NA
## FamilyDicksoniaceae:long                              NA         NA      NA
## FamilyDipterocarpaceae:long                           NA         NA      NA
## FamilyEbenaceae:long                                  NA         NA      NA
## FamilyElaeocarpaceae:long                             NA         NA      NA
## FamilyEricaceae:long                           6.409e-02  2.972e-02   2.156
## FamilyEuphorbiaceae:long                       4.371e-02  4.554e-02   0.960
## FamilyFabaceae - C:long                               NA         NA      NA
## FamilyFabaceae - M:long                               NA         NA      NA
## FamilyFabaceae - P:long                               NA         NA      NA
## FamilyFagaceae:long                                   NA         NA      NA
## FamilyGentianaceae:long                               NA         NA      NA
## FamilyHeliconiaceae:long                              NA         NA      NA
## FamilyJuglandaceae:long                               NA         NA      NA
## FamilyJuncaginaceae:long                              NA         NA      NA
## FamilyLamiaceae:long                                  NA         NA      NA
## FamilyLauraceae:long                                  NA         NA      NA
## FamilyMaesaceae:long                                  NA         NA      NA
## FamilyMalvaceae:long                           3.701e-02  4.705e-02   0.787
## FamilyMelastomataceae:long                            NA         NA      NA
## FamilyMoraceae:long                                   NA         NA      NA
## FamilyMyristicaceae:long                              NA         NA      NA
## FamilyMyrsinaceae:long                                NA         NA      NA
## FamilyMyrtaceae:long                           5.316e-02  3.216e-02   1.653
## FamilyOchnaceae:long                                  NA         NA      NA
## FamilyOnagraceae:long                                 NA         NA      NA
## FamilyOrobanchaceae:long                              NA         NA      NA
## FamilyPhyllanthaceae:long                             NA         NA      NA
## FamilyPicrodendraceae:long                            NA         NA      NA
## FamilyPinaceae:long                            1.006e-03  3.483e-02   0.029
## FamilyPoaceae:long                                    NA         NA      NA
## FamilyPolemoniaceae:long                              NA         NA      NA
## FamilyPolygonaceae:long                               NA         NA      NA
## FamilyProteaceae:long                          3.544e-02  2.338e-02   1.516
## FamilyRanunculaceae:long                              NA         NA      NA
## FamilyRhamnaceae:long                                 NA         NA      NA
## FamilyRosaceae:long                            6.083e-02  5.503e-02   1.105
## FamilyRubiaceae:long                                  NA         NA      NA
## FamilyRutaceae:long                                   NA         NA      NA
## FamilySalicaceae:long                                 NA         NA      NA
## FamilySapindaceae:long                         4.512e-02  2.639e-02   1.710
## FamilySapotaceae:long                                 NA         NA      NA
## FamilyScrophulariaceae:long                           NA         NA      NA
## FamilyThymelaeaceae:long                              NA         NA      NA
## FamilyUlmaceae:long                                   NA         NA      NA
## FamilyUrticaceae:long                                 NA         NA      NA
## FamilyViolaceae:long                                  NA         NA      NA
## FamilyXanthorrhoeaceae:long                           NA         NA      NA
## FamilyZygophyllaceae:long                             NA         NA      NA
## FamilyAsteraceae:alt                           2.977e-02  1.249e-02   2.384
## FamilyAtherospermataceae:alt                          NA         NA      NA
## FamilyBalsaminaceae:alt                               NA         NA      NA
## FamilyBetulaceae:alt                                  NA         NA      NA
## FamilyBrassicaceae:alt                                NA         NA      NA
## FamilyCactaceae:alt                           -4.211e-03  9.104e-03  -0.463
## FamilyCasuarinaceae:alt                               NA         NA      NA
## FamilyChloranthaceae:alt                              NA         NA      NA
## FamilyChrysobalanaceae:alt                            NA         NA      NA
## FamilyCistaceae:alt                                   NA         NA      NA
## FamilyCornaceae:alt                                   NA         NA      NA
## FamilyCrassulaceae:alt                                NA         NA      NA
## FamilyCunoniaceae:alt                                 NA         NA      NA
## FamilyCupressaceae:alt                                NA         NA      NA
## FamilyCyperaceae:alt                                  NA         NA      NA
## FamilyDennstaedtiaceae:alt                            NA         NA      NA
## FamilyDicksoniaceae:alt                               NA         NA      NA
## FamilyDipterocarpaceae:alt                            NA         NA      NA
## FamilyEbenaceae:alt                                   NA         NA      NA
## FamilyElaeocarpaceae:alt                              NA         NA      NA
## FamilyEricaceae:alt                            2.057e-03  9.530e-03   0.216
## FamilyEuphorbiaceae:alt                               NA         NA      NA
## FamilyFabaceae - C:alt                                NA         NA      NA
## FamilyFabaceae - M:alt                                NA         NA      NA
## FamilyFabaceae - P:alt                                NA         NA      NA
## FamilyFagaceae:alt                                    NA         NA      NA
## FamilyGentianaceae:alt                                NA         NA      NA
## FamilyHeliconiaceae:alt                               NA         NA      NA
## FamilyJuglandaceae:alt                                NA         NA      NA
## FamilyJuncaginaceae:alt                               NA         NA      NA
## FamilyLamiaceae:alt                                   NA         NA      NA
## FamilyLauraceae:alt                                   NA         NA      NA
## FamilyMaesaceae:alt                                   NA         NA      NA
## FamilyMalvaceae:alt                                   NA         NA      NA
## FamilyMelastomataceae:alt                             NA         NA      NA
## FamilyMoraceae:alt                                    NA         NA      NA
## FamilyMyristicaceae:alt                               NA         NA      NA
## FamilyMyrsinaceae:alt                                 NA         NA      NA
## FamilyMyrtaceae:alt                           -2.001e-03  8.486e-03  -0.236
## FamilyOchnaceae:alt                                   NA         NA      NA
## FamilyOnagraceae:alt                                  NA         NA      NA
## FamilyOrobanchaceae:alt                               NA         NA      NA
## FamilyPhyllanthaceae:alt                              NA         NA      NA
## FamilyPicrodendraceae:alt                             NA         NA      NA
## FamilyPinaceae:alt                            -7.229e-03  8.920e-03  -0.810
## FamilyPoaceae:alt                             -3.846e-03  8.857e-03  -0.434
## FamilyPolemoniaceae:alt                               NA         NA      NA
## FamilyPolygonaceae:alt                                NA         NA      NA
## FamilyProteaceae:alt                          -8.882e-03  8.984e-03  -0.989
## FamilyRanunculaceae:alt                               NA         NA      NA
## FamilyRhamnaceae:alt                                  NA         NA      NA
## FamilyRosaceae:alt                                    NA         NA      NA
## FamilyRubiaceae:alt                                   NA         NA      NA
## FamilyRutaceae:alt                                    NA         NA      NA
## FamilySalicaceae:alt                                  NA         NA      NA
## FamilySapindaceae:alt                                 NA         NA      NA
## FamilySapotaceae:alt                                  NA         NA      NA
## FamilyScrophulariaceae:alt                            NA         NA      NA
## FamilyThymelaeaceae:alt                               NA         NA      NA
## FamilyUlmaceae:alt                                    NA         NA      NA
## FamilyUrticaceae:alt                                  NA         NA      NA
## FamilyViolaceae:alt                                   NA         NA      NA
## FamilyXanthorrhoeaceae:alt                            NA         NA      NA
## FamilyZygophyllaceae:alt                              NA         NA      NA
## FamilyAsteraceae:temp                          9.741e+00  2.649e+00   3.678
## FamilyAtherospermataceae:temp                         NA         NA      NA
## FamilyBalsaminaceae:temp                              NA         NA      NA
## FamilyBetulaceae:temp                                 NA         NA      NA
## FamilyBrassicaceae:temp                               NA         NA      NA
## FamilyCactaceae:temp                                  NA         NA      NA
## FamilyCasuarinaceae:temp                              NA         NA      NA
## FamilyChloranthaceae:temp                             NA         NA      NA
## FamilyChrysobalanaceae:temp                           NA         NA      NA
## FamilyCistaceae:temp                                  NA         NA      NA
## FamilyCornaceae:temp                                  NA         NA      NA
## FamilyCrassulaceae:temp                               NA         NA      NA
## FamilyCunoniaceae:temp                                NA         NA      NA
## FamilyCupressaceae:temp                               NA         NA      NA
## FamilyCyperaceae:temp                                 NA         NA      NA
## FamilyDennstaedtiaceae:temp                           NA         NA      NA
## FamilyDicksoniaceae:temp                              NA         NA      NA
## FamilyDipterocarpaceae:temp                           NA         NA      NA
## FamilyEbenaceae:temp                                  NA         NA      NA
## FamilyElaeocarpaceae:temp                             NA         NA      NA
## FamilyEricaceae:temp                                  NA         NA      NA
## FamilyEuphorbiaceae:temp                              NA         NA      NA
## FamilyFabaceae - C:temp                               NA         NA      NA
## FamilyFabaceae - M:temp                               NA         NA      NA
## FamilyFabaceae - P:temp                               NA         NA      NA
## FamilyFagaceae:temp                                   NA         NA      NA
## FamilyGentianaceae:temp                               NA         NA      NA
## FamilyHeliconiaceae:temp                              NA         NA      NA
## FamilyJuglandaceae:temp                               NA         NA      NA
## FamilyJuncaginaceae:temp                              NA         NA      NA
## FamilyLamiaceae:temp                                  NA         NA      NA
## FamilyLauraceae:temp                                  NA         NA      NA
## FamilyMaesaceae:temp                                  NA         NA      NA
## FamilyMalvaceae:temp                                  NA         NA      NA
## FamilyMelastomataceae:temp                            NA         NA      NA
## FamilyMoraceae:temp                                   NA         NA      NA
## FamilyMyristicaceae:temp                              NA         NA      NA
## FamilyMyrsinaceae:temp                                NA         NA      NA
## FamilyMyrtaceae:temp                           1.113e+00  3.408e-01   3.266
## FamilyOchnaceae:temp                                  NA         NA      NA
## FamilyOnagraceae:temp                                 NA         NA      NA
## FamilyOrobanchaceae:temp                              NA         NA      NA
## FamilyPhyllanthaceae:temp                             NA         NA      NA
## FamilyPicrodendraceae:temp                            NA         NA      NA
## FamilyPinaceae:temp                                   NA         NA      NA
## FamilyPoaceae:temp                             6.425e-01  2.770e-01   2.320
## FamilyPolemoniaceae:temp                              NA         NA      NA
## FamilyPolygonaceae:temp                               NA         NA      NA
## FamilyProteaceae:temp                                 NA         NA      NA
## FamilyRanunculaceae:temp                              NA         NA      NA
## FamilyRhamnaceae:temp                                 NA         NA      NA
## FamilyRosaceae:temp                                   NA         NA      NA
## FamilyRubiaceae:temp                                  NA         NA      NA
## FamilyRutaceae:temp                                   NA         NA      NA
## FamilySalicaceae:temp                                 NA         NA      NA
## FamilySapindaceae:temp                                NA         NA      NA
## FamilySapotaceae:temp                                 NA         NA      NA
## FamilyScrophulariaceae:temp                           NA         NA      NA
## FamilyThymelaeaceae:temp                              NA         NA      NA
## FamilyUlmaceae:temp                                   NA         NA      NA
## FamilyUrticaceae:temp                                 NA         NA      NA
## FamilyViolaceae:temp                                  NA         NA      NA
## FamilyXanthorrhoeaceae:temp                           NA         NA      NA
## FamilyZygophyllaceae:temp                             NA         NA      NA
## lat:long                                      -4.850e-04  3.239e-04  -1.497
## lat:temp                                      -5.288e-03  1.011e-03  -5.231
## lat:NPP                                       -4.384e-05  4.734e-05  -0.926
## long:alt                                      -1.635e-06  2.683e-06  -0.609
## long:temp                                     -1.650e-03  5.755e-04  -2.868
## long:NPP                                       7.363e-06  2.523e-06   2.918
## alt:temp                                      -5.835e-05  1.708e-05  -3.416
## temp:NPP                                      -1.631e-04  8.021e-05  -2.034
##                                               Pr(>|t|)    
## (Intercept)                                   0.471210    
## growthformHerb                                0.401018    
## growthformHerb/Shrub                          0.959364    
## growthformShrub                               0.459650    
## growthformShrub/Tree                          0.576973    
## growthformTree                                0.209409    
## FamilyAsteraceae                              0.024915 *  
## FamilyAtherospermataceae                      0.240297    
## FamilyBalsaminaceae                           0.445806    
## FamilyBetulaceae                              0.781378    
## FamilyBrassicaceae                            0.520831    
## FamilyCactaceae                               0.374789    
## FamilyCasuarinaceae                           0.457542    
## FamilyChloranthaceae                          0.451076    
## FamilyChrysobalanaceae                        0.531209    
## FamilyCistaceae                               0.457466    
## FamilyCornaceae                               0.444251    
## FamilyCrassulaceae                            0.461491    
## FamilyCunoniaceae                             0.580889    
## FamilyCupressaceae                            0.306038    
## FamilyCyperaceae                              0.397873    
## FamilyDennstaedtiaceae                              NA    
## FamilyDicksoniaceae                           0.468095    
## FamilyDipterocarpaceae                        0.398966    
## FamilyEbenaceae                               0.397522    
## FamilyElaeocarpaceae                          0.323838    
## FamilyEricaceae                               0.451810    
## FamilyEuphorbiaceae                           0.438290    
## FamilyFabaceae - C                            0.453083    
## FamilyFabaceae - M                            0.244394    
## FamilyFabaceae - P                            0.340876    
## FamilyFagaceae                                0.799425    
## FamilyGentianaceae                            0.406314    
## FamilyHeliconiaceae                           0.418080    
## FamilyJuglandaceae                            0.307939    
## FamilyJuncaginaceae                           0.445752    
## FamilyLamiaceae                               0.474074    
## FamilyLauraceae                               0.548932    
## FamilyMaesaceae                               0.447295    
## FamilyMalvaceae                               0.537119    
## FamilyMelastomataceae                         0.447635    
## FamilyMoraceae                                0.380649    
## FamilyMyristicaceae                           0.349574    
## FamilyMyrsinaceae                             0.393440    
## FamilyMyrtaceae                               0.011832 *  
## FamilyOchnaceae                               0.448174    
## FamilyOnagraceae                                    NA    
## FamilyOrobanchaceae                           0.553364    
## FamilyPhyllanthaceae                          0.393311    
## FamilyPicrodendraceae                         0.359305    
## FamilyPinaceae                                0.087314 .  
## FamilyPoaceae                                 0.529061    
## FamilyPolemoniaceae                           0.415849    
## FamilyPolygonaceae                            0.444282    
## FamilyProteaceae                              0.106800    
## FamilyRanunculaceae                           0.447850    
## FamilyRhamnaceae                              0.466950    
## FamilyRosaceae                                0.016301 *  
## FamilyRubiaceae                               0.441496    
## FamilyRutaceae                                0.456598    
## FamilySalicaceae                              0.433258    
## FamilySapindaceae                             0.607885    
## FamilySapotaceae                              0.276654    
## FamilyScrophulariaceae                        0.448464    
## FamilyThymelaeaceae                           0.336091    
## FamilyUlmaceae                                0.231765    
## FamilyUrticaceae                              0.338942    
## FamilyViolaceae                               0.447326    
## FamilyXanthorrhoeaceae                        0.455414    
## FamilyZygophyllaceae                          0.467051    
## lat                                           0.325004    
## long                                          0.860913    
## alt                                           0.550965    
## temp                                          0.447922    
## NPP                                           0.298512    
## growthformHerb:FamilyAsteraceae               0.761650    
## growthformHerb/Shrub:FamilyAsteraceae               NA    
## growthformShrub:FamilyAsteraceae                    NA    
## growthformShrub/Tree:FamilyAsteraceae               NA    
## growthformTree:FamilyAsteraceae                     NA    
## growthformHerb:FamilyAtherospermataceae             NA    
## growthformHerb/Shrub:FamilyAtherospermataceae       NA    
## growthformShrub:FamilyAtherospermataceae            NA    
## growthformShrub/Tree:FamilyAtherospermataceae       NA    
## growthformTree:FamilyAtherospermataceae             NA    
## growthformHerb:FamilyBalsaminaceae                  NA    
## growthformHerb/Shrub:FamilyBalsaminaceae            NA    
## growthformShrub:FamilyBalsaminaceae                 NA    
## growthformShrub/Tree:FamilyBalsaminaceae            NA    
## growthformTree:FamilyBalsaminaceae                  NA    
## growthformHerb:FamilyBetulaceae                     NA    
## growthformHerb/Shrub:FamilyBetulaceae               NA    
## growthformShrub:FamilyBetulaceae              0.483172    
## growthformShrub/Tree:FamilyBetulaceae               NA    
## growthformTree:FamilyBetulaceae                     NA    
## growthformHerb:FamilyBrassicaceae                   NA    
## growthformHerb/Shrub:FamilyBrassicaceae             NA    
## growthformShrub:FamilyBrassicaceae                  NA    
## growthformShrub/Tree:FamilyBrassicaceae             NA    
## growthformTree:FamilyBrassicaceae                   NA    
## growthformHerb:FamilyCactaceae                      NA    
## growthformHerb/Shrub:FamilyCactaceae                NA    
## growthformShrub:FamilyCactaceae                     NA    
## growthformShrub/Tree:FamilyCactaceae                NA    
## growthformTree:FamilyCactaceae                      NA    
## growthformHerb:FamilyCasuarinaceae                  NA    
## growthformHerb/Shrub:FamilyCasuarinaceae            NA    
## growthformShrub:FamilyCasuarinaceae                 NA    
## growthformShrub/Tree:FamilyCasuarinaceae            NA    
## growthformTree:FamilyCasuarinaceae                  NA    
## growthformHerb:FamilyChloranthaceae                 NA    
## growthformHerb/Shrub:FamilyChloranthaceae           NA    
## growthformShrub:FamilyChloranthaceae                NA    
## growthformShrub/Tree:FamilyChloranthaceae           NA    
## growthformTree:FamilyChloranthaceae                 NA    
## growthformHerb:FamilyChrysobalanaceae               NA    
## growthformHerb/Shrub:FamilyChrysobalanaceae         NA    
## growthformShrub:FamilyChrysobalanaceae              NA    
## growthformShrub/Tree:FamilyChrysobalanaceae         NA    
## growthformTree:FamilyChrysobalanaceae               NA    
## growthformHerb:FamilyCistaceae                      NA    
## growthformHerb/Shrub:FamilyCistaceae                NA    
## growthformShrub:FamilyCistaceae                     NA    
## growthformShrub/Tree:FamilyCistaceae                NA    
## growthformTree:FamilyCistaceae                      NA    
## growthformHerb:FamilyCornaceae                      NA    
## growthformHerb/Shrub:FamilyCornaceae                NA    
## growthformShrub:FamilyCornaceae                     NA    
## growthformShrub/Tree:FamilyCornaceae                NA    
## growthformTree:FamilyCornaceae                      NA    
## growthformHerb:FamilyCrassulaceae                   NA    
## growthformHerb/Shrub:FamilyCrassulaceae             NA    
## growthformShrub:FamilyCrassulaceae                  NA    
## growthformShrub/Tree:FamilyCrassulaceae             NA    
## growthformTree:FamilyCrassulaceae                   NA    
## growthformHerb:FamilyCunoniaceae                    NA    
## growthformHerb/Shrub:FamilyCunoniaceae              NA    
## growthformShrub:FamilyCunoniaceae                   NA    
## growthformShrub/Tree:FamilyCunoniaceae        0.448590    
## growthformTree:FamilyCunoniaceae                    NA    
## growthformHerb:FamilyCupressaceae                   NA    
## growthformHerb/Shrub:FamilyCupressaceae             NA    
## growthformShrub:FamilyCupressaceae            0.447817    
## growthformShrub/Tree:FamilyCupressaceae             NA    
## growthformTree:FamilyCupressaceae                   NA    
## growthformHerb:FamilyCyperaceae                     NA    
## growthformHerb/Shrub:FamilyCyperaceae               NA    
## growthformShrub:FamilyCyperaceae                    NA    
## growthformShrub/Tree:FamilyCyperaceae               NA    
## growthformTree:FamilyCyperaceae                     NA    
## growthformHerb:FamilyDennstaedtiaceae               NA    
## growthformHerb/Shrub:FamilyDennstaedtiaceae         NA    
## growthformShrub:FamilyDennstaedtiaceae              NA    
## growthformShrub/Tree:FamilyDennstaedtiaceae         NA    
## growthformTree:FamilyDennstaedtiaceae               NA    
## growthformHerb:FamilyDicksoniaceae                  NA    
## growthformHerb/Shrub:FamilyDicksoniaceae            NA    
## growthformShrub:FamilyDicksoniaceae                 NA    
## growthformShrub/Tree:FamilyDicksoniaceae            NA    
## growthformTree:FamilyDicksoniaceae                  NA    
## growthformHerb:FamilyDipterocarpaceae               NA    
## growthformHerb/Shrub:FamilyDipterocarpaceae         NA    
## growthformShrub:FamilyDipterocarpaceae              NA    
## growthformShrub/Tree:FamilyDipterocarpaceae         NA    
## growthformTree:FamilyDipterocarpaceae               NA    
## growthformHerb:FamilyEbenaceae                      NA    
## growthformHerb/Shrub:FamilyEbenaceae                NA    
## growthformShrub:FamilyEbenaceae                     NA    
## growthformShrub/Tree:FamilyEbenaceae                NA    
## growthformTree:FamilyEbenaceae                      NA    
## growthformHerb:FamilyElaeocarpaceae                 NA    
## growthformHerb/Shrub:FamilyElaeocarpaceae           NA    
## growthformShrub:FamilyElaeocarpaceae                NA    
## growthformShrub/Tree:FamilyElaeocarpaceae           NA    
## growthformTree:FamilyElaeocarpaceae                 NA    
## growthformHerb:FamilyEricaceae                      NA    
## growthformHerb/Shrub:FamilyEricaceae                NA    
## growthformShrub:FamilyEricaceae               0.724294    
## growthformShrub/Tree:FamilyEricaceae                NA    
## growthformTree:FamilyEricaceae                      NA    
## growthformHerb:FamilyEuphorbiaceae                  NA    
## growthformHerb/Shrub:FamilyEuphorbiaceae            NA    
## growthformShrub:FamilyEuphorbiaceae           0.444658    
## growthformShrub/Tree:FamilyEuphorbiaceae            NA    
## growthformTree:FamilyEuphorbiaceae                  NA    
## growthformHerb:FamilyFabaceae - C                   NA    
## growthformHerb/Shrub:FamilyFabaceae - C             NA    
## growthformShrub:FamilyFabaceae - C                  NA    
## growthformShrub/Tree:FamilyFabaceae - C             NA    
## growthformTree:FamilyFabaceae - C                   NA    
## growthformHerb:FamilyFabaceae - M                   NA    
## growthformHerb/Shrub:FamilyFabaceae - M             NA    
## growthformShrub:FamilyFabaceae - M                  NA    
## growthformShrub/Tree:FamilyFabaceae - M       0.432247    
## growthformTree:FamilyFabaceae - M                   NA    
## growthformHerb:FamilyFabaceae - P             0.412368    
## growthformHerb/Shrub:FamilyFabaceae - P             NA    
## growthformShrub:FamilyFabaceae - P            0.434788    
## growthformShrub/Tree:FamilyFabaceae - P             NA    
## growthformTree:FamilyFabaceae - P                   NA    
## growthformHerb:FamilyFagaceae                       NA    
## growthformHerb/Shrub:FamilyFagaceae                 NA    
## growthformShrub:FamilyFagaceae                0.451214    
## growthformShrub/Tree:FamilyFagaceae                 NA    
## growthformTree:FamilyFagaceae                       NA    
## growthformHerb:FamilyGentianaceae                   NA    
## growthformHerb/Shrub:FamilyGentianaceae             NA    
## growthformShrub:FamilyGentianaceae                  NA    
## growthformShrub/Tree:FamilyGentianaceae             NA    
## growthformTree:FamilyGentianaceae                   NA    
## growthformHerb:FamilyHeliconiaceae                  NA    
## growthformHerb/Shrub:FamilyHeliconiaceae            NA    
## growthformShrub:FamilyHeliconiaceae                 NA    
## growthformShrub/Tree:FamilyHeliconiaceae            NA    
## growthformTree:FamilyHeliconiaceae                  NA    
## growthformHerb:FamilyJuglandaceae                   NA    
## growthformHerb/Shrub:FamilyJuglandaceae             NA    
## growthformShrub:FamilyJuglandaceae                  NA    
## growthformShrub/Tree:FamilyJuglandaceae             NA    
## growthformTree:FamilyJuglandaceae                   NA    
## growthformHerb:FamilyJuncaginaceae                  NA    
## growthformHerb/Shrub:FamilyJuncaginaceae            NA    
## growthformShrub:FamilyJuncaginaceae                 NA    
## growthformShrub/Tree:FamilyJuncaginaceae            NA    
## growthformTree:FamilyJuncaginaceae                  NA    
## growthformHerb:FamilyLamiaceae                      NA    
## growthformHerb/Shrub:FamilyLamiaceae                NA    
## growthformShrub:FamilyLamiaceae                     NA    
## growthformShrub/Tree:FamilyLamiaceae                NA    
## growthformTree:FamilyLamiaceae                      NA    
## growthformHerb:FamilyLauraceae                      NA    
## growthformHerb/Shrub:FamilyLauraceae                NA    
## growthformShrub:FamilyLauraceae                     NA    
## growthformShrub/Tree:FamilyLauraceae                NA    
## growthformTree:FamilyLauraceae                      NA    
## growthformHerb:FamilyMaesaceae                      NA    
## growthformHerb/Shrub:FamilyMaesaceae                NA    
## growthformShrub:FamilyMaesaceae                     NA    
## growthformShrub/Tree:FamilyMaesaceae                NA    
## growthformTree:FamilyMaesaceae                      NA    
## growthformHerb:FamilyMalvaceae                      NA    
## growthformHerb/Shrub:FamilyMalvaceae                NA    
## growthformShrub:FamilyMalvaceae                     NA    
## growthformShrub/Tree:FamilyMalvaceae                NA    
## growthformTree:FamilyMalvaceae                      NA    
## growthformHerb:FamilyMelastomataceae                NA    
## growthformHerb/Shrub:FamilyMelastomataceae          NA    
## growthformShrub:FamilyMelastomataceae               NA    
## growthformShrub/Tree:FamilyMelastomataceae          NA    
## growthformTree:FamilyMelastomataceae                NA    
## growthformHerb:FamilyMoraceae                       NA    
## growthformHerb/Shrub:FamilyMoraceae                 NA    
## growthformShrub:FamilyMoraceae                      NA    
## growthformShrub/Tree:FamilyMoraceae                 NA    
## growthformTree:FamilyMoraceae                       NA    
## growthformHerb:FamilyMyristicaceae                  NA    
## growthformHerb/Shrub:FamilyMyristicaceae            NA    
## growthformShrub:FamilyMyristicaceae                 NA    
## growthformShrub/Tree:FamilyMyristicaceae            NA    
## growthformTree:FamilyMyristicaceae                  NA    
## growthformHerb:FamilyMyrsinaceae                    NA    
## growthformHerb/Shrub:FamilyMyrsinaceae              NA    
## growthformShrub:FamilyMyrsinaceae                   NA    
## growthformShrub/Tree:FamilyMyrsinaceae        0.444600    
## growthformTree:FamilyMyrsinaceae                    NA    
## growthformHerb:FamilyMyrtaceae                      NA    
## growthformHerb/Shrub:FamilyMyrtaceae                NA    
## growthformShrub:FamilyMyrtaceae               0.441684    
## growthformShrub/Tree:FamilyMyrtaceae                NA    
## growthformTree:FamilyMyrtaceae                      NA    
## growthformHerb:FamilyOchnaceae                      NA    
## growthformHerb/Shrub:FamilyOchnaceae                NA    
## growthformShrub:FamilyOchnaceae                     NA    
## growthformShrub/Tree:FamilyOchnaceae                NA    
## growthformTree:FamilyOchnaceae                      NA    
## growthformHerb:FamilyOnagraceae                     NA    
## growthformHerb/Shrub:FamilyOnagraceae               NA    
## growthformShrub:FamilyOnagraceae                    NA    
## growthformShrub/Tree:FamilyOnagraceae               NA    
## growthformTree:FamilyOnagraceae                     NA    
## growthformHerb:FamilyOrobanchaceae                  NA    
## growthformHerb/Shrub:FamilyOrobanchaceae            NA    
## growthformShrub:FamilyOrobanchaceae                 NA    
## growthformShrub/Tree:FamilyOrobanchaceae            NA    
## growthformTree:FamilyOrobanchaceae                  NA    
## growthformHerb:FamilyPhyllanthaceae                 NA    
## growthformHerb/Shrub:FamilyPhyllanthaceae           NA    
## growthformShrub:FamilyPhyllanthaceae                NA    
## growthformShrub/Tree:FamilyPhyllanthaceae           NA    
## growthformTree:FamilyPhyllanthaceae                 NA    
## growthformHerb:FamilyPicrodendraceae                NA    
## growthformHerb/Shrub:FamilyPicrodendraceae          NA    
## growthformShrub:FamilyPicrodendraceae               NA    
## growthformShrub/Tree:FamilyPicrodendraceae          NA    
## growthformTree:FamilyPicrodendraceae                NA    
## growthformHerb:FamilyPinaceae                       NA    
## growthformHerb/Shrub:FamilyPinaceae                 NA    
## growthformShrub:FamilyPinaceae                      NA    
## growthformShrub/Tree:FamilyPinaceae                 NA    
## growthformTree:FamilyPinaceae                       NA    
## growthformHerb:FamilyPoaceae                  0.441587    
## growthformHerb/Shrub:FamilyPoaceae                  NA    
## growthformShrub:FamilyPoaceae                       NA    
## growthformShrub/Tree:FamilyPoaceae                  NA    
## growthformTree:FamilyPoaceae                        NA    
## growthformHerb:FamilyPolemoniaceae                  NA    
## growthformHerb/Shrub:FamilyPolemoniaceae            NA    
## growthformShrub:FamilyPolemoniaceae                 NA    
## growthformShrub/Tree:FamilyPolemoniaceae            NA    
## growthformTree:FamilyPolemoniaceae                  NA    
## growthformHerb:FamilyPolygonaceae                   NA    
## growthformHerb/Shrub:FamilyPolygonaceae             NA    
## growthformShrub:FamilyPolygonaceae                  NA    
## growthformShrub/Tree:FamilyPolygonaceae             NA    
## growthformTree:FamilyPolygonaceae                   NA    
## growthformHerb:FamilyProteaceae                     NA    
## growthformHerb/Shrub:FamilyProteaceae               NA    
## growthformShrub:FamilyProteaceae              0.441738    
## growthformShrub/Tree:FamilyProteaceae               NA    
## growthformTree:FamilyProteaceae                     NA    
## growthformHerb:FamilyRanunculaceae                  NA    
## growthformHerb/Shrub:FamilyRanunculaceae            NA    
## growthformShrub:FamilyRanunculaceae                 NA    
## growthformShrub/Tree:FamilyRanunculaceae            NA    
## growthformTree:FamilyRanunculaceae                  NA    
## growthformHerb:FamilyRhamnaceae                     NA    
## growthformHerb/Shrub:FamilyRhamnaceae               NA    
## growthformShrub:FamilyRhamnaceae                    NA    
## growthformShrub/Tree:FamilyRhamnaceae               NA    
## growthformTree:FamilyRhamnaceae                     NA    
## growthformHerb:FamilyRosaceae                 0.423206    
## growthformHerb/Shrub:FamilyRosaceae                 NA    
## growthformShrub:FamilyRosaceae                0.451569    
## growthformShrub/Tree:FamilyRosaceae                 NA    
## growthformTree:FamilyRosaceae                       NA    
## growthformHerb:FamilyRubiaceae                      NA    
## growthformHerb/Shrub:FamilyRubiaceae                NA    
## growthformShrub:FamilyRubiaceae               0.094423 .  
## growthformShrub/Tree:FamilyRubiaceae                NA    
## growthformTree:FamilyRubiaceae                      NA    
## growthformHerb:FamilyRutaceae                       NA    
## growthformHerb/Shrub:FamilyRutaceae                 NA    
## growthformShrub:FamilyRutaceae                      NA    
## growthformShrub/Tree:FamilyRutaceae                 NA    
## growthformTree:FamilyRutaceae                       NA    
## growthformHerb:FamilySalicaceae                     NA    
## growthformHerb/Shrub:FamilySalicaceae               NA    
## growthformShrub:FamilySalicaceae              0.195221    
## growthformShrub/Tree:FamilySalicaceae               NA    
## growthformTree:FamilySalicaceae                     NA    
## growthformHerb:FamilySapindaceae                    NA    
## growthformHerb/Shrub:FamilySapindaceae              NA    
## growthformShrub:FamilySapindaceae                   NA    
## growthformShrub/Tree:FamilySapindaceae              NA    
## growthformTree:FamilySapindaceae                    NA    
## growthformHerb:FamilySapotaceae                     NA    
## growthformHerb/Shrub:FamilySapotaceae               NA    
## growthformShrub:FamilySapotaceae                    NA    
## growthformShrub/Tree:FamilySapotaceae               NA    
## growthformTree:FamilySapotaceae                     NA    
## growthformHerb:FamilyScrophulariaceae               NA    
## growthformHerb/Shrub:FamilyScrophulariaceae         NA    
## growthformShrub:FamilyScrophulariaceae              NA    
## growthformShrub/Tree:FamilyScrophulariaceae         NA    
## growthformTree:FamilyScrophulariaceae               NA    
## growthformHerb:FamilyThymelaeaceae                  NA    
## growthformHerb/Shrub:FamilyThymelaeaceae            NA    
## growthformShrub:FamilyThymelaeaceae                 NA    
## growthformShrub/Tree:FamilyThymelaeaceae            NA    
## growthformTree:FamilyThymelaeaceae                  NA    
## growthformHerb:FamilyUlmaceae                       NA    
## growthformHerb/Shrub:FamilyUlmaceae                 NA    
## growthformShrub:FamilyUlmaceae                      NA    
## growthformShrub/Tree:FamilyUlmaceae                 NA    
## growthformTree:FamilyUlmaceae                       NA    
## growthformHerb:FamilyUrticaceae                     NA    
## growthformHerb/Shrub:FamilyUrticaceae               NA    
## growthformShrub:FamilyUrticaceae                    NA    
## growthformShrub/Tree:FamilyUrticaceae               NA    
## growthformTree:FamilyUrticaceae                     NA    
## growthformHerb:FamilyViolaceae                      NA    
## growthformHerb/Shrub:FamilyViolaceae                NA    
## growthformShrub:FamilyViolaceae                     NA    
## growthformShrub/Tree:FamilyViolaceae                NA    
## growthformTree:FamilyViolaceae                      NA    
## growthformHerb:FamilyXanthorrhoeaceae               NA    
## growthformHerb/Shrub:FamilyXanthorrhoeaceae         NA    
## growthformShrub:FamilyXanthorrhoeaceae              NA    
## growthformShrub/Tree:FamilyXanthorrhoeaceae         NA    
## growthformTree:FamilyXanthorrhoeaceae               NA    
## growthformHerb:FamilyZygophyllaceae                 NA    
## growthformHerb/Shrub:FamilyZygophyllaceae           NA    
## growthformShrub:FamilyZygophyllaceae                NA    
## growthformShrub/Tree:FamilyZygophyllaceae           NA    
## growthformTree:FamilyZygophyllaceae                 NA    
## growthformHerb:lat                            0.420417    
## growthformHerb/Shrub:lat                            NA    
## growthformShrub:lat                           0.045980 *  
## growthformShrub/Tree:lat                      0.360912    
## growthformTree:lat                                  NA    
## growthformHerb:long                           0.369731    
## growthformHerb/Shrub:long                           NA    
## growthformShrub:long                          0.665856    
## growthformShrub/Tree:long                           NA    
## growthformTree:long                                 NA    
## growthformHerb:NPP                            0.329320    
## growthformHerb/Shrub:NPP                            NA    
## growthformShrub:NPP                           0.123317    
## growthformShrub/Tree:NPP                            NA    
## growthformTree:NPP                                  NA    
## FamilyAsteraceae:lat                          0.283138    
## FamilyAtherospermataceae:lat                        NA    
## FamilyBalsaminaceae:lat                             NA    
## FamilyBetulaceae:lat                          0.785252    
## FamilyBrassicaceae:lat                        0.482497    
## FamilyCactaceae:lat                           0.046821 *  
## FamilyCasuarinaceae:lat                             NA    
## FamilyChloranthaceae:lat                            NA    
## FamilyChrysobalanaceae:lat                    0.570866    
## FamilyCistaceae:lat                                 NA    
## FamilyCornaceae:lat                                 NA    
## FamilyCrassulaceae:lat                              NA    
## FamilyCunoniaceae:lat                               NA    
## FamilyCupressaceae:lat                              NA    
## FamilyCyperaceae:lat                                NA    
## FamilyDennstaedtiaceae:lat                          NA    
## FamilyDicksoniaceae:lat                             NA    
## FamilyDipterocarpaceae:lat                          NA    
## FamilyEbenaceae:lat                           0.402213    
## FamilyElaeocarpaceae:lat                            NA    
## FamilyEricaceae:lat                           0.223129    
## FamilyEuphorbiaceae:lat                       0.344510    
## FamilyFabaceae - C:lat                        0.399668    
## FamilyFabaceae - M:lat                              NA    
## FamilyFabaceae - P:lat                        0.339143    
## FamilyFagaceae:lat                            0.367245    
## FamilyGentianaceae:lat                        0.426756    
## FamilyHeliconiaceae:lat                             NA    
## FamilyJuglandaceae:lat                              NA    
## FamilyJuncaginaceae:lat                             NA    
## FamilyLamiaceae:lat                                 NA    
## FamilyLauraceae:lat                                 NA    
## FamilyMaesaceae:lat                                 NA    
## FamilyMalvaceae:lat                           0.451583    
## FamilyMelastomataceae:lat                           NA    
## FamilyMoraceae:lat                                  NA    
## FamilyMyristicaceae:lat                             NA    
## FamilyMyrsinaceae:lat                               NA    
## FamilyMyrtaceae:lat                           0.675074    
## FamilyOchnaceae:lat                                 NA    
## FamilyOnagraceae:lat                                NA    
## FamilyOrobanchaceae:lat                       0.476749    
## FamilyPhyllanthaceae:lat                            NA    
## FamilyPicrodendraceae:lat                           NA    
## FamilyPinaceae:lat                            0.159861    
## FamilyPoaceae:lat                             0.420438    
## FamilyPolemoniaceae:lat                       0.433041    
## FamilyPolygonaceae:lat                        0.444087    
## FamilyProteaceae:lat                          0.149712    
## FamilyRanunculaceae:lat                             NA    
## FamilyRhamnaceae:lat                                NA    
## FamilyRosaceae:lat                            0.160271    
## FamilyRubiaceae:lat                           0.221123    
## FamilyRutaceae:lat                                  NA    
## FamilySalicaceae:lat                          0.240215    
## FamilySapindaceae:lat                         0.380059    
## FamilySapotaceae:lat                                NA    
## FamilyScrophulariaceae:lat                          NA    
## FamilyThymelaeaceae:lat                             NA    
## FamilyUlmaceae:lat                            0.282140    
## FamilyUrticaceae:lat                          0.358817    
## FamilyViolaceae:lat                                 NA    
## FamilyXanthorrhoeaceae:lat                          NA    
## FamilyZygophyllaceae:lat                            NA    
## FamilyAsteraceae:long                         0.001268 ** 
## FamilyAtherospermataceae:long                       NA    
## FamilyBalsaminaceae:long                            NA    
## FamilyBetulaceae:long                         0.374743    
## FamilyBrassicaceae:long                             NA    
## FamilyCactaceae:long                          0.155711    
## FamilyCasuarinaceae:long                            NA    
## FamilyChloranthaceae:long                           NA    
## FamilyChrysobalanaceae:long                         NA    
## FamilyCistaceae:long                                NA    
## FamilyCornaceae:long                                NA    
## FamilyCrassulaceae:long                             NA    
## FamilyCunoniaceae:long                              NA    
## FamilyCupressaceae:long                             NA    
## FamilyCyperaceae:long                               NA    
## FamilyDennstaedtiaceae:long                         NA    
## FamilyDicksoniaceae:long                            NA    
## FamilyDipterocarpaceae:long                         NA    
## FamilyEbenaceae:long                                NA    
## FamilyElaeocarpaceae:long                           NA    
## FamilyEricaceae:long                          0.056466 .  
## FamilyEuphorbiaceae:long                      0.359718    
## FamilyFabaceae - C:long                             NA    
## FamilyFabaceae - M:long                             NA    
## FamilyFabaceae - P:long                             NA    
## FamilyFagaceae:long                                 NA    
## FamilyGentianaceae:long                             NA    
## FamilyHeliconiaceae:long                            NA    
## FamilyJuglandaceae:long                             NA    
## FamilyJuncaginaceae:long                            NA    
## FamilyLamiaceae:long                                NA    
## FamilyLauraceae:long                                NA    
## FamilyMaesaceae:long                                NA    
## FamilyMalvaceae:long                          0.449768    
## FamilyMelastomataceae:long                          NA    
## FamilyMoraceae:long                                 NA    
## FamilyMyristicaceae:long                            NA    
## FamilyMyrsinaceae:long                              NA    
## FamilyMyrtaceae:long                          0.129339    
## FamilyOchnaceae:long                                NA    
## FamilyOnagraceae:long                               NA    
## FamilyOrobanchaceae:long                            NA    
## FamilyPhyllanthaceae:long                           NA    
## FamilyPicrodendraceae:long                          NA    
## FamilyPinaceae:long                           0.977523    
## FamilyPoaceae:long                                  NA    
## FamilyPolemoniaceae:long                            NA    
## FamilyPolygonaceae:long                             NA    
## FamilyProteaceae:long                         0.160569    
## FamilyRanunculaceae:long                            NA    
## FamilyRhamnaceae:long                               NA    
## FamilyRosaceae:long                           0.294900    
## FamilyRubiaceae:long                                NA    
## FamilyRutaceae:long                                 NA    
## FamilySalicaceae:long                               NA    
## FamilySapindaceae:long                        0.118118    
## FamilySapotaceae:long                               NA    
## FamilyScrophulariaceae:long                         NA    
## FamilyThymelaeaceae:long                            NA    
## FamilyUlmaceae:long                                 NA    
## FamilyUrticaceae:long                               NA    
## FamilyViolaceae:long                                NA    
## FamilyXanthorrhoeaceae:long                         NA    
## FamilyZygophyllaceae:long                           NA    
## FamilyAsteraceae:alt                          0.038360 *  
## FamilyAtherospermataceae:alt                        NA    
## FamilyBalsaminaceae:alt                             NA    
## FamilyBetulaceae:alt                                NA    
## FamilyBrassicaceae:alt                              NA    
## FamilyCactaceae:alt                           0.653610    
## FamilyCasuarinaceae:alt                             NA    
## FamilyChloranthaceae:alt                            NA    
## FamilyChrysobalanaceae:alt                          NA    
## FamilyCistaceae:alt                                 NA    
## FamilyCornaceae:alt                                 NA    
## FamilyCrassulaceae:alt                              NA    
## FamilyCunoniaceae:alt                               NA    
## FamilyCupressaceae:alt                              NA    
## FamilyCyperaceae:alt                                NA    
## FamilyDennstaedtiaceae:alt                          NA    
## FamilyDicksoniaceae:alt                             NA    
## FamilyDipterocarpaceae:alt                          NA    
## FamilyEbenaceae:alt                                 NA    
## FamilyElaeocarpaceae:alt                            NA    
## FamilyEricaceae:alt                           0.833476    
## FamilyEuphorbiaceae:alt                             NA    
## FamilyFabaceae - C:alt                              NA    
## FamilyFabaceae - M:alt                              NA    
## FamilyFabaceae - P:alt                              NA    
## FamilyFagaceae:alt                                  NA    
## FamilyGentianaceae:alt                              NA    
## FamilyHeliconiaceae:alt                             NA    
## FamilyJuglandaceae:alt                              NA    
## FamilyJuncaginaceae:alt                             NA    
## FamilyLamiaceae:alt                                 NA    
## FamilyLauraceae:alt                                 NA    
## FamilyMaesaceae:alt                                 NA    
## FamilyMalvaceae:alt                                 NA    
## FamilyMelastomataceae:alt                           NA    
## FamilyMoraceae:alt                                  NA    
## FamilyMyristicaceae:alt                             NA    
## FamilyMyrsinaceae:alt                               NA    
## FamilyMyrtaceae:alt                           0.818358    
## FamilyOchnaceae:alt                                 NA    
## FamilyOnagraceae:alt                                NA    
## FamilyOrobanchaceae:alt                             NA    
## FamilyPhyllanthaceae:alt                            NA    
## FamilyPicrodendraceae:alt                           NA    
## FamilyPinaceae:alt                            0.436576    
## FamilyPoaceae:alt                             0.673346    
## FamilyPolemoniaceae:alt                             NA    
## FamilyPolygonaceae:alt                              NA    
## FamilyProteaceae:alt                          0.346160    
## FamilyRanunculaceae:alt                             NA    
## FamilyRhamnaceae:alt                                NA    
## FamilyRosaceae:alt                                  NA    
## FamilyRubiaceae:alt                                 NA    
## FamilyRutaceae:alt                                  NA    
## FamilySalicaceae:alt                                NA    
## FamilySapindaceae:alt                               NA    
## FamilySapotaceae:alt                                NA    
## FamilyScrophulariaceae:alt                          NA    
## FamilyThymelaeaceae:alt                             NA    
## FamilyUlmaceae:alt                                  NA    
## FamilyUrticaceae:alt                                NA    
## FamilyViolaceae:alt                                 NA    
## FamilyXanthorrhoeaceae:alt                          NA    
## FamilyZygophyllaceae:alt                            NA    
## FamilyAsteraceae:temp                         0.004264 ** 
## FamilyAtherospermataceae:temp                       NA    
## FamilyBalsaminaceae:temp                            NA    
## FamilyBetulaceae:temp                               NA    
## FamilyBrassicaceae:temp                             NA    
## FamilyCactaceae:temp                                NA    
## FamilyCasuarinaceae:temp                            NA    
## FamilyChloranthaceae:temp                           NA    
## FamilyChrysobalanaceae:temp                         NA    
## FamilyCistaceae:temp                                NA    
## FamilyCornaceae:temp                                NA    
## FamilyCrassulaceae:temp                             NA    
## FamilyCunoniaceae:temp                              NA    
## FamilyCupressaceae:temp                             NA    
## FamilyCyperaceae:temp                               NA    
## FamilyDennstaedtiaceae:temp                         NA    
## FamilyDicksoniaceae:temp                            NA    
## FamilyDipterocarpaceae:temp                         NA    
## FamilyEbenaceae:temp                                NA    
## FamilyElaeocarpaceae:temp                           NA    
## FamilyEricaceae:temp                                NA    
## FamilyEuphorbiaceae:temp                            NA    
## FamilyFabaceae - C:temp                             NA    
## FamilyFabaceae - M:temp                             NA    
## FamilyFabaceae - P:temp                             NA    
## FamilyFagaceae:temp                                 NA    
## FamilyGentianaceae:temp                             NA    
## FamilyHeliconiaceae:temp                            NA    
## FamilyJuglandaceae:temp                             NA    
## FamilyJuncaginaceae:temp                            NA    
## FamilyLamiaceae:temp                                NA    
## FamilyLauraceae:temp                                NA    
## FamilyMaesaceae:temp                                NA    
## FamilyMalvaceae:temp                                NA    
## FamilyMelastomataceae:temp                          NA    
## FamilyMoraceae:temp                                 NA    
## FamilyMyristicaceae:temp                            NA    
## FamilyMyrsinaceae:temp                              NA    
## FamilyMyrtaceae:temp                          0.008481 ** 
## FamilyOchnaceae:temp                                NA    
## FamilyOnagraceae:temp                               NA    
## FamilyOrobanchaceae:temp                            NA    
## FamilyPhyllanthaceae:temp                           NA    
## FamilyPicrodendraceae:temp                          NA    
## FamilyPinaceae:temp                                 NA    
## FamilyPoaceae:temp                            0.042787 *  
## FamilyPolemoniaceae:temp                            NA    
## FamilyPolygonaceae:temp                             NA    
## FamilyProteaceae:temp                               NA    
## FamilyRanunculaceae:temp                            NA    
## FamilyRhamnaceae:temp                               NA    
## FamilyRosaceae:temp                                 NA    
## FamilyRubiaceae:temp                                NA    
## FamilyRutaceae:temp                                 NA    
## FamilySalicaceae:temp                               NA    
## FamilySapindaceae:temp                              NA    
## FamilySapotaceae:temp                               NA    
## FamilyScrophulariaceae:temp                         NA    
## FamilyThymelaeaceae:temp                            NA    
## FamilyUlmaceae:temp                                 NA    
## FamilyUrticaceae:temp                               NA    
## FamilyViolaceae:temp                                NA    
## FamilyXanthorrhoeaceae:temp                         NA    
## FamilyZygophyllaceae:temp                           NA    
## lat:long                                      0.165174    
## lat:temp                                      0.000384 ***
## lat:NPP                                       0.376311    
## long:alt                                      0.555973    
## long:temp                                     0.016734 *  
## long:NPP                                      0.015356 *  
## alt:temp                                      0.006587 ** 
## temp:NPP                                      0.069341 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.14 on 10 degrees of freedom
##   (15 observations deleted due to missingness)
## Multiple R-squared:  0.9981, Adjusted R-squared:  0.9695 
## F-statistic: 34.87 on 152 and 10 DF,  p-value: 5.06e-07

R2 = 0.998 is a warning sign: the selected model has about 150 parameters for 163 observations (10 residual degrees of freedom). AIC can overfit with a small n and a huge candidate set (consider AICc or a smaller candidate set), and R2 says nothing about prediction on new data: validate with cross-validation (Section 8.4).

8.6.2 Exercise — Life satisfaction #3

CautionExercise

Revisit our previous analysis on life satisfaction (case study in the chapter on causal inference, Section 7.7.2). Now, build a predictive model of life satisfaction. If you want, you can also tune this model in a validation or cross-validation split.

Here the goal is prediction, not causal inference, so - unlike the causal analysis of life satisfaction (Section 7.7.2) - we are free to throw in all available predictors and interactions and let a selection/validation procedure trim the model. Because we will not interpret the individual coefficients causally, the Table 2 fallacy (Section 7.4.1) does not apply.

library(EcoData); library(MASS)
dat = soep[complete.cases(soep[,c("lebensz_org","einkommenj1","syear","sex",
        "alter","anz_pers","bildung","erwerb","gesund_org")]), ]

# large candidate model incl. 2-way interactions
fullModel = lm(lebensz_org ~ (sqrt(einkommenj1) + syear + sex + alter +
        anz_pers + bildung + erwerb + gesund_org)^2, data = dat)

predModel = stepAIC(fullModel, trace = 0)   # AIC-optimal predictive model
summary(predModel)

Call:
lm(formula = lebensz_org ~ sqrt(einkommenj1) + syear + sex + 
    alter + anz_pers + bildung + erwerb + gesund_org + sqrt(einkommenj1):gesund_org + 
    syear:alter + syear:bildung + sex:alter + sex:anz_pers + 
    sex:bildung + sex:erwerb + alter:anz_pers + alter:bildung + 
    alter:erwerb + anz_pers:gesund_org + bildung:erwerb + erwerb:gesund_org, 
    data = dat)

Residuals:
    Min      1Q  Median      3Q     Max 
-9.1817 -0.7683  0.1289  0.9417  5.3018 

Coefficients:
                               Estimate Std. Error t value Pr(>|t|)    
(Intercept)                  -2.087e+01  7.699e+01  -0.271 0.786358    
sqrt(einkommenj1)             1.044e-04  4.721e-04   0.221 0.825000    
syear                         1.464e-02  3.818e-02   0.384 0.701337    
sex                           4.441e-01  1.375e-01   3.231 0.001235 ** 
alter                        -2.822e+00  8.246e-01  -3.422 0.000622 ***
anz_pers                      4.779e-02  2.866e-02   1.668 0.095418 .  
bildung                       1.061e+01  5.111e+00   2.075 0.037990 *  
erwerb                       -2.905e-01  4.085e-02  -7.112 1.18e-12 ***
gesund_org                   -8.126e-01  4.962e-02 -16.378  < 2e-16 ***
sqrt(einkommenj1):gesund_org  2.779e-04  1.842e-04   1.509 0.131342    
syear:alter                   1.404e-03  4.089e-04   3.434 0.000596 ***
syear:bildung                -5.238e-03  2.534e-03  -2.067 0.038783 *  
sex:alter                    -5.857e-03  1.400e-03  -4.184 2.87e-05 ***
sex:anz_pers                  2.274e-02  1.569e-02   1.450 0.147098    
sex:bildung                  -1.959e-02  7.571e-03  -2.588 0.009672 ** 
sex:erwerb                    4.077e-02  1.279e-02   3.189 0.001430 ** 
alter:anz_pers               -9.854e-04  5.395e-04  -1.826 0.067807 .  
alter:bildung                -7.471e-04  2.480e-04  -3.013 0.002592 ** 
alter:erwerb                  4.687e-03  4.226e-04  11.090  < 2e-16 ***
anz_pers:gesund_org           2.375e-02  7.601e-03   3.124 0.001787 ** 
bildung:erwerb                9.510e-03  2.323e-03   4.094 4.25e-05 ***
erwerb:gesund_org            -3.223e-02  9.458e-03  -3.408 0.000656 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.471 on 21598 degrees of freedom
Multiple R-squared:  0.2262,    Adjusted R-squared:  0.2255 
F-statistic: 300.7 on 21 and 21598 DF,  p-value: < 2.2e-16

Evaluate the predictive error with cross-validation. Note that this CV error is optimistic, because the model was selected on the same data, and rows of the same person are not independent:

library(boot)
set.seed(1)
cvFit = glm(formula(predModel), data = dat)
causalFit = glm(lebensz_org ~ sqrt(einkommenj1) + syear + sex + alter + anz_pers + bildung + erwerb + gesund_org, data = dat)
c(predModel = cv.glm(dat, cvFit, K = 10)$delta[1],       # 10-fold CV mean squared error
  causalModel = cv.glm(dat, causalFit, K = 10)$delta[1])  # the same for the causal model
  predModel causalModel 
   2.165679    2.187814

Note that we filtered to complete.cases on the modeled variables first, so that stepAIC compares all candidate models on the same rows. Compared with the simpler causal model (Section 7.7.2), the selected model is only slightly better on prediction (about 2.17 vs 2.19), a very small gain for 13 additional coefficients: a more complex model often wins slightly on prediction even where it would be the wrong choice for inference.