This is a simple example of how to use \textbf{'MixFishSim'} to generate simulations of the dynamics in a mixed fishery. We describe how to calibrate the habitat fields, the population models, the fishery model and implement a simple fixed spatial closure. \
First, load the packages and set a seed for reproducibility.
library(MixFishSim) library(knitr) opts_chunk$set(tidy = TRUE) set.seed(123)
This vignette is a paired down example of how to construct a simulation using MixFishSim. We include only a basic example and encourage users to explore the other features of the package. \
First we specify the basic parameters of the simulation. This includes the dimensions of the spatial domain, the number of years to simulate, the number of fleets and vessels per fleet and the number of species and how often (in weeks) the fish move.
The object returned is used internally by MixFishSim a list with two levels:
sim <- init_sim(nrows = 10, ncols = 10, n_years = 10, n_tows_day = 4, n_days_wk_fished = 5, n_fleets = 2, n_vessels = 20, n_species = 2, move_freq = 2) class(sim) sim$idx names(sim$brk.idx)
This function creates the spatial fields which support the fish populations and determine their spatial distributions. You define the parameters for the matern covariance function for each population and optionally the location of any spawning closure areas.
It returns a list of suitable habitat for each species (hab), the habitat as adjusted during the spawning period (spwn_hab) and the binary location of spawning areas (spwn_loc). It also returns the locations as x1,x2,y1,y2 and the multiplier of attractiveness to the spawning area during spawning periods (spwn_mult).
If plot.dist = TRUE, it returns the plots to a file.
hab <- create_hab(sim_init = sim, spp.ctrl = list( "spp.1" = list('nu' = 1/0.015, 'var' = 1, 'scale' = 1, 'Aniso' = matrix(nc = 2, c(1.5, 3, -3, 4))), "spp.2" = list('nu' = 1/0.05, 'var' = 2, 'scale' = 12, 'Aniso' = matrix(nc = 2, c(1, 2, -1, 2))) ), spawn_areas = list( "spp1" = list( 'area1' = c(2,5,2,5), 'area2' = c(6,8,6,8) ), "spp2" = list( 'area1' = c(5,6,6,6) )), spwn_mult = 10, plot.dist = FALSE) print(hab) ## Plot the unadjusted habitat fields plot_habitat(hab$hab) ## Plot the adjusted habitat fields plot_habitat(hab$spwn_hab)
Now we need to set up the population models for the simulations. We do this with the init_pop function. We set the initial population biomasses, movement rates, recruitment parameter and growth and natural mortality rates.
The object created stores all the starting conditions and containers for recording the changes in the populations during the simulations.
We can plot the starting distributions for each population as a check.
Pop <- init_pop(sim_init = sim, Bio = c("spp1" = 1e5, "spp2" = 1e5), hab = hab[["hab"]], start_cell = c(5,5), lambda = c("spp1" = 0.1, "spp2" = 0.1), init_move_steps = 20, rec_params = list("spp1" = c("model" = "BH", "a" = 54, "b" = 2, "cv" = 0.7), "spp2" = c("model" = "BH", "a" = 27, "b" = 4,"cv" = 0.3)), rec_wk = list("spp1" = 3:6, "spp2" = 4:8), spwn_wk = list("spp1" = 4:8, "spp2" = 4:8), M = c("spp1" = 0.2/365, "spp2" = 0.2/365), wt = c("spp1"= 1, "spp2" = 1), wtm1 = c("spp1" = 0.1, "spp2" = 0.1), K = c("spp1" = 0.3, "spp2" = 0.3) ) names(Pop) Pop$dem_params image(Pop$Start_pop[[1]], main = "spp1 starting biomass") image(Pop$Start_pop[[2]], main = "spp2 starting biomass")
Now we set up the population tolerance to different temperatures which determines how the populations move during the course of a year. We can then plot the combined spatiotemporal suitable habitat to examine how these interact.
moveCov <- init_moveCov(sim_init = sim, steps = 52, spp_tol = list("spp1" = list("mu" = 12, "va" = 8), "spp2" = list("mu" = 15, "va" = 7) ) ) plot(norm_fun(x = 0:25, mu = 12, va = 8)/max(norm_fun(0:25, 12, 8)), type = "l", xlab = "Temperature", ylab = "Tolerance", lwd = 2) lines(norm_fun(x = 0:25, mu = 15, va = 7)/ max(norm_fun(0:25, 15, 7)), type = "l", col = "blue", lwd = 2) legend(x = 2, y = 0.9, legend = c("spp1", "spp2"), lwd = 2, col = c("black", "blue")) plot_spatiotemp_hab(hab = hab, moveCov = moveCov, spwn_wk = list("spp1" = 4:8, "spp2" = 4:8))
Here we initialise the fleet with fish landings price per tonne, catchability coefficients per population, fuel cost, the coefficients for the step function and fleet behaviour.
We can plot the behaviour of the step function to check its suitable for our simulations. This determines the relationship between the monetary value gained from a fishing tow and the next move by the vessel when using the correlated random walk function.
fleets <- init_fleet(sim_init = sim, VPT = list("spp1" = 4, "spp2" = 3), Qs = list("fleet 1" = c("spp1" = 1e-5, "spp2" = 3e-5), "fleet 2" = c("spp1" = 5e-5, "spp2" = 1e-5) ), fuelC = list("fleet1" = 3, "fleet 2" = 8), step_params = list("fleet 1" = c("rate" = 3, "B1" = 1, "B2" = 2, "B3" = 3), "fleet 2" = c("rate" = 3, "B1" = 2, "B2" = 4, "B3" = 4) ), past_knowledge = TRUE, past_year_month = TRUE, past_trip = TRUE, threshold = 0.7 ) test_step(step_params = fleets$fleet_params[[1]]$step_params, rev.max = 1e2) test_step(step_params = fleets$fleet_params[[2]]$step_params, rev.max = 1e2)
We set up a spatial closure. There are multiple options in defining this, but we simply define a static fixed site closure for demonstration purposes.
closure <- init_closure(input_coords = data.frame(x = c(9,10), y = c(6,10)), spp1 = "spp1", year_start = 5)
Its also possible to define a survey design using the init_survey function, but we do not do so for this demonstration. Please refer to the function help file if this is required.
Finally we run the simulation. The output is a list of objects containing all the information on fisheries catches, the population dynamics and population distributions. These can be examined with some inbuilt plotting functions.
res <- run_sim(sim_init = sim, pop_init = Pop, move_cov = moveCov, fleets_init = fleets, hab_init = hab, save_pop_bio = TRUE, survey = NULL, closure = closure)
There are a series of input plotting functions to visualise the results of the simulation. For example, we can explore:
Users will wish to define their own plots, depending on the issues of interest and all the results are saved in the output from the run_sim function.
## Biological p1 <- plot_pop_summary(results = res, timestep = "annual", save = FALSE ) p1 p2 <- plot_daily_fdyn(res) p2 ## Fishery logs <- combine_logs(res[["fleets_catches"]]) p3 <- plot_vessel_move(sim_init = sim, logs = logs, fleet_no = 1, vessel_no = 5, year_trip = 5, trip_no = 10) p3 p4 <- plot_realised_stepF(logs = logs, fleet_no = 1, vessel_no = 1) p4
Note in our example how the fishing mortality rate for species 2 changes following the spatial closure, which was set to cover some of the core distribution of the population.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.