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 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