Ex 4: Killer Whale

Here we step through the Killer Whale Example using the script version of MixSIAR. For a demonstration using the GUI version, see the MixSIAR Manual. For a thorough walkthrough of how to use MixSIAR in a script, see the Wolves Example, which provides more commentary and explanation.

For a clean, runnable .R script, look at mixsiar_script_killerwhale.R in the example_scripts folder of the MixSIAR package install:

library(MixSIAR)
mixsiar.dir <- find.package("MixSIAR")
paste0(mixsiar.dir,"/example_scripts")

You can run the killer whale example script directly with:

source(paste0(mixsiar.dir,"/example_scripts/mixsiar_script_killerwhale.R"))

Killer Whale Example

The Killer Whale Example demonstrates the difference informative priors make, and illustrates how to construct them:

Load MixSIAR package

Load mixture data

See ?load_mix_data for details.

The killer whale consumer data has no covariates—we are estimating the diet of one population where all individuals are considered identical (factors=NULL, fac_random=NULL, fac_nested=NULL, cont_effects=NULL).

Load source data

See ?load_source_data for details.

We do not have any fixed/random/continuous effects or concentration dependence in the source data (source_factors=NULL, conc_dep=FALSE). We only have source means, SD, and sample size—not the original “raw” data (data_type="means").

Load discrimination data

See ?load_discr_data for details.

Plot data

This is your chance to check:

Notice that the high degree of correlation places the sources on a line, which makes it difficult for the model to determine if killer whales are eating Chinook vs. Coho, or Sockeye vs. Steelhead vs. Chum. We can give the model additional information to resolve this problem by using an informative prior.

Run the model with uninformative prior

Run the model with informative prior

From Semmens et al. (in prep):

One of the benefits to conducting mixture models in a Bayesian framework is that information from other data sources can be included via informative prior distributions Moore and Semmens 2008. Suppose we have collected 14 fecal samples from the killer whale population we are studying, and we find 10 chinook, 1 chum, 0 coho, 0 sockeye, and 3 steelhead. The sum of the Dirichlet hyperparameters roughly correspond to prior sample size, so one approach would be to construct a prior with \(\alpha\) = (10, 1, 0, 0, 3). A downside of this prior is that a sample size of 14 represents a very informative prior, with much of the parameter space given very little weight. Keeping the relative contributions the same, the hyperparameters can be rescaled to have the same mean, but different variance. An alternative is to scale the prior to have a weight of 5, equal to the same weight as the “uninformative” prior \(\alpha\) = (1, 1, 1, 1, 1). This prior could be represented as \(\alpha=\left(\frac{10*5}{14},\frac{1*5}{14},\frac{0*5}{14},\frac{0*5}{14},\frac{3*5}{14}\right)=\left(3.57,0.36,0.01,0.01,1.07\right)\).

You cannot have \(\alpha\) parameters = 0, so while we had no coho and sockeye in our fecal samples, we set their α parameters = 0.01. Remember that the sources are sorted alphabetically.

The plot_prior function shows your prior next to the uninformative prior. Compare the model results with the uninformative prior with those from the informative prior.