This notebook presents allometric analyses of Height on Girth and Bark thickness on Girth for eucalypt trees from 5 regions: Grampians Vic, South and east of the Hume Highway, the Victorian Mallee and Stirling Range, WA.

knitr::opts_chunk$set(cache = TRUE, dev = c("png", "pdf"))
library(eucs)
library(tidyverse)
library(lme4)

Plot data

Plot some scatter plots.

ggplot(traits_stirling, aes(x = girth_cm, y = height_m, colour = taxon)) +
geom_point(size = 2) +
theme_classic() +
theme(legend.position = "none") +
geom_smooth(method = "lm")

ggplot(data = traits_stirling, aes(girth_cm, height_m)) +
geom_point() +
facet_wrap(~taxon) +
xlab("Girth (cm)") +
ylab("Height (m)") +
ggtitle("Eucalypt height vs. girth")

ggplot(
  data = traits_stirling, aes(x = girth_cm, y = height_m, colour = taxon)
) +
geom_point(size = 2) +
xlab("Girth (cm)") +
ylab("Height (m)") +
scale_x_log10() +
scale_y_log10() +
theme_classic() +
theme(legend.position = "none") +
geom_smooth(method = "lm")

lm1 <- lm(height_m ~ girth_cm, data = traits_stirling)
summary(lm1)

Call:
lm(formula = height_m ~ girth_cm, data = traits_stirling)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.8648 -0.9588 -0.1770  0.5375  8.3571 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 2.883192   0.154696   18.64   <2e-16 ***
girth_cm    0.071951   0.002354   30.57   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.99 on 280 degrees of freedom
  (12 observations deleted due to missingness)
Multiple R-squared:  0.7694,    Adjusted R-squared:  0.7686 
F-statistic: 934.4 on 1 and 280 DF,  p-value: < 2.2e-16
confint(lm1)
                2.5 %     97.5 %
(Intercept) 2.5786768 3.18770728
girth_cm    0.0673172 0.07658404
lm2 <- lm(log(height_m) ~ log(girth_cm), data = traits_stirling)
summary(lm2)

Call:
lm(formula = log(height_m) ~ log(girth_cm), data = traits_stirling)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.75083 -0.14694  0.01071  0.17203  0.59775 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)   -0.16013    0.05592  -2.864  0.00451 ** 
log(girth_cm)  0.54291    0.01661  32.692  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.2672 on 280 degrees of freedom
  (12 observations deleted due to missingness)
Multiple R-squared:  0.7924,    Adjusted R-squared:  0.7917 
F-statistic:  1069 on 1 and 280 DF,  p-value: < 2.2e-16
confint(lm2)
                   2.5 %      97.5 %
(Intercept)   -0.2702127 -0.05005195
log(girth_cm)  0.5102201  0.57559965

Fit a varying intercepts model.

mod1 <- lmer(
  log(height_m) ~ log(girth_cm) + (1 | taxon), data = traits_stirling
)
summary(mod1)
Linear mixed model fit by REML ['lmerMod']
Formula: log(height_m) ~ log(girth_cm) + (1 | taxon)
   Data: traits_stirling

REML criterion at convergence: -105.4

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.79134 -0.63681 -0.03042  0.59425  2.44124 

Random effects:
 Groups   Name        Variance Std.Dev.
 taxon    (Intercept) 0.04233  0.2058  
 Residual             0.02896  0.1702  
Number of obs: 282, groups:  taxon, 33

Fixed effects:
              Estimate Std. Error t value
(Intercept)   -0.16790    0.08537  -1.967
log(girth_cm)  0.55197    0.02367  23.318

Correlation of Fixed Effects:
            (Intr)
lg(grth_cm) -0.899
plot(mod1)  # looks alright, no paterns evident

qqnorm(resid(mod1))
qqline(resid(mod1))  # points fall nicely onto the line - good!

Varying slopes & intecepts.

mod2 <- lmer(
  log(height_m) ~ log(girth_cm) + (1 + log(girth_cm) | taxon),
  data = traits_stirling
)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00310319 (tol = 0.002, component 1)
summary(mod2)
Linear mixed model fit by REML ['lmerMod']
Formula: log(height_m) ~ log(girth_cm) + (1 + log(girth_cm) | taxon)
   Data: traits_stirling

REML criterion at convergence: -111.9

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.78715 -0.48546 -0.00704  0.47148  2.54947 

Random effects:
 Groups   Name          Variance Std.Dev. Corr 
 taxon    (Intercept)   0.26663  0.5164        
          log(girth_cm) 0.02732  0.1653   -0.91
 Residual               0.02508  0.1584        
Number of obs: 282, groups:  taxon, 33

Fixed effects:
              Estimate Std. Error t value
(Intercept)   -0.11893    0.12803  -0.929
log(girth_cm)  0.53564    0.04105  13.048

Correlation of Fixed Effects:
            (Intr)
lg(grth_cm) -0.942
convergence code: 0
Model failed to converge with max|grad| = 0.00310319 (tol = 0.002, component 1)
plot(mod2)  # looks alright, no paterns evident

qqnorm(resid(mod2))
qqline(resid(mod2))

# not as good as varying intercept, use mod1. also too litle infor in the data
coef(mod1)
$taxon
                                             (Intercept) log(girth_cm)
Corymbia calophylla                         -0.360680312     0.5519723
Eucalyptus angulosa                         -0.176727452     0.5519723
Eucalyptus buprestium                       -0.109704089     0.5519723
Eucalyptus conglobata subsp. perata         -0.044896119     0.5519723
Eucalyptus cornuta                           0.009172194     0.5519723
Eucalyptus decipiens                        -0.378571463     0.5519723
Eucalyptus decipiens subsp. chalara         -0.506924623     0.5519723
Eucalyptus decurva                           0.018869050     0.5519723
Eucalyptus doratoxylon                      -0.322403914     0.5519723
Eucalyptus erectifolia                      -0.161827744     0.5519723
Eucalyptus falcata                          -0.121167598     0.5519723
Eucalyptus flocktoniae                       0.110336441     0.5519723
Eucalyptus hebetifolia                      -0.025705386     0.5519723
Eucalyptus incrassata                       -0.139095852     0.5519723
Eucalyptus lehmannii                        -0.250095696     0.5519723
Eucalyptus ligulata subsp. stirlingica      -0.158119203     0.5519723
Eucalyptus loxophleba subsp. loxophleba     -0.235886147     0.5519723
Eucalyptus marginata subsp. marginata       -0.206357975     0.5519723
Eucalyptus megacarpa                        -0.374027782     0.5519723
Eucalyptus occidentalis                     -0.007786424     0.5519723
Eucalyptus pachyloma                         0.036369700     0.5519723
Eucalyptus phaenophylla subsp. phaenophylla -0.223250922     0.5519723
Eucalyptus pleurocarpa                       0.032969541     0.5519723
Eucalyptus pluricaulis                      -0.145059354     0.5519723
Eucalyptus preissiana                       -0.264848126     0.5519723
Eucalyptus rudis                             0.296019031     0.5519723
Eucalyptus staeri                           -0.673095917     0.5519723
Eucalyptus talyuberlup                       0.057843693     0.5519723
Eucalyptus thamnoides                       -0.480826944     0.5519723
Eucalyptus uncinata                         -0.301658313     0.5519723
Eucalyptus vergrandis                       -0.060373977     0.5519723
Eucalyptus wandoo                           -0.153475963     0.5519723
Eucalyptus xanthonema                       -0.219693530     0.5519723

attr(,"class")
[1] "coef.mer"
confint(mod1)
Computing profile confidence intervals ...
                   2.5 %      97.5 %
.sig01         0.1569273  0.26627377
.sigma         0.1560433  0.18602399
(Intercept)   -0.3347535 -0.00103789
log(girth_cm)  0.5057015  0.59823878
rHt <- traits_stirling$height_m / traits_stirling$girth_cm
rBT <- traits_stirling$bark_thickness_mm / traits_stirling$girth_cm
ggplot(data = traits_stirling) +
aes(x = girth_cm, y = bark_thickness_mm, colour = taxon) +
geom_point(size = 2) +
xlab("Girth (cm)") +
ylab("Bark (mm)") +
scale_x_log10() +
scale_y_log10() +
theme_classic() +
theme(legend.position = "none") +
geom_smooth(method = "lm")

Fit a varying intercepts model for Bark Thickness on Girth.

mod3 <- lmer(
  log(bark_thickness_mm) ~ log(girth_cm) + (1 | taxon), data = traits_stirling
)
summary(mod3)
Linear mixed model fit by REML ['lmerMod']
Formula: log(bark_thickness_mm) ~ log(girth_cm) + (1 | taxon)
   Data: traits_stirling

REML criterion at convergence: 181.3

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.47503 -0.46307  0.09969  0.51555  2.45687 

Random effects:
 Groups   Name        Variance Std.Dev.
 taxon    (Intercept) 0.13030  0.361   
 Residual             0.08011  0.283   
Number of obs: 279, groups:  taxon, 33

Fixed effects:
              Estimate Std. Error t value
(Intercept)   -0.72974    0.14578  -5.006
log(girth_cm)  0.71179    0.04013  17.735

Correlation of Fixed Effects:
            (Intr)
lg(grth_cm) -0.894
IycgLS0tCiMnIHRpdGxlOiAiRXVjYWx5cHQgYWxsb21ldHJ5IgojJyBhdXRob3I6ICJXaWxsaWFtIEsuIE1vcnJpcyAmIFBldGVyIEEuIFZlc2siCiMnIGRhdGU6ICJgciBTeXMuRGF0ZSgpYCIKIycgb3V0cHV0OgojJyAgIHJtYXJrZG93bjo6aHRtbF9ub3RlYm9vazoKIycgICAgIGNvZGVfZm9sZGluZzogaGlkZQojJyAtLS0KCiMnIFRoaXMgbm90ZWJvb2sgcHJlc2VudHMgYWxsb21ldHJpYyBhbmFseXNlcyBvZiBIZWlnaHQgb24gR2lydGggYW5kIEJhcmsKIycgdGhpY2tuZXNzIG9uIEdpcnRoIGZvciBldWNhbHlwdCB0cmVlcyBmcm9tIDUgcmVnaW9uczogR3JhbXBpYW5zIFZpYywgU291dGgKIycgYW5kIGVhc3Qgb2YgdGhlIEh1bWUgSGlnaHdheSwgdGhlIFZpY3RvcmlhbiBNYWxsZWUgYW5kIFN0aXJsaW5nIFJhbmdlLCBXQS4KCiMrIHBrZ3MsIG1lc3NhZ2UgPSBGQUxTRQprbml0cjo6b3B0c19jaHVuayRzZXQoY2FjaGUgPSBUUlVFLCBkZXYgPSBjKCJwbmciLCAicGRmIikpCmxpYnJhcnkoZXVjcykKbGlicmFyeSh0aWR5dmVyc2UpCmxpYnJhcnkobG1lNCkKCiMnICMgUGxvdCBkYXRhCiMnIFBsb3Qgc29tZSBzY2F0dGVyIHBsb3RzLgojKyBwbG90LWhlaWdodC1vbi1naXJ0aCwgbWVzc2FnZSA9IEZBTFNFLCB3YXJuaW5nID0gRkFMU0UKZ2dwbG90KHRyYWl0c19zdGlybGluZywgYWVzKHggPSBnaXJ0aF9jbSwgeSA9IGhlaWdodF9tLCBjb2xvdXIgPSB0YXhvbikpICsKZ2VvbV9wb2ludChzaXplID0gMikgKwp0aGVtZV9jbGFzc2ljKCkgKwp0aGVtZShsZWdlbmQucG9zaXRpb24gPSAibm9uZSIpICsKZ2VvbV9zbW9vdGgobWV0aG9kID0gImxtIikKCiMrIHRyeS1wbG90LWZhY2V0LCBtZXNzYWdlID0gRkFMU0UsIHdhcm5pbmcgPSBGQUxTRQpnZ3Bsb3QoZGF0YSA9IHRyYWl0c19zdGlybGluZywgYWVzKGdpcnRoX2NtLCBoZWlnaHRfbSkpICsKZ2VvbV9wb2ludCgpICsKZmFjZXRfd3JhcCh+dGF4b24pICsKeGxhYigiR2lydGggKGNtKSIpICsKeWxhYigiSGVpZ2h0IChtKSIpICsKZ2d0aXRsZSgiRXVjYWx5cHQgaGVpZ2h0IHZzLiBnaXJ0aCIpCgojKyBsb2dfaHQtb24tbG9nX2dpcnRoLCBtZXNzYWdlID0gRkFMU0UsIHdhcm5pbmcgPSBGQUxTRQpnZ3Bsb3QoCiAgZGF0YSA9IHRyYWl0c19zdGlybGluZywgYWVzKHggPSBnaXJ0aF9jbSwgeSA9IGhlaWdodF9tLCBjb2xvdXIgPSB0YXhvbikKKSArCmdlb21fcG9pbnQoc2l6ZSA9IDIpICsKeGxhYigiR2lydGggKGNtKSIpICsKeWxhYigiSGVpZ2h0IChtKSIpICsKc2NhbGVfeF9sb2cxMCgpICsKc2NhbGVfeV9sb2cxMCgpICsKdGhlbWVfY2xhc3NpYygpICsKdGhlbWUobGVnZW5kLnBvc2l0aW9uID0gIm5vbmUiKSArCmdlb21fc21vb3RoKG1ldGhvZCA9ICJsbSIpCgojKyBsbTEtSC1vbi1HCmxtMSA8LSBsbShoZWlnaHRfbSB+IGdpcnRoX2NtLCBkYXRhID0gdHJhaXRzX3N0aXJsaW5nKQpzdW1tYXJ5KGxtMSkKY29uZmludChsbTEpCgojKyBsbTItbG9nSC1vbi1sb2dHCmxtMiA8LSBsbShsb2coaGVpZ2h0X20pIH4gbG9nKGdpcnRoX2NtKSwgZGF0YSA9IHRyYWl0c19zdGlybGluZykKc3VtbWFyeShsbTIpCmNvbmZpbnQobG0yKQoKIycgRml0IGEgdmFyeWluZyBpbnRlcmNlcHRzIG1vZGVsLgojKyBtb2QxLUgtb24tRwptb2QxIDwtIGxtZXIoCiAgbG9nKGhlaWdodF9tKSB+IGxvZyhnaXJ0aF9jbSkgKyAoMSB8IHRheG9uKSwgZGF0YSA9IHRyYWl0c19zdGlybGluZwopCnN1bW1hcnkobW9kMSkKCiMrIGRpYWcKcGxvdChtb2QxKSAgIyBsb29rcyBhbHJpZ2h0LCBubyBwYXRlcm5zIGV2aWRlbnQKcXFub3JtKHJlc2lkKG1vZDEpKQpxcWxpbmUocmVzaWQobW9kMSkpICAjIHBvaW50cyBmYWxsIG5pY2VseSBvbnRvIHRoZSBsaW5lIC0gZ29vZCEKCiMnIFZhcnlpbmcgc2xvcGVzICYgaW50ZWNlcHRzLgojKyBtb2QyLUgtb24tRy12YXJ5aW5nLXNsb3Blcy1hbmQtaW50ZXJjZXB0cwptb2QyIDwtIGxtZXIoCiAgbG9nKGhlaWdodF9tKSB+IGxvZyhnaXJ0aF9jbSkgKyAoMSArIGxvZyhnaXJ0aF9jbSkgfCB0YXhvbiksCiAgZGF0YSA9IHRyYWl0c19zdGlybGluZwopCnN1bW1hcnkobW9kMikKCiMrIGRpYWcxfQpwbG90KG1vZDIpICAjIGxvb2tzIGFscmlnaHQsIG5vIHBhdGVybnMgZXZpZGVudApxcW5vcm0ocmVzaWQobW9kMikpCnFxbGluZShyZXNpZChtb2QyKSkKIyBub3QgYXMgZ29vZCBhcyB2YXJ5aW5nIGludGVyY2VwdCwgdXNlIG1vZDEuIGFsc28gdG9vIGxpdGxlIGluZm9yIGluIHRoZSBkYXRhCgojKyBjb2Vmc01vZDEKY29lZihtb2QxKQpjb25maW50KG1vZDEpCgojKyBjYWxjLWh0LW9uLWdpcnRoCnJIdCA8LSB0cmFpdHNfc3RpcmxpbmckaGVpZ2h0X20gLyB0cmFpdHNfc3RpcmxpbmckZ2lydGhfY20KckJUIDwtIHRyYWl0c19zdGlybGluZyRiYXJrX3RoaWNrbmVzc19tbSAvIHRyYWl0c19zdGlybGluZyRnaXJ0aF9jbQoKIysgcGxvdC1CVC1vbi1HLCBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nID0gRkFMU0UKZ2dwbG90KGRhdGEgPSB0cmFpdHNfc3RpcmxpbmcpICsKYWVzKHggPSBnaXJ0aF9jbSwgeSA9IGJhcmtfdGhpY2tuZXNzX21tLCBjb2xvdXIgPSB0YXhvbikgKwpnZW9tX3BvaW50KHNpemUgPSAyKSArCnhsYWIoIkdpcnRoIChjbSkiKSArCnlsYWIoIkJhcmsgKG1tKSIpICsKc2NhbGVfeF9sb2cxMCgpICsKc2NhbGVfeV9sb2cxMCgpICsKdGhlbWVfY2xhc3NpYygpICsKdGhlbWUobGVnZW5kLnBvc2l0aW9uID0gIm5vbmUiKSArCmdlb21fc21vb3RoKG1ldGhvZCA9ICJsbSIpCgojJyBGaXQgYSB2YXJ5aW5nIGludGVyY2VwdHMgbW9kZWwgZm9yIEJhcmsgVGhpY2tuZXNzIG9uIEdpcnRoLgojKyBtb2QzLUItb24tRwptb2QzIDwtIGxtZXIoCiAgbG9nKGJhcmtfdGhpY2tuZXNzX21tKSB+IGxvZyhnaXJ0aF9jbSkgKyAoMSB8IHRheG9uKSwgZGF0YSA9IHRyYWl0c19zdGlybGluZwopCnN1bW1hcnkobW9kMykK