fitGP_fish: Fit a "fisheries" GP model

In tanyalrogers/GPEDM: Gaussian Process regression for Empirical Dynamic Modeling

Fit a "fisheries" GP model

Description

Fits an alternate version of the GP model with one additonal parameter (catchability), designed for use in fisheries applications.

Usage

fitGP_fish(
  data,
  y,
  m,
  h,
  z = NULL,
  pop = NULL,
  time = NULL,
  scaling = c("global", "local", "none"),
  initpars = NULL,
  modeprior = 1,
  fixedpars = NULL,
  rhofixed = NULL,
  rhomatrix = NULL,
  augdata = NULL,
  linprior = "none",
  predictmethod = NULL,
  newdata = NULL,
  xname = "escapement",
  ytrans = c("none", "log", "gr1", "gr2"),
  bfixed = NULL,
  bshared = TRUE
)

Arguments

data	A data frame, required.
y	The response variable (required). Typically CPUE.
m	Lags of the response variable (required). Typically lags of CPUE.
h	Lags of the variable to be multipled by b and subtracted from m (required). Typically harvest. The dimensions of m and h must match.
z	Other predictor variables, e.g. covariates, that are unmodified (optional).
pop	Identifies separate populations (optional, if not supplied, defaults to 1 population). Population values can be either numeric, character, or factor.
time	The time index (recommended). If not supplied, defaults to a numeric index.
scaling	How the variables should be scaled (see Details). Scaling can be "global" (default), "local" (only applicable if more than 1 pop), or "none". Scaling is applied to y and each x. For a mix of scalings, do it manually beforehand and set scaling to "none". All outputs are backtransformed to original scale.
initpars	Starting values for hyperparameters (see Details) in the order c(phi,ve,sigma2,rho), in constrained (non-transformed) space. Should be a numeric vector with length ncol(x)+3. Defaults used if not supplied: c(rep(0.1,ncol(x)), 0.001, 1-0.001, 0.5)
modeprior	This value is used by the phi prior and sets the expected number of modes over the unit interval. Defaults to 1.
fixedpars	Fixes values of the hyperparameters phi, ve, and sigma2 (if desired). Should be a numeric vector with length ncol(x)+2 in the order c(phi,ve,sigma2) in constrained (non-transformed) space. Hyperparameters to estimate should be input as NA. The value of rho is fixed separately using argument rhofixed.
rhofixed	Fixes the rho parameter, if desired (see Details).
rhomatrix	Symmetrical square matrix of fixed pairwise rho values to use, with 1's on the diagonal. The rows and columns should be named with the population identifier. The output of getrhomatrix can be used here. All populations in the dataset must appear in rhomatrix. Populations in rhomatrix that are not in the dataset are allowed and will not be used.
augdata	A data frame with augmentation data (see Details).
linprior	Fit GP model to the residuals of a linear relationship between y and the first variable of x (prior to scaling) and backtransforms as appropriate. If y is growth rate, acts like a Ricker prior. Defaults to "none". Option "local" fits a separate regression for each population, "global" fits a single regression.
predictmethod	Using the training data, loo does leave-one-out prediction, lto does leave-timepoint-out prediction, and sequential does sequential (leave-future-out) prediction.
newdata	Data frame containing the same columns supplied in the original model.
xname	What the composite variable m-bh should be called. Defaults to "escapement".
ytrans	Tranformation to apply to y before fitting (m remains untransformed). Either "none" (default), "log", "gr1", or "gr2". R2 is calculated on y in its original units.
bfixed	Fixes b and bypasses optimization. If there are multiple pops, can be a scalar (shared across pops), or a named vector with different b's for each pop.
bshared	If there are multiple pops, should they share the same b (TRUE, default) or have different values of b (FALSE)? Note that optimizing multiple b values can be slow. Ignored if bfixed is supplied. If bshared=FALSE and either pop is not supplied or there is only one pop, reverts back to a single b estimate.

Details

This fits a GP model of the form

y=f(m-bh,z)

where y is catch per unit effort (CPUE), m are lags of CPUE, h are lags of harvest (in numbers or biomass), b is catchability (scalar), and z are optional covariates. CPUE is assumed proportional to total biomass (or numbers), with proportionality constant b. The composite variable m-bh is the biomass or number of individuals remaining after harvesting (escapement).

Parameter b is found using optimize applied to the posterior likelihood. Alternatively, a fixed value for b can be provided under bfixed.

If fitting a hierarchical model, the default behavior is to estimate a single value ofb shared by all pops (bshared=TRUE). You can fix a single shared value of b by providing a single value under bfixed. Alternatively, you can estimate different values of b for each population by setting bshared=FALSE. This will use the Nelder-Mead method of optim (and will be quite a bit slower than the single parameter optimization). You can fix different values of b for each population by supplying a named vector under bfixed.

Parameter ytrans applies a transformation to y before fitting the model. Inputs y and m should be in untransformed CPUE. ytrans="none" (the default) will apply no tranformation, ytrans="log" with compute log(y), ytrans="gr1" will compute log(y_t/m_{t-1}), and ytrans="gr2" will compute log(y_t/(m_{t-1}-bh_{t-1}).

Using this method requires the use of data with pre-generated lags (option A1 in fitGP). For more details on fitting a fisheries model and an example see the vignette. For more elaboration on the inputs, (e.g. scaling, augdata) see fitGP.

Value

A list (class GP and GPpred) with the same elements as fitGP and with additonal element b, the names of m, h, z stored under inputs, and the composite variable (escapement) included in the insampresults table.

tanyalrogers/GPEDM documentation built on June 1, 2025, 8:10 p.m.