Newer
Older
##############################################
## Hide/Show : level 1
## Fatalities
output$hide_fatalities <- eventReactive({
input$button_fatalities
},{
if(input$button_fatalities%%2 == 1) TRUE else FALSE
}, ignoreInit = TRUE)
outputOptions(output, "hide_fatalities", suspendWhenHidden = FALSE)
## Population Size
output$hide_pop_size <- eventReactive({
input$button_pop_size
},{
if(input$button_pop_size%%2 == 1) TRUE else FALSE
}, ignoreInit = TRUE)
outputOptions(output, "hide_pop_size", suspendWhenHidden = FALSE)
## Population Growth
output$hide_pop_growth <- eventReactive({
input$button_pop_growth
},{
if(input$button_pop_growth%%2 == 1) TRUE else FALSE
}, ignoreInit = TRUE)
outputOptions(output, "hide_pop_growth", suspendWhenHidden = FALSE)
## Carrying capacity
output$hide_carrying_cap <- eventReactive({
input$button_carrying_cap
},{
if(input$button_carrying_cap%%2 == 1) TRUE else FALSE
}, ignoreInit = TRUE)
outputOptions(output, "hide_carrying_cap", suspendWhenHidden = FALSE)
# Display Carrying capacity Unit Info
output$carrying_cap_unit_info <- renderText({
if(input$pop_size_unit == "Npair"){
paste0("Nombre de couple")
} else {
paste0("Effectif total")
}
})
##############################################
## Hide/Show : level 2
#------------
# Hide all
#------------
shinyjs::hide("fatalities_mean")
shinyjs::hide("fatalities_se")
thierrychambert
committed
shinyjs::hide("fatalities_lower")
shinyjs::hide("fatalities_upper")
shinyjs::hide("farm_number_cumulated")
shinyjs::hide("fatalities_mat_cumulated")
shinyjs::hide("pop_size_lower")
shinyjs::hide("pop_size_upper")
shinyjs::hide("pop_size_mean")
shinyjs::hide("pop_size_se")
shinyjs::hide("pop_size_mat_expert")
shinyjs::hide("pop_growth_lower")
shinyjs::hide("pop_growth_upper")
shinyjs::hide("pop_growth_mean")
shinyjs::hide("pop_growth_se")
shinyjs::hide("pop_growth_mat_expert")
shinyjs::hide("pop_trend")
shinyjs::hide("pop_trend_strength")
shinyjs::hide("carrying_capacity")
shinyjs::hide("carrying_cap_mat_expert")
shinyjs::hide("carrying_cap_run_expert")
#------------
# Show some
#------------
# Show inputs for fatalities part
# Show inputs for none cumulated impacts scenario
if(input$analysis_choice == "scenario"){
shinyjs::show("fatalities_input_type")
if(input$fatalities_input_type == "itvl"){
thierrychambert
committed
shinyjs::show("fatalities_lower")
shinyjs::show("fatalities_upper")
if(input$fatalities_input_type == "val"){
shinyjs::show("fatalities_mean")
shinyjs::show("fatalities_se")
}
}
}
# Show inputs for cumulated scenario
if(input$analysis_choice == "cumulated"){
shinyjs::hide("fatalities_input_type")
shinyjs::show("farm_number_cumulated")
shinyjs::show("fatalities_mat_cumulated")
# Show inputs for population size part
if(input$button_pop_size%%2 == 1){
shinyjs::show("pop_size_input_type")
if(input$pop_size_input_type == "itvl"){
shinyjs::show("pop_size_lower")
shinyjs::show("pop_size_upper")
}
shinyjs::show("pop_size_mean")
shinyjs::show("pop_size_se")
}
# Show inputs for population trend/growth part
if(input$button_pop_growth%%2 == 1){
shinyjs::show("pop_growth_input_type")
if(input$pop_growth_input_type == "itvl"){
shinyjs::show("pop_growth_lower")
shinyjs::show("pop_growth_upper")
}
shinyjs::show("pop_growth_mean")
shinyjs::show("pop_growth_se")
}
if(input$pop_trend != "stable"){
shinyjs::show("pop_trend_strength")
}
# Show inputs for carrying capacity part
if(input$button_carrying_cap%%2 == 1){
shinyjs::show("carrying_cap_input_type")
if(input$carrying_cap_input_type == "val"){
shinyjs::show("carrying_capacity")
}
if(input$carrying_cap_input_type == "eli_exp"){
shinyjs::show("carrying_cap_mat_expert")
shinyjs::show("carrying_cap_run_expert")
}
}
# Show inputs vital rates part
}) # en observe show/hide
###~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~###
##############################################
## Reactive values
thierrychambert
committed
out <- reactiveValues(run = NULL, msg = NULL)
ready <- reactiveValues(fatalities = TRUE, pop_size = TRUE, pop_growth = TRUE, carrying_capacity = TRUE)
nsim = NULL,
cumulated_impacts = FALSE,
fatalities_mean_nb = NULL,
fatalities_se = NULL,
fatalities_se_nb = NULL,
onset_time = NULL,
onset_year = NULL,
pop_size_mean = NULL,
pop_size_se = NULL,
pop_growth_mean = NULL,
pop_growth_se = NULL,
fecundities = NULL,
survivals = NULL,
s_calibrated = NULL,
f_calibrated = NULL,
vr_calibrated = NULL,
theta = NULL,
model_demo = NULL,
time_horzion = NULL,
coeff_var_environ = NULL,
fatal_constant = NULL,
fatalities_eli_result = NULL,
pop_size_eli_result = NULL,
pop_growth_eli_result = NULL,
carrying_cap_eli_result = NULL
###~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~###
##############################################
## Update matrix cumulated impact
##-------------------------------------------
observeEvent({
input$farm_number_cumulated
}, {
nfarm <- input$farm_number_cumulated
init_cumul_new <- init_cumul[1:nfarm,]
updateMatrixInput(session, inputId = "fatalities_mat_cumulated",
value = matrix(init_cumul_new, nrow = nfarm, ncol = 3, byrow = FALSE,
dimnames = list(paste("Parc", c(1:nfarm)),
c("Moyenne",
"Erreur-type",
#####
##--------------------------------------------
##--------------------------------------------
# Function to run the elication analysis
func_eli <- function(mat_expert){
t_mat_expert <- t(mat_expert)
vals <- t_mat_expert[2:4,]
Cp <- t_mat_expert[5,]
weights <- t_mat_expert[1,]
out <- elicitation(vals, Cp, weights)
return(list(out = out, mean = out$mean_smooth, SE = sqrt(out$var_smooth)))
}
# Function to plot the elication analysis output
plot_expert <- function(out, show_se = TRUE, ...){
plot_elicitation(out, ylab = "", xlab = "Valeur du paramtre", cex.lab = 1.2, yaxt = "n")
mtext(text = "Densit de probabilit", side = 2, line = 2, cex = 1.2)
y2 <- dgamma(x = out$mean_smooth, shape = out$shape_smooth, rate = out$rate_smooth)
xx <- qgamma(p = c(0.01,0.99), shape = out$shape_smooth, rate = out$rate_smooth)
clip(xx[1], xx[2], -100, y2)
abline(v = out$mean_smooth, lwd = 3, col = "darkblue")
mtext(text = paste("Moyenne = ", round(out$mean_smooth,2)), side = 3, line = 2.5, cex = 1.2, adj = 0)
if(show_se) mtext(text = paste("Erreur-type = ", round(sqrt(out$var_smooth), 2)), side = 3, line = 1, cex = 1.2, adj = 0)
}
## Fatalities
##----------------------
observeEvent({
input$fatalities_run_expert
}, {
if( all(!is.na(input$fatalities_mat_expert)) ) {
param$fatalities_eli_result <- func_eli(input$fatalities_mat_expert)
output$title_distri_plot <- renderText({ "Mortalits annuelles" })
output$distri_plot <- renderPlot({ plot_expert(param$fatalities_eli_result$out) })
} else {
print("missing value")
} # end if
}) # end observeEvent
## Population size
##----------------------
observeEvent({
input$pop_size_run_expert
}, {
if(all(!is.na(input$pop_size_mat_expert))) {
## run elicitation analysis
param$pop_size_eli_result <- func_eli(input$pop_size_mat_expert)
output$title_distri_plot <- renderText({ "Taille de population" })
output$distri_plot <- renderPlot({ plot_expert(param$pop_size_eli_result$out) })
} else {
print("missing value")
} # end if
}) # end observeEvent
## Population growth
##----------------------
observeEvent({
input$pop_growth_run_expert
},{
if(all(!is.na(input$pop_growth_mat_expert))) {
## run elicitation analysis
param$pop_growth_eli_result <- func_eli(input$pop_growth_mat_expert)
output$title_distri_plot <- renderText({ "Taux de croissance de la population" })
output$distri_plot <- renderPlot({ plot_expert(param$pop_growth_eli_result$out) })
} else {
print("missing value")
} # end if
}) # end observeEvent
## Carrying capacity
##----------------------
observeEvent({
input$carrying_cap_run_expert
},{
if(all(!is.na(input$carrying_cap_mat_expert))) {
## run elicitation analysis
param$carrying_cap_eli_result <- func_eli(input$carrying_cap_mat_expert)
output$title_distri_plot <- renderText({ "Capacit de charge" })
output$distri_plot <- renderPlot({ plot_expert(param$carrying_cap_eli_result$out, show_se = FALSE) })
} else {
print("missing value")
} # end if
}) # end observeEvent
###~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~###
#####
#####
##--------------------------------------------
## Display parameter distribution
##--------------------------------------------
# Function to plot a gamma distribution
plot_gamma <- function(mu, se, show_mode = TRUE, show_mean = TRUE, show_se = TRUE, ...){
## Define shape and scale parameter of gamma distribution
shape = (mu/se)^2
scale = se^2/mu
par(mar = c(5, 4, 6, 2))
curve(dgamma(x, shape=shape, scale=scale), from = max(0,mu-3*se), to = mu+4*se, lwd = 3, col = "darkblue", yaxt = "n",
ylab = "", xlab = "Valeur du paramtre", cex.lab = 1.2)
mtext(text = "Densit de probabilit", side = 2, line = 2, cex = 1.2)
# show mode
MU <- (shape-1)*scale
y_MU <- dgamma(x = MU, shape = shape, scale = scale)
xx <- qgamma(p = c(0.01,0.99), shape = shape, scale = scale)
clip(xx[1], xx[2], -100, y_MU)
abline(v = MU, lwd = 3, col = "darkblue")
# show mean
y_mu <- dgamma(x = mu, shape = shape, scale = scale)
clip(xx[1], xx[2], -100, y_mu)
abline(v = mu, lwd = 2, col = "darkblue", lty = 2)
if(show_mode) mtext(text = paste("Mode = ", round(MU, 2)), side = 3, line = 4, cex = 1.2, adj = 0)
if(show_mean) mtext(text = paste("Moyenne = ", round(mu, 2)), side = 3, line = 2.5, cex = 1.2, adj = 0)
if(show_se) mtext(text = paste("Erreur-type = ", round(se, 2)), side = 3, line = 1, cex = 1.2, adj = 0)
thierrychambert
committed
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
} # end function plot_gamma
plot_gamma_cumulated_impacts <- function(mu, se, nparc, ...){
## Define shape and scale parameter of gamma distribution
shape = (mu/se)^2
scale = se^2/mu
## Define x and y lim
xx = yy = list()
for(j in 1:nparc){
xx[[j]] = seq(from = max(0,mu[j]-4*se[j]), to = mu[j]+4*se[j], length.out = 1e3)
yy[[j]] = dgamma(xx[[j]], shape=shape[j], scale=scale[j])
}
ylim = c(min(unlist(yy)), max(unlist(yy))*1.4)
xlim = c(min(unlist(xx)), max(unlist(xx)))
## Plot
j=1
curve(dgamma(x, shape=shape[j], scale=scale[j]),
from = max(0,mu[j]-4*se[j]), to = mu[j]+4*se[j], n = 1e4,
xlim = xlim, ylim = ylim,
lwd = 3, col = j, yaxt = "n", xaxt = "n",
#xaxp = c(round(xlim[1]), round(xlim[2]), n = 10),
ylab = "", xlab = "Valeur du paramtre", cex.lab = 1.2)
axis(side = 1, at = seq(round(xlim[1]), round(xlim[2]),
by = max(round((round(xlim[2])-round(xlim[1]))/10),1) ))
mtext(text = "Densit de probabilit", side = 2, line = 2, cex = 1.2)
y1 <- dgamma(x = mu[j], shape = shape[j], scale = scale[j])
segments(x0 = mu[j], y0 = 0, y1 = y1, lty = 2, lwd = 3, col = j)
points(x = mu[j], y = y1, pch = 19, cex = 1.5, col = j)
for(j in 2:nparc){
curve(dgamma(x, shape=shape[j], scale=scale[j]),
from = max(0,mu[j]-4*se[j]), to = mu[j]+4*se[j], n = 1e4,
lwd = 3, col = j, yaxt = "n",
ylab = "", xlab = "Valeur du paramtre", cex.lab = 1.2, add = TRUE)
y1 <- dgamma(x = mu[j], shape = shape[j], scale = scale[j])
segments(x0 = mu[j], y0 = 0, y1 = y1, lty = 2, lwd = 3, col = j)
points(x = mu[j], y = y1, pch = 19, cex = 1.5, col = j)
}
legend(x = xlim[1], y = ylim[2], legend = paste("Parc", 1:nparc),
lwd = 3, col = 1:nparc, text.col = 1:nparc, cex = 1.5,
bty = "n", horiz = TRUE)
} # end function plot_gamma_cumulated_impacts
########################
## Fatalities
##----------------------
observeEvent({
input$analysis_choice
input$button_fatalities
input$fatalities_input_type
input$fatalities_run_expert
input$farm_number_cumulated
input$fatalities_mat_cumulated
},{
if(input$analysis_choice != "cumulated"){
# Show from input values: if button is ON and input_type is set on "value" or "itvl" (thus not "eli_exp")
if(input$button_fatalities%%2 == 1 & input$fatalities_input_type != "eli_exp"){
output$title_distri_plot <- renderText({ "Mortalits annuelles" })
output$distri_plot <- renderPlot({
if(input$fatalities_input_type == "itvl"){
req(input$fatalities_lower, input$fatalities_upper)
plot_gamma(mu = tail(param$fatalities_mean, -1), se = tail(param$fatalities_se, -1))
}else{
req(input$fatalities_mean, input$fatalities_se)
plot_gamma(mu = tail(param$fatalities_mean, -1), se = tail(param$fatalities_se, -1))
}
})
} else {
# Show from elicitation expert: if button is ON and input_type is set on "expert elicitation"
if(input$button_fatalities%%2 == 1 & input$fatalities_input_type == "eli_exp"){
if(!is.null(param$fatalities_eli_result)){
output$title_distri_plot <- renderText({ "Mortalits annuelles" })
output$distri_plot <- renderPlot({ plot_expert(param$fatalities_eli_result$out) })
} else {
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
# Hide otherwise (when button is OFF)
}else{
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
}
thierrychambert
committed
# When analysis = cumulated impacts
}else{
output$title_distri_plot <- renderText({ "Mortalits annuelles par parc (impacts cumuls)" })
# Plot: note we use the "NULL + delay" sequence only to avoid error message in R console
output$distri_plot <- NULL
delay(5,
output$distri_plot <- renderPlot({
plot_gamma_cumulated_impacts(mu = input$fatalities_mat_cumulated[,1],
thierrychambert
committed
se = input$fatalities_mat_cumulated[,2],
nparc = input$farm_number_cumulated)
thierrychambert
committed
} # end "if"
########################
## Population size
##----------------------
observeEvent({
input$button_pop_size
},{
# Show from input values: if button is ON and input_type is set on "value"
if(input$button_pop_size%%2 == 1 & input$pop_size_input_type != "eli_exp"){
output$title_distri_plot <- renderText({ "Taille initiale de la population" })
output$distri_plot <- renderPlot({
req(param$pop_size_mean, param$pop_size_se)
plot_gamma(mu = param$pop_size_mean, se = param$pop_size_se)
})
} else {
# Show from elicitation expert: if button is ON and input_type is set on "expert elicitation"
if(input$button_pop_size%%2 == 1 & input$pop_size_input_type == "eli_exp"){
if(!is.null(param$pop_size_eli_result)){
output$title_distri_plot <- renderText({ "Taille initiale de la population" })
output$distri_plot <- renderPlot({ plot_expert(param$pop_size_eli_result$out) })
} else {
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
# Hide otherwise (when button is OFF)
}else{
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
########################
## Population growth
##----------------------
observeEvent({
input$pop_growth_input_type
input$button_pop_growth
},{
# Show from input values: if button is ON and input_type is set on "value" or "interval"
if(input$button_pop_growth%%2 == 1 & input$pop_growth_input_type != "eli_exp" & input$pop_growth_input_type != "trend"){
output$title_distri_plot <- renderText({ "Taux de croissance de la population" })
output$distri_plot <- renderPlot({
req(param$pop_growth_mean, param$pop_growth_se > 0)
plot_gamma(mu = param$pop_growth_mean, se = param$pop_growth_se)
})
} else {
# Show from elicitation expert: if button is ON and input_type is set on "expert elicitation"
if(input$button_pop_growth%%2 == 1 & input$pop_growth_input_type == "eli_exp"){
if(!is.null(param$pop_growth_eli_result)){
output$title_distri_plot <- renderText({ "Taux de croissance de la population" })
output$distri_plot <- renderPlot({ plot_expert(param$pop_growth_eli_result$out) })
} else {
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
# Hide otherwise (when button is OFF)
}else{
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
}
}, ignoreInit = FALSE)
########################
## Carrying capacity
##----------------------
observeEvent({
input$carrying_cap_input_type
input$button_carrying_cap
},{
# Show from elicitation expert: if button is ON and input_type is set on "expert elicitation"
if(input$button_carrying_cap%%2 == 1 & input$carrying_cap_input_type == "eli_exp"){
if(!is.null(param$carrying_cap_eli_result)){
output$title_distri_plot <- renderText({ "Capacit de charge" })
output$distri_plot <- renderPlot({ plot_expert(param$carrying_cap_eli_result$out) })
} else {
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
# Hide otherwise (when button is OFF)
}else{
output$title_distri_plot <- NULL
output$distri_plot <- NULL
}
}, ignoreInit = FALSE)
#####
#####
##-------------------------------------------------
## Display parameter values (on the side panel)
##-------------------------------------------------
#################################
##-------------------------------
## UNIT
output$fatalities_unit_info <- renderText({
if(!is.null(input$fatalities_unit)){
if(input$fatalities_unit == "h"){
paste0("Taux de mortalit")
} else {
paste0("Nombre de mortalits")
}
}
})
## Values
output$fatalities_mean_info <- renderText({
if(input$fatalities_unit == "h") add_perc <- "%" else add_perc <- ""
paste0(c("Moyenne : ",
paste0(tail(param$fatalities_mean, -1), add_perc, collapse = ", ")
), collapse = "")
})
output$fatalities_se_info <- renderText({
if(input$fatalities_unit == "h") add_perc <- "%" else add_perc <- ""
paste0(c("Erreur-type : ",
paste0(tail(param$fatalities_se, -1), add_perc, collapse = ", ")
), collapse = "")
})
thierrychambert
committed
#################################
## Poplutation size
##-------------------------------
## UNIT
thierrychambert
committed
output$pop_size_unit_info <- renderText({
if(!is.null(param$pop_size_unit)){
if(param$pop_size_unit == "Npair"){
paste0("Nombre de couple")
thierrychambert
committed
paste0("Effectif total")
thierrychambert
committed
output$pop_size_mean_info <- renderText({ paste0("Moyenne : ", param$pop_size_mean) })
output$pop_size_se_info <- renderText({ paste0("Erreur-type : ", param$pop_size_se) })
## Show Popsize by age (table)
# Function to create the table
make_mat_popsizes <- function(data_sf, species, pop_size, pop_size_unit, survivals, fecundities){
nam <- data_sf %>%
filter(Nom_espece == species) %>%
select(classes_age) %>%
unlist %>%
as.vector
matrix(round(pop_vector(pop_size = pop_size, pop_size_type = pop_size_unit, s = survivals, f = fecundities)),
nrow = 1,
dimnames = list("Effectifs", nam)
)
}
# Display the table (Note the delay : piece is just there to avoid an error message - time for parameters to be "loaded in")
output$pop_size_by_age <- renderTable({
if(any(is.na(param$survivals)) | any(is.na(param$fecundities))){
matrix("Valeurs de survies et/ ou de fcondits manquantes",
nrow = 1, dimnames = list(NULL, "Erreur"))
}else{
make_mat_popsizes(data_sf = data_sf, species = input$species_choice, pop_size = param$pop_size_mean,
pop_size_unit = input$pop_size_unit, s = param$survivals, f = param$fecundities)
} # end if
},
width = "500px",
rownames = FALSE,
digits = 0)
)
#################################
## Population growth
##-------------------------------
thierrychambert
committed
output$pop_growth_mean_info <- renderText({ paste0("Moyenne : ", param$pop_growth_mean) })
output$pop_growth_se_info <- renderText({ paste0("Erreur-type : ", param$pop_growth_se) })
#################################
thierrychambert
committed
## Carrying capacity
##-------------------------------
# UNIT (like pop size)
thierrychambert
committed
# Source info "unit"
if(is.null(param$pop_size_unit)){
unit1 <- input$pop_size_unit
}else{
unit1 <- param$pop_size_unit
}
thierrychambert
committed
thierrychambert
committed
if(unit1 == "Npair"){
info1 <- paste0("Nombre de couple")
thierrychambert
committed
info1 <- paste0("Effectif total")
thierrychambert
committed
paste0(info1, " : ", param$carrying_capacity)
#################################
##-------------------------------
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
# Function to create the matrix
make_mat_vr <- function(data_sf, species){
out_mat <- data_sf %>%
filter(Nom_espece == species) %>%
select(classes_age, survie, fecondite)
return(out_mat)
}
# Update the vital rate matrix (mat_fill_vr) when changing species in the list
observeEvent({
input$species_choice
}, {
if(input$species_choice == "Espce gnrique") {} else {
tab_species <- make_mat_vr(data_sf = data_sf, species = input$species_choice)
if(all(is.na(tab_species))) {
updateMatrixInput(session, inputId = "mat_fill_vr",
value = matrix(data = NA,
nrow = 4,
ncol = 2,
dimnames = list(c("Juv 0", "Sub 1", "Sub 2", "Adulte"), c("Survie", "Fcondit"))))
} else {
number_age_class <- nrow(tab_species)
ages <- tab_species$classes_age
survivals <- tab_species$survie
fecundities <- tab_species$fecondite
updateMatrixInput(session, inputId = "mat_fill_vr",
value = matrix(data = c(survivals, fecundities),
nrow = number_age_class,
ncol = 2,
dimnames = list(ages, c("Survie", "Fcondit"))))
} # end if 2
} # end if 1
}) # end observeEvent species_list
# Display vital rates output table
# Display intrinsic lambda (based solely on Leslie matrix)
delay(ms = 300,
output$lambda0_info <- renderUI({
lam <- lambda(build_Leslie(s = input$mat_fill_vr[,1], f = input$mat_fill_vr[,2]))
withMathJax(sprintf("$$\\lambda = %.02f$$", lam))
})
)
#####
##--------------------------------------------
## Select parameter values for simulations
##--------------------------------------------
thierrychambert
committed
# Functions to calculate mean and SD from lower & upper values
get_mu <- function(lower, upper) (lower + upper)/2
get_sd <- function(lower, upper, coverage) ((abs(upper - lower)/2))/qnorm(1-((1-coverage)/2))
#################################
## Cumulated impacts or not ?
##-------------------------------
observeEvent({
input$run
}, {
if(input$analysis_choice == "scenario"){
param$cumulated_impacts = FALSE
} else {
param$cumulated_impacts = TRUE
} # end if
}) # end observeEvent
#################################
## Fatalities
##-------------------------------
observe({
# Case 1 : Not cumulated effects (if1)
if(input$analysis_choice == "scenario"){
# Case 1.1 : Values from expert elicitation (if2)
if(input$fatalities_input_type == "eli_exp"){
if(!(is.null(param$fatalities_eli_result))){
param$fatalities_mean <- c(0, round(param$fatalities_eli_result$mean, 2))
param$onset_time <- NULL
thierrychambert
committed
param$fatalities_se <- c(0, round(param$fatalities_eli_result$SE, 3))
ready$fatalities <- TRUE
} else {
ready$fatalities <- FALSE
}
} else {
if(input$fatalities_input_type == "val"){
# Case 1.2 : Values directly provided as mean & SE
param$fatalities_mean <- c(0, input$fatalities_mean)
param$onset_time <- NULL
param$fatalities_se <- c(0, input$fatalities_se)
ready$fatalities <- TRUE
}else{
# Case 1.3 : Values directly provided as lower/upper interval
thierrychambert
committed
param$fatalities_mean <- c(0, round(get_mu(lower = input$fatalities_lower, upper = input$fatalities_upper), 2))
thierrychambert
committed
param$fatalities_se <- c(0, round(get_sd(lower = input$fatalities_lower, upper = input$fatalities_upper, coverage = CP), 3))
} # end (if2)
# Case 2 : Cumulated effects (if-else 1)
} else {
ready$fatalities <- TRUE
param$fatalities_mean <- c(0, input$fatalities_mat_cumulated[,1])
param$fatalities_se <- c(0, input$fatalities_mat_cumulated[,2])
param$onset_year <- c(min(input$fatalities_mat_cumulated[,3]), input$fatalities_mat_cumulated[,3])
param$onset_time <- param$onset_year - min(param$onset_year) + 1
} # end (if1)
}) # end observe
###~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~###
# Make sure fatalities are expressed as "number" (not rate) for the run_simul function
se_prod2 <- function(mu1, se1, mu2, se2) sqrt((se1^2 * se2^2) + (se1^2 * mu2^2) + (mu1^2 * se2^2))
observeEvent({
input$run
},{
if(input$fatalities_unit == "h"){
pop_size_tot <- sum(pop_vector(pop_size = param$pop_size_mean, pop_size_type = param$pop_size_type, s = param$survivals, f = param$fecundities)[-1])
param$fatalities_mean_nb <- (param$fatalities_mean/100) * pop_size_tot
param$fatalities_se_nb <- se_prod2(mu1 = param$fatalities_mean/100,
se1 = param$fatalities_se/100,
mu2 = pop_size_tot,
se2 = (pop_size_tot/param$pop_size_mean) * param$pop_size_se)
}else{
param$fatalities_mean_nb <- param$fatalities_mean
param$fatalities_se_nb <- param$fatalities_se
}
})
#################################
## Population size
##-------------------------------
observe({
# Case 1 : Values from expert elicitation
if(input$pop_size_input_type == "eli_exp"){
if(!(is.null(param$pop_size_eli_result))){
param$pop_size_mean <- round(param$pop_size_eli_result$mean)
param$pop_size_se <- round(param$pop_size_eli_result$SE)
ready$pop_size <- TRUE
} else {
ready$pop_size <- FALSE
}
} else {
if(input$pop_size_input_type == "val"){
# Case 2 : Values directly provided as mean & SE
ready$pop_size <- TRUE
param$pop_size_mean <- input$pop_size_mean
param$pop_size_se <- input$pop_size_se
}else{
# Case 3 : Values directly provided as lower/upper interval
ready$pop_size <- TRUE
param$pop_size_mean <- round(get_mu(lower = input$pop_size_lower, upper = input$pop_size_upper), 2)
param$pop_size_se <- round(get_sd(lower = input$pop_size_lower, upper = input$pop_size_upper, coverage = CP), 3)
} # end (if3)
}
param$pop_size_unit <- input$pop_size_unit
})
#################################
## Population growth
##-------------------------------
observe({
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
# Case 1 : Values from expert elicitation
if(input$pop_growth_input_type == "eli_exp"){
if(!(is.null(param$pop_growth_eli_result))){
param$pop_growth_mean <- round(min(1 + param$rMAX_species, round(param$pop_growth_eli_result$mean, 2)), 2)
param$pop_growth_se <- round(param$pop_growth_eli_result$SE, 2)
ready$pop_growth <- TRUE
} else {
ready$pop_growth <- FALSE
}
} else {
# Case 2 : Trend information
if(input$pop_growth_input_type == "trend"){
ready$pop_growth <- TRUE
if(input$pop_trend == "growth") {
if(input$pop_trend_strength == "weak") {
param$pop_growth_mean <- 1.01
} else if(input$pop_trend_strength == "average"){
param$pop_growth_mean <- 1.03
} else {
param$pop_growth_mean <- 1.06
}
} else if(input$pop_trend == "decline"){
if(input$pop_trend_strength == "weak") {
param$pop_growth_mean <- 0.99
} else if(input$pop_trend_strength == "average"){
param$pop_growth_mean <- 0.97
} else {
param$pop_growth_mean <- 0.94
}
} else {
param$pop_growth_mean <- 1
}
param$pop_growth_se <- 0
# Case 3 : Values directly provided (i.e., not from expert elicitation)
} else {
if(input$pop_growth_input_type == "val"){
# Case 2 : Values directly provided as mean & SE
ready$pop_growth <- TRUE
param$pop_growth_mean <- round(min(1 + param$rMAX_species, input$pop_growth_mean), 3)
param$pop_growth_se <- input$pop_growth_se
}else{
# Case 3 : Values directly provided as lower/upper interval
ready$pop_growth <- TRUE
param$pop_growth_mean <- round(min(1 + param$rMAX_species,
round(get_mu(lower = input$pop_growth_lower, upper = input$pop_growth_upper), 2)
), 3)
param$pop_growth_se <- round(get_sd(lower = input$pop_growth_lower, upper = input$pop_growth_upper, coverage = CP), 3)
} # end (if3)
}
}
})
#################################
## Carrying capacity
##------------------------------
observe({
if(input$carrying_cap_input_type == "eli_exp"){
if(!(is.null(param$carrying_cap_eli_result))){
param$carrying_capacity <- round(param$carrying_cap_eli_result$mean)
ready$carrying_capacity <- TRUE
} else {
ready$carrying_capacity <- FALSE
}
} else {
ready$carrying_capacity <- TRUE
param$carrying_capacity <- input$carrying_capacity
}