It is not a common practice to fertilize trees, but some situation arises that demands the fertilization of trees. Nutrient deficiency is a big one, but, isn’t that the reason why we fertilize anything–to make up for deficient nutrient. Well, fertilizing trees destabilizes them since they are evolved to thrive and survive in certain soil fertility spectrum. The detriments to fertilizing trees could include:
That’s why the advice of undertaking a soil test before planting should be adhered to.
This all is not to say you shouldn’t fertilize trees, just do so after a soil test to see if its necessary. In the northern part of Sweden fertilization is common, it is however prohibited in the South. Swedish forest fertilization mainly involves the application of Nitrogen which is normally the limiting nutrient for high stand growth.
Different fertilization experiments have been set see the effect of fertilization application frequency on the stand development.
The data used in this post is from an experiment set to see the effect of fertilization frequency on the growth of a stand. The experiment is a young Norway Spruce stands which was established with 5 blocks having randomly distributed treatments in 0.1 ha plots. The treatments are with 3 different intensities in fertilization
F1: Fertilized every year
F2: Fertilized every second year
F3: Fertilized every third year
C: Control without fertilzation.The experiment was measured initially in year 1972, first measurement revision (revision 1) and there after there were several revisions, but the important revisions are the focus here which is an interval of 5 years period (rev 1, 4, 5, and 6). This means that the difference between revision 1 and 4 is 15 years, then the addition of the usual 5 years interval.
For the post the main thing is to check if there’s any difference between the different treatments. Let’s start by importing the necessary packages to be used in this project.
library(pacman)
p_load(dplyr, ggplot2, readr, ggthemes, ggtext, tidyr)
my_theme <- theme_clean() +
theme(
plot.title = element_markdown(hjust = .75, size = 14, face = "plain"),
plot.subtitle = element_text(size = 9),
legend.title = element_text(size = 9),
legend.background = element_blank()
)
theme_set(my_theme)The data is available here
fert_raw <- read_delim(
"https://raw.githubusercontent.com/xrander/SLU-Plantation-Experimentation/master/Data/Lab6/expfert.txt",
delim = "\t"
)
fert_raw |>
head()# A tibble: 6 × 7
block revision treatment volume CAI domheight age
<dbl> <dbl> <chr> <dbl> <dbl> <dbl> <dbl>
1 1523 1 C 11.3 0 5.77 20
2 1523 1 F1 8.07 0 5.37 20
3 1523 1 F2 8.13 0 5 20
4 1523 1 F3 7.53 0 5.2 20
5 1523 4 C 28.8 6.43 7.83 25
6 1523 4 F1 43.1 13.7 8.53 25Let’s see what we can understand from the data.
str(fert_raw)spc_tbl_ [80 × 7] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
$ block : num [1:80] 1523 1523 1523 1523 1523 ...
$ revision : num [1:80] 1 1 1 1 4 4 4 4 5 5 ...
$ treatment: chr [1:80] "C" "F1" "F2" "F3" ...
$ volume : num [1:80] 11.27 8.07 8.13 7.53 28.83 ...
$ CAI : num [1:80] 0 0 0 0 6.43 ...
$ domheight: num [1:80] 5.77 5.37 5 5.2 7.83 ...
$ age : num [1:80] 20 20 20 20 25 25 25 25 30 30 ...
- attr(*, "spec")=
.. cols(
.. block = col_double(),
.. revision = col_double(),
.. treatment = col_character(),
.. volume = col_double(),
.. CAI = col_double(),
.. domheight = col_double(),
.. age = col_double()
.. )
- attr(*, "problems")=<externalptr> The variables block, revision, and treatment should be factors as they are categorical.
fert_tbl <- fert_raw |>
mutate(
across(block:treatment, as.factor),
age = as.integer(age)
)summary(fert_tbl) block revision treatment volume CAI
1523:16 1:20 C :20 Min. : 1.767 Min. : 0.00
1524:16 4:20 F1:20 1st Qu.: 14.667 1st Qu.: 1.75
1525:16 5:20 F2:20 Median : 54.050 Median :11.68
1526:16 6:20 F3:20 Mean : 86.636 Mean :10.76
1527:16 3rd Qu.:134.154 3rd Qu.:17.71
Max. :280.400 Max. :28.57
domheight age
Min. : 0.000 Min. :13.0
1st Qu.: 5.667 1st Qu.:20.0
Median : 9.200 Median :25.0
Mean : 8.926 Mean :26.1
3rd Qu.:12.858 3rd Qu.:30.0
Max. :16.867 Max. :35.0 From the data, there are no missing data, and the number of replications as seen with revision and treatment seems to match. This will be better understood if presented in another way.
with(
fert_tbl,
ftable(block, revision, treatment)
) treatment C F1 F2 F3
block revision
1523 1 1 1 1 1
4 1 1 1 1
5 1 1 1 1
6 1 1 1 1
1524 1 1 1 1 1
4 1 1 1 1
5 1 1 1 1
6 1 1 1 1
1525 1 1 1 1 1
4 1 1 1 1
5 1 1 1 1
6 1 1 1 1
1526 1 1 1 1 1
4 1 1 1 1
5 1 1 1 1
6 1 1 1 1
1527 1 1 1 1 1
4 1 1 1 1
5 1 1 1 1
6 1 1 1 1The experiment setup shows that we have 5 blocks number 1523 to 1526, with each having the 4 revisions on 4 treatments. This is a classic randomized block design.
Let’s get a sense of the volume distribution. How height and volume developed overtime is presented in Figure 2.
fert_tbl |>
ggplot(aes(volume)) +
geom_density() 
It should be noted that this are repeated measurements on the same experimental unit. Let’s investigate to see if the repeated measurements of volume are correlated
fert_tbl |>
pivot_wider(
id_cols = c(block, treatment),
names_from = revision,
values_from = volume,
names_glue = "revision_{revision}_vol"
) |>
select(where(is.double)) |>
GGally::ggpairs(
title = "Correlation between Repeated Volume Measurements",
xlab = "Volume"
)
Figure 1 shows volume to correlated to each other with the exception of revision_1 which is weak (rev 5 and 6) to moderate (rev 4).
vol_plot <- fert_tbl |>
ggplot(aes(age, volume, colour = treatment)) +
geom_line() +
scale_color_colorblind() +
labs(
x = "Age",
y = expression(paste("Volume (", m^{3}/ha, ")")),
title = "Growth m^3^/ha of *Picea Abies* across Four Fertilization Treatments",
subtitle = "Growth varies across blocks for fertilized treatments, control treatments produces less volume"
) +
facet_wrap(~block)
height_plot <- fert_tbl |>
ggplot(aes(age, domheight, colour = treatment)) +
geom_line() +
geom_point() +
scale_color_colorblind() +
labs(
x = "Age",
y = "Height (m)",
title = "Height Development of *P. Abies* for Four Fertilization Treatments"
) +
facet_wrap(~block)
vol_plot
height_plot
We have an insight of how the treatments influence growth of the trees across the blocks, Figure 2. Let’s investigate volume and height distribution across the treatments.
fert_tbl |>
ggplot(aes(treatment, volume)) +
geom_boxplot(
fill = "whitesmoke",
col = "maroon4"
) +
labs(
x = "Treatment",
y = expression(paste("Volume (", m^{3}/ha, ")")),
title = "Volume distribution of *P. Abies* for the different treatments"
)
fert_tbl |>
ggplot(aes(treatment, domheight)) +
geom_boxplot(
fill = "whitesmoke",
col = "seagreen"
) +
labs(
x = "Treatment",
y = "Height (m)",
title = "Height distribution of *P. Abies* for the different treatments"
)
Figure 3 shows there’s a marked difference between the treatment. This will be tested.
Given the setup of the research, and the result from Figure 1, simple ANOVA isn’t the right choose to answer the question posed. ANOVA can only be used if a particular revisions is investigated, instead for this, we will use a Linear Mixed Effect Model with the effects:
Fixed effects: treatment, revision, and their interaction (treatment × revision will answer the question do fertilization regimes differ in volume development?).
Random effects: block (to handle repeated measures across revisions).
We will use the lme4 package.
p_load(lme4)fert_mod <- lmer(volume ~ treatment * revision + (1 | block), data = fert_tbl)
summary(fert_mod)Linear mixed model fit by REML ['lmerMod']
Formula: volume ~ treatment * revision + (1 | block)
Data: fert_tbl
REML criterion at convergence: 570.8
Scaled residuals:
Min 1Q Median 3Q Max
-2.50043 -0.55433 -0.06999 0.38196 3.00524
Random effects:
Groups Name Variance Std.Dev.
block (Intercept) 204.7 14.31
Residual 247.8 15.74
Number of obs: 80, groups: block, 5
Fixed effects:
Estimate Std. Error t value
(Intercept) 5.9467 9.5133 0.625
treatmentF1 -0.9833 9.9565 -0.099
treatmentF2 0.0200 9.9565 0.002
treatmentF3 -1.4733 9.9565 -0.148
revision4 17.5333 9.9565 1.761
revision5 52.9200 9.9565 5.315
revision6 103.0933 9.9565 10.354
treatmentF1:revision4 20.7233 14.0807 1.472
treatmentF2:revision4 22.2667 14.0807 1.581
treatmentF3:revision4 12.2800 14.0807 0.872
treatmentF1:revision5 75.4767 14.0807 5.360
treatmentF2:revision5 72.8900 14.0807 5.177
treatmentF3:revision5 44.9067 14.0807 3.189
treatmentF1:revision6 133.3233 14.0807 9.469
treatmentF2:revision6 130.2433 14.0807 9.250
treatmentF3:revision6 94.4867 14.0807 6.710
Correlation matrix not shown by default, as p = 16 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need itThe result shows that there are differences in volume between the blocks, Std. Dev. of 14.31. We can estimate the ICC of the blocks by measuring the proportion of the volume that is due to the differences within groups vs between groups.
icc <- 14.31^2 / ( (14.31^2) + (15.74^2) )
icc[1] 0.4525202About 45% of the variance is due to block effect. The effect of the treatments at revision one is mostly negative (-0.98 for F1 and -1.47 for F3) or almost similar to control treatment (0.02 for F1) These are small differences and changes was more pronounced overtime than in the first 15 years. revision4 to revision6 shows the growth of control overtime.
The most interesting findings is at the interaction terms, it shows how the treatments changes growth relative to control at each revision. F3 consistently gave the lowest addition of growth over control when compared to the other treatments at each revisions.
F1 added the highest growth ever at revision6 with about 133 m3/ha added, while F2 added a similar growth with 130 m3/ha. This shows that up to 30 years fertilization matters.
Fertilization significantly boosts growth and the effect gets stronger as years pass.