7.8 Set up the hypotheses as different models
Only the \(\mathbf{Z}\) matrices change for our model. We will set up a base model list used for all models.
<- list(
mod.list B = "identity",
U = "unequal",
Q = "equalvarcov",
Z = "placeholder",
A = "scaling",
R = "diagonal and equal",
x0 = "unequal",
tinitx = 0
)
Then we set up the \(\mathbf{Z}\) matrices using the factor short-cut.
<- list(
Z.models H1 = factor(c("pnw", "pnw", rep("ps", 4), "ca", "ca", "pnw", "pnw", "ps")),
H2 = factor(c(rep("coast", 2), rep("ps", 4), rep("coast", 4), "ps")),
H3 = factor(c(rep("N", 6), "S", "S", "N", "S", "N")),
H4 = factor(c("nc", "nc", "is", "is", "ps", "ps", "sc", "sc", "nc", "sc", "is")),
H5 = factor(rep("pan", 11)),
H6 = factor(1:11) # site
)names(Z.models) <-
c("stock", "coast+PS", "N+S", "NC+strait+PS+SC", "panmictic", "site")
7.8.1 Fit the models
We loop through the models, fit and store the results:
<- NULL
out.tab <- list()
fits for (i in 1:length(Z.models)) {
$Z <- Z.models[[i]]
mod.list<- MARSS::MARSS(sealData, model = mod.list, silent = TRUE,
fit control = list(maxit = 1000))
<- data.frame(H = names(Z.models)[i], logLik = fit$logLik,
out AICc = fit$AICc, num.param = fit$num.params, m = length(unique(Z.models[[i]])),
num.iter = fit$numIter, converged = !fit$convergence)
<- rbind(out.tab, out)
out.tab <- c(fits, list(fit))
fits }
We will use AICc and AIC weights to summarize the data support for the different hypotheses. First we will sort the fits based on AICc:
<- order(out.tab$AICc)
min.AICc .1 <- out.tab[min.AICc, ] out.tab
Next we add the \(\Delta\)AICc values by subtracting the lowest AICc:
.1 <- cbind(out.tab.1, delta.AICc = out.tab.1$AICc - out.tab.1$AICc[1]) out.tab
Relative likelihood is defined as \(\,\text{exp}(-\Delta \mathrm{AICc}/2)\).
.1 <- cbind(out.tab.1, rel.like = exp(-1 * out.tab.1$delta.AICc/2)) out.tab
The AIC weight for a model is its relative likelihood divided by the sum of all the relative likelihoods.
.1 <- cbind(out.tab.1, AIC.weight = out.tab.1$rel.like/sum(out.tab.1$rel.like)) out.tab
Let’s look at the model weights (out.tab.1
):
H delta.AICc AIC.weight converged
NC+strait+PS+SC 0.00 0.979 TRUE
site 7.65 0.021 TRUE
N+S 36.97 0.000 TRUE
stock 47.02 0.000 TRUE
coast+PS 48.78 0.000 TRUE
panmictic 71.67 0.000 TRUE