MAR 580: Models for Marine Ecosystem-Based Management

TMB workshop

Session III

Programming bits, State-space models

​

08 September 2022

Review

• Mixed effects models
• Poisson GLMM

Programming bits�(The “if” statement)

​

• The if statement:

if (condition)

Statements-1;

else

Statements-2;

• Combining statements:

{ statement; statement; }

​

​

Constructing conditions

• RELATIONAL OPERATORS:

== Equal to (NOT =)

<= Less than or equal to

>= Greater than or equal to

< Less than

> Greater than

!= Not equal to

​

• LOGICAL OPERATORS:

&& Condition is true if both sub-conditions are true

|| Condition is true if one of the sub-conditions is true

! Condition is true if sub-condition is false

​

• Note: TMB may not work if a parameter or derived variable appears in a condition.

​

The “for” statement

for (initial state; terminal condition; increment)

{ statements; }

​

• Initial state: Usually involves setting a “loop control variable” to some initial value.
• Terminal condition: The loop will run until this condition is TRUE (note it need not involve the “loop control variable”).
• Increment: This usually involves updating the “loop control variable” (e.g. X++; X--).
• Parameters and derived variables should not appear in the terminal condition.

​

biomass(1) = 1;

for (year=0; year<nyear; year++)

biomass(year+1)=biomass(year)*r;

Vector & matrix functions

http://kaskr.github.io/adcomp/matrix_arrays_8cpp-example.html

​

Extracting from arrays and matrices:

vector = matrix.col(index)

vector = matrix.row(index)

vector = matrix.diagonal();

Subvectors: extract a subset of a vector:

vector2 = vector1(index1,index2)

NB: arrays start at 0, so index1=2 is element 3!

​

Process error modeling in fisheries

• Measurement OR Process error estimators
• Total error estimators (errors in variables models)
• e.g. linear model: error in both the X’s and the Y’s
• Penalized likelihood models
• Penalty added to likelihood function that describes the deviation of a parameter from distribution, estimated like a normal parameter.
• e.g. recruitment deviations around a stock-recruit curve.
• Many statistical stock assessment models use this approach to obtain MLEs (integration over the random effects is achieved when doing Bayesian estimation).
• State-space models
• Formally account for both process and observation error.

​

• Reading: Comparing among error models, Punt (2003)

State-space models

• Special case of mixed-effects models
• As with mixed-effects models, the aim is to partition variance.
• Two discrete time processes running in parallel.
• State process model
• Observation model

​

State-space modeling, hidden process

Observations

Fish population

Time

t = 0

t = 3

t = 2

t = 1

?

State-space modeling, hidden process

Observations

Fish population

Time

t = 0

t = 3

t = 2

t = 1

?

State-space modeling, hidden process

Observations

Fish population

Time

t = 0

t = 3

t = 2

t = 1

?

We are interested in separating variability associated with population process from that due to the observations.

State-space model: exponential growth

• State process model

​

​

​

​

• Observation model

State-space model: exponential growth

• State process model

​

​

​

​

• Observation model

Population growth rate

Population size at time t

Annual deviation in population growth, process error

Process error variance

State-space model: exponential growth

• State process model

​

​

​

​

• Observation model

Population growth rate

Population size at time t

Annual deviation in population growth, process error

Observationat time t

Constant relating population size to observed values

Observation error residual

Process error variance

Observation error variance

State-space models: exponential growth

• Lamda = 1.05, tau=2, sigma=4, n0=10
• Deterministic dynamics

State-space models: exponential growth

• Lamda = 1.05, tau=2, sigma=4, n0=10
• Process error, one realization

State-space models: exponential growth

• Lamda = 1.05, tau=2, sigma=4, n0=10
• Process error, 10 realizations

State-space models: exponential growth

• Lamda = 1.05, tau=2, sigma=4, n0=10
• Observations around the process

​

We might be interested in estimating the values for these parameters!

State-space models, more generally

• Process and Observation modeled as conditional probability density functions.
• State process model describes change in state vector from time t-1 to t.

​

• Observation model relates the state to the observed data, given values for the states.

​

• are the model parameters.
• We also need the initial state,

Relationship between models and parameters

Thorson & Minto (2015)

Inference for state-space models

• Joint distribution for the observations and states:

​

​

Probability of the initial state.

Product over time steps t=1 to t=T.

Probability of state at time t conditional on state at time t-1 and the parameters of the state process model.

Probability of the observation at time t conditional on the state at time t and the parameters of the observation model.

Inference for state-space models

• Joint distribution for the observations and states:

​

​

• We are interested in the marginal distribution of the observations,

​

​

Maximum likelihood estimates for the model parameters

Likelihood of the parameters given the data from time t=1 to time=T.

Inference for state-space models

• Joint distribution for the observations and states:

​

​

• We are interested in the marginal distribution of the observations,

​

​

• As before with mixed-effects models, we achieve this by integrating out the states in the joint distribution.

Estimating state-space models

• Various methods exist for obtaining the (log) likelihood, and then to maximize it.
• Kalman filter (Kalman, 1960)

Seminal work, for linear state and observation models when process and observation errors are Gaussian.

• Additional filters relaxing some assumptions

(Extended Kalman filter, unscented Kalman filter, etc.)

• Expectation-Maximization (EM) algorithm
• Laplace approximation (e.g. TMB)
• Importance sampling
• Estimation can also be achieved in a Bayesian framework. (e.g. with Markov Chain Monte Carlo)

​

Estimating state-space models in TMB

• Take advantage of the conditional distribution for the states:
• i.e. make the state variables the random effects, e.g.

​

• For a Schaefer biomass dynamics model:

​

• Conditioned on the biomass at the previous time step, this model is separable, i.e. ln(B_t) does not depend on ln(B_t-2).

​

Example: �SE Alaska Steller sea lion pup production

• Data: Counts at Forrester Island Complex
• exponential growth model

Exercise: Linear state-space production model

• The Gompertz population dynamics model:

​

​

• Similar to other models, except density dependence is on ln(B) rather than B.
• Application to �lingcod survey data�(Schnute 1994)

​

​

​

​

​

​

​

​

​

​

Schnute, J. T. (1994). A general framework for developing sequential fisheries models. Canadian Journal of Fisheries and Aquatic Sciences, 51(8), 1676-1688.

Initial biomass

• Bt = Bt exp(r + beta*ln(Bt))

LN(Bt) = LN(Bt) + r + beta*LN(Bt)

LN(Bt) = -r/beta

​

Pseudocode for cpp

• Read in data
• Length of data
• Declare parameters
• R, beta, Obs variance, proc variance, random_effects
• Transformation of beta
• Declare Vector of logB
• Populate logB(1) with –r/beta
• Loop over time (Population update equation)
• Likelihood function
• Observations. [normal – dnorm()]
• Process errors. [normal – dnorm()

​

Take homes

• State-space approach provides robust framework for simultaneously considering uncertainty associated with population processes and observations.
• Increasingly being used for management advice.(e.g. SAM, age-structured state-space fisheries model (Nielsen & Berg 2014) ).
• Readily handles missing data.
• Easy to extend the univariate models considered here:
• Multiple state variables for population stages
• Parameters or environmental covariates modeled as state variables, and include process model for these.
• Hierarchical approach to share information among time series (e.g. common trends, Dynamic Factor analysis)

Recommended Reading

de Valpine, P., & Hastings, A. (2002). Fitting population models incorporating process noise and observation error. Ecological Monographs, 72(1), 57-76.

Harvey, A. C. (1990). Forecasting, structural time series models and the Kalman filter. Cambridge university press.

Holmes, E. E., Ward, E. J., & Wills, K. (2012). Marss: Multivariate autoregressive state-space models for analyzing time-series data. The R Journal, 4(1), 11-19.

Kalman, R. E. 1960. A new approach to linear filtering and prediction problems. Transactions of the ASME–Journal of Basic Engineering 82(Series D): 35–45.

Meinhold, R. J., & Singpurwalla, N. D. (1983). Understanding the Kalman filter. The American Statistician, 37(2), 123-127.

Nielsen, A. and Berg, C. W. 2014. Estimation of time-varying selectivity in stock assessments using state–space models. Fisheries Research 158(1):96–101.

Newman, K. B., Buckland, S. T., Morgan, B. J., King, R., Borchers, D. L., Cole, D. J., ... & Thomas, L. (2014). Modelling Population Dynamics. Springer.

Pedersen, M.W., Berg, C.W., Thygesen,U. H., Nielsen, A., and Madsen, H. 2011. Estimation methods for nonlinear state-space models in ecology. Ecological Modelling, 222: 1394–1400.

Punt, A. E. (2003). Extending production models to include process error in the population dynamics. Canadian Journal of Fisheries and Aquatic Sciences, 60(10), 1217-1228.

Schnute, J. T. (1994). A general framework for developing sequential fisheries models. Canadian Journal of Fisheries and Aquatic Sciences, 51(8), 1676-1688.

Thorson, J. T., & Minto, C. (2015). Mixed effects: a unifying framework for statistical modelling in fisheries biology. ICES Journal of Marine Science: Journal du Conseil, fsu213.