Compare trophic and classical SDMs: difference and equalities
Here, we present the differences between a classic SDM and our trophic SDM, starting from the simplest case when these two approaches coincide to more realistic cases where they do not.
The simplest case
Let’s consider a very simple case with one environmental covariate \(x\), and two species \(y_A\) and \(y_B\). We assume species \(y_A\) to be basal and species \(y_B\) to feed on \(y_A\). This leads to the following directed acyclic graph that resumes the dependencies across variables, assuming bottom-up control.
library(igraph)
G=graph_from_edgelist(as.matrix(data.frame(from = c("x","x","A"),
to = c("A","B","B"))))
plot(G)
We now assume that data on species distributions are continuous
observable entities (e.g. basal area, biomass), that we can model with a
linear regression. We name \(y_{iA}\)
and \(y_{iB}\) the observations of
species \(y_A\) and \(y_B\) at site \(i\), respectively, and we assume that the
environment has a quadratic effect on both species. Therefore, using
classical SDMs, model equations are:
\[
\begin{aligned}
y_{iA} &= \beta^{SDM}_{0A} + \beta^{SDM}_{1A} x_i +
\beta^{SDM}_{2A}x_i^2 + \epsilon_{iA}, \ \ \epsilon_{iA} \sim N(0,
\sigma_A^{SDM}) \\
y_{iB} &= \beta^{SDM}_{0B} + \beta^{SDM}_{1B} x_i +
\beta^{SDM}_{2B}x_i^2 + \epsilon_{iB}, \ \ \epsilon_{iB} \sim N(0,
\sigma_B^{SDM})
\end{aligned}
\]
Instead, our trophic model includes a dependence of \(y_B\) on its prey \(y_A\):
\[
\begin{aligned}
y_{iA} &= \beta^{tSDM}_{0A} + \beta^{tSDM}_{1A} x_i +
\beta^{tSDM}_{2A}x_i^2 + \epsilon_{iA}, \ \ \epsilon_{iA} \sim N(0,
\sigma_A^{tSDM}) \\
y_{iB} &= \beta^{tSDM}_{0B} + \beta^{tSDM}_{1B} x_i +
\beta^{tSDM}_{2B}x_i^2 + \alpha y_{iA} + \epsilon_{iB}, \ \
\epsilon_{iB} \sim N(0, \sigma_B^{tSDM})
\end{aligned}
\]
We hereafter focus on the expectation of the models (corresponding to
the predicted values), focusing in particular on species \(y_B\) for which the two models are
different. For SDMs:
\[E[y_{iB} | \ x_i] = \beta^{SDM}_{0B} +
\beta^{SDM}_{1B} x_i + \beta^{SDM}_{2B}x_i^2 \] For trophic SDM,
we need to condition on a value of \(y_{iA}\). We here choose to condition on
the expected value of \(y_{iA}\), which
corresponds to using predictions from \(y_A\) to predict \(y_B\). In a linear model, with Gaussian
data and a linear interaction between prey and predator, we have that
\(E_A[E_B[y_B | \ y_A]] = E_B[y_B | \
E[y_A]]\) (which is not necessarily true in other cases).
We thus have:
\[E[y_{iB} | \ x_i, E[y_{iA}] ] =
\beta^{tSDM}_{0B} + \beta^{tSDM}_{1B} x_i + \beta^{tSDM}_{2B}x_i^2 +
\alpha E[y_{iA}] \] We can rewrite the expected value of \(y_B\) only as a function of the environment
\(x\) only. \[
\begin{aligned}
E[y_{iB} | \ x_i, E[y_{iA}] ] &= \beta^{tSDM}_{0B} +
\beta^{tSDM}_{1B} x_i + \beta^{tSDM}_{2B}x_i^2 + \alpha
(\beta^{tSDM}_{0A} + \beta^{tSDM}_{1A} x_i + \beta^{tSDM}_{2A}x_i^2 ) \\
& = \beta^{tSDM}_{0B} + \alpha \beta^{tSDM}_{0A} + (
\beta^{tSDM}_{1B} + \alpha \beta^{tSDM}_{1A}) x_i + (\beta^{tSDM}_{2B} +
\alpha \beta^{tSDM}_{2A}) x_i^2 \\
&= \beta^{SDM}_{0B} + \beta^{SDM}_{1B} x_i + \beta^{SDM}_{2B}x_i^2
\end{aligned}
\]
In this very simple case of a linear regression model with a linear
interactions between the predator and the preys, the estimation of the
regression coefficients is such that \[\beta^{SDM}_{.B} = \beta^{tSDM}_{.B} + \alpha
\beta^{tSDM}_{.A} \] As a consequence the models estimate the
same realised niche, and predictions from trophic and classical SDMs
coincide.
However, our trophic model allows a separation of the direct effect of
the environment on \(y_B\) (i.e., the
potential niche of \(y_B\)) from its
indirect effect due to the effect of \(y_A\) on \(y_B\). In other words, using the jargon of
structural equation models, \(\beta^{tSDM}_{1B}\) and \(\beta^{tSDM}_{2B}\) are the direct effects
of the environment on species \(y_B\)
(i.e., its potential niche), the terms \(\alpha \beta^{tSDM}_{1A}\) and \(\alpha \beta^{tSDM}_{2A}\) are the indirect
effect of the environment on species \(y_B\). The sum of these two effects, \(\beta^{tSDM}_{1B} + \alpha
\beta^{tSDM}_{1A}\) and \(\beta^{tSDM}_{2B} + \alpha
\beta^{tSDM}_{2A}\) corresponds to the total effect of the
environment on \(y_B\), i.e., the
realised niche of \(y_B\). In this very
particular case, SDM can correctly estimate the total effect (the
realised niche), but cannot separate the direct and indirect effect of
the environment on \(y_B\) and thus
they do not correctly estimate the potential niche.
Simulation
We hereafter present a short example confirming the analytical derivation of these results.
First, we simulate species distribution data accordingly to the directed acyclic graph presented above. Therefore, the environment has a quadratic effect on species \(y_A\) and species \(y_B\), with a linear effect of \(y_A\) on \(y_B\).
set.seed(1712)
# Create the environmental gradient
X = seq(0,5,length.out=100)
# Simulate data for the basal species A
A = 3*X + 2*X^2 + rnorm(100, sd = 5)
# Simulate data for the predator species B
B = 2*X -3*X^2 + 3*A + rnorm(100, sd = 5)
oldpar = par()
par(mfrow=c(1,3))
plot(X,A, main = "Species A realised \n and potential niche")
plot(X,B, main = "Species B realised niche")
plot(X, 2*X -3*X^2, main = "Species B potential niche")
We now run both a classic SDM and a trophic SDM
# For species A, the two models coincide
m_A = lm(A ~ X + I(X^2))
# For species B, the two models are different:
## trophic model
m_trophic = lm(B ~ X + I(X^2) + A)
## classic SDM
m_SDM = lm(B ~ X + I(X^2))We can now test the relationships between the coefficients of SDM and trophic SDM.
# Show the equality of coefficients
all.equal(coef(m_SDM), coef(m_trophic)[-4] + coef(m_trophic)[4]*coef(m_A))
#> [1] TRUE
# Prediction coincide
pred_SDM = m_SDM$fitted.values
pred_trophic = predict(m_trophic,
newdata = data.frame(X = X, A = m_A$fitted.values))
all.equal(pred_SDM, pred_trophic)
#> [1] TRUEThe coefficient of classic SDM for species \(y_B\) are equal to the environmental coefficients of trophic SDM for species \(y_B\) plus the indirect effect of the environment on species \(y_B\) (given by the product of the biotic term and the environmental coefficients of species \(y_A\)). As a consequence, the predictions of both models coincide when the predicted value of species \(y_A\) is used to predict \(y_B\). However, unlike the traditional SDM, our trophic model allows to retrieve the correct potential niche of species \(y_B\)
# Potential niche infered by SDM
plot(X, m_SDM$fitted.values, ylim = c(-65, 120), ylab = "B",
main = "Potential niche of yB")
points(X, cbind(1, X, X^2) %*% coef(m_trophic)[-4], col ="red")
points(X, 2*X -3*X^2, col ="blue")
legend(-0.1, 120, legend=c("Potential niche inferred by SDM",
"Potential niche inferred by trophic SDM",
"True fundamenta niche"), pch =1,
col=c("black","red", "blue"), cex = 0.5)
