In jae0/snowcrab: Snow crab assessment suite for aegis* and associated projects.
require(knitr)
knitr::opts_chunk$set(
root.dir = data_root,
echo = FALSE,
out.width="6.2in",
# dev.args = list(type = "cairo"),
fig.retina = 2,
dpi=192
)
# inits and data loading (front load all required data)
require(aegis)
year.assessment = params$year.assessment
year_previous = year.assessment - 1
p = bio.snowcrab::load.environment( year.assessment=year.assessment )
SCD = project.datadirectory("bio.snowcrab")
media_loc = params$media_loc
# fishery_model_results = file.path( "/home", "jae", "projects", "dynamical_model", "snowcrab", "outputs" )
fishery_model_results = file.path( SCD, "fishery_model" )
sn_env = snowcrab_load_key_results_to_memory( year.assessment, debugging=params$debugging, loc_dde=params$loc_dde, return_as_list=TRUE )
attach(sn_env)
# predator diet data
diet_data_dir = file.path( SCD, "data", "diets" )
require(data.table) # for speed
require(lubridate)
require(stringr)
require(gt) # table formatting
library(janitor)
require(ggplot2)
require(aegis) # map-related
require(bio.taxonomy) # handle species codes
# assimilate the CSV data tables:
# diet = get_feeding_data( diet_data_dir, redo=TRUE ) # if there is a data update
diet = get_feeding_data( diet_data_dir, redo=FALSE )
tx = taxa_to_code("snow crab")
# matching codes are
# spec tsn tx vern tx_index
#1 528 172379 BENTHODESMUS BENTHODESMUS 1659
#2 2522 98427 CHIONOECETES SPIDER QUEEN SNOW UNID 728
#3 2526 98428 CHIONOECETES OPILIO SNOW CRAB QUEEN 729
# 2 and 3 are correct
snowcrab_predators = diet[ preyspeccd %in% c(2522, 2526), ] # n=159 oservations out of a total of 58287 observations in db (=0.28% of all data)
snowcrab_predators$Species = code_to_taxa(snowcrab_predators$spec)$vern
snowcrab_predators$Predator = factor(snowcrab_predators$Species)
counts = snowcrab_predators[ , .(Frequency=.N), by=.(Species)]
setorderv(counts, "Frequency", order=-1)
# species composition
psp = speciescomposition_parameters( yrs=p$yrs, runlabel="1999_present" )
pca = speciescomposition_db( DS="pca", p=psp )
pcadata = as.data.frame( pca$loadings )
pcadata$vern = stringr::str_to_title( taxonomy.recode( from="spec", to="taxa", tolookup=rownames( pcadata ) )$vern )
# bycatch summaries
o_cfaall = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfaall" )
o_cfanorth = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfanorth" )
o_cfasouth = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfasouth" )
o_cfa4x = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfa4x" )
Life history {.c}
Photos {.c}
::: columns
:::: column
fn1=file.path( media_loc, "snowcrab_zoea.png" )
knitr::include_graphics( c(fn1 ) )
# \@ref(fig:photos)
::::
:::: column
fn2=file.path( media_loc, "snowcrab_male.png" )
fn3=file.path( media_loc, "snowcrab_male_and_female.png" )
knitr::include_graphics( c( fn2, fn3) )
# \@ref(fig:photos)
::::
:::
Life history stages{.c}
loc = file.path( Sys.getenv("HOME"), "projects", "dynamical_model", "snowcrab", "media" )
fn1=file.path( loc, "life_history.png" )
knitr::include_graphics( fn1 )
# \@ref(fig:lifehistory)
Growth stanzas {.c}
loc = file.path( Sys.getenv("HOME"), "projects", "dynamical_model", "snowcrab", "media" )
fn1=file.path( loc, "life_history_male.png" )
knitr::include_graphics( fn1 )
# \@ref(fig:lifehistory_male)
Growth modes{.c}
loc = file.path( Sys.getenv("HOME"), "bio.data", "bio.snowcrab", "output" )
fn1=file.path( loc, "size_structure", "growth_summary1.png" )
knitr::include_graphics( c(fn1) )
Growth modes 2 {.c}
loc = file.path( Sys.getenv("HOME"), "bio.data", "bio.snowcrab", "output" )
fn1=file.path( loc, "size_structure", "growth_summary2.png" )
knitr::include_graphics( c(fn1) )
Notable traits
- Sexual dimorphism
- Pelagic in larval stages; benthic in pre-adolescent and adult stages
- Biannual, annual molts depending upon size/age and environmental conditions
- Terminal molt to maturity and survives for up to 5 years
- Total life span: up to 15 years.
- Stenothermic (narrow) temperature requirements < 6 C
- Ontogenetic shifts in habitat preferences: warmer, complex substrates to mud.
- Cannibalism of immature crab by mature female snow crab are known
- Primiparous female (57.4 mm CW) produces between 35,000 to 46,000 eggs
- Multiparous females more fecund (>100,000 eggs)
- Eggs extruded February and April and brooded for up to two years (> 80% in SSE follow an annual cycle) and hatched/released from April to June.
Life history: Spatially and temporally complex
- 14+ instars/life stages with different preferences;
- Multiply prey and predators (ontogenetic shift); many ecosystem changes
- Growth: biannual (up to instar 5), annual, biennial (skip moulters)
- Reproduction: annual and biennial; strong sexual selection (dimorphism) and parental investment for up to 2 years; lags of up to 10 years
- Mortality: pulsed due to moult vulnerability, cannibalism
- Movement: seasonal, annual, ontogenetic
- Evolution: sexual dimorphism, parental investment of brood, stenothermic
- Population structure integrates all these effects ("Cumulative Effects")
Population structure integrates everything
- Spatial structures: scales > 100's km; local crab holes ( 10's to 100's m )
- Confluence of major oceanic currents (temperature, salinity, nutrients, primary and secondary production)
- Temporal processes: scales > 10 yrs; daily, tidal currents,
- 10-15 y life cycle, seasonal, annual, multiyear, decadal, evolutionary, etc.
- Biological (individual, age classes, population, etc.): cross multiple space, time and organismal scales
- Ontogeny (developmental change and habitat requirements)
- Life history (habitat preferences)
- Inter-specific interactions (feeding relationships change with biological stage)
- Interactions with humans
Multiple processes at multiple scales
- Populations integrate these processes
- Induces auto-correlations at multiple spatial-temporal-organisational scales
- We perceive this as "autocorrelation"
- First law of geography: "everything is related to everything else, but near things are more related than distant things" (Tobler 1970)
- "Nonstationary" processes (mean and variance not stable in space and time)
- Ignoring correlated errors causes inappropriate parameter estimates
Observations error (precision, accuracy)
- Pragmatic compromise between costs vs information gain
- SSE snow crab domain with ~400 stations and ~109,120 $\text{km}^{2}$
- Each station represents 273 $\text{km}^{2}$ , but actually samples ~0.0039 $\text{km}^{2}$
- A factor of 1:70,000 (spatially)
- Each sample is about 5 minutes but represents a whole year (525,600 minutes)
- A factor of 1:105,120
- In space-time:
- A factor of ~ 1:10 billion (note: 1 Angstrom = 10$^{-10}$m) .. the same scale as an atom!
- A few samples will do if it is spatially homogenous ("well mixed" = constant forcings, no ecosystem variability, no spatiotemporal autocorrelation)
- Spatial bias/aliasing
- Temporal bias/aliasing ("Spring" before 2004)
Solutions: Experimental design
- Naive (Ignore everything) "Random" Sampling / Fixed Station Sampling (IID)
- Stratified (Hope this works) "Random" Sampling -- depth is all that matters (IID within, IID between)
- Ignore spatiotemporal structure
- Ignore Ecosystem variability
- Assume no autocorrelated errors
- Geostatistical design (gridded, pseudo-random)
- Address spatial autocorrelation as a stationary process
- Abundance is NOT first and second order stationary and non-Gaussian (distributional bias)
- Ignores time (temporal discretization/aliasing)
- Ignores ecosystem variability
- Variogram solutions unstable when population size declines
Solutions: Experimental design 2
- Hierarchical process temporal and covariate as an external forcing and spatial as a local process (UKED)
- Variogram instability across time/stage when abundance declined and distribution was spotty
- NOT second-order stationary
- Mosaic process:
- Hierarchical first and second-order non-stationary model (modular)
- Global Hurdle/INLA/GAM/GLM
- Local Random spatial effects (residuals) Matern (SPDE; FFT)
- Independent space and time process
- Slow
- Model-based inference: Bayesian hierarchical mixed-effects model
- Conditional Autoregressive Space-Time Models (CARSTM, separable: Kronecker product space BYM2, time AR1)
- Hurdle/Poisson and Binomial number process, with Gaussian weight process
Latent ecological process
-
Real (latent, unobserved) ecological processes are generaly spatiotemporal processes
-
Observations from samples/surveys are used to infer the real (latent, unobserved) state
-
Stock assessments, almost always, focus only upon a purely temporal process by integrating the latent spatiotemporal process
-
Experimental design assumes samples are random or random-stratified, almost always, as a purely spatial process
-
This is a problem:
- temporal aliasing is likely (as it is usually ignored)
- spatial aliasing is likely when sampling environments that cannot be randomly sampled
- many locations are not directly accessible to trawls (rocks and bedrock)
- easier to sample locations (softer, gravel or mud substrates) coincide with preferred snow habitats
Snow crab survey locations {.c}
::: columns
:::: column
loc = file.path( SCD, "output", "maps", "survey.locations" )
yrsplot = setdiff( year.assessment + c(0:-9), 2020)
fn6 = file.path( loc, paste( "survey.locations", yrsplot[6], "png", sep=".") )
fn5 = file.path( loc, paste( "survey.locations", yrsplot[5], "png", sep=".") )
fn4 = file.path( loc, paste( "survey.locations", yrsplot[4], "png", sep=".") )
fn3 = file.path( loc, paste( "survey.locations", yrsplot[3], "png", sep=".") )
fn2 = file.path( loc, paste( "survey.locations", yrsplot[2], "png", sep=".") )
fn1 = file.path( loc, paste( "survey.locations", yrsplot[1], "png", sep=".") )
knitr::include_graphics( c( fn2, fn1) )
# \@ref(fig:survey-locations-map)
::::
:::: column
fn1 = file.path( media_loc, "carstm_prediction_domain.png" )
knitr::include_graphics( c( fn1) )
::::
:::
Clustering {.c}
```r aggregation for moulting and migration and Alaska red king crab \emph{Paralithodes camtschaticus} aggregation in Alaska for egg release, migrations.' }
fn1=file.path( media_loc, "australian_leptomithrax_gaimardii.png" )
fn2=file.path( media_loc, "kingcrab_aggregation.png" )
knitr::include_graphics( c(fn1, fn2 ) )
\@ref(fig:aggregation)
- *Other* crab species show "Mounding" for protection from predation (larval, moulting and females)
- Narrow habitat preferences force them to move and cluster when environment is poor
## Clustering2 {.c}
```r
fn = file.path(p$project.outputdir, "maps", "map_highdensity_locations.png" )
knitr::include_graphics( fn )
# \@ref(fig:aggregation)
if (0) {
# high density locations directly from databases
M = snowcrab.db( DS="set.complete", p=p )
setDT(M)
i = which(M$totno.all > 2.5*10^5)
H = M[i, .( plon, plat, towquality, dist, distance, surfacearea, vessel, yr, z, julian, no.male.all, no.female.all, cw.mean, totno.all, totno.male.imm, totno.male.mat, totno.female.imm, totno.female.mat, totno.female.primiparous, totno.female.multiparous, totno.female.berried)]
H$log10density = log10(H$totno.all)
library(ggplot2)
cst = coastline_db( p=p, project_to=st_crs(pg) )
isodepths = c(100, 200, 300)
isob = isobath_db( DS="isobath", depths=isodepths, project_to=st_crs(pg))
isob$level = as.factor( isob$level)
plt = ggplot() +
geom_sf( data=cst, show.legend=FALSE ) +
geom_sf( data=isob, aes( alpha=0.1, fill=level), lwd=0.1, show.legend=FALSE) +
geom_point(data=H, aes(x=plon, y=plat, colour=log10density), size=5) +
coord_sf(xlim = c(270, 940 ), ylim = c(4780, 5200 )) +
theme(legend.position=c(0.08, 0.8))
png(filename=fn, width=1000,height=600, res=144)
(plt)
dev.off()
}
Sampling bias: depth
loc = file.path( homedir, "projects", "snowcrabframework" )
knitr::include_graphics( file.path( loc, "bias_depth.png" ) )
Sampling bias: substrate grain size
loc = file.path( homedir, "projects", "snowcrabframework" )
knitr::include_graphics( file.path( loc, "bias_substrate.png" ) )
Sampling bias: bottom temperature
loc = file.path( homedir, "projects", "snowcrabframework" )
knitr::include_graphics( file.path( loc, "bias_temp.png" ) )
Sampling bias: species composition 1 (temperature related)
loc = file.path( homedir, "projects", "snowcrabframework" )
knitr::include_graphics( file.path( loc, "bias_pca1.png" ) )
Sampling bias: species composition 2 (depth related)
loc = file.path( homedir, "projects", "snowcrabframework" )
knitr::include_graphics( file.path( loc, "bias_pca2.png" ) )
Sampling bias: summary
Indications of sampling bias relative to spatial domain of snow crab:
- deeper and more favourable locations are being sampled
- statistically overlapping but biologically important
Generalized linear models (GLM)
\small
::: columns
:::: column
$S={S_{1},...,S_{K}}$ is a set of $k=1,...,K$ non-overlapping (areal) units.
$\boldsymbol{Y}=(y_{1},...,y_{K})$ are observations on $S$, then:
$$\begin{aligned}
Y & \sim f(y|\Omega)
g(\mu) &=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}.
\end{aligned}$$
-
$\mu = \text{E}(Y)$ is the expected value of $\boldsymbol{Y}$
-
$f(\cdot)$ indicates an exponential function
-
$\Omega$ is the set of the parameters of the function $f(\cdot)$
::::
:::: column
-
$g(\cdot) = f^{-1}(\cdot)$ is the linearizing link function
-
$\boldsymbol{x=}(x_{kv})\in\Re^{K\times V}$ is the matrix of covariates
-
$\boldsymbol{\beta}$ are the $V$ covariate parameters with MVN prior, with mean $\mu_{\beta}$ and diagonal variance matrix $\Sigma_{\beta}$
-
$\boldsymbol{O=}(o_{1},...,o_{K}\boldsymbol{)}$ are offsets, if any
-
$\boldsymbol{\varepsilon}=(\varepsilon_{1},...,\varepsilon_{K})$ are residual errors, if any
::::
:::
\normalsize
Generalized linear models (GLM) ...
\small
::: columns
:::: column
For each distributional family:
\vspace{2mm}
$Y\sim\text{Normal}(\mu,\sigma^{2})$
- $\mu=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}$
$Y\sim\text{Poisson}(\mu)$, and
- $\text{ln}(\mu)=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}$.
$Y\sim\text{Binomial}(\eta,\theta)$
- $\text{ln}(\theta/(1-\theta))=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}$,
- $\eta$ is the vector of number of trials
- $\theta$ the vector of probabilities of success in each trial
::::
:::: column
Stratified random sampling assumes:
$$\begin{aligned}
\varepsilon_{s} & \sim\text{N}(0,{}^{\varepsilon}!\sigma_{s}^{2})
\end{aligned}$$
i.e., IID errors in space (time is usually ignored)
::::
:::
\normalsize
Autocorrelated spatial errors
\small
Kriging extends $\varepsilon$ to spatial contraints ("variograms").
$$\begin{aligned}
Y_{t} & \sim\text{MVN}(\boldsymbol{\mu}{t},\boldsymbol{\Sigma}{t})\
g(\mu_{t}) & =\boldsymbol{x_{t}}^{T}\boldsymbol{\beta_{t}}+\boldsymbol{\omega}{t}+\boldsymbol{\varepsilon}{t}\
\varepsilon_{t} & \sim N(0,{}^{\varepsilon}!\sigma_{t}^{2})\
\omega_{t} & \sim\text{GP}(\boldsymbol{0},C(s_{t},s_{t}';^{\omega}!\theta_{t}))\
\boldsymbol{\Sigma_{t}} & =\left[C(\text{s}{it},\text{s}{jt};^{\omega}!\theta_{t})\right]{i,j=1}^{K}+{}^{\varepsilon}!\sigma{t}^{2}I_{S} \
C(h){\text{Mat\'{e}rn}} & = ;^{\omega}\sigma^{2}\frac{1}{2^{\nu-1}\Gamma(\nu)}(\sqrt{2\nu}h/\phi)^{\nu}\ K{\nu}(\sqrt{2\nu}h/\phi).
\end{aligned}$$
\footnotesize
- ignores temporal structure, focus upon spatial-only process
-
assumes first and second order stationarity (incorrect)
-
problems when abundance declines and spatial autocorrelation changes
across time or becomes unstable or impossible to estimable
\normalsize
Reference Points ... {.c}
# odir = file.path( fishery_model_results, year.assessment, "logistic_discrete_historical" )
odir = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" )
fn1 = file.path( odir, 'plot_hcr_cfanorth.png' )
fn2 = file.path( odir, 'plot_hcr_cfasouth.png' )
fn3 = file.path( odir, 'plot_hcr_cfa4x.png' )
include_graphics(c(fn1, fn2, fn3) )
# \@ref(fig:logistic-hcr)
Prior and posterior comparisons: K
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' )
loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" )
fn1 = file.path(loc, 'plot_prior_K_cfanorth.png')
fn2 = file.path(loc, 'plot_prior_K_cfasouth.png')
fn3 = file.path(loc, 'plot_prior_K_cfa4x.png')
# knitr::include_graphics( c(fn1) )
knitr::include_graphics( c(fn1, fn2, fn3) )
Prior and posterior comparisons: r
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' )
loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" )
fn1 = file.path(loc, 'plot_prior_r_cfanorth.png')
fn2 = file.path(loc, 'plot_prior_r_cfasouth.png')
fn3 = file.path(loc, 'plot_prior_r_cfa4x.png')
knitr::include_graphics( c(fn1, fn2, fn3) )
Prior and posterior comparisons: q
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' )
loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" )
fn1 = file.path(loc, 'plot_prior_q1_cfanorth.png')
fn2 = file.path(loc, 'plot_prior_q1_cfasouth.png')
fn3 = file.path(loc, 'plot_prior_q1_cfa4x.png')
knitr::include_graphics( c(fn1, fn2, fn3) )
Prior and posterior comparisons: observaton error
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' )
loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" )
fn1 = file.path(loc, 'plot_prior_bosd_cfanorth.png')
fn2 = file.path(loc, 'plot_prior_bosd_cfasouth.png')
fn3 = file.path(loc, 'plot_prior_bosd_cfa4x.png')
knitr::include_graphics( c( fn1, fn2, fn3) )
Prior and posterior comparisons: process error
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' )
loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" )
fn1 = file.path(loc, 'plot_prior_bpsd_cfanorth.png')
fn2 = file.path(loc, 'plot_prior_bpsd_cfasouth.png')
fn3 = file.path(loc, 'plot_prior_bpsd_cfa4x.png')
knitr::include_graphics( c(fn1, fn2, fn3) )
jae0/snowcrab documentation built on Feb. 27, 2024, 2:42 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
require(knitr) knitr::opts_chunk$set( root.dir = data_root, echo = FALSE, out.width="6.2in", # dev.args = list(type = "cairo"), fig.retina = 2, dpi=192 ) # inits and data loading (front load all required data) require(aegis) year.assessment = params$year.assessment year_previous = year.assessment - 1 p = bio.snowcrab::load.environment( year.assessment=year.assessment ) SCD = project.datadirectory("bio.snowcrab") media_loc = params$media_loc # fishery_model_results = file.path( "/home", "jae", "projects", "dynamical_model", "snowcrab", "outputs" ) fishery_model_results = file.path( SCD, "fishery_model" ) sn_env = snowcrab_load_key_results_to_memory( year.assessment, debugging=params$debugging, loc_dde=params$loc_dde, return_as_list=TRUE ) attach(sn_env) # predator diet data diet_data_dir = file.path( SCD, "data", "diets" ) require(data.table) # for speed require(lubridate) require(stringr) require(gt) # table formatting library(janitor) require(ggplot2) require(aegis) # map-related require(bio.taxonomy) # handle species codes # assimilate the CSV data tables: # diet = get_feeding_data( diet_data_dir, redo=TRUE ) # if there is a data update diet = get_feeding_data( diet_data_dir, redo=FALSE ) tx = taxa_to_code("snow crab") # matching codes are # spec tsn tx vern tx_index #1 528 172379 BENTHODESMUS BENTHODESMUS 1659 #2 2522 98427 CHIONOECETES SPIDER QUEEN SNOW UNID 728 #3 2526 98428 CHIONOECETES OPILIO SNOW CRAB QUEEN 729 # 2 and 3 are correct snowcrab_predators = diet[ preyspeccd %in% c(2522, 2526), ] # n=159 oservations out of a total of 58287 observations in db (=0.28% of all data) snowcrab_predators$Species = code_to_taxa(snowcrab_predators$spec)$vern snowcrab_predators$Predator = factor(snowcrab_predators$Species) counts = snowcrab_predators[ , .(Frequency=.N), by=.(Species)] setorderv(counts, "Frequency", order=-1) # species composition psp = speciescomposition_parameters( yrs=p$yrs, runlabel="1999_present" ) pca = speciescomposition_db( DS="pca", p=psp ) pcadata = as.data.frame( pca$loadings ) pcadata$vern = stringr::str_to_title( taxonomy.recode( from="spec", to="taxa", tolookup=rownames( pcadata ) )$vern ) # bycatch summaries o_cfaall = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfaall" ) o_cfanorth = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfanorth" ) o_cfasouth = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfasouth" ) o_cfa4x = observer.db( DS="bycatch_summary", p=p, yrs=p$yrs, region="cfa4x" )
Life history {.c}
Photos {.c}
::: columns :::: column
fn1=file.path( media_loc, "snowcrab_zoea.png" ) knitr::include_graphics( c(fn1 ) ) # \@ref(fig:photos)
:::: :::: column
fn2=file.path( media_loc, "snowcrab_male.png" ) fn3=file.path( media_loc, "snowcrab_male_and_female.png" ) knitr::include_graphics( c( fn2, fn3) ) # \@ref(fig:photos)
:::: :::
Life history stages{.c}
loc = file.path( Sys.getenv("HOME"), "projects", "dynamical_model", "snowcrab", "media" ) fn1=file.path( loc, "life_history.png" ) knitr::include_graphics( fn1 ) # \@ref(fig:lifehistory)
Growth stanzas {.c}
loc = file.path( Sys.getenv("HOME"), "projects", "dynamical_model", "snowcrab", "media" ) fn1=file.path( loc, "life_history_male.png" ) knitr::include_graphics( fn1 ) # \@ref(fig:lifehistory_male)
Growth modes{.c}
loc = file.path( Sys.getenv("HOME"), "bio.data", "bio.snowcrab", "output" ) fn1=file.path( loc, "size_structure", "growth_summary1.png" ) knitr::include_graphics( c(fn1) )
Growth modes 2 {.c}
loc = file.path( Sys.getenv("HOME"), "bio.data", "bio.snowcrab", "output" ) fn1=file.path( loc, "size_structure", "growth_summary2.png" ) knitr::include_graphics( c(fn1) )
Notable traits
- Sexual dimorphism
- Pelagic in larval stages; benthic in pre-adolescent and adult stages
- Biannual, annual molts depending upon size/age and environmental conditions
- Terminal molt to maturity and survives for up to 5 years
- Total life span: up to 15 years.
- Stenothermic (narrow) temperature requirements < 6 C
- Ontogenetic shifts in habitat preferences: warmer, complex substrates to mud.
- Cannibalism of immature crab by mature female snow crab are known
- Primiparous female (57.4 mm CW) produces between 35,000 to 46,000 eggs
- Multiparous females more fecund (>100,000 eggs)
- Eggs extruded February and April and brooded for up to two years (> 80% in SSE follow an annual cycle) and hatched/released from April to June.
Life history: Spatially and temporally complex
- 14+ instars/life stages with different preferences;
- Multiply prey and predators (ontogenetic shift); many ecosystem changes
- Growth: biannual (up to instar 5), annual, biennial (skip moulters)
- Reproduction: annual and biennial; strong sexual selection (dimorphism) and parental investment for up to 2 years; lags of up to 10 years
- Mortality: pulsed due to moult vulnerability, cannibalism
- Movement: seasonal, annual, ontogenetic
- Evolution: sexual dimorphism, parental investment of brood, stenothermic
- Population structure integrates all these effects ("Cumulative Effects")
Population structure integrates everything
- Spatial structures: scales > 100's km; local crab holes ( 10's to 100's m )
- Confluence of major oceanic currents (temperature, salinity, nutrients, primary and secondary production)
- Temporal processes: scales > 10 yrs; daily, tidal currents,
- 10-15 y life cycle, seasonal, annual, multiyear, decadal, evolutionary, etc.
- Biological (individual, age classes, population, etc.): cross multiple space, time and organismal scales
- Ontogeny (developmental change and habitat requirements)
- Life history (habitat preferences)
- Inter-specific interactions (feeding relationships change with biological stage)
- Interactions with humans
Multiple processes at multiple scales
- Populations integrate these processes
- Induces auto-correlations at multiple spatial-temporal-organisational scales
- We perceive this as "autocorrelation"
- First law of geography: "everything is related to everything else, but near things are more related than distant things" (Tobler 1970)
- "Nonstationary" processes (mean and variance not stable in space and time)
- Ignoring correlated errors causes inappropriate parameter estimates
Observations error (precision, accuracy)
- Pragmatic compromise between costs vs information gain
- SSE snow crab domain with ~400 stations and ~109,120 $\text{km}^{2}$
- Each station represents 273 $\text{km}^{2}$ , but actually samples ~0.0039 $\text{km}^{2}$
- A factor of 1:70,000 (spatially)
- Each sample is about 5 minutes but represents a whole year (525,600 minutes)
- A factor of 1:105,120
- In space-time:
- A factor of ~ 1:10 billion (note: 1 Angstrom = 10$^{-10}$m) .. the same scale as an atom!
- A few samples will do if it is spatially homogenous ("well mixed" = constant forcings, no ecosystem variability, no spatiotemporal autocorrelation)
- Spatial bias/aliasing
- Temporal bias/aliasing ("Spring" before 2004)
Solutions: Experimental design
- Naive (Ignore everything) "Random" Sampling / Fixed Station Sampling (IID)
- Stratified (Hope this works) "Random" Sampling -- depth is all that matters (IID within, IID between)
- Ignore spatiotemporal structure
- Ignore Ecosystem variability
- Assume no autocorrelated errors
- Geostatistical design (gridded, pseudo-random)
- Address spatial autocorrelation as a stationary process
- Abundance is NOT first and second order stationary and non-Gaussian (distributional bias)
- Ignores time (temporal discretization/aliasing)
- Ignores ecosystem variability
- Variogram solutions unstable when population size declines
Solutions: Experimental design 2
- Hierarchical process temporal and covariate as an external forcing and spatial as a local process (UKED)
- Variogram instability across time/stage when abundance declined and distribution was spotty
- NOT second-order stationary
- Mosaic process:
- Hierarchical first and second-order non-stationary model (modular)
- Global Hurdle/INLA/GAM/GLM
- Local Random spatial effects (residuals) Matern (SPDE; FFT)
- Independent space and time process
- Slow
- Model-based inference: Bayesian hierarchical mixed-effects model
- Conditional Autoregressive Space-Time Models (CARSTM, separable: Kronecker product space BYM2, time AR1)
- Hurdle/Poisson and Binomial number process, with Gaussian weight process
Latent ecological process
-
Real (latent, unobserved) ecological processes are generaly spatiotemporal processes
-
Observations from samples/surveys are used to infer the real (latent, unobserved) state
-
Stock assessments, almost always, focus only upon a purely temporal process by integrating the latent spatiotemporal process
-
Experimental design assumes samples are random or random-stratified, almost always, as a purely spatial process
-
This is a problem:
- temporal aliasing is likely (as it is usually ignored)
- spatial aliasing is likely when sampling environments that cannot be randomly sampled
- many locations are not directly accessible to trawls (rocks and bedrock)
- easier to sample locations (softer, gravel or mud substrates) coincide with preferred snow habitats
Snow crab survey locations {.c}
::: columns :::: column
loc = file.path( SCD, "output", "maps", "survey.locations" ) yrsplot = setdiff( year.assessment + c(0:-9), 2020) fn6 = file.path( loc, paste( "survey.locations", yrsplot[6], "png", sep=".") ) fn5 = file.path( loc, paste( "survey.locations", yrsplot[5], "png", sep=".") ) fn4 = file.path( loc, paste( "survey.locations", yrsplot[4], "png", sep=".") ) fn3 = file.path( loc, paste( "survey.locations", yrsplot[3], "png", sep=".") ) fn2 = file.path( loc, paste( "survey.locations", yrsplot[2], "png", sep=".") ) fn1 = file.path( loc, paste( "survey.locations", yrsplot[1], "png", sep=".") ) knitr::include_graphics( c( fn2, fn1) ) # \@ref(fig:survey-locations-map)
:::: :::: column
fn1 = file.path( media_loc, "carstm_prediction_domain.png" ) knitr::include_graphics( c( fn1) )
:::: :::
Clustering {.c}
```r aggregation for moulting and migration and Alaska red king crab \emph{Paralithodes camtschaticus} aggregation in Alaska for egg release, migrations.' } fn1=file.path( media_loc, "australian_leptomithrax_gaimardii.png" ) fn2=file.path( media_loc, "kingcrab_aggregation.png" ) knitr::include_graphics( c(fn1, fn2 ) )
\@ref(fig:aggregation)
- *Other* crab species show "Mounding" for protection from predation (larval, moulting and females)
- Narrow habitat preferences force them to move and cluster when environment is poor
## Clustering2 {.c}
```r
fn = file.path(p$project.outputdir, "maps", "map_highdensity_locations.png" )
knitr::include_graphics( fn )
# \@ref(fig:aggregation)
if (0) {
# high density locations directly from databases
M = snowcrab.db( DS="set.complete", p=p )
setDT(M)
i = which(M$totno.all > 2.5*10^5)
H = M[i, .( plon, plat, towquality, dist, distance, surfacearea, vessel, yr, z, julian, no.male.all, no.female.all, cw.mean, totno.all, totno.male.imm, totno.male.mat, totno.female.imm, totno.female.mat, totno.female.primiparous, totno.female.multiparous, totno.female.berried)]
H$log10density = log10(H$totno.all)
library(ggplot2)
cst = coastline_db( p=p, project_to=st_crs(pg) )
isodepths = c(100, 200, 300)
isob = isobath_db( DS="isobath", depths=isodepths, project_to=st_crs(pg))
isob$level = as.factor( isob$level)
plt = ggplot() +
geom_sf( data=cst, show.legend=FALSE ) +
geom_sf( data=isob, aes( alpha=0.1, fill=level), lwd=0.1, show.legend=FALSE) +
geom_point(data=H, aes(x=plon, y=plat, colour=log10density), size=5) +
coord_sf(xlim = c(270, 940 ), ylim = c(4780, 5200 )) +
theme(legend.position=c(0.08, 0.8))
png(filename=fn, width=1000,height=600, res=144)
(plt)
dev.off()
}
Sampling bias: depth
loc = file.path( homedir, "projects", "snowcrabframework" ) knitr::include_graphics( file.path( loc, "bias_depth.png" ) )
Sampling bias: substrate grain size
loc = file.path( homedir, "projects", "snowcrabframework" ) knitr::include_graphics( file.path( loc, "bias_substrate.png" ) )
Sampling bias: bottom temperature
loc = file.path( homedir, "projects", "snowcrabframework" ) knitr::include_graphics( file.path( loc, "bias_temp.png" ) )
Sampling bias: species composition 1 (temperature related)
loc = file.path( homedir, "projects", "snowcrabframework" ) knitr::include_graphics( file.path( loc, "bias_pca1.png" ) )
Sampling bias: species composition 2 (depth related)
loc = file.path( homedir, "projects", "snowcrabframework" ) knitr::include_graphics( file.path( loc, "bias_pca2.png" ) )
Sampling bias: summary
Indications of sampling bias relative to spatial domain of snow crab:
- deeper and more favourable locations are being sampled
- statistically overlapping but biologically important
Generalized linear models (GLM)
\small
::: columns :::: column
$S={S_{1},...,S_{K}}$ is a set of $k=1,...,K$ non-overlapping (areal) units.
$\boldsymbol{Y}=(y_{1},...,y_{K})$ are observations on $S$, then:
$$\begin{aligned} Y & \sim f(y|\Omega) g(\mu) &=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}. \end{aligned}$$
-
$\mu = \text{E}(Y)$ is the expected value of $\boldsymbol{Y}$
-
$f(\cdot)$ indicates an exponential function
-
$\Omega$ is the set of the parameters of the function $f(\cdot)$
::::
:::: column
-
$g(\cdot) = f^{-1}(\cdot)$ is the linearizing link function
-
$\boldsymbol{x=}(x_{kv})\in\Re^{K\times V}$ is the matrix of covariates
-
$\boldsymbol{\beta}$ are the $V$ covariate parameters with MVN prior, with mean $\mu_{\beta}$ and diagonal variance matrix $\Sigma_{\beta}$
-
$\boldsymbol{O=}(o_{1},...,o_{K}\boldsymbol{)}$ are offsets, if any
-
$\boldsymbol{\varepsilon}=(\varepsilon_{1},...,\varepsilon_{K})$ are residual errors, if any
:::: :::
\normalsize
Generalized linear models (GLM) ...
\small ::: columns :::: column For each distributional family: \vspace{2mm}
$Y\sim\text{Normal}(\mu,\sigma^{2})$
- $\mu=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}$
$Y\sim\text{Poisson}(\mu)$, and
- $\text{ln}(\mu)=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}$.
$Y\sim\text{Binomial}(\eta,\theta)$
- $\text{ln}(\theta/(1-\theta))=\boldsymbol{x}^{T}\boldsymbol{\beta}+\boldsymbol{O}+\boldsymbol{\varepsilon}$,
- $\eta$ is the vector of number of trials
- $\theta$ the vector of probabilities of success in each trial ::::
:::: column
Stratified random sampling assumes:
$$\begin{aligned} \varepsilon_{s} & \sim\text{N}(0,{}^{\varepsilon}!\sigma_{s}^{2}) \end{aligned}$$
i.e., IID errors in space (time is usually ignored)
:::: :::
\normalsize
Autocorrelated spatial errors
\small
Kriging extends $\varepsilon$ to spatial contraints ("variograms").
$$\begin{aligned} Y_{t} & \sim\text{MVN}(\boldsymbol{\mu}{t},\boldsymbol{\Sigma}{t})\ g(\mu_{t}) & =\boldsymbol{x_{t}}^{T}\boldsymbol{\beta_{t}}+\boldsymbol{\omega}{t}+\boldsymbol{\varepsilon}{t}\ \varepsilon_{t} & \sim N(0,{}^{\varepsilon}!\sigma_{t}^{2})\ \omega_{t} & \sim\text{GP}(\boldsymbol{0},C(s_{t},s_{t}';^{\omega}!\theta_{t}))\ \boldsymbol{\Sigma_{t}} & =\left[C(\text{s}{it},\text{s}{jt};^{\omega}!\theta_{t})\right]{i,j=1}^{K}+{}^{\varepsilon}!\sigma{t}^{2}I_{S} \ C(h){\text{Mat\'{e}rn}} & = ;^{\omega}\sigma^{2}\frac{1}{2^{\nu-1}\Gamma(\nu)}(\sqrt{2\nu}h/\phi)^{\nu}\ K{\nu}(\sqrt{2\nu}h/\phi). \end{aligned}$$
\footnotesize - ignores temporal structure, focus upon spatial-only process
-
assumes first and second order stationarity (incorrect)
-
problems when abundance declines and spatial autocorrelation changes across time or becomes unstable or impossible to estimable
\normalsize
Reference Points ... {.c}
# odir = file.path( fishery_model_results, year.assessment, "logistic_discrete_historical" ) odir = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" ) fn1 = file.path( odir, 'plot_hcr_cfanorth.png' ) fn2 = file.path( odir, 'plot_hcr_cfasouth.png' ) fn3 = file.path( odir, 'plot_hcr_cfa4x.png' ) include_graphics(c(fn1, fn2, fn3) ) # \@ref(fig:logistic-hcr)
Prior and posterior comparisons: K
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' ) loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" ) fn1 = file.path(loc, 'plot_prior_K_cfanorth.png') fn2 = file.path(loc, 'plot_prior_K_cfasouth.png') fn3 = file.path(loc, 'plot_prior_K_cfa4x.png') # knitr::include_graphics( c(fn1) ) knitr::include_graphics( c(fn1, fn2, fn3) )
Prior and posterior comparisons: r
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' ) loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" ) fn1 = file.path(loc, 'plot_prior_r_cfanorth.png') fn2 = file.path(loc, 'plot_prior_r_cfasouth.png') fn3 = file.path(loc, 'plot_prior_r_cfa4x.png') knitr::include_graphics( c(fn1, fn2, fn3) )
Prior and posterior comparisons: q
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' ) loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" ) fn1 = file.path(loc, 'plot_prior_q1_cfanorth.png') fn2 = file.path(loc, 'plot_prior_q1_cfasouth.png') fn3 = file.path(loc, 'plot_prior_q1_cfa4x.png') knitr::include_graphics( c(fn1, fn2, fn3) )
Prior and posterior comparisons: observaton error
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' ) loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" ) fn1 = file.path(loc, 'plot_prior_bosd_cfanorth.png') fn2 = file.path(loc, 'plot_prior_bosd_cfasouth.png') fn3 = file.path(loc, 'plot_prior_bosd_cfa4x.png') knitr::include_graphics( c( fn1, fn2, fn3) )
Prior and posterior comparisons: process error
# loc = file.path(data_root, 'bio.snowcrab', 'fishery_model', '2023', 'logistic_discrete_historical' ) loc = file.path( homedir, "projects", "dynamical_model", "snowcrab", "outputs_alt", year.assessment, "logistic_discrete_historical" ) fn1 = file.path(loc, 'plot_prior_bpsd_cfanorth.png') fn2 = file.path(loc, 'plot_prior_bpsd_cfasouth.png') fn3 = file.path(loc, 'plot_prior_bpsd_cfa4x.png') knitr::include_graphics( c(fn1, fn2, fn3) )
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.