Skip to content

Commit

Permalink
Changed some font sizes, fixed plot bug, edited text
Browse files Browse the repository at this point in the history
  • Loading branch information
drfinlayscott committed Sep 16, 2020
1 parent 18fc1b3 commit 075dc0a
Show file tree
Hide file tree
Showing 5 changed files with 105 additions and 27 deletions.
6 changes: 3 additions & 3 deletions AMPLE/R/modules.R
Original file line number Diff line number Diff line change
Expand Up @@ -39,15 +39,15 @@ stoch_params_setterUI <- function(id, show_biol_prod_sigma = TRUE, show_biol_est
options <- list()
if (show_biol_prod_sigma){
options[[length(options)+1]] <-
tags$span(title="Natural variability in the productivity of the stock through changes in growth and natural mortality", numericInput(ns("biol_prod_sigma"), label = "Biological productivity variability", value = init_prod_sigma, min=0, max=1, step=0.05))
tags$span(title="Natural variability in the stock biological processes (e.g. recruitment and growth)", numericInput(ns("biol_prod_sigma"), label = "Biological variability", value = init_prod_sigma, min=0, max=1, step=0.05))
}
if (show_biol_est_sigma){
options[[length(options)+1]] <-
tags$span(title="Simulating the difference between the 'true' biomass and the 'estimated' biomass used by the HCR by applying randomly generated noise", numericInput(ns("biol_est_sigma"), label = "Estimation error variability", value = init_est_sigma, min=0, max=1, step=0.05))
tags$span(title="Simulating the difference between the 'true' biomass and the 'estimated' biomass used by the HCR by applying randomly generated noise", numericInput(ns("biol_est_sigma"), label = "Estimation variability", value = init_est_sigma, min=0, max=1, step=0.05))
}
if (show_biol_est_bias){
options[[length(options)+1]] <-
tags$span(title="Simulating the difference between the 'true' biomass and the 'estimated' biomass used by the HCR by applying a continuous bias (positive or negative)", numericInput(ns("biol_est_bias"), label = "Estimation error bias", value = init_est_bias, min=-0.5, max=0.5, step=0.05))
tags$span(title="Simulating the difference between the 'true' biomass and the 'estimated' biomass used by the HCR by applying a continuous bias (positive or negative)", numericInput(ns("biol_est_bias"), label = "Estimation bias", value = init_est_bias, min=-0.5, max=0.5, step=0.05))
}

vars <- conditionalPanel(condition="input.show_var == true", ns=ns, options)
Expand Down
86 changes: 79 additions & 7 deletions AMPLE/R/plots.R
Original file line number Diff line number Diff line change
Expand Up @@ -184,7 +184,7 @@ sideways_histogram <- function(dat, range, lhist=20, num.dnorm=5*lhist, dcol="bl
#' # Yield curve
#' plot_yieldcurve_projections(stock=stock2, stock_params=stock_params, app_params=app_params)
#' @export
plot_hcr <- function(stock, stock_params, mp_params, app_params, timestep=NULL, show_last=TRUE){
plot_hcr <- function(stock, stock_params, mp_params, app_params, timestep=NULL, show_last=TRUE, ...){

if (mp_params$output_type=="catch"){
ymax <- get_catch_ymax(stock$catch, mp_params)
Expand All @@ -202,7 +202,7 @@ plot_hcr <- function(stock, stock_params, mp_params, app_params, timestep=NULL,
#xlab <- "Estimated biomass / K"
xlab <- "Estimated SB/SBF=0"
# Plot empty axes
plot(x=xrange,y=yrange,type="n",xlab=xlab, ylab=ylab, xaxs="i", yaxs="i", main="The HCR")
plot(x=xrange,y=yrange,type="n",xlab=xlab, ylab=ylab, xaxs="i", yaxs="i", main="The HCR", ...)

# If model based MP (or NA) add the B/K based reference points
if (mp_params$mp_type %in% c(as.character(NA), "model")){
Expand Down Expand Up @@ -362,13 +362,14 @@ plot_biomass <- function(stock, stock_params, mp_params, timestep=NULL, show_las
#' @rdname front_page_plots
#' @name Front page plots
#' @export
# Used for IntroHCR and IntroIndicators and maybe others
plot_catch <- function(stock, stock_params, mp_params, app_params=NULL, timestep=NULL, show_last=TRUE, max_spaghetti_iters=50, quantiles, nspaghetti=5, add_grid=TRUE, xlab="Year", ghost_col="grey", true_col="black", ...){
years <- as.numeric(dimnames(stock$biomass)$year)
# Set Ylim - use same as HCR plot
ymax <- get_catch_ymax(stock$catch, mp_params)
yrange <- c(0, ymax)
# Empty axis
plot(x=years, y=years, type="n", ylim=c(0, ymax), ylab="Catch", xlab=xlab, xaxs="i", yaxs="i",...)
plot(x=years, y=years, type="n", ylim=c(0, ymax), ylab="Catch", xlab=xlab, xaxs="i", yaxs="i", ...)
if (add_grid){
grid()
}
Expand Down Expand Up @@ -412,6 +413,77 @@ plot_catch <- function(stock, stock_params, mp_params, app_params=NULL, timestep
}
}

# Generic timeseries plot
plot_indiv_timeseries_base <- function(data, stock, stock_params, mp_params, app_params=NULL, show_last=TRUE, max_spaghetti_iters=50, quantiles, nspaghetti=5, yrange, ylab, add_grid=TRUE, xlab="Year", ghost_col="grey", true_col="black", ...){

years <- as.numeric(dimnames(stock$biomass)$year)
# Plot empty axis
plot(x=years, y=years, type="n", ylim=c(yrange[1], yrange[2]), ylab=ylab, xlab=xlab, xaxs="i", yaxs="i",...)
if (add_grid){
grid()
}
# Get last iteration
last_iter <- dim(data)[1]
# If we have more than X iters, draw envelope of iters
if(last_iter > max_spaghetti_iters){
# Draw ribbon
draw_ribbon(x=years, y=data, quantiles=quantiles)
# Add spaghetti
for (iter in 1:nspaghetti){
lines(x=years, y=data[iter,], lty=spaghetti_lty, lwd=spaghetti_lwd, col=spaghetti_col)
}
}
# Else plot individual iters
else{
# Plot all iters as ghosts
if (last_iter > 1){
for (i in 1:last_iter){
lines(x=years, y=data[i,], col=ghost_col, lwd=ghost_lwd, lty=ghost_lty)
}
}
}
# Current iteration
if(show_last){
lines(x=years, y=data[last_iter,], col=true_col, lwd=last_lwd, lty=last_lty)
}

}

# Relative CPUE
plot_relcpue <- function(stock, stock_params, mp_params, app_params=NULL, show_last=TRUE, max_spaghetti_iters=50, quantiles, nspaghetti=5, add_grid=TRUE, ymax=NA, ...){

years <- as.numeric(dimnames(stock$biomass)$year)
cpue <- stock$catch / stock$effort
rel_cpue <- sweep(cpue, 1, cpue[,app_params$last_historical_timestep], "/")
if(is.na(ymax)){
ymax <- max(c(rel_cpue * 1.1, 1.0), na.rm=TRUE)
}
yrange <- c(0, ymax)

# Plot it
plot_indiv_timeseries_base(data=rel_cpue, stock=stock, stock_params=stock_params, mp_params=mp_params, app_params=app_params, show_last=show_last, max_spaghetti_iters=max_spaghetti_iters, quantiles=quantiles, nspaghetti=nspaghetti, yrange=yrange, ylab="Relative CPUE", add_grid=add_grid, ...)

# Add 1 line
lines(x=years,y=rep(1,length(years)), lty=2)

}


plot_releffort <- function(stock, stock_params, mp_params, app_params=NULL, show_last=TRUE, max_spaghetti_iters=50, quantiles, nspaghetti=5, add_grid=TRUE, ...){
years <- as.numeric(dimnames(stock$biomass)$year)
rel_effort <- sweep(stock$effort, 1, stock$effort[,app_params$last_historical_timestep], "/")
# Set Ylim - use same as HCR plot
ymax <- max(c(rel_effort * 1.1, 1.0), na.rm=TRUE)
ymax <- min(10, ymax, na.rm=TRUE)
yrange <- c(0, ymax)

# Plot it
plot_indiv_timeseries_base(data=rel_effort, stock=stock, stock_params=stock_params, mp_params=mp_params, app_params=app_params, show_last=show_last, max_spaghetti_iters=max_spaghetti_iters, quantiles=quantiles, nspaghetti=nspaghetti, yrange=yrange, ylab="Relative effort", add_grid=add_grid, ...)

# Add 1 line
lines(x=years,y=rep(1,length(years)), lty=2)

}



Expand All @@ -427,7 +499,7 @@ plot_catch <- function(stock, stock_params, mp_params, app_params=NULL, timestep
#' @rdname front_page_plots
#' @name Front page plots
#' @export
plot_metric_with_histo <- function(stock, stock_params, mp_params, metric, app_params=NULL, show_last=TRUE, percentile_range = c(20,80)){
plot_metric_with_histo <- function(stock, stock_params, mp_params, metric, app_params=NULL, show_last=TRUE, percentile_range = c(20,80), ...){
def.par <- par(no.readonly = TRUE) # as seen in layout doc
on.exit(par(def.par))
# Plot the metric with an extra sideways histogram
Expand All @@ -438,15 +510,15 @@ plot_metric_with_histo <- function(stock, stock_params, mp_params, metric, app_p
par(mar=c(pext, pext, bspc, bspc), oma=rep(ospc, 4)) # plot parameters
if (metric == "biomass"){
# The timeseries of biomass in the big window
plot_biomass(stock=stock, stock_params=stock_params, mp_params=mp_params, show_last=show_last, quantiles=percentile_range / 100)
plot_biomass(stock=stock, stock_params=stock_params, mp_params=mp_params, show_last=show_last, quantiles=percentile_range / 100, ...)
# The histogram should be the true B/K
final_yr <- dim(stock$biomass)[2]
dat <- stock$biomass[,final_yr] / stock_params$k
range <- c(0,1)
}
else if (metric == "catch"){
# The timeseries of biomass in the big window
plot_catch(stock=stock, stock_params=stock_params, mp_params=mp_params, show_last=show_last, quantiles=percentile_range / 100)
plot_catch(stock=stock, stock_params=stock_params, mp_params=mp_params, show_last=show_last, quantiles=percentile_range / 100, ...)
final_yr <- dim(stock$catch)[2]
dat <- stock$catch[,final_yr]
range <- c(0, get_catch_ymax(stock$catch, mp_params))
Expand All @@ -458,7 +530,7 @@ plot_metric_with_histo <- function(stock, stock_params, mp_params, metric, app_p
max_rel_cpue <- 2
rel_cpue[rel_cpue > max_rel_cpue] <- max_rel_cpue
ymax <- max(c(rel_cpue * 1.1, 1.0), na.rm=TRUE)
plot_relcpue(stock=stock, stock_params=stock_params, mp_params=mp_params, app_params=app_params, show_last=show_last, quantiles=percentile_range / 100, ymax=ymax)
plot_relcpue(stock=stock, stock_params=stock_params, mp_params=mp_params, app_params=app_params, show_last=show_last, quantiles=percentile_range / 100, ymax=ymax, ...)
# data for the histogram
final_yr <- dim(rel_cpue)[2]
dat <- rel_cpue[,final_yr]
Expand Down
10 changes: 5 additions & 5 deletions ComparingPerformance/app.R
Original file line number Diff line number Diff line change
Expand Up @@ -154,7 +154,7 @@ ui <- navbarPage(
# Total number of years (including historical years)
numericInput("nyears", "Number of years", value = 30, min=20, max=100, step=1),
# Number of iteration
numericInput("niters", "Number of iterations", value = 1000, min=10, max=1000, step=10),
numericInput("niters", "Number of iterations", value = 500, min=10, max=1000, step=10),
br(),
actionButton("dump", "Dump data")
)
Expand Down Expand Up @@ -366,19 +366,19 @@ server <- function(input, output,session) {

# Call the HCR plot function
output$plothcr <- renderPlot({
plot_hcr(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params, show_last = FALSE)
plot_hcr(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params, show_last = FALSE, cex.axis=1.1, cex.lab=1.3)
})

output$plotbiomasshisto <- renderPlot({
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="biomass", show_last = FALSE, percentile_range = pi_percentiles)
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="biomass", show_last = FALSE, percentile_range = pi_percentiles, cex.axis=1.1, cex.lab=1.3)
})

output$plotcatchhisto <- renderPlot({
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="catch", show_last = FALSE, percentile_range = pi_percentiles)
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="catch", show_last = FALSE, percentile_range = pi_percentiles, cex.axis=1.1, cex.lab=1.3)
})

output$plotrelcpuehisto <- renderPlot({
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="relcpue", app_params=app_params, show_last = FALSE, percentile_range = pi_percentiles)
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="relcpue", app_params=app_params, show_last = FALSE, percentile_range = pi_percentiles, cex.axis=1.1, cex.lab=1.3)
})


Expand Down
22 changes: 14 additions & 8 deletions IntroHCR/app.R
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ library(AMPLE)

#------------------------------------------------------------
# Plot stuff
plot_hcr_intro_arrow <- function(stock, timestep){
plot_hcr_intro_arrow <- function(stock, timestep, guide_lwd=3, guide_col="blue"){
# Arrow from B/K to HCR
btimestep <- min(timestep, dim(stock$hcr_ip)[2])
last_hcr_ip <- stock$hcr_ip[1, btimestep]
Expand Down Expand Up @@ -112,10 +112,16 @@ ui <- navbarPage(
p("The ghosts of catch limits from the past are shown as grey dashed lines on the catch plot and as grey dots on the HCR plot. These allow you to see which parts of the HCR have been active."),
h2("Variability"),
p("Variability can be included in the projection in two ways: through variability in the stock productivity and through the estimated level of stock biomass being different to the true level of the stock biomass. These options are initially turned off. The options can be seen by clicking on the ", strong("Show variability options"), "box."),
p("Biological productivity variability represents the variability of the natural procesess of the stock, for example growth and natural mortality. Increasing the variability will increase the 'bumpiness' of the stock trajectory. As biological variability is always encountered in fisheries it is essential that a selected HCR is robust to the variability."),
p("Estimation error simulates the difference between the true level of the stock biomass and the estimated level. Unfortunately, the true abundance of a fish stock is never known. Instead, estimates of abundance are made, for example using stock assessment models. The HCR uses the estimated biomass, not the true biomass. This means that the catch limit that is set by the HCR is based on estimated biomass. If the biomass is estimated poorly the resulting catch limit set by the HCR may not be appropriate."),
p("Here, estimation error is modelled using two different processes: random error and consistent bias (positive or negative). The bias represents situations where the biomass is consistently over or under estimated."),
p("When estimation error is active the biomass plot shows two lines. The black line shows the true biomass, the blue line shwows the estimated biomass. It is the blue line that feeds the HCR. Increasing the estiamation bias and variability will increase the difference between these lines."),
p("Biological variability is the natural variability in the biological processes of the stock (e.g. recruitment and growth). Increasing the variability will increase the 'bumpiness' of the stock trajectory. As biological variability is always encountered in fisheries it is important that a selected HCR is robust to the variability."),
p("The true stock status (e.g. abundance) of a fish stock is never known.
Instead, the stock status must be estimated. For example, stock assessment models can be used to estimate stock abundance.
HCRs use the estimated value of stock status, not the true value.
In this app the HCR uses an estimated value of biomass (SB/SBF=0) to set the catches.
This means that if the biomass is estimated poorly, the resulting catch limit set by the HCR may not be appropriate."),
p("Here, estimation variability and estimation bias can be used to simulate the estimation process.
They can be adjusted to change the difference between the true level of the biomass and the estimated level. The bias represents situations where the biomass is consistently over or under estimated."),
p("When estimation variability or bias is active the biomass plot shows two lines. The black line shows the true biomass, the blue line shows the estimated biomass. It is the blue line that feeds the HCR. Increasing the estimation bias and variability will increase the difference between these lines.
If you change the values for variability and bias you will see changes in the performance of the HCR. This is why HCRs are designed to be used with a specific estimation process."),
h1("Tutorial"),
p("A more detailed tutorial can be found at these links:"),
a("Tutorial (pdf)",target="_blank",href= "introHCR.pdf"),
Expand Down Expand Up @@ -203,15 +209,15 @@ server <- function(input, output,session) {
})

output$plotcatch <- renderPlot({
plot_catch(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params, timestep=timestep(), main="Catch", add_grid=TRUE)
plot_catch(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params, timestep=timestep(), main="Catch", add_grid=TRUE, cex.axis=1.1, cex.lab=1.3)
})

output$plothcr <- renderPlot({
plot_hcr(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params, timestep=timestep()+1)
plot_hcr(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params, timestep=timestep()+1, cex.axis=1.1, cex.lab=1.3)
})

output$plotbiomass <- renderPlot({
plot_biomass(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), timestep=timestep()+1, main="SB / SBF=0", add_grid=TRUE)
plot_biomass(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), timestep=timestep()+1, main="SB / SBF=0", add_grid=TRUE, cex.axis=1.1, cex.lab=1.3)
})

output$plotarrow <- renderPlot({
Expand Down
8 changes: 4 additions & 4 deletions IntroIndicators/app.R
Original file line number Diff line number Diff line change
Expand Up @@ -216,19 +216,19 @@ server <- function(input, output,session) {

# Call the HCR plot function
output$plothcr <- renderPlot({
plot_hcr(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params)
plot_hcr(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), app_params=app_params, cex.axis=1.1, cex.lab=1.3)
})

output$plotbiomasshisto <- renderPlot({
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="biomass", percentile_range = pi_percentiles)
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="biomass", percentile_range = pi_percentiles, cex.axis=1.1, cex.lab=1.3)
})

output$plotcatchhisto <- renderPlot({
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="catch", percentile_range = pi_percentiles)
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="catch", percentile_range = pi_percentiles, cex.axis=1.1, cex.lab=1.3)
})

output$plotrelcpuehisto <- renderPlot({
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="relcpue", app_params=app_params, percentile_range = pi_percentiles)
plot_metric_with_histo(stock=stock, stock_params=get_stock_params(), mp_params=get_mp_params(), metric="relcpue", app_params=app_params, percentile_range = pi_percentiles, cex.axis=1.1, cex.lab=1.3)
})

output$itercount <- renderText({
Expand Down

0 comments on commit 075dc0a

Please sign in to comment.