Diversity-interaction model
To disentangle plant guilds effects on SOC, we built a
Diversity-Interaction Model (DIM) (Kirwan et al. 2009) in two
steps. Preliminary normality tests (distribution histogram and normal
Q-Q plot, and the Shapiro Wilk W test (Shapiro et al. 1968))
showed that the variable SOC was not normaly distributed (W = 0.948;
p-value < 0.001). We log-transformed SOC to prevent a breach
of the normality assumption by the residuals of the models (W = 0.99;
p-value = 0.18). From now on, it should be understood that the mention
of SOC in the context of modelling refers always to log(SOC).
As a first step, to control the wide climatic and topographical
variability of the Pyrenees (Garcia-Pausas et al. 2007, 2017) we
fitted a linear model (Legendre & Legendre 1998) for SOC considering
only environmental drivers other than plant guilds. We included both
main effects and pair-wise interactions. We selected several candidate
terms by a semi-automatic procedure by a genetic algorithm included in
the R package glmulti (Calcagno & Mazancourt 2010). We selected the
best model including only main effects of the predictors, using SOC as
the response variable. Then we modelled the residuals of that model with
the same predictors to look for significant candidate level-two
interaction terms. Once we obtained a set of candidate model terms (with
both main effects and interactions), those were included into a single
model on which we performed backward-forward selection. We used several
methods to compare and determine the final model, including the
corrected Akaike information criterion (AICc; (Symonds & Moussalli
2011)), the adjusted determination coefficient
(Radj2) and model comparison
techniques with the “anova()” function in R, using Chi-square tests to
test whether the reduction in the residual sum of squares was
statistically significant.
At the second step, to complete the DIM, we added to our null model the
guild proportion variables and the guild pair-wise interactions. Because
all guild proportions sum 100%, we built a model without intercept at
this step (Kirwan et al. 2009). We discarded curved relationships
between the model terms and SOC by comparing this final DIM model with
the corresponding generalized interaction model (Connolly et al.2013) using a likelihood test (Chi-squared = 1.27 x
10-3; p-value = 0.97; Hothorn et al., 2019).
In order to show the effects of plant guilds we calculated the predicted
effects of the model for all possible combinations of grasses, forbs and
legumes with the emmeans package (Lenth et al. 2019), fixing the
rest of the variables at their mean value. Sedges were fixed at 0%
because of their low effects on SOC (Table 1), their low mean proportion
(2%) and to obtain clearer plots. We built a ternary plot with the
ggtern package (Hamilton & Ferry 2018) to show predicted SOC variation
across plant guild proportions. We also plotted SOC variation across a
legume proportion gradient under: 1) a balanced scenario where grasses
and forbs accounted equal proportions; 2) a forb-dominance scenario
where forb:grass proportions were 80:20; 3) a grass-dominance scenario
were forb:grass proportion was 20:80; and 4) a co-dominance scenario
where grasses and forbs accounted equal proportions (50:50), of the
remaining proportion unaccounted by legumes.