8.9 Problems

Read Section 8.8 for the data and tips on answering the questions and setting up your models. Note the questions asking about the effects on growth rate are asking about the C matrix in \[\mathbf{x}_t=\mathbf{B}\mathbf{x}_{t-1}+\mathbf{C}\mathbf{c}_t+\mathbf{w}_t\] The \(\mathbf{C}\mathbf{c}_t+\mathbf{w}_t\) are the process errors and represent the growth rates (growth above or below what you would expect given \(\mathbf{x}_{t-1}\)). Use your raw data in the MARSS model. You do not need to difference the data to get at the growth rates since the process model is modeling that.

  1. How does month affect the mean phytoplankton population growth rates? Show a plot of the estimated mean growth rate versus month for each taxon using three approaches to estimate the month effect (factor, polynomial, Fourier series). Estimate seasonal effects without any covariate (Temp, TP) effects.

  2. It is likely that both temperature and total phosphorus (TP) affect phytoplankton population growth rates. Using MARSS models, estimate the effect of Temp and TP on growth rates of each taxon. Leave out the seasonal covariates from question 1, i.e. only use Temp and TP as covariates. Make a plot of the point estimates of the Temp and TP effects with the 95% CIs added to the plot. tidy() is an easy way to get the parameters CIs.

  3. Estimate the Temp and TP effects using B="unconstrained".

    1. Compare the \(\mathbf{B}\) matrix for the fit from question 2 and from question 3. Describe the species interactions modeled by the \(\mathbf{B}\) matrix when B="unconstrained". How is it different than the \(\mathbf{B}\) matrix from question 2? Note, you can retrieve the matrix using coef(fit, type="matrix")$B.

    2. Do the Temp and TP effects change when you use B="unconstrained"? Make sure to look at the CIs also.

  4. Using MARSS models, evaluate which (Temp or TP) is the more important driver or if both are important. Again, leave out the seasonal covariates from question 1, i.e. only use Temp and TP as covariates. Compare two approaches: comparison of effect sizes in a model with both Temp and TP and model selection using a set of models.

  5. Evaluate whether the effect of temperature (Temp) on the taxa manifests itself via their underlying physiology (by affecting growth rates and thus abundance) or because physical changes in the water stratification makes them easier/harder to sample in some months. Leave out the seasonal covariates from question 1, i.e. only use Temp and TP as the covariates. For TP, assume it always affects growth rates, never the observation errors.

  6. Is there support for temperature or TP affecting all functional groups’ growth rates the same, or are the effects on one taxon different from another? Make sure to test all possibilities: the Temp and TP effects are the same for all taxa, and one covariate effect is the same across taxa while the other’s effects are unique across taxa.

  7. Compare your results for question 2 using linear regression, by using the lm() function. You’ll need to look at the response of each taxon separately, i.e. one response variable. You can have a multivariate response variable with lm() but the functions will be doing 6 independent linear regressions. In your lm() model, use only Temp and TP (and intercept) as covariates. Compare the estimated effects to those from question 2. How are they different? How is this model different from the model you fit in question 2?

  1. Temp and TP are negatively correlated (cor = -0.66). A common threshold for collinearity in regression models is 0.7. Temp and TP fall below that but are close. One approach to collinearity is sequential regression (Dormann et al. 2013). The first (most influential) covariate is included ‘as is’ and the second covariate appears as the residuals of a regression of the second against the first. The covariates are now orthogonal however the second covariate is conditioned on the first. If we see an effect of the residuals covariate, it is the effect of TP additional to the contribution it already made through its relationship with temperature. Rerun question 2 using sequential regression (see code below).

Make your Temp and TP covariates orthogonal using sequential regression. Do your conclusions about the effects of Temperature and TP change?

Below is code to construct your orthogonal covariates for sequential regression.

fit <- lm(covars[1, ] ~ covars[2, ])
seqcovs <- rbind(covars[1, ], residuals(fit))
avg <- apply(seqcovs, 1, mean)
sd <- sqrt(apply(seqcovs, 1, var, na.rm = TRUE))
seqcovs <- (seqcovs - avg)/sd
rownames(seqcovs) <- c("Temp", "TPresids")
  1. Compare the AICc’s of the 3 seasonal models from question 1 and the 4 Temp/TP models from question 5. What does this tell you about the Temp and TP only models?

  2. We cannot just fit a model with season and Temp plus TP since Temp and TP are highly seasonal. That will cause problems if we have something that explain season (a polynomial) and a covariate that has seasonality. Instead, use sequential regression to fit a model with seasonality, Temp and TP. Use a 3rd order polynomial with the poly() function to create orthogonal season covariates and then use sequential regression (code in problem 8) to create Temp and TP covariates that are orthogonal to your season covariates. Fit the model and compare a model with only season to a model with season and Temp plus TP.

  3. Another approach to looking at effects of covariates which have season cycles is to examine if the seasonal anomalies of the independent variable can be explained by the seasonal anomalies of the dependent variables. In other words, can an unusually high February abundance (higher than expected) be explained by an unusually high or low February temperature? In this approach, you remove season so you do not need to model it (with factor, polynomial, etc). The stl() function can be used to decompose a time series using LOESS. We’ll use stl() since it can handle missing values.

    1. Decompose the Diatom time series using stl() and plot. Use na.action=zoo::na.approx to deal with the NAs. Use s.window="periodic". Other than that you can use the defaults.
    i <- "Diatoms"
    dati <- ts(dat[i, ], frequency = 12)
    a <- stl(dati, "periodic", na.action = zoo::na.approx)
    1. Create dependent variables and covariates that are anomolies by modifying the following code. For the anomaly, you will use the remainder plus the trend. You will need to adapt this code to create the anomalies for Temp and TP and for dat (your data).
    i <- "Diatoms"
    a <- stl(ts(dat[i, ], frequency = 12), "periodic", na.action = zoo::na.approx)
    anom <- a$time.series[, "remainder"] + a$time.series[, "trend"]
    1. Notice that you have simply removed the seasonal cycle from the data. Using the seasonal anomalies (from part b), estimate the effect of Temp and TP on each taxon’s growth rate. You will use the same model as in question 2, but use the seasonal anomalies as data and covariates.
anoms <- matrix(NA, dim(dat)[1] + dim(covars)[1], dim(dat)[2])
rownames(anoms) <- c(rownames(dat), rownames(covars))
for (i in 1:dim(dat)[1]) {
    a <- stl(ts(dat[i, ], frequency = 12), "periodic", na.action = zoo::na.approx)
    anoms[i, ] <- a$time.series[, "remainder"] + a$time.series[, 
for (i in 1:dim(covars)[1]) {
    a <- stl(ts(covars[i, ], frequency = 12), "periodic", na.action = zoo::na.approx)
    anoms[i + dim(dat)[1], ] <- a$time.series[, "remainder"] + 
        a$time.series[, "trend"]