| Title: | Statistical Catch-at-Age Plotting Environment |
|---|---|
| Description: | Import, plot, and diagnose results from statistical catch-at-age models, used in fisheries stock assessment. |
| Authors: | Arni Magnusson [aut, cre] |
| Maintainer: | Arni Magnusson <[email protected]> |
| License: | GPL-3 |
| Version: | 2.3.5 |
| Built: | 2024-11-21 05:37:04 UTC |
| Source: | https://github.com/arni-magnusson/scape |
Import and plot results from statistical catch-at-age models, used in fisheries stock assessments.
Import model results:
importADCAM |
ADCAM model results |
importCol |
Coleraine model results |
Plot model fit to data:
plotCA |
catch at age |
plotCL |
catch at length |
plotIndex |
abundance index |
plotLA |
length at age |
Plot derived quantities:
plotB |
biomass, recruitment, and landings |
plotN |
numbers at age |
plotSel |
selectivity and maturity |
Sigmas and sample sizes:
getN, getSigmaI, getSigmaR
|
extract |
estN, estSigmaI, estSigmaR
|
estimate |
iterate |
combine all get* and
est*
|
Import MCMC results:
importMCMC |
traces of likelihoods, parameters, biomass and recruitment |
importProj |
medium-term projections of biomass and catch |
Examples:
x.cod, x.ling, x.oreo,
x.saithe, x.sbw |
cod, ling, oreo, and whiting assessments |
xmcmc, xproj |
MCMC results |
browseVignettes() shows the vignettes, found in the
‘scape/doc’ directory.
Arni Magnusson.
Magnusson, A. (2005). R goes fishing: Analyzing fisheries data using AD Model Builder and R. Proceedings of the 4th International Workshop on Distributed Statistical Computing.
Magnusson, A. and Hilborn, R. (2007). What makes fisheries data informative? Fish and Fisheries, 8, 337–358. doi:10.1111/j.1467-2979.2007.00258.x
Magnusson, A., Punt, A.E. and Hilborn, R. (2013). Measuring uncertainty in fisheries stock assessment: the delta method, bootstrap, and MCMC. Fish and Fisheries, 14, 325–342. doi:10.1111/j.1467-2979.2012.00473.x
Most scape graphics are trellis plots, rendered with the
lattice package.
The functions Args and ll (package gdata) can be
useful for browsing unwieldy functions and objects.
The plotMCMC package can be used to diagnose the results of MCMC analysis.
Estimate the effective sample size for catch-at-age or catch-at-length data, based on the multinomial distribution.
estN(model, what="CAc", series=NULL, init=NULL, FUN=mean, ceiling=Inf, digits=0) estN.int(P, Phat) # internal functionestN(model, what="CAc", series=NULL, init=NULL, FUN=mean, ceiling=Inf, digits=0) estN.int(P, Phat) # internal function
model |
fitted |
what |
name of model element: |
series |
vector of strings indicating which gears or surveys to analyze (all by default). |
init |
initial sample size, determining the relative pattern of the effective sample size between years. |
FUN |
function to standardize the effective sample size. |
ceiling |
largest possible sample size in one year. |
digits |
number of decimal places to use when rounding, or
|
P |
observed catch-at-age or catch-at-length matrix. |
Phat |
fitted catch-at-age or catch-at-length matrix. |
The init sample sizes set a fixed pattern for the relative
sample sizes between years. For example, if there are two years of
catch-at-age data and the initial sample sizes are 100 in year 1 and
200 in year 2, the effective sample size will be two times greater in
year 2 than in year 1, although both will be scaled up or down
depending on how closely the model fits the catch-at-age data. The
value of init can be one of the following:
NULLmeans read the initial sample sizes from the
existing SS column (default).
means read the initial sample sizes from the SS
column in that model (object of class scape).
means those are the initial sample sizes (same length as the number of years).
FALSEmeans ignore the initial sample sizes and use
the empirical multinomial sample size () in each
year.
1means calculate one effective sample size to use
across all years, e.g. the mean or median of .
The idea behind FUN=mean is to guarantee that regardless of the
value of init, the mean effective sample size will always be
the same. Other functions can be used to a similar effect, such as
FUN=median.
The estN function is implemented for basic single-sex datasets.
If the data are sex-specific, estN pools (averages) the sexes
before estimating effective sample sizes. The general function
estN.int, on the other hand, is suitable for analyzing any
datasets in matrix format. The ‘int’ in estN.int
stands for internal (not integer), analogous to rep.int,
seq.int, sort.int, and similar functions.
Numeric vector of effective sample sizes (one value if init=1),
or a list of such vectors when analyzing multiple series.
This function uses the empirical multinomial sample size to estimate an effective sample size, which may be appropriate as likelihood weights for catch-at-age and catch-at-length data. The better the model fits the data, the larger the effective sample size. See McAllister and Ianelli (1997), Gavaris and Ianelli (2002), and Magnusson et al. (2013).
estN can be used iteratively, along with
estSigmaI and estSigmaR to assign
likelihood weights that are indicated by the model fit to the data.
Sigmas and sample sizes are then adjusted between model runs, until
they converge. The iterate function facilitates this procedure.
If is the observed proportion of fish at age (or
length bin) in year , and
is the fitted proportion, then the estimated sample size in that year
is:
Due to the non-random and non-independent nature of sampling fish, the effective sample size, for statistical purposes, is much less than the number of fish sampled. Common starting points include using the number of tows as the sample size in each year, or using the empirical multinomial sample sizes. Those “initial” sample sizes can then be scaled up or down. Sample sizes between 20 and 200 are common in the stock assessment literature.
Gavaris, S. and Ianelli, J.N. (2002). Statistical issues in fisheries' stock assessments. Scandinavian Journal of Statistics, 29, 245–271. doi:10.1111/1467-9469.00282
Magnusson, A., Punt, A.E., and Hilborn, R. (2013). Measuring uncertainty in fisheries stock assessment: the delta method, bootstrap, and MCMC. Fish and Fisheries, 14, 325–342. doi:10.1111/j.1467-2979.2012.00473.x
McAllister, M.K. and Ianelli, J.N. (1997). Bayesian stock assessment using catch-age data and the sampling-importance resampling algorithm. Canadian Journal of Fisheries and Aquaticic Sciences, 54, 284–300. doi:10.1139/f96-285
getN, getSigmaI, getSigmaR,
estN, estSigmaI, and estSigmaR
extract and estimate sample sizes and sigmas.
iterate combines all the get* and est*
functions in one call.
plotCA and plotCL show what is behind the
sample-size estimation.
scape-package gives an overview of the package.
## Exploring candidate sample sizes: getN(x.sbw) # sample sizes used in assessment: number of tows estN(x.sbw) # effective sample size, given data (tows) and model fit estN(x.sbw, ceiling=200) # could use this estN(x.sbw, init=FALSE) # from model fit, disregarding tows plotCA(x.sbw) # years with good fit => large sample size estN(x.sbw, init=1) # one sample size across all years estN(x.sbw, init=c(rep(1,14),rep(2,9))) # two sampling periods ## Same mean, regardless of init: mean(estN(x.sbw, digits=NULL)) mean(estN(x.sbw, digits=NULL, init=FALSE)) mean(estN(x.sbw, digits=NULL, init=1)) mean(estN(x.sbw, digits=NULL, init=c(rep(1,14),rep(2,9)))) ## Same median, regardless of init: median(estN(x.sbw, FUN=median, digits=NULL)) median(estN(x.sbw, FUN=median, digits=NULL, init=FALSE)) median(estN(x.sbw, FUN=median, digits=NULL, init=1)) median(estN(x.sbw, FUN=median, digits=NULL, init=c(rep(1,14),rep(2,9)))) ## Multiple series: getN(x.ling, "CLc") # sample size used in assessment getN(x.ling, "CLc", digits=0) # rounded estN(x.ling, "CLc") # model fit implies larger sample sizes getN(x.ling, "CLc", series="1", digits=0) # get one series estN(x.ling, "CLc", series="1") # estimate one series## Exploring candidate sample sizes: getN(x.sbw) # sample sizes used in assessment: number of tows estN(x.sbw) # effective sample size, given data (tows) and model fit estN(x.sbw, ceiling=200) # could use this estN(x.sbw, init=FALSE) # from model fit, disregarding tows plotCA(x.sbw) # years with good fit => large sample size estN(x.sbw, init=1) # one sample size across all years estN(x.sbw, init=c(rep(1,14),rep(2,9))) # two sampling periods ## Same mean, regardless of init: mean(estN(x.sbw, digits=NULL)) mean(estN(x.sbw, digits=NULL, init=FALSE)) mean(estN(x.sbw, digits=NULL, init=1)) mean(estN(x.sbw, digits=NULL, init=c(rep(1,14),rep(2,9)))) ## Same median, regardless of init: median(estN(x.sbw, FUN=median, digits=NULL)) median(estN(x.sbw, FUN=median, digits=NULL, init=FALSE)) median(estN(x.sbw, FUN=median, digits=NULL, init=1)) median(estN(x.sbw, FUN=median, digits=NULL, init=c(rep(1,14),rep(2,9)))) ## Multiple series: getN(x.ling, "CLc") # sample size used in assessment getN(x.ling, "CLc", digits=0) # rounded estN(x.ling, "CLc") # model fit implies larger sample sizes getN(x.ling, "CLc", series="1", digits=0) # get one series estN(x.ling, "CLc", series="1") # estimate one series
Estimate the effective sigma (magnitude of observation noise) for a survey or commercial abundance index, based on the empirical standard deviation.
estSigmaI(model, what="s", series=NULL, init=NULL, FUN=mean, p=1, digits=2)estSigmaI(model, what="s", series=NULL, init=NULL, FUN=mean, p=1, digits=2)
model |
fitted |
what |
which effective sigma to estimate: |
series |
vector of strings indicating which gears or surveys to analyze (all by default). |
init |
initial sigma, determining the relative pattern of the effective sigmas between years. |
FUN |
function to use when scaling a vector of sigmas. |
p |
effective number of parameters estimated in the model. |
digits |
number of decimal places to use when rounding, or
|
The init sigmas set a fixed pattern for the relative sigmas
between years. For example, if there are two years of abundance index
data and the initial sigmas are 0.1 in year 1 and 0.2 in year 2, the
effective sigma will be two times greater in year 2 than in year 1,
although both will be scaled up or down depending on how closely the
model fits the abundance index. The value of init can be one of
the following:
NULLmeans read the initial sigmas from the existing
CV column (default).
means read the initial sigmas from the CV column
in that model (object of class scape).
means those are the initial sigmas (same length as the number of years).
FALSE or 1
means use one effective sigma
() across all years.
The idea behind FUN=mean is to guarantee that regardless of the
value of init, the mean effective sigma will always be the
same. Other functions can be used to a similar effect, such as
FUN=median.
Numeric vector of effective sigmas (one value if init=1), or a
list of such vectors when analyzing multiple series.
This function uses the empirical standard deviation to estimate an effective sigma, which may be appropriate as likelihood weights for abundance index data. The better the model fits the data, the smaller the effective sigma.
estSigmaI can be used iteratively, along with
estN and estSigmaR to assign likelihood
weights that are indicated by the model fit to the data. Sigmas and
sample sizes are then adjusted between model runs, until they
converge. The iterate function facilitates this procedure.
If is the residual sum of squares in log space, is
the number of abundance index data points, and is the
effective number of parameters estimated in the model, then the
estimated effective sigma is:
There is no simple way to calculate for statistical
catch-at-age models. The default value of 1 is likely to underestimate
the true magnitude of observation noise.
getN, getSigmaI, getSigmaR,
estN, estSigmaI, and estSigmaR
extract and estimate sample sizes and sigmas.
iterate combines all the get* and est*
functions in one call.
plotIndex shows what is behind the sigma estimation.
scape-package gives an overview of the package.
## Exploring candidate sigmas: getSigmaI(x.cod) # sigma used in assessment 0.20 estSigmaI(x.cod) # model fit implies 0.17 plotIndex(x.cod) # model fit estSigmaI(x.cod, p=8) # eight estimated parameters implies 0.22 getSigmaI(x.sbw) # sigma used in assessment estSigmaI(x.sbw) # model fit implies smaller sigma estSigmaI(x.sbw, init=1) # could use 0.17 in all years ## Same mean, regardless of init: mean(estSigmaI(x.sbw, digits=NULL)) mean(estSigmaI(x.sbw, digits=NULL, init=1)) ## Same median, regardless of init: median(estSigmaI(x.sbw, FUN=median, digits=NULL)) median(estSigmaI(x.sbw, FUN=median, digits=NULL, init=1)) ## Multiple series: getSigmaI(x.oreo, "c") # sigma used in assessment getSigmaI(x.oreo, "c", digits=2) # rounded estSigmaI(x.oreo, "c") # model fit implies smaller sigma estSigmaI(x.oreo, "c", init=1) # could use 0.19 in all years estSigmaI(x.oreo, "c", init=1, digits=3) # series 2 slightly worse fit # estSigmaI(x.oreo, "c", init=1, p=11) # more parameters than datapoints getSigmaI(x.oreo, "c", series="Series 2-1") # get one series estSigmaI(x.oreo, "c", series="Series 2-1") # estimate one series## Exploring candidate sigmas: getSigmaI(x.cod) # sigma used in assessment 0.20 estSigmaI(x.cod) # model fit implies 0.17 plotIndex(x.cod) # model fit estSigmaI(x.cod, p=8) # eight estimated parameters implies 0.22 getSigmaI(x.sbw) # sigma used in assessment estSigmaI(x.sbw) # model fit implies smaller sigma estSigmaI(x.sbw, init=1) # could use 0.17 in all years ## Same mean, regardless of init: mean(estSigmaI(x.sbw, digits=NULL)) mean(estSigmaI(x.sbw, digits=NULL, init=1)) ## Same median, regardless of init: median(estSigmaI(x.sbw, FUN=median, digits=NULL)) median(estSigmaI(x.sbw, FUN=median, digits=NULL, init=1)) ## Multiple series: getSigmaI(x.oreo, "c") # sigma used in assessment getSigmaI(x.oreo, "c", digits=2) # rounded estSigmaI(x.oreo, "c") # model fit implies smaller sigma estSigmaI(x.oreo, "c", init=1) # could use 0.19 in all years estSigmaI(x.oreo, "c", init=1, digits=3) # series 2 slightly worse fit # estSigmaI(x.oreo, "c", init=1, p=11) # more parameters than datapoints getSigmaI(x.oreo, "c", series="Series 2-1") # get one series estSigmaI(x.oreo, "c", series="Series 2-1") # estimate one series
Estimate sigma R (recruitment variability), based on the empirical standard deviation of recruitment deviates in log space.
estSigmaR(model, digits=2)estSigmaR(model, digits=2)
model |
fitted |
digits |
number of decimal places to use when rounding, or
|
Vector of two numbers, estimating recruitment variability based on (1) the estimated age composition in the first year, and (2) subsequent annual recruitment.
This function uses the empirical standard deviation to estimate sigma R, which may be appropriate as likelihood penalty (or Bayesian prior distribution) for recruitment deviates from the stock-recruitment curve. The smaller the estimated recruitment deviates, the smaller the estimated sigma R.
estSigmaR can be used iteratively, along with
estN and estSigmaI to assign likelihood
weights that are indicated by the model fit to the data. Sigmas and
sample sizes are then adjusted between model runs, until they
converge. The iterate function facilitates this procedure.
If is the sum of squared recruitment deviates in log space
and is the number of estimated recruitment deviates, then the
estimated sigma R is:
The denominator is neither 1 nor , since is
based on deviates from zero and not the mean, and the deviates do not
converge to zero as the number of model parameters increases.
getN, getSigmaI, getSigmaR,
estN, estSigmaI, and estSigmaR
extract and estimate sample sizes and sigmas.
iterate combines all the get* and est*
functions in one call.
plotN and plotB(..., what="s") show what
is behind the sigma R estimation.
scape-package gives an overview of the package.
getSigmaR(x.cod) # sigmaR used in assessment 0.5 and 1.0 estSigmaR(x.cod) # model estimates imply 0.20 and 0.52 getSigmaR(x.ling) # 0.6, deterministic age distribution in first year estSigmaR(x.ling) # model estimates imply 0.36 getSigmaR(x.sbw) estSigmaR(x.sbw) # large deviates in first year plotN(x.sbw) # enormous plus group and 1991 cohort # x.oreo assessment had deterministic recruitment, so no deviatesgetSigmaR(x.cod) # sigmaR used in assessment 0.5 and 1.0 estSigmaR(x.cod) # model estimates imply 0.20 and 0.52 getSigmaR(x.ling) # 0.6, deterministic age distribution in first year estSigmaR(x.ling) # model estimates imply 0.36 getSigmaR(x.sbw) estSigmaR(x.sbw) # large deviates in first year plotN(x.sbw) # enormous plus group and 1991 cohort # x.oreo assessment had deterministic recruitment, so no deviates
Extract the sample size that was used in a model, from catch-at-age or catch-at-length data.
getN(model, what="CAc", series=NULL, digits=NULL)getN(model, what="CAc", series=NULL, digits=NULL)
model |
fitted |
what |
name of model element: |
series |
vector of strings indicating which gears or surveys to analyze (all by default). |
digits |
number of decimal places to use when rounding, or
|
Numeric vector of year-specific sample sizes, or a list of such vectors when analyzing multiple series.
Thin wrapper to access model$element$SS, providing a uniform
interface with other get* and est* functions.
See discussion in the estN documentation.
getN, getSigmaI, getSigmaR,
estN, estSigmaI, and
estSigmaR extract and estimate sample sizes and sigmas.
scape-package gives an overview of the package.
## Exploring candidate sample sizes: getN(x.sbw) # sample sizes used in assessment: number of tows estN(x.sbw) # effective sample size, given data (tows) and model fit ## Multiple series: getN(x.ling, "CLc") # sample size used in assessment getN(x.ling, "CLc", digits=0) # rounded estN(x.ling, "CLc") # model fit implies larger sample sizes getN(x.ling, "CLc", series="1", digits=0) # get one series estN(x.ling, "CLc", series="1") # estimate one series## Exploring candidate sample sizes: getN(x.sbw) # sample sizes used in assessment: number of tows estN(x.sbw) # effective sample size, given data (tows) and model fit ## Multiple series: getN(x.ling, "CLc") # sample size used in assessment getN(x.ling, "CLc", digits=0) # rounded estN(x.ling, "CLc") # model fit implies larger sample sizes getN(x.ling, "CLc", series="1", digits=0) # get one series estN(x.ling, "CLc", series="1") # estimate one series
Extract the sigma (magnitude of observation noise) that was used in a model, from survey or commercial abundance index data.
getSigmaI(model, what="s", series=NULL, digits=NULL)getSigmaI(model, what="s", series=NULL, digits=NULL)
model |
fitted |
what |
which sigma to extract: |
series |
vector of strings indicating which gears or surveys to analyze (all by default). |
digits |
number of decimal places to use when rounding, or
|
Numeric vector of year-specific sigmas, or a list of such vectors when analyzing multiple series.
Thin wrapper to access model$element$CV, providing a uniform
interface with other get* and est* functions.
See discussion in the estSigmaI documentation.
getN, getSigmaI, getSigmaR,
estN, estSigmaI, and
estSigmaR extract and estimate sample sizes and sigmas.
scape-package gives an overview of the package.
## Exploring candidate sigmas: getSigmaI(x.cod) # sigma used in assessment 0.20 estSigmaI(x.cod) # model fit implies 0.17 ## Multiple series: getSigmaI(x.oreo, "c") # sigma used in assessment getSigmaI(x.oreo, "c", digits=2) # rounded estSigmaI(x.oreo, "c") # model fit implies smaller sigma getSigmaI(x.oreo, "c", series="Series 2-1") # get one series estSigmaI(x.oreo, "c", series="Series 2-1") # estimate one series## Exploring candidate sigmas: getSigmaI(x.cod) # sigma used in assessment 0.20 estSigmaI(x.cod) # model fit implies 0.17 ## Multiple series: getSigmaI(x.oreo, "c") # sigma used in assessment getSigmaI(x.oreo, "c", digits=2) # rounded estSigmaI(x.oreo, "c") # model fit implies smaller sigma getSigmaI(x.oreo, "c", series="Series 2-1") # get one series estSigmaI(x.oreo, "c", series="Series 2-1") # estimate one series
Extract sigma R (recruitment variability) that was used in a model, as
indicated in the Dev$sigmaR model component.
getSigmaR(model, digits=NULL)getSigmaR(model, digits=NULL)
model |
fitted |
digits |
number of decimal places to use when rounding, or
|
Vector of two numbers, representing recruitment variability in (1) the estimated age composition in the first year, and (2) subsequent annual recruitment.
Thin wrapper to access model$Dev$sigmaR, providing a uniform
interface with other get* and est* functions.
See discussion in the estSigmaR documentation.
getN, getSigmaI, getSigmaR,
estN, estSigmaI, and
estSigmaR extract and estimate sample sizes and sigmas.
getSigmaR(x.cod) # sigmaR used in assessment 0.5 and 1.0 estSigmaR(x.cod) # model estimates imply 0.20 and 0.52 getSigmaR(x.ling) # 0.6, deterministic age distribution in first year estSigmaR(x.ling) # model estimates imply 0.36 getSigmaR(x.sbw) estSigmaR(x.sbw) # large deviates in first year plotN(x.sbw) # enormous plus group and 1991 cohort # x.oreo assessment had deterministic recruitment, so no deviatesgetSigmaR(x.cod) # sigmaR used in assessment 0.5 and 1.0 estSigmaR(x.cod) # model estimates imply 0.20 and 0.52 getSigmaR(x.ling) # 0.6, deterministic age distribution in first year estSigmaR(x.ling) # model estimates imply 0.36 getSigmaR(x.sbw) estSigmaR(x.sbw) # large deviates in first year plotN(x.sbw) # enormous plus group and 1991 cohort # x.oreo assessment had deterministic recruitment, so no deviates
Import ADCAM model results into a standard format suitable for plotting.
importADCAM(dir)importADCAM(dir)
dir |
directory containing the results files (‘resultsbyyearandage’, ‘resultsbyyear’, ‘resultsbyage’), as well as the input file ‘catchparameters’ specifying selectivity periods. |
A list of class scape containing the following elements:
N |
predicted numbers at age |
B |
predicted biomass, recruitment, and observed landings (year things) |
Sel |
predicted selectivity and observed maturity (age things) |
Dev |
predicted recruitment deviates from the stock-recruitment curve |
Survey |
survey abundance index and fit |
CAc |
commercial C@A (catch at age) and fit |
CAs |
survey C@A (catch at age) and fit |
B$VB is reference biomass.
Survey$CV is 0 to give plotIndex sensible y-axis
limits.
CAc$SS and CAs$SS are the
log(obs+const) - log(fit+const) residuals used in the ADCAM
likelihood.
Bjornsson, H. and A. Magnusson. (2009) ADCAM user manual (draft version). ICES CM 2009/ACOM:56, Annex 6.
read.table, readLines, and
scan import any data.
x.saithe was created using importADCAM.
importCol imports Coleraine model results.
scape-package gives an overview of the package.
path <- system.file("example/saithe", package="scape") saithe <- importADCAM(path)path <- system.file("example/saithe", package="scape") saithe <- importADCAM(path)
Import Coleraine model results from ‘.res’ file, and rearrange into a standard format suitable for plotting.
importCol(res.file, Dev=FALSE, CPUE=FALSE, Survey=FALSE, CAc=FALSE, CAs=FALSE, CLc=FALSE, CLs=FALSE, LA=FALSE, quiet=TRUE)importCol(res.file, Dev=FALSE, CPUE=FALSE, Survey=FALSE, CAc=FALSE, CAs=FALSE, CLc=FALSE, CLs=FALSE, LA=FALSE, quiet=TRUE)
res.file |
name of Coleraine model results file to import. |
Dev |
whether recruitment deviates were estimated in model. |
CPUE |
whether model was fitted to catch-per-unit-effort data. |
Survey |
whether model was fitted to survey abundance index data. |
CAc |
whether model was fitted to commercial catch-at-age data. |
CAs |
whether model was fitted to survey catch-at-age data. |
CLc |
whether model was fitted to commercial catch-at-length data. |
CLs |
whether model was fitted to survey catch-at-length data. |
LA |
whether model was fitted to length-at-age data. |
quiet |
whether to report progress while parsing file. |
A list of class scape containing at least N, B,
and Sel. The other elements may or may not be included in the
list, depending on how importCol was called:
N |
predicted numbers at age |
B |
predicted biomass, recruitment, and observed landings (year things) |
Sel |
predicted selectivity and observed maturity (age things) |
Dev |
predicted recruitment deviates from the stock-recruitment curve |
CPUE, Survey
|
commercial and survey abundance index and fit |
CAc, CAs
|
commercial and survey C@A (catch at age) and fit |
CLc, CLs
|
commercial and survey C@L (catch at length) and fit |
LA |
observed L@A and fit |
This import function is implemented for the Coleraine statistical catch-at-age software, and can serve as a template for scape users who would like to implement import functions for specific stock assessment software.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(x.cod),
ll(x.cod$N), and head(x.cod$N).
Hilborn, R., Maunder, M., Parma, A., Ernst, B., Payne, J., and Starr, P. (2003). Coleraine: A generalized age-structured stock assessment model. User's manual version 2.0. University of Washington Report SAFS-UW-0116.
read.table, readLines, and
scan import any data.
x.cod, x.ling, x.oreo, and
x.sbw were created using importCol.
importADCAM imports ADCAM model results.
scape-package gives an overview of the package.
path <- system.file("example/cod.res", package="scape") cod <- importCol(path, Dev=TRUE, Survey=TRUE, CAc=TRUE, CAs=TRUE)path <- system.file("example/cod.res", package="scape") cod <- importCol(path, Dev=TRUE, Survey=TRUE, CAc=TRUE, CAs=TRUE)
Import Coleraine MCMC traces for likelihoods, parameters, spawning biomass, and recruitment.
importMCMC(dir, coda=FALSE, quiet=TRUE, pretty.labels=FALSE, l.choose=NULL, p.choose=NULL)importMCMC(dir, coda=FALSE, quiet=TRUE, pretty.labels=FALSE, l.choose=NULL, p.choose=NULL)
dir |
directory containing the files ‘mcmclike.out’, ‘params.pst’, ‘spawbiom.pst’ and ‘recruits.pst’. |
coda |
whether data frames should be coerced to class |
quiet |
whether to report progress while parsing files in directory. |
pretty.labels |
whether likelihood and parameter columns should be renamed |
l.choose |
vector of strings, indicating which likelihood
components to import, or |
p.choose |
vector of strings, indicating which parameters to
import, or |
A list containing:
L |
likelihoods |
P |
parameters |
B |
biomass by year |
R |
recruitment by year |
as data frames, or mcmc objects if coda=TRUE.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(xmcmc),
ll(xmcmc$P), and head(xmcmc$P).
The plotMCMC package is recommended for plotting MCMC diagnostics and posteriors.
Hilborn, R., Maunder, M., Parma, A., Ernst, B., Payne, J., and Starr, P. (2003). Coleraine: A generalized age-structured stock assessment model. User's manual version 2.0. University of Washington Report SAFS-UW-0116.
read.table, readLines, and
scan import any data.
importMCMC and importProj import Coleraine MCMC
results.
xmcmc was created using importMCMC.
scape-package gives an overview of the package.
path <- system.file("example/mcmc", package="scape") mcmc <- importMCMC(path) # or rename and select particular elements: mcmc <- importMCMC(path, pretty.labels=TRUE, l.choose=c("CAc","CAs","Survey","Prior","Total"), p.choose=c("R0","Rinit","uinit","cSleft","cSfull", "sSleft","sSfull","logq"))path <- system.file("example/mcmc", package="scape") mcmc <- importMCMC(path) # or rename and select particular elements: mcmc <- importMCMC(path, pretty.labels=TRUE, l.choose=c("CAc","CAs","Survey","Prior","Total"), p.choose=c("R0","Rinit","uinit","cSleft","cSfull", "sSleft","sSfull","logq"))
Import Coleraine MCMC traces for spawning biomass and catch, projected into the near future.
importProj(dir, coda=FALSE, quiet=TRUE)importProj(dir, coda=FALSE, quiet=TRUE)
dir |
directory containing the files ‘strategy.out’, ‘projspbm.pst’ and ‘procatch.pst’. |
coda |
whether data frames should be coerced to class |
quiet |
whether to report progress while parsing files in directory. |
A list containing:
B |
biomass by catch policy and year |
Y |
catch by catch policy and year |
as lists of data frames, or mcmc objects if coda=TRUE.
MCMC projections can be used to evaluate the short-term outcome of harvest policies (constant catch or constant harvest rate), given the uncertainty about parameter values and random future recruitment.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(xproj),
ll(xproj$B), ll(xproj$B$"0.2"), and
head(xproj$B$"0.2").
The plotMCMC package is recommended for plotting MCMC diagnostics and posteriors.
Hilborn, R., Maunder, M., Parma, A., Ernst, B., Payne, J., and Starr, P. (2003). Coleraine: A generalized age-structured stock assessment model. User's manual version 2.0. University of Washington Report SAFS-UW-0116.
read.table, readLines, and
scan import any data.
importMCMC and importProj import Coleraine MCMC
results.
xproj was created using importProj.
scape-package gives an overview of the package.
path <- system.file("example/mcmc", package="scape") proj <- importProj(path)path <- system.file("example/mcmc", package="scape") proj <- importProj(path)
Compare current sigmas and sample sizes with candidate values, by
running variations of estSigmaR, estN, and
estSigmaI on all model components.
iterate(model, ceiling=Inf, p=1, digits.n=0, digits.sigma=2)iterate(model, ceiling=Inf, p=1, digits.n=0, digits.sigma=2)
model |
fitted |
ceiling |
largest possible sample size in one year, passed to
|
p |
effective number of parameters estimated in the model, passed
to |
digits.n |
number of decimal places to use when rounding sample
sizes, or |
digits.sigma |
number of decimal places to use when rounding
sigmas, or |
List containing data frames summarizing current sigmas and sample sizes, as well as candidate values. The following abbreviations are used in column names:
candidate sigma, the empirical standard deviation.
candidate sample sizes, the empirical multinomial sample sizes.
vector of candidate values, whose mean equals
sigmahat or nhat.
vector of candidate values, whose median equals
sigmahat or nhat.
vector of identical candidate values, the mean of
nhat.
vector of identical candidate values, the median of
nhat.
getN, getSigmaI, getSigmaR,
estN, estSigmaI, and
estSigmaR extract and estimate sample sizes and sigmas.
iterate combines all the get* and est*
functions in one call.
scape-package gives an overview of the package.
iterate(x.cod) iterate(x.ling) iterate(x.oreo) iterate(x.sbw)iterate(x.cod) iterate(x.ling) iterate(x.oreo) iterate(x.sbw)
Plot scape model predicted biomass, stock recruitment, and
landings.
plotB(model, what="d", series=NULL, years=NULL, axes=TRUE, div=1, legend="bottom", main="", xlab="", ylab="", cex.main=1.2, cex.legend=1, cex.lab=1, cex.axis=0.8, las=1, tck=c(1,what=="d")/2, tick.number=5, lty.grid=3, col.grid="white", pch=16, cex.points=0.8, col.points="black", lty.lines=1:3, lwd.lines=2, col.lines="black", ratio.bars=3, col.bars="gray", plot=TRUE, ...)plotB(model, what="d", series=NULL, years=NULL, axes=TRUE, div=1, legend="bottom", main="", xlab="", ylab="", cex.main=1.2, cex.legend=1, cex.lab=1, cex.axis=0.8, las=1, tck=c(1,what=="d")/2, tick.number=5, lty.grid=3, col.grid="white", pch=16, cex.points=0.8, col.points="black", lty.lines=1:3, lwd.lines=2, col.lines="black", ratio.bars=3, col.bars="gray", plot=TRUE, ...)
model |
fitted |
what |
what to plot: |
series |
vector of strings indicating which column names in
|
years |
vector of numbers indicating which years to include (all by default). |
axes |
whether to plot axis values. |
div |
denominator to shorten values on the y axis, or a vector with two elements referring to x and y axis. |
legend |
legend location: |
main |
main title. |
xlab |
x-axis label. |
ylab |
y-axis label. |
cex.main |
size of main title. |
cex.legend |
size of legend text. |
cex.lab |
size of axis labels. |
cex.axis |
size of tick labels. |
las |
orientation of tick labels: 0=parallel, 1=horizontal, 2=perpendicular, 3=vertical. |
tck |
tick mark length. |
tick.number |
number of tick marks. |
lty.grid |
line type of gridlines. |
col.grid |
color of gridlines. |
pch |
symbol for points. |
cex.points |
size of points. |
col.points |
color of points. |
lty.lines |
line type of main lines. |
lwd.lines |
line width of main lines. |
col.lines |
color of main lines. |
ratio.bars |
width of bars. |
col.bars |
color of bars. |
plot |
whether to draw plot. |
... |
passed to |
The "d"[efault] plot shows spawning biomass and vulnerable
biomass as lines, and landings as bars, on the same scale.
When plot=TRUE, a trellis plot is drawn and a data frame is
returned, containing the data used for plotting. When
plot=FALSE, a trellis object is returned.
The Args function from the gdata package is recommended
for reviewing the arguments, instead of args.
xyplot, panel.barchart,
and panel.superpose are the underlying drawing
functions.
plotCA, plotCL, plotIndex,
and plotLA plot model fit and data.
plotB, plotN, and plotSel plot
derived quantities.
scape-package gives an overview of the package.
plotB(x.ling, series=c("VB.1","VB.2","Y"), div=1000, xlab="Year\n", ylab="Biomass and landings (kt)") plotB(x.ling, "s", div=1000, xlab="Biomass age 4+ (kt)", ylab="Recruitment (million one-year-olds)")plotB(x.ling, series=c("VB.1","VB.2","Y"), div=1000, xlab="Year\n", ylab="Biomass and landings (kt)") plotB(x.ling, "s", div=1000, xlab="Biomass age 4+ (kt)", ylab="Recruitment (million one-year-olds)")
Plot scape model fit to catch-at-age data.
plotCA(model, what="c", fit=TRUE, swap=FALSE, series=NULL, sex=NULL, years=NULL, ages=NULL, axes=TRUE, same.limits=TRUE, log=FALSE, base=10, eps.log=1e-5, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=!fit, tck=c(1,fit)/2, tick.number=5, lty.grid=3, col.grid="gray", pch=16, cex.points=0.5, col.points="black", lty.lines=1, lwd.lines=2, col.lines=c("red","blue"), plot=TRUE, ...)plotCA(model, what="c", fit=TRUE, swap=FALSE, series=NULL, sex=NULL, years=NULL, ages=NULL, axes=TRUE, same.limits=TRUE, log=FALSE, base=10, eps.log=1e-5, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=!fit, tck=c(1,fit)/2, tick.number=5, lty.grid=3, col.grid="gray", pch=16, cex.points=0.5, col.points="black", lty.lines=1, lwd.lines=2, col.lines=c("red","blue"), plot=TRUE, ...)
model |
fitted |
what |
what to plot: |
fit |
whether to overlay fitted values on observed data. |
swap |
whether to swap ages and years between axes or panels. |
series |
string indicating which gear or survey to plot (first by default). |
sex |
string indicating which sex to plot (both by default). |
years |
vector of numbers indicating which years to plot (all by default). |
ages |
vector of numbers indicating which ages to plot (all by default). |
axes |
whether to plot axis values. |
same.limits |
whether panels should have same y-axis limits. |
log |
whether to log-transform values. |
base |
logarithm base. |
eps.log |
small number to add before log-transforming to avoid log 0. |
main |
main title. |
xlab |
x-axis label. |
ylab |
y-axis label. |
cex.main |
size of main title. |
cex.lab |
size of axis labels. |
cex.axis |
size of tick labels. |
cex.strip |
size of strip labels. |
col.strip |
color of strip labels. |
strip |
logical flag (whether to plot strip labels), or a
function passed to |
las |
orientation of tick labels: 0=parallel, 1=horizontal, 2=perpendicular, 3=vertical. |
tck |
tick mark length. |
tick.number |
number of tick marks. |
lty.grid |
line type of gridlines. |
col.grid |
color of gridlines. |
pch |
symbol for points. |
cex.points |
size of points. |
col.points |
color of points. |
lty.lines |
line type of main lines, possibly a vector where element 2 refers to males. |
lwd.lines |
line width of main lines, possibly a vector where element 2 refers to males. |
col.lines |
color of main lines, possibly a vector where element 2 refers to males. |
plot |
whether to draw plot. |
... |
passed to |
When plot=TRUE, a trellis plot is drawn and a data frame is
returned, containing the data used for plotting. When
plot=FALSE, a trellis object is returned.
The Args function from the gdata package is recommended
for reviewing the arguments, instead of args.
xyplot, panel.xyplot,
and panel.superpose are the underlying drawing
functions.
plotCA, plotCL, plotIndex, and
plotLA plot model fit and data.
plotB, plotN, and plotSel
plot derived quantities.
scape-package gives an overview of the package.
plotCA(x.sbw, fit=FALSE, strip=FALSE, xlab="Age", ylab="Year", tick.number=10) plotCA(x.cod, xlab="Age", ylab="Proportion in catch", cex.strip=0.7, cex.axis=0.7, col.lines="brown", layout=c(8,4)) plotCA(x.cod, xlab="Age", ylab="Proportion in catch", cex.strip=0.7, cex.axis=0.7, col.lines="brown", layout=c(2,4), swap=TRUE, ages=3:10, same.limits=FALSE) plotCA(x.ling, "s", col.points=c("red","blue"), lty.lines=0, xlab="Age", ylab="Observed proportion in survey", tck=0.5, cex.strip=0.7, cex.axis=0.7) plotCA(x.ling, "s", xlab="Age", ylab="Observed proportion in survey", fit=FALSE, cex.strip=0.7, cex.axis=0.7, tck=0.5, layout=c(5,2)) plotCA(x.ling, "s", xlab="Age", ylab="Observed proportion in survey", fit=FALSE, cex.strip=0.7, cex.axis=0.7, tck=0.5, layout=c(5,6), swap=TRUE)plotCA(x.sbw, fit=FALSE, strip=FALSE, xlab="Age", ylab="Year", tick.number=10) plotCA(x.cod, xlab="Age", ylab="Proportion in catch", cex.strip=0.7, cex.axis=0.7, col.lines="brown", layout=c(8,4)) plotCA(x.cod, xlab="Age", ylab="Proportion in catch", cex.strip=0.7, cex.axis=0.7, col.lines="brown", layout=c(2,4), swap=TRUE, ages=3:10, same.limits=FALSE) plotCA(x.ling, "s", col.points=c("red","blue"), lty.lines=0, xlab="Age", ylab="Observed proportion in survey", tck=0.5, cex.strip=0.7, cex.axis=0.7) plotCA(x.ling, "s", xlab="Age", ylab="Observed proportion in survey", fit=FALSE, cex.strip=0.7, cex.axis=0.7, tck=0.5, layout=c(5,2)) plotCA(x.ling, "s", xlab="Age", ylab="Observed proportion in survey", fit=FALSE, cex.strip=0.7, cex.axis=0.7, tck=0.5, layout=c(5,6), swap=TRUE)
Plot scape model fit to catch-at-length data.
plotCL(model, what="c", fit=TRUE, swap=FALSE, series=NULL, sex=NULL, years=NULL, lengths=NULL, axes=TRUE, same.limits=TRUE, log=FALSE, base=10, eps.log=1e-5, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=!fit, tck=c(1,fit)/2, tick.number=5, lty.grid=3, col.grid="gray", pch=16, cex.points=0.5, col.points="black", lty.lines=1, lwd.lines=2, col.lines=c("red","blue"), plot=TRUE, ...)plotCL(model, what="c", fit=TRUE, swap=FALSE, series=NULL, sex=NULL, years=NULL, lengths=NULL, axes=TRUE, same.limits=TRUE, log=FALSE, base=10, eps.log=1e-5, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=!fit, tck=c(1,fit)/2, tick.number=5, lty.grid=3, col.grid="gray", pch=16, cex.points=0.5, col.points="black", lty.lines=1, lwd.lines=2, col.lines=c("red","blue"), plot=TRUE, ...)
model |
fitted |
what |
what to plot: |
fit |
whether to overlay fitted values on observed data. |
swap |
whether to swap lengths and years between axes or panels. |
series |
string indicating which gear or survey to plot (first by default). |
sex |
string indicating which sex to plot (both by default). |
years |
vector of numbers indicating which years to plot plot (all by default). |
lengths |
vector of numbers indicating which lengths to plot (all by default). |
axes |
whether to plot axis values. |
same.limits |
whether panels should have same y-axis limits. |
log |
whether to log-transform values. |
base |
logarithm base. |
eps.log |
small number to add before log-transforming to avoid log 0. |
main |
main title. |
xlab |
x-axis label. |
ylab |
y-axis label. |
cex.main |
size of main title. |
cex.lab |
size of axis labels. |
cex.axis |
size of tick labels. |
cex.strip |
size of strip labels. |
col.strip |
color of strip labels. |
strip |
logical flag (whether to plot strip labels), or a
function passed to |
las |
orientation of tick labels: 0=parallel, 1=horizontal, 2=perpendicular, 3=vertical. |
tck |
tick mark length. |
tick.number |
number of tick marks. |
lty.grid |
line type of gridlines. |
col.grid |
color of gridlines. |
pch |
symbol for points. |
cex.points |
size of points. |
col.points |
color of points. |
lty.lines |
line type of main lines, possibly a vector where element 2 refers to males. |
lwd.lines |
line width of main lines, possibly a vector where element 2 refers to males. |
col.lines |
color of main lines, possibly a vector where element 2 refers to males. |
plot |
whether to draw plot. |
... |
passed to |
When plot=TRUE, a trellis plot is drawn and a data frame is
returned, containing the data used for plotting. When
plot=FALSE, a trellis object is returned.
The Args function from the gdata package is recommended
for reviewing the arguments, instead of args.
xyplot, panel.xyplot,
and panel.superpose are the underlying drawing
functions.
plotCA, plotCL, plotIndex, and
plotLA plot model fit and data.
plotB, plotN, and plotSel
plot derived quantities.
scape-package gives an overview of the package.
plotCL(x.ling, fit=FALSE, strip=FALSE, series="1", sex="Female", xlab="Length (cm)", ylab="Year") plotCL(x.oreo, xlab="Length (cm)", ylab="Proportion in catch") plotCL(x.oreo, "s", layout=c(2,1), xlab="Length (cm)", ylab="Observed proportion in survey", cex.points=0.8, col.points=c("red","blue"), lty.lines=0) plotCL(x.ling, fit=FALSE, series="2", xlab="Length (cm)", ylab="Observed proportion in trawl catch", tck=0.5) plotCL(x.ling, series="2", swap=TRUE, lengths=70:150, lty.grid=0)plotCL(x.ling, fit=FALSE, strip=FALSE, series="1", sex="Female", xlab="Length (cm)", ylab="Year") plotCL(x.oreo, xlab="Length (cm)", ylab="Proportion in catch") plotCL(x.oreo, "s", layout=c(2,1), xlab="Length (cm)", ylab="Observed proportion in survey", cex.points=0.8, col.points=c("red","blue"), lty.lines=0) plotCL(x.ling, fit=FALSE, series="2", xlab="Length (cm)", ylab="Observed proportion in trawl catch", tck=0.5) plotCL(x.ling, series="2", swap=TRUE, lengths=70:150, lty.grid=0)
Plot scape model fit to abundance index data.
plotIndex(model, what="s", series=NULL, axes=TRUE, same.limits=FALSE, between=list(x=axes,y=axes), ylim=NULL, q=1, bar=1, log=FALSE, base=10, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=1, tck=c(1,0)/2, tick.number=5, lty.grid=3, col.grid="white", pch=16, cex.points=1.2, col.points="black", lty.lines=1, lwd.lines=4, col.lines="dimgray", lty.bar=1, plot=TRUE, ...)plotIndex(model, what="s", series=NULL, axes=TRUE, same.limits=FALSE, between=list(x=axes,y=axes), ylim=NULL, q=1, bar=1, log=FALSE, base=10, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=1, tck=c(1,0)/2, tick.number=5, lty.grid=3, col.grid="white", pch=16, cex.points=1.2, col.points="black", lty.lines=1, lwd.lines=4, col.lines="dimgray", lty.bar=1, plot=TRUE, ...)
model |
fitted |
what |
what to plot: |
series |
vector of strings indicating which gears or surveys to plot (all by default). |
axes |
whether to plot axis values. |
same.limits |
whether panels should have same y-axis limits. |
between |
list with |
ylim |
vector with lower and upper y-axis limits. |
q |
denominator to scale the y axis, e.g. to vulnerable biomass.
Similar to the |
bar |
extent of error bars relative to standard error. |
log |
whether to log-transform values. |
base |
logarithm base. |
main |
main title. |
xlab |
x-axis label. |
ylab |
y-axis label. |
cex.main |
size of main title. |
cex.lab |
size of axis labels. |
cex.axis |
size of tick labels. |
cex.strip |
size of strip labels. |
col.strip |
color of strip labels. |
strip |
logical flag (whether to plot strip labels), or a
function passed to |
las |
orientation of tick labels: 0=parallel, 1=horizontal, 2=perpendicular, 3=vertical. |
tck |
tick mark length. |
tick.number |
number of tick marks. |
lty.grid |
line type of gridlines. |
col.grid |
color of gridlines. |
pch |
symbol for points. |
cex.points |
size of points. |
col.points |
color of points and error bars. |
lty.lines |
line type of main lines. |
lwd.lines |
line width of main lines. |
col.lines |
color of main lines. |
lty.bar |
line type of error bars. |
plot |
whether to draw plot. |
... |
passed to |
When plot=TRUE, a trellis plot is drawn and a data frame is
returned, containing the data used for plotting. When
plot=FALSE, a trellis object is returned.
The Args function from the gdata package is recommended
for reviewing the arguments, instead of args.
xyplot, panel.xyplot,
and panel.xYplot are the underlying drawing
functions.
plotCA, plotCL, plotIndex, and
plotLA plot model fit and data.
plotB, plotN, and plotSel
plot derived quantities.
scape-package gives an overview of the package.
plotIndex(x.cod, xlab="Year", ylab="Survey abundance index", strip=FALSE) plotIndex(x.oreo, "c", series="Series 1-1", xlim=c(1981,1990)) plotIndex(x.oreo, "c", xlim=list(c(1981,1990),c(1992,2002)), xlab="Year", ylab="Observed CPUE", col.points=c("salmon","seagreen"), lty.lines=0)plotIndex(x.cod, xlab="Year", ylab="Survey abundance index", strip=FALSE) plotIndex(x.oreo, "c", series="Series 1-1", xlim=c(1981,1990)) plotIndex(x.oreo, "c", xlim=list(c(1981,1990),c(1992,2002)), xlab="Year", ylab="Observed CPUE", col.points=c("salmon","seagreen"), lty.lines=0)
Plot scape model fit to length-at-age data.
plotLA(model, together=FALSE, sex=NULL, axes=TRUE, same.limits=TRUE, between=list(x=axes,y=axes), ylim=NULL, bands=1, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=1, tck=0, tick.number=5, lty.grid=3, col.grid="gray", pch=16, cex.points=0.5, col.points="black", lty.lines=1, lwd.lines=4, col.lines=c("red","blue"), lty.bands=2*(!together), lwd.bands=1, col.bands="black", plot=TRUE, ...)plotLA(model, together=FALSE, sex=NULL, axes=TRUE, same.limits=TRUE, between=list(x=axes,y=axes), ylim=NULL, bands=1, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=1, tck=0, tick.number=5, lty.grid=3, col.grid="gray", pch=16, cex.points=0.5, col.points="black", lty.lines=1, lwd.lines=4, col.lines=c("red","blue"), lty.bands=2*(!together), lwd.bands=1, col.bands="black", plot=TRUE, ...)
model |
fitted |
together |
whether to plot both sexes in one panel. |
sex |
string indicating which sex to plot (both by default). |
axes |
whether to plot axis values. |
same.limits |
whether panels should have same y-axis limits. |
between |
list with |
ylim |
vector with lower and upper y-axis limits. |
bands |
extent of error bands relative to standard error. |
main |
main title. |
xlab |
x-axis label. |
ylab |
y-axis label. |
cex.main |
size of main title. |
cex.lab |
size of axis labels. |
cex.axis |
size of tick labels. |
cex.strip |
size of strip labels. |
col.strip |
color of strip labels. |
strip |
logical flag (whether to plot strip labels), or a
function passed to |
las |
orientation of tick labels: 0=parallel, 1=horizontal, 2=perpendicular, 3=vertical. |
tck |
tick mark length. |
tick.number |
number of tick marks. |
lty.grid |
line type of gridlines. |
col.grid |
color of gridlines. |
pch |
symbol for points, possibly a vector where element 2 refers to males. |
cex.points |
size of points, possibly a vector where element 2 refers to males. |
col.points |
color of points, possibly a vector where element 2 refers to males. |
lty.lines |
line type of main lines, possibly a vector where element 2 refers to males. |
lwd.lines |
line width of main lines, possibly a vector where element 2 refers to males. |
col.lines |
color of main lines, possibly a vector where element 2 refers to males. |
lty.bands |
line type of error bands. |
lwd.bands |
line width of error bands. |
col.bands |
color of error bands, possibly a vector where element 2 refers to males. |
plot |
whether to draw plot. |
... |
passed to |
When plot=TRUE, a trellis plot is drawn and a data frame is
returned, containing the data used for plotting. When
plot=FALSE, a trellis object is returned.
The Args function from the gdata package is recommended
for reviewing the arguments, instead of args.
xyplot, panel.xyplot,
and panel.superpose are the underlying drawing
functions.
plotCA, plotCL, plotIndex,
and plotLA plot model fit and data.
plotB, plotN, and plotSel
plot derived quantities.
scape-package gives an overview of the package.
plotLA(x.oreo, xlab="Age", ylab="Length (cm)") mykey <- list(text=list(lab=c("Female","Male")), space="right", lines=list(lwd=4,col=c("red","blue"))) plotLA(x.oreo, together=TRUE, xlab="Age", ylab="Length (cm)", pch=NA, key=mykey) mykey <- list(text=list(lab=c("Female","Male")), space="right", points=list(pch=16,cex=0.5,col=c("red","blue"))) plotLA(x.oreo, together=TRUE, xlab="Age", ylab="Length (cm)", col.points=c("red","blue"), lty.lines=0, key=mykey)plotLA(x.oreo, xlab="Age", ylab="Length (cm)") mykey <- list(text=list(lab=c("Female","Male")), space="right", lines=list(lwd=4,col=c("red","blue"))) plotLA(x.oreo, together=TRUE, xlab="Age", ylab="Length (cm)", pch=NA, key=mykey) mykey <- list(text=list(lab=c("Female","Male")), space="right", points=list(pch=16,cex=0.5,col=c("red","blue"))) plotLA(x.oreo, together=TRUE, xlab="Age", ylab="Length (cm)", col.points=c("red","blue"), lty.lines=0, key=mykey)
Plot scape model predicted numbers at age.
plotN(model, what="d", swap=FALSE, years=NULL, ages=NULL, axes=TRUE, same.limits=TRUE, div=1, log=FALSE, base=10, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=(what=="b"), tck=c(1,what=="b")/2, tick.number=10, lty.grid=3, col.grid="white", pch=16, cex.points=1, col.points="black", ratio.bars=3, col.bars="gray", plot=TRUE, ...)plotN(model, what="d", swap=FALSE, years=NULL, ages=NULL, axes=TRUE, same.limits=TRUE, div=1, log=FALSE, base=10, main="", xlab="", ylab="", cex.main=1.2, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=(what=="b"), tck=c(1,what=="b")/2, tick.number=10, lty.grid=3, col.grid="white", pch=16, cex.points=1, col.points="black", ratio.bars=3, col.bars="gray", plot=TRUE, ...)
model |
fitted |
what |
what to plot: |
swap |
whether to swap ages and years between axes or panels. |
years |
vector of numbers indicating which years to plot (all by default). |
ages |
vector of numbers indicating which ages to plot (all by default). |
axes |
whether to plot axis values. |
same.limits |
whether panels should have same y-axis limits. |
div |
denominator to shorten values on the y axis. |
log |
whether to log-transform values. |
base |
logarithm base. |
main |
main title. |
xlab |
x-axis label. |
ylab |
y-axis label. |
cex.main |
size of main title. |
cex.lab |
size of axis labels. |
cex.axis |
size of tick labels. |
cex.strip |
size of strip labels. |
col.strip |
color of strip labels. |
strip |
logical flag (whether to plot strip labels), or a
function passed to |
las |
orientation of tick labels: 0=parallel, 1=horizontal, 2=perpendicular, 3=vertical. |
tck |
tick mark length. |
tick.number |
number of tick marks. |
lty.grid |
line type of gridlines. |
col.grid |
color of gridlines. |
pch |
symbol for points. |
cex.points |
size of points. |
col.points |
color of points. |
ratio.bars |
width of bars. |
col.bars |
color of bars. |
plot |
whether to draw plot. |
... |
passed to |
The "d"[efault] plot is a combination of
"i"[nitial year] and "r"[ecruitment].
When plot=TRUE, a trellis plot is drawn and a data frame is
returned, containing the data used for plotting. When
plot=FALSE, a trellis object is returned.
The Args function from the gdata package is recommended
for reviewing the arguments, instead of args.
xyplot, panel.xyplot,
and panel.barchart are the underlying drawing
functions.
plotCA, plotCL, plotIndex,
and plotLA plot model fit and data.
plotB, plotN, and plotSel plot
derived quantities.
scape-package gives an overview of the package.
plotN(x.cod, div=1000, xlab=c("Age (years)","Year"), ylab="Individuals (million)") plotN(x.cod, "l", div=1000, xlab="Age", ylab="Individuals (million)") plotN(x.cod, "r", age=3, div=1000, xlim=c(1967,2002)) plotN(x.cod, "p", div=1000, ages=3:10, xlim=c(2,11), xlab="Age", ylab="Individuals (million)", cex.strip=0.7, cex.axis=0.7, tck=0.5) plotN(x.cod, "b", xlab="Age (years)", ylab="Year", cex.points=0.7)plotN(x.cod, div=1000, xlab=c("Age (years)","Year"), ylab="Individuals (million)") plotN(x.cod, "l", div=1000, xlab="Age", ylab="Individuals (million)") plotN(x.cod, "r", age=3, div=1000, xlim=c(1967,2002)) plotN(x.cod, "p", div=1000, ages=3:10, xlim=c(2,11), xlab="Age", ylab="Individuals (million)", cex.strip=0.7, cex.axis=0.7, tck=0.5) plotN(x.cod, "b", xlab="Age (years)", ylab="Year", cex.points=0.7)
Plot scape model predicted selectivity and observed maturity.
plotSel(model, together=FALSE, series=NULL, sex=NULL, axes=TRUE, legend="bottom", main="", xlab="", ylab="", cex.main=1.2, cex.legend=1, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=1, tck=0, tick.number=5, lty.grid=3, col.grid="gray", pch="m", cex.points=1, col.points="black", lty.lines=1, lwd.lines=4, col.lines=c("red","blue"), plot=TRUE, ...)plotSel(model, together=FALSE, series=NULL, sex=NULL, axes=TRUE, legend="bottom", main="", xlab="", ylab="", cex.main=1.2, cex.legend=1, cex.lab=1, cex.axis=0.8, cex.strip=0.8, col.strip="gray95", strip=strip.custom(bg=col.strip), las=1, tck=0, tick.number=5, lty.grid=3, col.grid="gray", pch="m", cex.points=1, col.points="black", lty.lines=1, lwd.lines=4, col.lines=c("red","blue"), plot=TRUE, ...)
model |
fitted |
together |
whether to plot gears in one panel. |
series |
vector of strings indicating which gears or surveys to plot (all by default). |
sex |
string indicating which sex to plot (both by default). |
axes |
whether to plot axis values. |
legend |
legend location: |
main |
main title. |
xlab |
x-axis label. |
ylab |
y-axis label. |
cex.main |
size of main title. |
cex.legend |
size of legend text. |
cex.lab |
size of axis labels. |
cex.axis |
size of tick labels. |
cex.strip |
size of strip labels. |
strip |
logical flag (whether to plot strip labels), or a
function passed to |
col.strip |
color of strip labels. |
las |
orientation of tick labels: 0=parallel, 1=horizontal, 2=perpendicular, 3=vertical. |
tck |
tick mark length. |
tick.number |
number of tick marks. |
lty.grid |
line type of gridlines. |
col.grid |
color of gridlines. |
pch |
symbol for points. |
cex.points |
size of points. |
col.points |
color of points. |
lty.lines |
line type of main lines. |
lwd.lines |
line width of main lines. |
col.lines |
color of main lines. |
plot |
whether to draw plot. |
... |
passed to |
When plot=TRUE, a trellis plot is drawn and a data frame is
returned, containing the data used for plotting. When
plot=FALSE, a trellis object is returned.
The Args function from the gdata package is recommended
for reviewing the arguments, instead of args.
xyplot, panel.points,
panel.lines, and
panel.superpose are the underlying drawing
functions.
plotCA, plotCL, plotIndex,
and plotLA plot model fit and data.
plotB, plotN, and plotSel plot
derived quantities.
scape-package gives an overview of the package.
plotSel(x.ling, xlab="Age", ylab="Selectivity and maturity") plotSel(x.cod, together=TRUE, xlab="Age\n", ylab="Selectivity", pch=NA, col.lines=c("coral","navyblue"), strip=FALSE)plotSel(x.ling, xlab="Age", ylab="Selectivity and maturity") plotSel(x.cod, together=TRUE, xlab="Age\n", ylab="Selectivity", pch=NA, col.lines=c("coral","navyblue"), strip=FALSE)
Stock assessment data and model fit for cod (Gadus morhua) in Icelandic waters, using a Coleraine statistical catch-at-age model.
This is a single-sex model with 10 age classes, the catch data starting in 1971 and ending in 2003. The model was fitted to three data components: survey abundance index, commercial catch at age, and survey catch at age.
x.codx.cod
List of class scape containing:
N
|
predicted numbers at age |
B
|
predicted biomass, recruitment, and observed landings (year things) |
Sel
|
predicted selectivity and observed maturity (age things) |
Dev
|
predicted recruitment deviates from the stock-recruitment curve |
Survey
|
survey abundance index and fit |
CAc
|
commercial C@A (catch at age) and fit |
CAs
|
survey C@A and fit |
Hilborn et al. (2003) give a general description of the Coleraine generalized model. Magnusson and Hilborn (2007) describe the model used in this assessment.
A maturity vector of zeros and ones was used to predict the biomass of age 4 and older, the quantity of main interest for the management of this stock.
Estimated parameters: R0, Rinit, uinit, Sleft[commercial], Sfull[c], Sleft[survey], Sfull[s], q, and 41 recruitment deviates.
The list was imported from the file ‘scape/example/cod.res’ using
importCol.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(x.cod),
ll(x.cod$N), and head(x.cod$N).
Magnusson, A. (2003) Coleraine assessment of the Icelandic cod stock. ICES North Western Working Group (NWWG) WD 31.
Marine Research Institute (2003) State of marine stocks in Icelandic waters 2002/2003. Marine Research in Iceland 97.
Hilborn, R., Maunder, M., Parma, A., Ernst, B., Payne, J., and Starr, P. (2003). Coleraine: A generalized age-structured stock assessment model. User's manual version 2.0. University of Washington Report SAFS-UW-0116.
Magnusson, A. and Hilborn, R. (2007). What makes fisheries data informative? Fish and Fisheries, 8, 337–358. doi:10.1111/j.1467-2979.2007.00258.x
importCol was used to import the fitted model.
x.cod, x.ling, x.oreo,
x.saithe, and x.sbw are fitted
scape models to explore.
xmcmc and xproj are MCMC results and
projections for the x.cod model run.
scape-package gives an overview of the package.
plotB(x.cod) plotCA(x.cod, "c") plotCA(x.cod, "s") plotIndex(x.cod, "s") plotN(x.cod) plotSel(x.cod)plotB(x.cod) plotCA(x.cod, "c") plotCA(x.cod, "s") plotIndex(x.cod, "s") plotN(x.cod) plotSel(x.cod)
Stock assessment data and model fit for ling (Genypterus blacodes) in New Zealand waters, using a Coleraine statistical catch-at-age model.
This is a two-sex model with 30 age classes and 29 length classes, the catch data starting in 1973 and ending in 2000. The model was fitted to five data components: longline abundance index, survey abundance index, survey catch at age, longline catch at length, and trawl catch at length.
x.lingx.ling
List of class scape containing:
N
|
predicted numbers at age |
B
|
predicted biomass, recruitment, and observed landings (year things) |
Sel
|
predicted selectivity and observed maturity (age things) |
Dev
|
predicted recruitment deviates from the stock-recruitment curve |
CPUE
|
commercial abundance index and fit |
Survey
|
survey abundance index and fit |
CAs
|
survey C@A (catch at age) and fit |
CLc
|
commercial C@L (catch at length) and fit |
Hilborn et al. (2003) give a general description of the Coleraine generalized model.
Estimated parameters: R0, Rinit, Sleft[trawl], Sfemale[t], Smale[t], Sright[t], Sleft[longline], Sfemale[l], Smale[l], Sright[l], Sleft[survey], Sfemale[s], Smale[s], Sright[s], q[l], q[s], and 29 recruitment deviates.
The list was imported from the file ‘scape/example/ling.res’
using importCol.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(x.ling),
ll(x.ling$N), and head(x.ling$N).
Annala, J.H., Sullivan, K.J., O'Brien, C.J., and Smith, N. W. M., eds. (2001). Report from the Fishery Assessment Plenary: Stock assessments and yield estimates. NIWA, Wellington.
Magnusson, A. (2001). SeaFIC assessment of Chatham Rise ling (LIN 3 and 4). Middle Depths Working Group Doc. 11. Report for the New Zealand Ministry of Fisheries.
Hilborn, R., Maunder, M., Parma, A., Ernst, B., Payne, J. and Starr, P. (2003) Coleraine: A generalized age-structured stock assessment model. User's manual version 2.0. University of Washington Report SAFS-UW-0116.
importCol was used to import the fitted model.
x.cod, x.ling, x.oreo,
x.saithe, and x.sbw are fitted
scape models to explore.
scape-package gives an overview of the package.
plotB(x.ling) plotCA(x.ling, "s") plotCL(x.ling, "c", series="1") plotCL(x.ling, "c", series="2") plotIndex(x.ling, "c") plotIndex(x.ling, "s") plotN(x.ling) plotSel(x.ling)plotB(x.ling) plotCA(x.ling, "s") plotCL(x.ling, "c", series="1") plotCL(x.ling, "c", series="2") plotIndex(x.ling, "c") plotIndex(x.ling, "s") plotN(x.ling) plotSel(x.ling)
Stock assessment data and model fit for smooth oreo (Pseudocyttus maculatus) in New Zealand waters, using a Coleraine statistical catch-at-age model.
This is a two-sex model with 80 age classes and 44 length classes, the catch data starting in 1979 and ending in 2001. The model was fitted to seven data components: pre-GPS commercial abundance index, post-GPS commercial abundance index, survey abundance index, commercial catch at length, survey catch at length, female length at age, and male length at age.
x.oreox.oreo
List of class scape containing:
N
|
predicted numbers at age |
B
|
predicted biomass, recruitment, and observed landings (year things) |
Sel
|
predicted selectivity and observed maturity (age things) |
CPUE
|
commercial abundance index and fit |
Survey
|
survey abundance index and fit |
CLc
|
commercial C@L (catch at length) and fit |
CLs
|
survey C@L and fit |
LA
|
observed L@A (length at age) |
Hilborn et al. (2003) give a general description of the Coleraine generalized model.
Since relatively few smooth oreo individuals have been aged, this assessment admits uncertainty about the von Bertalanffy growth curve, which is estimated for each sex. The acoustic survey abundance estimate is considered absolute, so q[survey] is fixed at 1.
Estimated parameters: R0, Sfull[commercial], Sfull[survey], q[pre-GPS], q[post-GPS], L80female, L80male, Kfemale, Kmale, CVfemale, and CVmale.
The list was imported from the files ‘scape/example/oreo.res’,
‘oreo.txt’, and ‘l_at_age.dat’ using importCol.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(x.oreo),
ll(x.oreo$N), and head(x.oreo$N).
Annala, J.H., Sullivan, K.J., O'Brien, C.J., Smith, N.W.M., and Grayling S.M., eds. (2003). Report from the Fishery Assessment Plenary: Stock assessments and yield estimates. Ministry of Fisheries, Wellington.
Hilborn, R., Maunder, M., Parma, A., Ernst, B., Payne, J., and Starr, P. (2003). Coleraine: A generalized age-structured stock assessment model. User's manual version 2.0. University of Washington Report SAFS-UW-0116.
Magnusson, A. (2003). Stock assessment of Chatham Rise smooth oreo (SSO4). Deepwater Working Group Doc. 16. Report for the New Zealand Ministry of Fisheries.
importCol was used to import the fitted model.
x.cod, x.ling, x.oreo,
x.saithe, and x.sbw are fitted
scape models to explore.
scape-package gives an overview of the package.
plotB(x.oreo) plotCL(x.oreo, "c") plotCL(x.oreo, "s") plotIndex(x.oreo, "c") plotIndex(x.oreo, "s") plotLA(x.oreo) plotN(x.oreo) plotSel(x.oreo)plotB(x.oreo) plotCL(x.oreo, "c") plotCL(x.oreo, "s") plotIndex(x.oreo, "c") plotIndex(x.oreo, "s") plotLA(x.oreo) plotN(x.oreo) plotSel(x.oreo)
Stock assessment data and model fit for saithe (Pollachius virens) in Icelandic waters, using an ADCAM-type (separable) catch-at-age model.
This is a single-sex model with 14 age classes, the catch data starting in 1980 and ending in 2014. The model was fitted to two data components: survey abundance index and commercial catch at age.
x.saithex.saithe
List of class scape containing:
N
|
predicted numbers at age |
B
|
predicted biomass, recruitment, and observed landings (year things) |
Sel
|
predicted selectivity and observed maturity (age things) |
Dev
|
predicted recruitment deviates from the stock-recruitment curve |
Survey
|
survey abundance index and fit |
CAc
|
commercial C@A (catch at age) and fit |
CAs
|
survey C@A and fit |
Bjornsson and Magnusson (2009) give a general description of the ADCAM modelling platform. The specific model used in this assessment is described in the ICES (2015) assessment report.
Fleet selectivity is modelled with three periods: 1980-1996, 1997-2003, and 2004 onwards.
The list was imported from results files in the directory
‘scape/example/saithe’ using importADCAM.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(x.saithe),
ll(x.saithe$N), and head(x.saithe$N).
ICES (2015) Report of the North-Western Working Group (NWWG). ICES CM 2015/ACOM:07.
Bjornsson, H. and A. Magnusson. (2009) ADCAM user manual (draft version). ICES CM 2009/ACOM:56, Annex 6.
importADCAM was used to import the fitted model.
x.cod, x.ling, x.oreo,
x.saithe, and x.sbw are fitted scape
models to explore.
scape-package gives an overview of the package.
plotB(x.saithe) plotCA(x.saithe, "c") plotCA(x.saithe, "s") plotIndex(x.saithe, "s") plotN(x.saithe) plotSel(x.saithe)plotB(x.saithe) plotCA(x.saithe, "c") plotCA(x.saithe, "s") plotIndex(x.saithe, "s") plotN(x.saithe) plotSel(x.saithe)
Stock assessment data and model fit for southern blue whiting (Micromesistius australis) in New Zealand waters, using a Coleraine statistical catch-at-age model.
This is a single-sex model with 11 age classes, the catch data starting in 1979 and ending in 2002. The model was fitted to two data components: survey abundance index and commercial catch at age.
x.sbwx.sbw
List of class scape containing:
N
|
predicted numbers at age |
B
|
predicted biomass, recruitment, and observed landings (year things) |
Sel
|
predicted selectivity and observed maturity (age things) |
Dev
|
predicted recruitment deviates from the stock-recruitment curve |
Survey
|
survey abundance index and fit |
CAc
|
commercial C@A (catch at age) and fit |
Hilborn et al. (2003) give a general description of the Coleraine generalized model.
The survey abundance index was preprocessed so that it contains only age 4 and older.
Estimated parameters: R0, Rinit, Rplus, Sleft[commercial], Sfull[c], q, and 33 recruitment deviates.
The list was imported from the file ‘scape/example/sbw.res’ using
importCol.
The functions ll (package gdata) and head are
recommended for browsing nested objects, e.g. ll(x.sbw),
ll(x.sbw$N), and head(x.sbw$N).
Annala, J. H., Sullivan, K. J., O'Brien, C. J., Smith, N. W. M. and Grayling S. M., eds. (2003) Report from the Fishery Assessment Plenary: Stock assessments and yield estimates. Ministry of Fisheries, Wellington.
Branch, T. A., Magnusson, A., Hilborn, R. and Starr, P. J. (2002) Stock assessment of the Campbell Island Rise population of southern blue whiting (Micromesistius australis) for the 2000–01 fishing season. University of Washington Report SAFS-UW-0107.
Hilborn, R., Maunder, M., Parma, A., Ernst, B., Payne, J. and Starr, P. (2003) Coleraine: A generalized age-structured stock assessment model. User's manual version 2.0. University of Washington Report SAFS-UW-0116.
Magnusson, A. and Hilborn, R. 2004. What is it in fisheries data that tells us about population abundance? Poster presented at the 4th World Fisheries Congress, Vancouver, BC.
importCol was used to import the fitted model.
x.cod, x.ling, x.oreo,
x.saithe, and x.sbw are fitted scape
models to explore.
scape-package gives an overview of the package.
plotB(x.sbw) plotCA(x.sbw, "c") plotIndex(x.sbw, "s") plotN(x.sbw) plotSel(x.sbw)plotB(x.sbw) plotCA(x.sbw, "c") plotIndex(x.sbw, "s") plotN(x.sbw) plotSel(x.sbw)
Markov chain Monte Carlo results from stock assessment of cod (Gadus morhua) in Icelandic waters.
xmcmcxmcmc
List containing four data frames:
L
|
likelihood components: CAc (commercial catch at
age), CAs (survey catch at age), Survey (survey
abundance index), Prior (penalty for priors), Total
(sum of all likelihood components; the objective function). |
P
|
estimated parameters: R0 (average virgin
recruitment), Rinit (initial recruitment scaler),
uinit (initial harvest rate), cSleft (left-side slope
of commercial selectivity curve), cSfull (age at full
commercial selectivity), sSleft (left-side slope of survey
selectivity curve), sSfull (age at full survey selectivity),
logq (log-transformed survey catchability). |
B
|
predicted biomass (age 4+) by year. |
R
|
predicted recruitment by year. |
See the x.cod help page for details about the data and
model.
The list was imported from the files ‘mcmclike.out’,
‘params.pst’, ‘spawbiom.pst’ and ‘recruits.pst’, using
the importMCMC() function. These files can be found in the
‘scape/example/mcmc’ directory.
The functions ll (package gdata) and head
are recommended for browsing nested objects, e.g. ll(xmcmc),
ll(xmcmc$P), and head(xmcmc$P).
The plotMCMC package is recommended for plotting MCMC diagnostics and posteriors.
Magnusson, A., Punt, A.E. and Hilborn, R. (2013). Measuring uncertainty in fisheries stock assessment: the delta method, bootstrap, and MCMC. Fish and Fisheries, 14, 325–342. doi:10.1111/j.1467-2979.2012.00473.x
importMCMC was used to import the MCMC results.
xmcmc and xproj are MCMC results and projections
for the x.cod model run.
scape-package gives an overview of the package.
# See examples in package 'plotMCMC'# See examples in package 'plotMCMC'
Markov chain Monte Carlo projections from stock assessment of cod (Gadus morhua) in Icelandic waters.
xprojxproj
List containing two lists:
B
|
projected biomass by year, given a constant harvest
rate policy: "0", "0.05", ..., "0.50". |
Y
|
projected catch by year, given a constant harvest rate
policy: "0", "0.05", ..., "0.50". |
MCMC projections can be used to evaluate the medium-term outcome of harvest policies (in this case, a certain fixed harvest rate), given the uncertainty about parameter values and random future recruitment.
See the x.cod help page for details about the data and
model.
The list was imported from the files ‘strategy.out’,
‘projspbm.out’ and ‘procatch.out’, using the
importProj() function. These files can be found in the
‘scape/example/mcmc’ directory.
The functions ll (package gdata) and head
are recommended for browsing nested objects, e.g. ll(xproj),
ll(xproj$B), ll(xproj$B$"0.25"), and
head(xproj$B$"0.25").
The plotMCMC package is recommended for plotting MCMC diagnostics and posteriors.
Magnusson, A., Punt, A. E. and Hilborn, R. (2013) Measuring uncertainty in fisheries stock assessment: the delta method, bootstrap, and MCMC. Fish and Fisheries, 14, 325–342.
importProj was used to import the MCMC projections.
xmcmc and xproj are MCMC results and projections
for the x.cod model run.
scape-package gives an overview of the package.
# See examples in package 'plotMCMC'# See examples in package 'plotMCMC'