R/SimOP.R
In paulhegedus/OFPE: On-Farm Precision Experiments (OFPE) Data Management and Analysis Tools
#' @title R6 Class for plotting outputs of the OFPE simulation.
#'
#' @description R6 class for plotting outputs of the OFPE simulation such as
#' figures comparing management scenarios, probability tables, and maps.
#'
#' This class should be initialized and executed after the user has initialized
#' and executed all other classes and has ran the SimClass$executeSim() method.
#' These figures rely on outputs from the simulation, so will cause errors if
#' used prior to execution of the simulation.
#'
#' See the SimClass$executeSim() method documentation for more detail on the management
#' scenarios reported. NR is net-return. The abbreviations used are; NRmin: uniform minimum
#' rate ('AAmin'), NRfs: uniform farmer selected rate, NRact: if experiment applied in
#' simulated year, NRssopt: site-specific optimum rates, NRffopt: uniform full-
#' field optimum rate, NRopp: if the crop received the opposite system crop price (e.g.
#' 'conventional' vs. 'organic' or reverse, see SimClass$executeSim() method
#' documentation for specifics on the difference in the calculation of this net-return
#' between systems).
#'
#' Error checking inputs to these functions is scant because the user should only use
#' these after executing the simulation method in the 'SimClass' R6 class to generate
#' the necessary data for the output generation methods.
#' @seealso \code{\link{SimClass}} for the class that executes the Monte Carlo
#' simulation and provides the data for saving and plotting in this class.
#' @export
SimOP <- R6::R6Class(
"SimOP",
public = list(
#' @field simClass 'SimClass' R6 object that has been initialized and the OFPE
#' simulation has been executed. This class holds the information necessary
#' to produce the output figures, maps, and tables from the Monte Carlo
#' simulation.
simClass = NULL,
#' @field input_list Optional, list of inputs required for generating the
#' output figures. Recreation of relevant info from a 'SimClass' object. Only
#' needed if a 'SimClass' object is not passed into class.
#'
#' 'input_list' is a named list with the following fields;
#'
#' dat_path: The path to where the data from the simulation was exported.
#' Must contain 'Bp.var', 'NRffmax', and 'NRopt' tables exported from the
#' 'SimClass'. This class requires file names to be in the same format as
#' exported from the 'SimClass' to import.
#'
#' unique_fieldname: Unique name for the field(s) used in the simulation.
#' This is used for saving and labeling of outputs. e.g. "sec35middle" or
#' "sec1east & sec1west". This must match the 'unique_fieldname' used to
#' save the simulation outputs.
#'
#' unique_fxn: The function name generated by the 'SimClass' object that
#' was used for saving the simulation outputs. e.g. "yldGAM-proGAM".
#'
#' sim_years: Years that the simulation was run. The output tables from
#' the simulation must be labeled with these years as exported from the
#' 'SimClass' object.
#'
#' opt: The optimization method used in the simulation. Used to identify
#' and label data. Must be in the same format as from the 'SimClass', e.g.
#' "deriv".
#'
#' fieldsize: Size of the field in acres.
#'
#' fs: The uniform rate of the experimental input that the farmer would
#' have applied in the absence of variable rate application.
#'
#' expvar: The OFPE experimental input code for the experimental variable,
#' i.e. "aa_n" for as-applied nitrogen fertilizer or "aa_sr" for as-applied
#' seed rate.
#'
#' farmername: The name of the farmer, corresponding to the format in the
#' database, e.g. lowercase.
#'
#' respvar: Character vector of OFPE response variable codes that were
#' used in the model fitting and simulation processes, e.g. "yld" or "pro"
#' or 'c("yld", "pro")'.
#'
#' db: Connection to an OFPE formatted database.
#'
#' utm_fieldname: The name of the field (corresponding to the name in the
#' OFPE database) that should be used for identifying the UTM zone. In the
#' case of multiple fields (i.e. "sec1east & sec1west") choose one because
#' their UTM code will be the same.
#'
#' out_path: Path to the folder to create an 'Outputs' folder and save
#' figures and tables.
input_list = NULL,
#' @field out_path Output path to create folder for outputs. Taken from either
#' a 'SimClass' object or an optional 'input_list' field upon initialization.
out_path = NULL,
#' @field SAVE Whether to save simulation outputs If NULL uses the user selected
#' choice. If not NULL and is logical, argument replaces previously set SAVE
#' options for the entire class.
SAVE = NULL,
#' @field SI Logical, whether to plot in SI units. Defaults to FALSE where imperial
#' units (lbs, ac, bu) are used. If TRUE kg and ha are used.
SI = NULL,
#' @field dat_path The path to where the data from the simulation was exported.
#' Must contain 'Bp.var', 'NRffmax', and 'NRopt' tables exported from the
#' 'SimClass'. This class requires file names to be in the same format as
#' exported from the 'SimClass' to import.
dat_path = NULL,
#' @field unique_fieldname Unique name for the field(s) used in the simulation.
#' This is used for saving and labeling of outputs. e.g. "sec35middle" or
#' "sec1east & sec1west". This must match the 'unique_fieldname' used to
#' save the simulation outputs.
unique_fieldname = NULL,
#' @field unique_fxn The function name generated by the 'SimClass' object that
#' was used for saving the simulation outputs. e.g. "yldGAM-proGAM".
unique_fxn = NULL,
#' @field sim_years Years that the simulation was run. The output tables from
#' the simulation must be labeled with these years as exported from the
#' 'SimClass' object.
sim_years = NULL,
#' @field opt The optimization method used in the simulation. Used to identify
#' and label data. Must be in the same format as from the 'SimClass', e.g. "deriv".
opt = NULL,
#' @field fieldsize Size of the field. Note: this must be
#' in the same units as those specified in DatClass. i.e. if DatClass$SI == TRUE
#' than fieldsize should be in hectares, else in acres.
fieldsize = NULL,
#' @field fs The uniform rate of the experimental input that the farmer would
#' have applied in the absence of variable rate application. Note: this rate must be
#' in the same units as those specified in DatClass. i.e. if DatClass$SI == TRUE
#' than fs should be in kg/ha, else in lbs/ac.
fs = NULL,
#' @field expvar The OFPE experimental input code for the experimental variable,
#' i.e. "aa_n" for as-applied nitrogen fertilizer or "aa_sr" for as-applied
#' seed rate.
expvar = NULL,
#' @field farmername The name of the farmer, corresponding to the format in the
#' database, e.g. lowercase.
farmername = NULL,
#' @field respvar Character vector of OFPE response variable codes that were
#' used in the model fitting and simulation processes, e.g. "yld" or "pro"
#' or 'c("yld", "pro")'.
respvar = NULL,
#' @field db Connection to an OFPE formatted database.
db = NULL,
#' @field utm_fieldname The name of the field (corresponding to the name in the
#' OFPE database) that should be used for identifying the UTM zone. In the
#' case of multiple fields (i.e. "sec1east & sec1west") choose one because
#' their UTM code will be the same.
utm_fieldname = NULL,
#' @field AAmin The minimum as-applied rate to simulate management
#' outcomues from (i.e. 0 lbs N per acre or 25 lbs seed per acre).
AAmin = NULL,
#' @description
#' Use this class to save figures and tables from simulations performed using
#' the 'SimClass' R6 object. This class can be initialized with an executed
#' 'SimClass' object, or a list of inputs. See 'input_list' in this class'
#' description. This class can also be initialized with neither and methods
#' for generating individual figures and tables can be executed by passing
#' in arguments. When used in this manner this may be more susceptible to
#' errors as the user is required to specify arguments in the same format
#' as found in the 'SimClass' object.
#'
#' Also pass in a logical to the 'create' argument to select whether to save
#' figures. Depending on this choice, the output folder will be generated.
#' This is done by a private method that sets up the output location for
#' the figures that the model produces. This will not overwrite any previously
#' generated diagnostic or validation plots from the ModClass. Functions for
#' plots and maps all have save options which must be accompanied by a folder
#' path. The folder created is named 'Outputs'. This folder contains a
#' folder called 'Maps' which contains maps of estimated responses, net-return,
#' and site-specific optimized rates. Another folder called 'EXP' is created
#' that holds figures related to the amounts and distributions for as-applied
#' rates of the experimental variable. Within this folder, subfolders called
#' 'ffoptEXP', 'ssoptEXP', and 'EXPapplied' are created. The acronyms SSOPT and
#' FFOPT correspond to site-specific optimized and full-field optimized rates.
#' Another main folder named 'NR' is created to hold net-return figures comparing
#' management outcomes. Subfolders within 'NR' are 'NRboxplots', 'NRbarplots'
#' and 'NRprobabilities'. The two former contain figures, and the latter contains
#' a table with the probability that the SSOPT strategy yields a higher net-return
#' than each of the other management strategies.
#' @param simClass Optional, 'SimClass' R6 object that has been initialized and the OFPE
#' simulation has been executed. This class holds the information necessary
#' to produce the output figures, maps, and tables from the Monte Carlo
#' simulation.
#' @param input_list Optional, list of inputs required for generating the
#' output figures. Recreation of relevant info from a 'SimClass' object. Only
#' needed if a 'SimClass' object is not passed into class.
#' @param create create Logical, whether to create the 'Output' folder and save
#' figures, maps, and tables. Default is TRUE, pass FALSE to skip output folder
#' step and prevent any output generation.
#' @param SI Logical, whether to plot in SI units. Defaults to FALSE where imperial
#' units (lbs, ac, bu) are used. If TRUE kg and ha are used.
#' @return A folder created in the path for model output figures.
initialize = function(simClass = NULL,
input_list = NULL,
create = TRUE,
SI = FALSE) {
stopifnot(
is.logical(create)
)
if (!is.null(simClass)) {
self$simClass <- simClass
self$out_path <- self$simClass$out_path
self$SAVE <- self$simClass$SAVE
}
if (!is.null(input_list)) {
stopifnot(
!is.null(input_list$dat_path),
!is.null(input_list$unique_fieldname),
!is.null(input_list$unique_fxn),
!is.null(input_list$sim_years),
!is.null(input_list$opt),
!is.null(input_list$fieldsize),
!is.null(input_list$fs),
!is.null(input_list$EXPvec),
!is.null(input_list$expvar),
!is.null(input_list$farmername),
!is.null(input_list$respvar),
!is.null(input_list$db),
!is.null(input_list$utm_fieldname)
)
self$input_list <- input_list
stopifnot(!is.null(input_list$out_path))
self$out_path <- input_list$out_path
}
if (!is.null(self$simClass) | !is.null(self$input_list)) {
self$dat_path <- ifelse(!is.null(self$simClass),
self$simClass$dat_path,
self$input_list$dat_path)
self$unique_fieldname <- ifelse(!is.null(self$simClass),
self$simClass$unique_fieldname,
self$input_list$unique_fieldname)
self$unique_fxn <- ifelse(!is.null(self$simClass),
self$simClass$unique_fxn,
self$input_list$unique_fxn)
if (!is.null(self$simClass)) {
self$sim_years <- self$simClass$sim_years
} else {
self$sim_years <- self$input_list$sim_years
}
self$opt <- ifelse(!is.null(self$simClass),
self$simClass$opt,
self$input_list$opt)
self$fieldsize <- ifelse(!is.null(self$simClass),
self$simClass$fieldsize,
self$input_list$fieldsize)
self$fs <- ifelse(!is.null(self$simClass),
self$simClass$fs,
self$input_list$fs)
self$AAmin <- ifelse(!is.null(self$simClass),
self$simClass$EXPvec[1],
self$input_list$EXPvec[1])
self$expvar <- ifelse(!is.null(self$simClass),
self$simClass$datClass$expvar,
self$input_list$expvar)
self$farmername <- ifelse(!is.null(self$simClass),
self$simClass$datClass$farmername,
self$input_list$farmername)
if (!is.null(self$simClass)) {
self$respvar <- self$simClass$datClass$respvar
} else {
self$respvar <- self$input_list$respvar
}
if (!is.null(self$simClass)) {
self$db <- self$simClass$dbCon$db
} else {
self$db <- self$input_list$db
}
self$utm_fieldname <- ifelse(!is.null(self$simClass),
self$simClass$datClass$fieldname[1],
self$input_list$utm_fieldname)
}
if (create) {
private$.setupOP()
} else {
self$SAVE <- FALSE
}
stopifnot(is.logical(SI))
self$SI <- SI
},
#' @description
#' Method for saving simulation output plots (see below). This method is used
#' to call the public methods for each figure, map, and table generated. This
#' applies over the user selected simulation years. The user can pass a
#' logical argument to this function which will disable all plotting and
#' not save anything to your computer It will override any previous
#' setting you have set for 'SAVE'. If calling plots individually, a
#' SAVE option for each can be selected.
#'
#' This method will only work if the 'SimOP' class has been initialized with
#' a 'SimClass' or 'input_list' object. If the 'SimOP' class is initialized with
#' neither, than individual plotting functions must be executed with user
#' supplied inputs. This class uses internal fields in the 'SimOP' class that
#' are created on instantiation with a 'SimClass' or 'input_list' object.
#' class from these objects.
#'
#' The output methods and output created by this method are; 'plotNRbox',
#' 'plotNRbar', 'mgmtNRprobTable', 'plotFFOPThist', 'plotTotExpAppl',
#' and 'plotSSOPThist'. Maps are created with the 'plotSimMaps' methods.
#' The maps saved are net-return maps for the SSOPT, FFOPT, FS, and Act
#' management strategies, estimated yield and protein maps for the SSOPT
#' and FS management strategies, and the SSOPT experimental rate map.
#'
#' @param SAVE Whether to save simulation outputs If NULL uses the user selected
#' choice. If not NULL and is logical, argument replaces previously set SAVE
#' options for the entire class.
#' @return Simulation outputs in the 'Outputs' folder.
savePlots = function(SAVE = NULL) {
stopifnot(
!is.null(self$simClass) | !is.null(self$input_list)
)
if (is.null(SAVE)) {
stopifnot(!is.null(self$SAVE))
} else {
stopifnot(is.logical(SAVE))
self$SAVE <- SAVE
}
if (self$SI) {
exp_lab <- "kg/ha"
nr_lab <- "$/ha"
yld_lab <- "kg/ha"
} else {
exp_lab <- "lbs/ac"
nr_lab <- "$/ac"
yld_lab <- "bu/ac"
}
if (self$SAVE) {
for (i in 1:length(self$sim_years)) {
#### Bp.var plots ####
Bp.var <- self$loadBpVar(sim_year = self$sim_years[i])
## net-return boxplots
temp_plot <- self$plotNRbox(
Bp.var = Bp.var,
sim_year = self$sim_years[i]
)
## difference in net-return boxplots
temp_plot <- self$plotDiffNRbox(
Bp.var = Bp.var,
sim_year = self$sim_years[i]
)
## probability of higher avg net return
temp_plot <- self$mgmtNRprobTable(
mgmt = "NR.ssopt",
Bp.var = Bp.var,
rt = 1000,
sim_year = self$sim_years[i]
)
## net-return bar plot for all strategies
temp_plot <- self$plotNRbar(
Bp.var = Bp.var,
sim_year = self$sim_years[i]
)
#### Get NRopt & NRffmax & make TF4 ####
## gets a subset of all points for all rates and iterations
NRopt <- suppressWarnings(self$loadNRopt(sim_year = self$sim_years[i]))
NRffmax <- self$loadNRffmax(sim_year = self$sim_years[i])
TF4 <- self$makeTF4(NRopt = NRopt,
NRffmax = NRffmax)
rm(NRffmax)
#### NRopt & NRffmax plots ####
## hist of ffopt rate across iterations (price scens)
temp_plot <- self$plotFFOPThist(
Bp.var = Bp.var,
TF4 = TF4,
sim_year = self$sim_years[i]
)
rm(Bp.var)
#### End Bp.var Plots ####
## total n fert applied across strategies
temp_plot <- self$plotTotExpAppl(
TF4 = TF4,
sim_year = self$sim_years[i]
)
## hist of ssopt rates for mean price scen.
temp_plot <- self$plotSSOPThist(
NRopt = NRopt,
TF4 = TF4,
sim_year = self$sim_years[i]
)
#### Maps (avg all sim) ####
## SSOPT N map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = paste0("SSOPT_", ifelse(self$expvar == "aa_n", "N", "SR")),
var_col_name = "EXP.rate.ssopt",
var_label = paste0(ifelse(self$expvar == "aa_n", "N", "Seed Rate"), " (", exp_lab, ")"),
var_main_label = paste0("SS.opt ", ifelse(self$expvar == "aa_n", "N", "Seed"), " rates for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
## NR maps
# SSOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SSOPTNR",
var_col_name = "NR.ssopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("SS.Opt NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "FFOPTNR",
var_col_name = "NR.ffopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FF.Opt NR (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "FSNR",
var_col_name = "NR.fs",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FS NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "ActNR",
var_col_name = "NR.act",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR if exp in ", self$sim_years[i], " conditions (NR.act)"),
sim_year = self$sim_years[i]
)
# Min NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "MinNR",
var_col_name = "NR.min",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with uniform minimum rate in ", self$sim_years[i], " conditions (NR.min)"),
sim_year = self$sim_years[i]
)
# Opp NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "OppNR",
var_col_name = "NR.opp",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with opposite system pricing in ", self$sim_years[i], " conditions (NR.opp)"),
sim_year = self$sim_years[i]
)
## Yld maps
# SSOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldOpt",
var_col_name = "yld.opt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("SS.Opt. predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldFS",
var_col_name = "yld.fs",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FS predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldAct",
var_col_name = "yld.act",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldFFOPT",
var_col_name = "yld.ffopt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FF.Opt predicted yield (", NRopt$EXP.rate.ffopt[1], " ", exp_lab,") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldMin",
var_col_name = "yld.min",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
if (any(grepl("pro", self$respvar))) {
# SSOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProOpt",
var_col_name = "pro.opt",
var_label = "Protein (%)",
var_main_label = paste0("SS.Opt. predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProFS",
var_col_name = "pro.fs",
var_label = "Protein (%)",
var_main_label = paste0("FS predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProFFOPT",
var_col_name = "pro.ffopt",
var_label = "Protein (%)",
var_main_label = paste0("FF.Opt predicted protein (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProAct",
var_col_name = "pro.act",
var_label = "Protein (%)",
var_main_label = paste0("Protein if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProMin",
var_col_name = "pro.min",
var_label = "Protein (%)",
var_main_label = paste0("Protein with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
}
#### END Maps (avg all sim) ####
#### Maps (actual econ conditions of sim year) ####
NRopt <- data.table::fread(paste0(self$dat_path, self$sim_years[i], "/",
self$unique_fieldname, "_NRopt_",
self$unique_fxn, "_SimYr",
self$sim_years[i], "EconCond_",
self$opt, ".csv")) %>%
as.data.frame()
for (k in 1:ncol(NRopt)) {
if (!grepl("^x$|^y$", names(NRopt)[k])) {
if (all(is.na(NRopt[, k]))) {
NRopt[, k] <- 0
}
}
}
NRopt <- data.table::as.data.table(NRopt)
## SSOPT N map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = paste0("SimYrEconCond_SSOPT_", ifelse(self$expvar == "aa_n", "N", "SR")),
var_col_name = "EXP.rate.ssopt",
var_label = paste0(ifelse(self$expvar == "aa_n", "N", "Seed Rate"), " (", exp_lab, ")"),
var_main_label = paste0("SS.opt ", ifelse(self$expvar == "aa_n", "N", "Seed"), " rates for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
## NR maps
# SSOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_SSOPTNR",
var_col_name = "NR.ssopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("SS.Opt NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_FFOPTNR",
var_col_name = "NR.ffopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FF.Opt NR (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_FSNR",
var_col_name = "NR.fs",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FS NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_ActNR",
var_col_name = "NR.act",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR if exp in ", self$sim_years[i], " conditions (NR.act)"),
sim_year = self$sim_years[i]
)
# Min NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_MinNR",
var_col_name = "NR.min",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with uniform minimum rate in ", self$sim_years[i], " conditions (NR.min)"),
sim_year = self$sim_years[i]
)
# Opp NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_OppNR",
var_col_name = "NR.opp",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with opposite system pricing in ", self$sim_years[i], " conditions (NR.opp)"),
sim_year = self$sim_years[i]
)
## Yld maps
# SSOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldOpt",
var_col_name = "yld.opt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("SS.Opt. predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldFS",
var_col_name = "yld.fs",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FS predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldAct",
var_col_name = "yld.act",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldFFOPT",
var_col_name = "yld.ffopt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FF.Opt predicted yield (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldMin",
var_col_name = "yld.min",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
if (any(grepl("pro", self$respvar))) {
# SSOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProOpt",
var_col_name = "pro.opt",
var_label = "Protein (%)",
var_main_label = paste0("SS.Opt. predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProFS",
var_col_name = "pro.fs",
var_label = "Protein (%)",
var_main_label = paste0("FS predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProFFOPT",
var_col_name = "pro.ffopt",
var_label = "Protein (%)",
var_main_label = paste0("FF.Opt predicted protein (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProAct",
var_col_name = "pro.act",
var_label = "Protein (%)",
var_main_label = paste0("Protein if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProMin",
var_col_name = "pro.min",
var_label = "Protein (%)",
var_main_label = paste0("Protein with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
}
#### END Maps (actual econ conditions of sim year) ####
} # End sim year loop
rm(NRopt, TF4, temp_plot)
} else{
if (!is.null(self$simClass)) {
self$simClass$cleanUp()
}
}
},
#' @description
#' Plot the variation in net returns over different price years with a boxplot.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/NR/NRboxplots' if selected.
plotNRbox = function(Bp.var,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
nr_lab <- "$/ha"
} else {
nr_lab <- "$/ac"
}
Bp.plot <- as.data.frame(Bp.var)
Bp.plot <- tidyr::pivot_longer(Bp.plot,
-c("BaseP", "ffopt.EXPrate", "EXP.cost", "median.ssopt.EXPrate", "sim"),
names_to = "Method",
values_to = "NR")
Bp.plot$Method <- factor(Bp.plot$Method)
level_order <- c(grep("NR.min", levels(Bp.plot$Method)),
grep("NR.fs", levels(Bp.plot$Method)),
grep("NR.ssopt", levels(Bp.plot$Method)),
grep("NR.ffopt", levels(Bp.plot$Method)),
grep("NR.act", levels(Bp.plot$Method)),
grep("NR.opp", levels(Bp.plot$Method)))
Bp.plot$Method <- factor(Bp.plot$Method,
levels(Bp.plot$Method)[level_order])
mgmt_labels <- c("Min. Rate",
"FS",
"SS.Opt",
"FF.Opt",
"Actual",
"Alt. Price")
yMIN <- DescTools::RoundTo(min(Bp.plot$NR,na.rm=T),5,floor)
yMAX <- DescTools::RoundTo(max(Bp.plot$NR,na.rm=T),5,ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <-
ggplot2::ggplot(Bp.plot) +
ggplot2::geom_boxplot(ggplot2::aes(x = Method, y = NR),
fill = "green3",
notch = FALSE) +
ggplot2::scale_y_continuous(name = paste0("Average Net Return (", nr_lab, ")"),
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP)) +
ggplot2::scale_x_discrete(name = "Management Strategy", labels = mgmt_labels) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
if (any(Bp.plot$NR < 0, na.rm = TRUE)) {
p <- p +
ggplot2::geom_hline(yintercept = 0, color = "red", linetype = 2)
}
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/NR/NRboxPlots/",
fieldname, "_avgNR_box_",
fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Plot the variation in the difference in net returns over different price years.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/NR/NRboxplots' if selected.
plotDiffNRbox = function(Bp.var,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
nr_lab <- "$/ha"
} else {
nr_lab <- "$/ac"
}
Bp.plot <- as.data.frame(Bp.var)
Bp.plot$`SS.Opt - Min. Rate` <- Bp.plot$NR.ssopt - Bp.plot$NR.min
Bp.plot$`SS.Opt - FS` <- Bp.plot$NR.ssopt - Bp.plot$NR.fs
Bp.plot$`SS.Opt - FF.Opt` <- Bp.plot$NR.ssopt - Bp.plot$NR.ffopt
Bp.plot$`SS.Opt - Actual` <- Bp.plot$NR.ssopt - Bp.plot$NR.act
Bp.plot$`SS.Opt - Alt. Price` <- Bp.plot$NR.ssopt - Bp.plot$NR.opp
Bp.plot <- Bp.plot[, c("SS.Opt - Min. Rate", "SS.Opt - FS", "SS.Opt - FF.Opt", "SS.Opt - Actual", "SS.Opt - Alt. Price")]
Bp.plot <- tidyr::pivot_longer(Bp.plot,
c("SS.Opt - Min. Rate", "SS.Opt - FS", "SS.Opt - FF.Opt", "SS.Opt - Actual", "SS.Opt - Alt. Price"),
names_to = "Comparison",
values_to = "Diff.NR")
Bp.plot$Comparison <- factor(Bp.plot$Comparison,
levels = c("SS.Opt - Min. Rate", "SS.Opt - FS", "SS.Opt - FF.Opt", "SS.Opt - Actual", "SS.Opt - Alt. Price"))
yMIN <- DescTools::RoundTo(min(Bp.plot$Diff.NR,na.rm=T),5,floor)
yMAX <- DescTools::RoundTo(max(Bp.plot$Diff.NR,na.rm=T),5,ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <-
ggplot2::ggplot(Bp.plot) +
ggplot2::geom_boxplot(ggplot2::aes(x = Comparison, y = Diff.NR),
fill = "green3",
notch = FALSE) +
ggplot2::scale_y_continuous(name = paste0("Difference in Net Return (", nr_lab, ")"),
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP)) +
ggplot2::scale_x_discrete(name = "Strategy Comparison", labels = levels(Bp.plot$Comparison)) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
if (any(Bp.plot$Diff.NR < 0, na.rm = TRUE)) {
p <- p +
ggplot2::geom_hline(yintercept = 0, color = "red", linetype = 2)
}
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/NR/NRboxPlots/",
fieldname, "_diffAvgNR_box_",
fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png", scale = 1,
width = 10, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Calculate the probability that a given management method is more profitable
#' than the other management strategies. For example, compare the probability that
#' the site-specific optimization strategy yielded a higher net-retrun than the
#' other management strategies. This probability is generated from across all
#' iterations of the Monte Carlo simulation, so reflects the variation in economic
#' conditions.
#' @param mgmt Select the management strategy to compare, select from; 'NR.ssopt',
#' 'NR.min', 'NR.fs', 'NR.ffopt', 'NR.act', 'NR.opp'. See 'SimClass' documentation
#' or this class' description for acronym definitions.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param rt Number of times to sample the data with replacement to generate the
#' probability that the net-return of the selected strategy was higher than each
#' other strategy.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A data.frame and data saved in 'Outputs/NR/NRProbabilities' if selected.
mgmtNRprobTable = function(mgmt,
Bp.var,
rt,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
Bp.var <- as.data.frame(Bp.var)
mgmt_opts <- c("NR.min", "NR.fs", "NR.ssopt", "NR.ffopt", "NR.act", "NR.opp")
mgmt_opts <- mgmt_opts[-grep(mgmt, mgmt_opts)]
mgmt_prob <- rep(list(0), length(mgmt_opts)) %>%
`names<-`(mgmt_opts)
mgmt_col <- grep(mgmt, names(Bp.var))
for (nn in 1:rt) {
rp <- as.integer(runif(1, 1, nrow(Bp.var)))
for (j in 1:length(mgmt_prob)) {
alt_col <- grep(names(mgmt_prob)[j], names(Bp.var))
mgmt_prob[[j]] <- ifelse(Bp.var[rp, mgmt_col] > Bp.var[rp, alt_col],
1 + mgmt_prob[[j]],
0 + mgmt_prob[[j]])
}
}
mgmt_prob <- lapply(mgmt_prob, function(x) x / rt)
Table.p <- do.call(rbind, mgmt_prob) %>%
t() %>%
as.data.frame()
colnames(Table.p) <- paste0(mgmt, ">", names(Table.p))
row.names(Table.p) <- c("Probability")
if (SAVE) {
data.table::fwrite(Table.p,
paste0(out_path, "/Outputs/NR/NRProbabilities/",
fieldname, "_mgmtProbs_", fxn,
"_", sim_year, "_", opt, ".txt"),
sep="\t")
}
return(Table.p)
},
#' @description
#' Plot the mean net return for each management strategy as a
#' bar plot and a quarter standard deviation around each mean.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/NR/NRbarplots' if selected.
plotNRbar = function(Bp.var,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
nr_lab <- "$/ha"
} else {
nr_lab <- "$/ac"
}
Bp.var <- as.data.frame(Bp.var)
#Average net return over the field for the different N application methods without organic
TF3 <- data.frame(NR.ssopt = mean(Bp.var[, 'NR.ssopt'], na.rm = TRUE),
NR.min = mean(Bp.var[, 'NR.min'], na.rm = TRUE),
NR.fs = mean(Bp.var[, 'NR.fs'], na.rm = TRUE),
NR.ffopt = mean(Bp.var[, 'NR.ffopt'], na.rm = TRUE),
NR.act = mean(Bp.var[, 'NR.act'], na.rm = TRUE),
NR.opp = mean(Bp.var[, 'NR.opp'], na.rm = TRUE))
#If plotting metric values don't forget to relabel y-axis in barplot
#Put std. error bars on
TF3 <- tidyr::pivot_longer(TF3,
dplyr::everything(),
names_to = "Method",
values_to = "MeanNR")
TF3$Method <- factor(TF3$Method)
level_order <- c(grep("NR.min", levels(TF3$Method)),
grep("NR.fs", levels(TF3$Method)),
grep("NR.ssopt", levels(TF3$Method)),
grep("NR.ffopt", levels(TF3$Method)),
grep("NR.act", levels(TF3$Method)),
grep("NR.opp", levels(TF3$Method)))
TF3$Method <- factor(TF3$Method,
levels(TF3$Method)[level_order])
TF3$SDNR <- c(sd(Bp.var[, 'NR.ssopt'], na.rm = TRUE),
sd(Bp.var[, 'NR.min'], na.rm = TRUE),
sd(Bp.var[, 'NR.fs'], na.rm = TRUE),
sd(Bp.var[, 'NR.ffopt'], na.rm = TRUE),
sd(Bp.var[, 'NR.act'], na.rm = TRUE),
sd(Bp.var[, 'NR.opp'], na.rm = TRUE))
TF3$SDNR <- ifelse(is.na(TF3$SDNR), 0, TF3$SDNR)
TF3$max <- TF3$MeanNR + (1.96 * TF3$SDNR / 10)
TF3$min <- TF3$MeanNR - (1.96 * TF3$SDNR / 10)
if (any(any(na.omit(TF3$min) < 0) | any(na.omit(TF3$max) > 0))) { # if means above and below zero
MIN <- DescTools::RoundTo((min(TF3$MeanNR, na.rm = T) -
(1.96 * min(TF3$SDNR, na.rm = T) / 10)),
5, floor)
MAX <- DescTools::RoundTo((max(TF3$MeanNR, na.rm = T) +
(1.96 * max(TF3$SDNR, na.rm = T) / 10)),
5, ceiling)
STEP <- (MAX - MIN) / 10
}
if (all(na.omit(TF3$min) < 0)) { # if means below zero
MIN <- DescTools::RoundTo((min(TF3$MeanNR, na.rm = T) -
(1.96 * min(TF3$SDNR, na.rm = T) / 10)),
5, floor)
MAX <- 0
STEP <- (MAX - MIN) / 10
}
if (all(na.omit(TF3$min > 0))) { # if means above zero
MIN <- 0
MAX <- DescTools::RoundTo((max(TF3$MeanNR, na.rm = T) +
(1.96 * max(TF3$SDNR, na.rm = T) / 10)),
5, ceiling)
STEP <- (MAX - MIN) / 10
}
mgmt_labels <- c("Min. Rate",
"FS",
"SS.Opt",
"FF.Opt",
"Actual",
"Alt. Price")
p <-
ggplot2::ggplot(TF3, ggplot2::aes(x = Method, y = MeanNR)) +
ggplot2::geom_bar(stat = "identity", fill = "green3") +
ggplot2::geom_errorbar(ggplot2::aes(ymin = MeanNR - (1.96 * SDNR / 10),
ymax = MeanNR + (1.96 * SDNR / 10)),
width = .2) + #
ggplot2::scale_y_continuous(name = paste0("Average net return (", nr_lab, ")"),
limits = c(MIN, MAX),
breaks = seq(MIN, MAX, STEP)) +
ggplot2::scale_x_discrete(name = "Management Strategy",
labels = mgmt_labels) +
ggplot2::geom_text(data = TF3,
ggplot2::aes(x = Method, y = MeanNR, label = round(MeanNR, 2)),
hjust = 1.15,
vjust = 1.25) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/NR/NRbarPlots/",
fieldname, "_avgNR_bar_",
fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Create a bar plot of the total experimental input applied with
#' each strategy. The amount of the site-specific input applied
#' is averaged across the simulation for each point and summed to
#' get a total amount of experimental input applied with the site-
#' specific method. The total amount of experimental input applied
#' with all other strategies is calculated the same way.
#' @param TF4 Data.frame with columns for the amount of experimental rates applied
#' for each strategy.
#' @param expvar Experimental variable optimized, select/input
#' 'As-Applied Nitrogen' or 'As-Applied Seed Rate'. This is the type of
#' input that was experimentally varied across the field as part of the
#' on-farm experimentation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/EXP/EXPapplied' if selected.
plotTotExpAppl = function(TF4,
expvar = self$expvar,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
exp_lab <- "kg"
} else {
exp_lab <- "lbs"
}
if (TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] == 0) {
txt <- "FF.Opt Rate = 0"
} else {
if (TF4[which(TF4$Method == "EXP.ssopt"),"EXP"] >
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"]) {
txt <- paste0(as.integer(
(TF4[which(TF4$Method == "EXP.ssopt"), "EXP"] -
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"]) /
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] * 100),
"% more ",
ifelse(expvar == "aa_n", "N", "Seed"), " used in SS.Opt than FF.Opt")
}else{
txt <- paste0(as.integer(
(TF4[which(TF4$Method == "EXP.ssopt"), "EXP"] -
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"]) /
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] * -100),
"% less ",
ifelse(expvar == "aa_n", "N", "Seed"), " used in SS.Opt than FF.Opt")
}
}
y_lab <- paste0(ifelse(expvar == "aa_n", "Nitrogen", "Seed"),
" used on field (", exp_lab, "/field)")
yMIN <- 0
yMAX <- DescTools::RoundTo(max(TF4$EXP, na.rm = TRUE), 5, ceiling)
ySTEP <- (yMAX - yMIN) / 10
gg_title <- paste0(ifelse(expvar == "aa_n", "Nitrogen", "Seed"),
" used with each application strategy")
mgmt_labels <- c("Min. Rate", "FS", "SS.Opt", "FF.Opt")
p <-
ggplot2::ggplot(TF4, ggplot2::aes(x = Method, y = EXP)) +
ggplot2::geom_bar(stat = "identity", fill = "blue") +
ggplot2::scale_y_continuous(name = y_lab,
limits=c(yMIN, yMAX),
breaks=seq(yMIN, yMAX, ySTEP)) +
ggplot2::scale_x_discrete(name = "Management Strategy",
labels = mgmt_labels) +
ggplot2::geom_text(data = TF4,
ggplot2::aes(x = Method, y = EXP, label = round(EXP, 2)),
vjust = -.25) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14)) +
ggplot2::ggtitle(gg_title, subtitle = txt)
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/EXP/EXPapplied/",
fieldname, "_tot",
ifelse(expvar == "aa_n", "N", "SR"),
"applied_bar_", fxn, "_", sim_year, "_",
opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Generate a histogram of site specific experimental inptut rates within
#' a field. Uses the average site-specific rate for each point across the
#' simulation, indicating the mean optimum rate at each point across all
#' economic conditions.
#' @param NRopt Data frame with the net-returns and experimental optimums
#' for every point for every simulation iteration.
#' @param TF4 Data.frame with columns for the amount of experimental rates applied
#' for each strategy.
#' @param expvar Experimental variable optimized, select/input
#' 'As-Applied Nitrogen' or 'As-Applied Seed Rate'. This is the type of
#' input that was experimentally varied across the field as part of the
#' on-farm experimentation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/EXP/ssoptEXP' if selected.
plotSSOPThist = function(NRopt,
TF4,
expvar = self$expvar,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
exp_lab <- "kg/ha"
} else {
exp_lab <- "lbs/acre"
}
NRopt <- as.data.frame(NRopt)
gg_title <- paste0("Site-specific profit maximizing ",
ifelse(expvar=="aa_n", "nitrogen", "seed"), " rates")
sub_title <- paste0("Mean SS.Opt ",
ifelse(expvar == "aa_n", "N", "Seed"),
" Rate = ",
round(mean(NRopt$EXP.rate.ssopt, na.rm = TRUE), 2),
" ", exp_lab, ", SD = ",
round(sd(NRopt$EXP.rate.ssopt, na.rm = TRUE), 2),
" ", exp_lab)
if (TF4[which(TF4$Method == "EXP.ssopt"), "EXP"] == 0) {
p <-
ggplot2::ggplot(NRopt) +
ggplot2::geom_histogram(ggplot2::aes(x = EXP.rate.ssopt + 1),
bins = 1,
col = "white",
fill = "blue",
na.rm = TRUE) +
ggplot2::labs(x = exp_lab) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle = sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
} else {
bin_num <- DescTools::RoundTo(max(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 5, ceiling) / 5
xMIN <- DescTools::RoundTo(min(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 5, floor)
xMAX <- DescTools::RoundTo(max(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 5, ceiling)
xSTEP <- (xMAX - xMIN) / 10
p <-
ggplot2::ggplot(NRopt) +
ggplot2::geom_histogram(ggplot2::aes(x = EXP.rate.ssopt + 1),
bins = bin_num,
col = "white",
fill = "blue",
na.rm = TRUE) +
ggplot2::scale_x_continuous(name = exp_lab,
limits = c(xMIN, xMAX),
breaks = seq(xMIN, xMAX, xSTEP),
oob = scales::squish) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle = sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
}
yMIN <- 0
yMAX <- DescTools::RoundTo(max(ggplot2::ggplot_build(p)$data[[1]]$count), 5, ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <- p +
ggplot2::geom_vline(ggplot2::aes(xintercept = round(mean(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 0)),
col = "red",
linetype = 2) +
ggplot2::scale_y_continuous(name = "Frequency",
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP))
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/EXP/ssoptEXP/",
fieldname, "_ssopt",
ifelse(expvar == "aa_n", "N", "SR"),
"_hist_", fxn, "_", sim_year, "_",
opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Generate a histogram of full-field experimental input rates identified
#' in the simulation. Shows the distribution of full-field optimum rates
#' across all economic scenarios in the simulation.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param TF4 Data.frame with columns for the amount of experimental rates applied
#' for each strategy.
#' @param expvar Experimental variable optimized, select/input
#' 'As-Applied Nitrogen' or 'As-Applied Seed Rate'. This is the type of
#' input that was experimentally varied across the field as part of the
#' on-farm experimentation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/EXP/ffoptEXP' if selected.
plotFFOPThist = function(Bp.var,
TF4,
expvar = self$expvar,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
exp_lab <- "kg/ha"
} else {
exp_lab <- "lbs/acre"
}
Bp.var <- as.data.frame(Bp.var)
x_lab <- paste0("Profit Maximizing Top-Dress ",
ifelse(expvar == "aa_n", "N", "Seed"),
" Rate (", exp_lab, ")")
gg_title <- paste0("Variation in best full-field uniform ",
ifelse(expvar == "aa_n", "N", "seed"),
" rate over years")
sub_title <- paste0("Mean FF.Opt ",
ifelse(expvar == "aa_n", "N", "Seed"),
" Rate = ",
round(mean(Bp.var[, "ffopt.EXPrate"], na.rm = TRUE), 2),
" ", exp_lab, ", SD = ",
round(sd(Bp.var[, "ffopt.EXPrate"], na.rm = TRUE), 2),
" ", exp_lab)
if (TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] == 0) {
p <-
ggplot2::ggplot(Bp.var) +
ggplot2::geom_histogram(ggplot2::aes(x = ffopt.EXPrate + 1),
bins=1,
col="white",
fill="blue",
na.rm=TRUE) +
ggplot2::labs(x=x_lab) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle=sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
} else {
bin_num <- DescTools::RoundTo(max(Bp.var$ffopt.EXPrate,
na.rm = TRUE), 5, ceiling) / 5
xMIN <- DescTools::RoundTo(min(Bp.var$ffopt.EXPrate,
na.rm = TRUE), 5, floor)
xMAX <- DescTools::RoundTo(max(Bp.var$ffopt.EXPrate,
na.rm = TRUE), 5, ceiling)
xSTEP <- (xMAX - xMIN) / 10
p <-
ggplot2::ggplot(Bp.var) +
ggplot2::geom_histogram(ggplot2::aes(x = ffopt.EXPrate),
bins = bin_num,
col = "white",
fill = "blue",
na.rm = TRUE) +
ggplot2::scale_x_continuous(name = x_lab,
limits = c(xMIN, xMAX),
breaks = seq(xMIN, xMAX, xSTEP),
oob = scales::squish) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle = sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
}
yMIN <- DescTools::RoundTo(min(ggplot2::ggplot_build(p)$data[[1]]$count), 5, floor)
yMAX <- DescTools::RoundTo(max(ggplot2::ggplot_build(p)$data[[1]]$count), 5, ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <- p +
ggplot2::geom_vline(ggplot2::aes(xintercept = round(mean(Bp.var[, "ffopt.EXPrate"],
na.rm = TRUE), 0)),
col = "red",
linetype = 2) +
ggplot2::scale_y_continuous(name = "Frequency",
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP))
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/EXP/ffoptEXP/",
fieldname, "_ffopt",
ifelse(expvar == "aa_n", "N", "SR"),
"_hist_", fxn, "_", sim_year, "_",
opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' This method is for plotting maps of simulation outcomes. These
#' include the net-returns from the management strategies, predicted
#' responses, and optimized rates. Is a wrapper to OFPE::plotMaps().
#' @param dat Data frame with the net-returns and experimental optimums
#' for every point for every simulation iteration.
#' @param var The label of the variable to map. Used in figure name.
#' @param var_col_name The name of the column of the variable in the
#' supplied data ('dat').
#' @param var_label The label to be applied to the legend of the map
#' corresponding to the variable mapped.
#' @param var_main_label The main label to apply to the map.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param fxn The functional form of the models used for analysis.
#' @param opt The optimization method used in the simulation. For labeling.
#' @param SAVE Logical, whether to save figure.
#' @param farmername The name of the farmer that manages the field.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @param db Connection to the OFPE database to identify UTM zone.
#' @param utm_fieldname Name of the field for identifying the UTM zone.
#' @return A 'ggmap' object and maps saved in 'Outputs/Maps' folder if selected.
plotSimMaps = function(dat,
var,
var_col_name,
var_label,
var_main_label,
fieldname = self$unique_fieldname,
sim_year,
fxn = self$unique_fxn,
opt = self$opt,
SAVE = self$SAVE,
farmername = self$farmername,
out_path = self$out_path,
db = self$db,
utm_fieldname = self$utm_fieldname) {
utm_zone <- OFPE::findUTMzone(db,
fieldname = utm_fieldname)
p <- OFPE::plotMaps(dat,
var_col_name,
var_label,
var_main_label,
fieldname,
farmername,
utm_zone) %>%
suppressMessages() %>%
suppressWarnings()
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(
file = paste0(out_path, "/Outputs/Maps/", sim_year, "/",
fieldname, "_", tolower(var),
"_map_", fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png",
width = 7.5, height = 7.5, units = "in"
)
}
return(p)
},
#' @description
#' Gets the Bp.var output table from the 'SimClass' simulation. The Bp.var table contains means
#' across each field for each iteration of the simulation. This records base
#' price received, the cost of the input, the mean net-return per acre for
#' each management strategy, and the full-field optimum rate. Arguments are
#' used to fill in the file name in the same format as how the file was
#' saved from the simulation.
#' @param dat_path The path to the folder where the data is saved.
#' @param unique_fieldname Unique field name corresponding to the fieldname
#' in the file name.
#' @param unique_fxn The functional form of the models used for analysis
#' found in the file name.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used. Found in the file name.
#' @param opt The optimization method used in the simulation. Found in file
#' name.
#' @return Data frame with 10 columns and rows equal to the number of
#' iterations in the simulation.
loadBpVar = function(dat_path = self$dat_path,
unique_fieldname = self$unique_fieldname,
unique_fxn = self$unique_fxn,
sim_year,
opt = self$opt) {
Bp.var <- data.table::fread(
paste0(dat_path, sim_year, "/",
unique_fieldname, "_BpVar_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv")
)
return(Bp.var)
},
#' @description
#' Gets the NRffmax table that is an output from the 'SimClass' simulation with rows for
#' each iteration and the full-field optimum experimental input rate and the
#' total net-return from the full-field. Arguments are used to fill in the file
#' name in the same format as how the file was saved from the simulation.
#' @param dat_path The path to the folder where the data is saved.
#' @param unique_fieldname Unique field name corresponding to the fieldname
#' in the file name.
#' @param unique_fxn The functional form of the models used for analysis
#' found in the file name.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used. Found in the file name.
#' @param opt The optimization method used in the simulation. Found in file
#' name.
#' @return Data frame with 3 columns and rows equal to the number of
#' iterations in the simulation.
loadNRffmax = function(dat_path = self$dat_path,
unique_fieldname = self$unique_fieldname,
unique_fxn = self$unique_fxn,
sim_year,
opt = self$opt) {
NRffmax <- data.table::fread(
paste0(dat_path, sim_year, "/",
unique_fieldname, "_NRffMaxData_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv")
)
return(NRffmax)
},
#' @description
#' Create the table TF4, which is used for the output functions 'plotFFOPThist',
#' 'plotTotExpAppl', and 'plotSSOPThist'. This contains mean values for each
#' of the strategies to create bar charts of experimental inputs used.
#' @param NRopt 'SimClass' output table that contains every point with
#' the net-return and experimenal rates for each management strategy, and
#' the estimated crop responses for each location for the simulation year that
#' are the means across all simulations.
#' @param NRffmax 'SimClass' output table that contains the optimum full-field
#' experimental rate, the total net-return for each iteration of the simulation.
#' @param fieldsize Size of the field in acres.
#' @param fs The farmer selected experimental rate used in the simulation.
#' @return Data frame with SSOPT and FFOPT experimental input results for
#' each management strategy.
makeTF4 = function(NRopt,
NRffmax,
fieldsize = self$fieldsize,
fs = self$fs) {
TF4 <- data.frame(EXP.ssopt = (sum(NRopt[,'EXP.rate.ssopt'], na.rm = T) *
fieldsize) / nrow(NRopt),
EXP.min = self$AAmin * fieldsize,
EXP.fs = fs * fieldsize,
EXP.ffopt = mean(NRffmax$EXP.rate, na.rm = T) * fieldsize)
TF4 <- tidyr::pivot_longer(TF4,
dplyr::everything(),
names_to = "Method",
values_to = "EXP")
TF4$Method <- factor(TF4$Method)
mgmt_order <- c(grep("EXP.min", levels(TF4$Method)),
grep("EXP.fs", levels(TF4$Method)),
grep("EXP.ssopt", levels(TF4$Method)),
grep("EXP.ffopt", levels(TF4$Method)))
TF4$Method <- factor(TF4$Method, levels(TF4$Method)[mgmt_order])
return(TF4)
},
#' @description
#' Gets the NRopt table that is an output from the 'SimClass' simulation with rows for
#' each in the field. Columns contain the means across all iterations of the
#' simulation. Columns in this file are the price received, cost of nitrogen,
#' SSOPT and FFOPT experimental rates and the estimated net-return, yield, and
#' protein for each management strategy at each location in the field.
#' Arguments are used to fill in the file name in the same format as how the file was saved from the simulation.
#' @param dat_path The path to the folder where the data is saved.
#' @param unique_fieldname Unique field name corresponding to the fieldname
#' in the file name.
#' @param unique_fxn The functional form of the models used for analysis
#' found in the file name.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used. Found in the file name.
#' @param opt The optimization method used in the simulation. Found in file
#' name.
#' @return Data frame with 22 columns and rows equal to the number of
#' observations in the dataset.
loadNRopt = function(dat_path = self$dat_path,
unique_fieldname = self$unique_fieldname,
unique_fxn = self$unique_fxn,
sim_year,
opt = self$opt) {
on.exit(closeAllConnections())
## get col names
NRopt_names <- read.csv(file = paste0(dat_path, sim_year, "/",
unique_fieldname, "_NRopt_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv"),
nrows = 1,
header = TRUE) %>%
names()
## average all cols by x and y
NRopt_names <- NRopt_names[!1:length(NRopt_names) %in%
grep("^x$|^y$|^field$", NRopt_names)]
avg_query <- paste0("AVG([", NRopt_names, "])") %>%
paste(collapse = ", ")
NRopt <- sqldf::read.csv.sql(
paste0(dat_path, sim_year, "/",
unique_fieldname, "_NRopt_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv"),
sql = paste0("SELECT x, y, field, ", avg_query,
"FROM file GROUP BY x, y")
) %>%
`names<-`(c("x", "y", "field", NRopt_names))
NRcols <- grep("NR", names(NRopt))
for (i in 1:length(NRcols)) {
NRopt <- NRopt[NRopt[, NRcols[i]] > -50000 &
NRopt[, NRcols[i]] < 50000, ]
}
return(NRopt)
}
),
private = list(
.setupOP = function() {
cwd <- paste0(self$out_path, "/Outputs") # outputs working directory
if (!file.exists(cwd)) {
dir.create(cwd)
#dir.create(paste0(cwd, "/", "Exploratory"))
dir.create(paste0(cwd, "/", "Maps"))
for (i in 1:length(self$sim_years)) {
dir.create(paste0(cwd, "/", "Maps/", self$sim_years[i]))
}
dir.create(paste0(cwd, "/", "EXP"))
dir.create(paste0(cwd, "/", "EXP/ffoptEXP"))
dir.create(paste0(cwd, "/", "EXP/ssoptEXP"))
dir.create(paste0(cwd, "/", "EXP/EXPapplied"))
dir.create(paste0(cwd, "/", "NR"))
dir.create(paste0(cwd, "/", "NR/NRbarPlots"))
dir.create(paste0(cwd, "/", "NR/NRboxPlots"))
dir.create(paste0(cwd, "/", "NR/NRprobabilities"))
} else {
# if (!file.exists(paste0(cwd, "/", "Exploratory"))) {
# dir.create(paste0(cwd, "/", "Exploratory"))
# }
if (!file.exists(paste0(cwd, "/", "Maps"))) {
dir.create(paste0(cwd, "/", "Maps"))
for (i in 1:length(self$sim_years)) {
dir.create(paste0(cwd, "/", "Maps/", self$sim_years[i]))
}
} else {
for (i in 1:length(self$sim_years)) {
if (!file.exists(paste0(cwd, "/", "Maps/", self$sim_years[i]))) {
dir.create(paste0(cwd, "/", "Maps/", self$sim_years[i]))
}
}
}
if (!file.exists(paste0(cwd, "/", "EXP"))) {
dir.create(paste0(cwd, "/", "EXP"))
}
if (!file.exists(paste0(cwd, "/", "EXP/ffoptEXP"))) {
dir.create(paste0(cwd, "/", "EXP/ffoptEXP"))
}
if (!file.exists(paste0(cwd, "/", "EXP/ssoptEXP"))) {
dir.create(paste0(cwd, "/", "EXP/ssoptEXP"))
}
if (!file.exists(paste0(cwd, "/", "EXP/EXPapplied"))) {
dir.create(paste0(cwd, "/", "EXP/EXPapplied"))
}
if (!file.exists(paste0(cwd, "/", "NR"))) {
dir.create(paste0(cwd, "/", "NR"))
}
if (!file.exists(paste0(cwd, "/", "NR/NRbarPlots"))) {
dir.create(paste0(cwd, "/", "NR/NRbarPlots"))
}
if (!file.exists(paste0(cwd, "/", "NR/NRboxPlots"))) {
dir.create(paste0(cwd, "/", "NR/NRboxPlots"))
}
if (!file.exists(paste0(cwd, "/", "NR/NRprobabilities"))) {
dir.create(paste0(cwd, "/", "NR/NRprobabilities"))
}
}
}
)
)
paulhegedus/OFPE documentation built on Nov. 23, 2022, 5:09 a.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.
#' @title R6 Class for plotting outputs of the OFPE simulation.
#'
#' @description R6 class for plotting outputs of the OFPE simulation such as
#' figures comparing management scenarios, probability tables, and maps.
#'
#' This class should be initialized and executed after the user has initialized
#' and executed all other classes and has ran the SimClass$executeSim() method.
#' These figures rely on outputs from the simulation, so will cause errors if
#' used prior to execution of the simulation.
#'
#' See the SimClass$executeSim() method documentation for more detail on the management
#' scenarios reported. NR is net-return. The abbreviations used are; NRmin: uniform minimum
#' rate ('AAmin'), NRfs: uniform farmer selected rate, NRact: if experiment applied in
#' simulated year, NRssopt: site-specific optimum rates, NRffopt: uniform full-
#' field optimum rate, NRopp: if the crop received the opposite system crop price (e.g.
#' 'conventional' vs. 'organic' or reverse, see SimClass$executeSim() method
#' documentation for specifics on the difference in the calculation of this net-return
#' between systems).
#'
#' Error checking inputs to these functions is scant because the user should only use
#' these after executing the simulation method in the 'SimClass' R6 class to generate
#' the necessary data for the output generation methods.
#' @seealso \code{\link{SimClass}} for the class that executes the Monte Carlo
#' simulation and provides the data for saving and plotting in this class.
#' @export
SimOP <- R6::R6Class(
"SimOP",
public = list(
#' @field simClass 'SimClass' R6 object that has been initialized and the OFPE
#' simulation has been executed. This class holds the information necessary
#' to produce the output figures, maps, and tables from the Monte Carlo
#' simulation.
simClass = NULL,
#' @field input_list Optional, list of inputs required for generating the
#' output figures. Recreation of relevant info from a 'SimClass' object. Only
#' needed if a 'SimClass' object is not passed into class.
#'
#' 'input_list' is a named list with the following fields;
#'
#' dat_path: The path to where the data from the simulation was exported.
#' Must contain 'Bp.var', 'NRffmax', and 'NRopt' tables exported from the
#' 'SimClass'. This class requires file names to be in the same format as
#' exported from the 'SimClass' to import.
#'
#' unique_fieldname: Unique name for the field(s) used in the simulation.
#' This is used for saving and labeling of outputs. e.g. "sec35middle" or
#' "sec1east & sec1west". This must match the 'unique_fieldname' used to
#' save the simulation outputs.
#'
#' unique_fxn: The function name generated by the 'SimClass' object that
#' was used for saving the simulation outputs. e.g. "yldGAM-proGAM".
#'
#' sim_years: Years that the simulation was run. The output tables from
#' the simulation must be labeled with these years as exported from the
#' 'SimClass' object.
#'
#' opt: The optimization method used in the simulation. Used to identify
#' and label data. Must be in the same format as from the 'SimClass', e.g.
#' "deriv".
#'
#' fieldsize: Size of the field in acres.
#'
#' fs: The uniform rate of the experimental input that the farmer would
#' have applied in the absence of variable rate application.
#'
#' expvar: The OFPE experimental input code for the experimental variable,
#' i.e. "aa_n" for as-applied nitrogen fertilizer or "aa_sr" for as-applied
#' seed rate.
#'
#' farmername: The name of the farmer, corresponding to the format in the
#' database, e.g. lowercase.
#'
#' respvar: Character vector of OFPE response variable codes that were
#' used in the model fitting and simulation processes, e.g. "yld" or "pro"
#' or 'c("yld", "pro")'.
#'
#' db: Connection to an OFPE formatted database.
#'
#' utm_fieldname: The name of the field (corresponding to the name in the
#' OFPE database) that should be used for identifying the UTM zone. In the
#' case of multiple fields (i.e. "sec1east & sec1west") choose one because
#' their UTM code will be the same.
#'
#' out_path: Path to the folder to create an 'Outputs' folder and save
#' figures and tables.
input_list = NULL,
#' @field out_path Output path to create folder for outputs. Taken from either
#' a 'SimClass' object or an optional 'input_list' field upon initialization.
out_path = NULL,
#' @field SAVE Whether to save simulation outputs If NULL uses the user selected
#' choice. If not NULL and is logical, argument replaces previously set SAVE
#' options for the entire class.
SAVE = NULL,
#' @field SI Logical, whether to plot in SI units. Defaults to FALSE where imperial
#' units (lbs, ac, bu) are used. If TRUE kg and ha are used.
SI = NULL,
#' @field dat_path The path to where the data from the simulation was exported.
#' Must contain 'Bp.var', 'NRffmax', and 'NRopt' tables exported from the
#' 'SimClass'. This class requires file names to be in the same format as
#' exported from the 'SimClass' to import.
dat_path = NULL,
#' @field unique_fieldname Unique name for the field(s) used in the simulation.
#' This is used for saving and labeling of outputs. e.g. "sec35middle" or
#' "sec1east & sec1west". This must match the 'unique_fieldname' used to
#' save the simulation outputs.
unique_fieldname = NULL,
#' @field unique_fxn The function name generated by the 'SimClass' object that
#' was used for saving the simulation outputs. e.g. "yldGAM-proGAM".
unique_fxn = NULL,
#' @field sim_years Years that the simulation was run. The output tables from
#' the simulation must be labeled with these years as exported from the
#' 'SimClass' object.
sim_years = NULL,
#' @field opt The optimization method used in the simulation. Used to identify
#' and label data. Must be in the same format as from the 'SimClass', e.g. "deriv".
opt = NULL,
#' @field fieldsize Size of the field. Note: this must be
#' in the same units as those specified in DatClass. i.e. if DatClass$SI == TRUE
#' than fieldsize should be in hectares, else in acres.
fieldsize = NULL,
#' @field fs The uniform rate of the experimental input that the farmer would
#' have applied in the absence of variable rate application. Note: this rate must be
#' in the same units as those specified in DatClass. i.e. if DatClass$SI == TRUE
#' than fs should be in kg/ha, else in lbs/ac.
fs = NULL,
#' @field expvar The OFPE experimental input code for the experimental variable,
#' i.e. "aa_n" for as-applied nitrogen fertilizer or "aa_sr" for as-applied
#' seed rate.
expvar = NULL,
#' @field farmername The name of the farmer, corresponding to the format in the
#' database, e.g. lowercase.
farmername = NULL,
#' @field respvar Character vector of OFPE response variable codes that were
#' used in the model fitting and simulation processes, e.g. "yld" or "pro"
#' or 'c("yld", "pro")'.
respvar = NULL,
#' @field db Connection to an OFPE formatted database.
db = NULL,
#' @field utm_fieldname The name of the field (corresponding to the name in the
#' OFPE database) that should be used for identifying the UTM zone. In the
#' case of multiple fields (i.e. "sec1east & sec1west") choose one because
#' their UTM code will be the same.
utm_fieldname = NULL,
#' @field AAmin The minimum as-applied rate to simulate management
#' outcomues from (i.e. 0 lbs N per acre or 25 lbs seed per acre).
AAmin = NULL,
#' @description
#' Use this class to save figures and tables from simulations performed using
#' the 'SimClass' R6 object. This class can be initialized with an executed
#' 'SimClass' object, or a list of inputs. See 'input_list' in this class'
#' description. This class can also be initialized with neither and methods
#' for generating individual figures and tables can be executed by passing
#' in arguments. When used in this manner this may be more susceptible to
#' errors as the user is required to specify arguments in the same format
#' as found in the 'SimClass' object.
#'
#' Also pass in a logical to the 'create' argument to select whether to save
#' figures. Depending on this choice, the output folder will be generated.
#' This is done by a private method that sets up the output location for
#' the figures that the model produces. This will not overwrite any previously
#' generated diagnostic or validation plots from the ModClass. Functions for
#' plots and maps all have save options which must be accompanied by a folder
#' path. The folder created is named 'Outputs'. This folder contains a
#' folder called 'Maps' which contains maps of estimated responses, net-return,
#' and site-specific optimized rates. Another folder called 'EXP' is created
#' that holds figures related to the amounts and distributions for as-applied
#' rates of the experimental variable. Within this folder, subfolders called
#' 'ffoptEXP', 'ssoptEXP', and 'EXPapplied' are created. The acronyms SSOPT and
#' FFOPT correspond to site-specific optimized and full-field optimized rates.
#' Another main folder named 'NR' is created to hold net-return figures comparing
#' management outcomes. Subfolders within 'NR' are 'NRboxplots', 'NRbarplots'
#' and 'NRprobabilities'. The two former contain figures, and the latter contains
#' a table with the probability that the SSOPT strategy yields a higher net-return
#' than each of the other management strategies.
#' @param simClass Optional, 'SimClass' R6 object that has been initialized and the OFPE
#' simulation has been executed. This class holds the information necessary
#' to produce the output figures, maps, and tables from the Monte Carlo
#' simulation.
#' @param input_list Optional, list of inputs required for generating the
#' output figures. Recreation of relevant info from a 'SimClass' object. Only
#' needed if a 'SimClass' object is not passed into class.
#' @param create create Logical, whether to create the 'Output' folder and save
#' figures, maps, and tables. Default is TRUE, pass FALSE to skip output folder
#' step and prevent any output generation.
#' @param SI Logical, whether to plot in SI units. Defaults to FALSE where imperial
#' units (lbs, ac, bu) are used. If TRUE kg and ha are used.
#' @return A folder created in the path for model output figures.
initialize = function(simClass = NULL,
input_list = NULL,
create = TRUE,
SI = FALSE) {
stopifnot(
is.logical(create)
)
if (!is.null(simClass)) {
self$simClass <- simClass
self$out_path <- self$simClass$out_path
self$SAVE <- self$simClass$SAVE
}
if (!is.null(input_list)) {
stopifnot(
!is.null(input_list$dat_path),
!is.null(input_list$unique_fieldname),
!is.null(input_list$unique_fxn),
!is.null(input_list$sim_years),
!is.null(input_list$opt),
!is.null(input_list$fieldsize),
!is.null(input_list$fs),
!is.null(input_list$EXPvec),
!is.null(input_list$expvar),
!is.null(input_list$farmername),
!is.null(input_list$respvar),
!is.null(input_list$db),
!is.null(input_list$utm_fieldname)
)
self$input_list <- input_list
stopifnot(!is.null(input_list$out_path))
self$out_path <- input_list$out_path
}
if (!is.null(self$simClass) | !is.null(self$input_list)) {
self$dat_path <- ifelse(!is.null(self$simClass),
self$simClass$dat_path,
self$input_list$dat_path)
self$unique_fieldname <- ifelse(!is.null(self$simClass),
self$simClass$unique_fieldname,
self$input_list$unique_fieldname)
self$unique_fxn <- ifelse(!is.null(self$simClass),
self$simClass$unique_fxn,
self$input_list$unique_fxn)
if (!is.null(self$simClass)) {
self$sim_years <- self$simClass$sim_years
} else {
self$sim_years <- self$input_list$sim_years
}
self$opt <- ifelse(!is.null(self$simClass),
self$simClass$opt,
self$input_list$opt)
self$fieldsize <- ifelse(!is.null(self$simClass),
self$simClass$fieldsize,
self$input_list$fieldsize)
self$fs <- ifelse(!is.null(self$simClass),
self$simClass$fs,
self$input_list$fs)
self$AAmin <- ifelse(!is.null(self$simClass),
self$simClass$EXPvec[1],
self$input_list$EXPvec[1])
self$expvar <- ifelse(!is.null(self$simClass),
self$simClass$datClass$expvar,
self$input_list$expvar)
self$farmername <- ifelse(!is.null(self$simClass),
self$simClass$datClass$farmername,
self$input_list$farmername)
if (!is.null(self$simClass)) {
self$respvar <- self$simClass$datClass$respvar
} else {
self$respvar <- self$input_list$respvar
}
if (!is.null(self$simClass)) {
self$db <- self$simClass$dbCon$db
} else {
self$db <- self$input_list$db
}
self$utm_fieldname <- ifelse(!is.null(self$simClass),
self$simClass$datClass$fieldname[1],
self$input_list$utm_fieldname)
}
if (create) {
private$.setupOP()
} else {
self$SAVE <- FALSE
}
stopifnot(is.logical(SI))
self$SI <- SI
},
#' @description
#' Method for saving simulation output plots (see below). This method is used
#' to call the public methods for each figure, map, and table generated. This
#' applies over the user selected simulation years. The user can pass a
#' logical argument to this function which will disable all plotting and
#' not save anything to your computer It will override any previous
#' setting you have set for 'SAVE'. If calling plots individually, a
#' SAVE option for each can be selected.
#'
#' This method will only work if the 'SimOP' class has been initialized with
#' a 'SimClass' or 'input_list' object. If the 'SimOP' class is initialized with
#' neither, than individual plotting functions must be executed with user
#' supplied inputs. This class uses internal fields in the 'SimOP' class that
#' are created on instantiation with a 'SimClass' or 'input_list' object.
#' class from these objects.
#'
#' The output methods and output created by this method are; 'plotNRbox',
#' 'plotNRbar', 'mgmtNRprobTable', 'plotFFOPThist', 'plotTotExpAppl',
#' and 'plotSSOPThist'. Maps are created with the 'plotSimMaps' methods.
#' The maps saved are net-return maps for the SSOPT, FFOPT, FS, and Act
#' management strategies, estimated yield and protein maps for the SSOPT
#' and FS management strategies, and the SSOPT experimental rate map.
#'
#' @param SAVE Whether to save simulation outputs If NULL uses the user selected
#' choice. If not NULL and is logical, argument replaces previously set SAVE
#' options for the entire class.
#' @return Simulation outputs in the 'Outputs' folder.
savePlots = function(SAVE = NULL) {
stopifnot(
!is.null(self$simClass) | !is.null(self$input_list)
)
if (is.null(SAVE)) {
stopifnot(!is.null(self$SAVE))
} else {
stopifnot(is.logical(SAVE))
self$SAVE <- SAVE
}
if (self$SI) {
exp_lab <- "kg/ha"
nr_lab <- "$/ha"
yld_lab <- "kg/ha"
} else {
exp_lab <- "lbs/ac"
nr_lab <- "$/ac"
yld_lab <- "bu/ac"
}
if (self$SAVE) {
for (i in 1:length(self$sim_years)) {
#### Bp.var plots ####
Bp.var <- self$loadBpVar(sim_year = self$sim_years[i])
## net-return boxplots
temp_plot <- self$plotNRbox(
Bp.var = Bp.var,
sim_year = self$sim_years[i]
)
## difference in net-return boxplots
temp_plot <- self$plotDiffNRbox(
Bp.var = Bp.var,
sim_year = self$sim_years[i]
)
## probability of higher avg net return
temp_plot <- self$mgmtNRprobTable(
mgmt = "NR.ssopt",
Bp.var = Bp.var,
rt = 1000,
sim_year = self$sim_years[i]
)
## net-return bar plot for all strategies
temp_plot <- self$plotNRbar(
Bp.var = Bp.var,
sim_year = self$sim_years[i]
)
#### Get NRopt & NRffmax & make TF4 ####
## gets a subset of all points for all rates and iterations
NRopt <- suppressWarnings(self$loadNRopt(sim_year = self$sim_years[i]))
NRffmax <- self$loadNRffmax(sim_year = self$sim_years[i])
TF4 <- self$makeTF4(NRopt = NRopt,
NRffmax = NRffmax)
rm(NRffmax)
#### NRopt & NRffmax plots ####
## hist of ffopt rate across iterations (price scens)
temp_plot <- self$plotFFOPThist(
Bp.var = Bp.var,
TF4 = TF4,
sim_year = self$sim_years[i]
)
rm(Bp.var)
#### End Bp.var Plots ####
## total n fert applied across strategies
temp_plot <- self$plotTotExpAppl(
TF4 = TF4,
sim_year = self$sim_years[i]
)
## hist of ssopt rates for mean price scen.
temp_plot <- self$plotSSOPThist(
NRopt = NRopt,
TF4 = TF4,
sim_year = self$sim_years[i]
)
#### Maps (avg all sim) ####
## SSOPT N map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = paste0("SSOPT_", ifelse(self$expvar == "aa_n", "N", "SR")),
var_col_name = "EXP.rate.ssopt",
var_label = paste0(ifelse(self$expvar == "aa_n", "N", "Seed Rate"), " (", exp_lab, ")"),
var_main_label = paste0("SS.opt ", ifelse(self$expvar == "aa_n", "N", "Seed"), " rates for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
## NR maps
# SSOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SSOPTNR",
var_col_name = "NR.ssopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("SS.Opt NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "FFOPTNR",
var_col_name = "NR.ffopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FF.Opt NR (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "FSNR",
var_col_name = "NR.fs",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FS NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "ActNR",
var_col_name = "NR.act",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR if exp in ", self$sim_years[i], " conditions (NR.act)"),
sim_year = self$sim_years[i]
)
# Min NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "MinNR",
var_col_name = "NR.min",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with uniform minimum rate in ", self$sim_years[i], " conditions (NR.min)"),
sim_year = self$sim_years[i]
)
# Opp NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "OppNR",
var_col_name = "NR.opp",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with opposite system pricing in ", self$sim_years[i], " conditions (NR.opp)"),
sim_year = self$sim_years[i]
)
## Yld maps
# SSOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldOpt",
var_col_name = "yld.opt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("SS.Opt. predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldFS",
var_col_name = "yld.fs",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FS predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldAct",
var_col_name = "yld.act",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldFFOPT",
var_col_name = "yld.ffopt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FF.Opt predicted yield (", NRopt$EXP.rate.ffopt[1], " ", exp_lab,") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estYldMin",
var_col_name = "yld.min",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
if (any(grepl("pro", self$respvar))) {
# SSOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProOpt",
var_col_name = "pro.opt",
var_label = "Protein (%)",
var_main_label = paste0("SS.Opt. predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProFS",
var_col_name = "pro.fs",
var_label = "Protein (%)",
var_main_label = paste0("FS predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProFFOPT",
var_col_name = "pro.ffopt",
var_label = "Protein (%)",
var_main_label = paste0("FF.Opt predicted protein (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProAct",
var_col_name = "pro.act",
var_label = "Protein (%)",
var_main_label = paste0("Protein if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "estProMin",
var_col_name = "pro.min",
var_label = "Protein (%)",
var_main_label = paste0("Protein with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
}
#### END Maps (avg all sim) ####
#### Maps (actual econ conditions of sim year) ####
NRopt <- data.table::fread(paste0(self$dat_path, self$sim_years[i], "/",
self$unique_fieldname, "_NRopt_",
self$unique_fxn, "_SimYr",
self$sim_years[i], "EconCond_",
self$opt, ".csv")) %>%
as.data.frame()
for (k in 1:ncol(NRopt)) {
if (!grepl("^x$|^y$", names(NRopt)[k])) {
if (all(is.na(NRopt[, k]))) {
NRopt[, k] <- 0
}
}
}
NRopt <- data.table::as.data.table(NRopt)
## SSOPT N map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = paste0("SimYrEconCond_SSOPT_", ifelse(self$expvar == "aa_n", "N", "SR")),
var_col_name = "EXP.rate.ssopt",
var_label = paste0(ifelse(self$expvar == "aa_n", "N", "Seed Rate"), " (", exp_lab, ")"),
var_main_label = paste0("SS.opt ", ifelse(self$expvar == "aa_n", "N", "Seed"), " rates for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
## NR maps
# SSOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_SSOPTNR",
var_col_name = "NR.ssopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("SS.Opt NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_FFOPTNR",
var_col_name = "NR.ffopt",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FF.Opt NR (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_FSNR",
var_col_name = "NR.fs",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("FS NR for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_ActNR",
var_col_name = "NR.act",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR if exp in ", self$sim_years[i], " conditions (NR.act)"),
sim_year = self$sim_years[i]
)
# Min NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_MinNR",
var_col_name = "NR.min",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with uniform minimum rate in ", self$sim_years[i], " conditions (NR.min)"),
sim_year = self$sim_years[i]
)
# Opp NR map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_OppNR",
var_col_name = "NR.opp",
var_label = paste0("Net Return (", nr_lab, ")"),
var_main_label = paste0("NR with opposite system pricing in ", self$sim_years[i], " conditions (NR.opp)"),
sim_year = self$sim_years[i]
)
## Yld maps
# SSOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldOpt",
var_col_name = "yld.opt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("SS.Opt. predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldFS",
var_col_name = "yld.fs",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FS predicted yield for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldAct",
var_col_name = "yld.act",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldFFOPT",
var_col_name = "yld.ffopt",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("FF.Opt predicted yield (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Yld map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estYldMin",
var_col_name = "yld.min",
var_label = paste0("Yield (", yld_lab, ")"),
var_main_label = paste0("Yield with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
if (any(grepl("pro", self$respvar))) {
# SSOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProOpt",
var_col_name = "pro.opt",
var_label = "Protein (%)",
var_main_label = paste0("SS.Opt. predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FS Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProFS",
var_col_name = "pro.fs",
var_label = "Protein (%)",
var_main_label = paste0("FS predicted protein for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# FFOPT Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProFFOPT",
var_col_name = "pro.ffopt",
var_label = "Protein (%)",
var_main_label = paste0("FF.Opt predicted protein (", NRopt$EXP.rate.ffopt[1], " ", exp_lab, ") for ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Act Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProAct",
var_col_name = "pro.act",
var_label = "Protein (%)",
var_main_label = paste0("Protein if exp in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
# Min Est. Pro map
temp_plot <- self$plotSimMaps(
dat = NRopt,
var = "SimYrEconCond_estProMin",
var_col_name = "pro.min",
var_label = "Protein (%)",
var_main_label = paste0("Protein with uniform minimum rate in ", self$sim_years[i], " conditions"),
sim_year = self$sim_years[i]
)
}
#### END Maps (actual econ conditions of sim year) ####
} # End sim year loop
rm(NRopt, TF4, temp_plot)
} else{
if (!is.null(self$simClass)) {
self$simClass$cleanUp()
}
}
},
#' @description
#' Plot the variation in net returns over different price years with a boxplot.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/NR/NRboxplots' if selected.
plotNRbox = function(Bp.var,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
nr_lab <- "$/ha"
} else {
nr_lab <- "$/ac"
}
Bp.plot <- as.data.frame(Bp.var)
Bp.plot <- tidyr::pivot_longer(Bp.plot,
-c("BaseP", "ffopt.EXPrate", "EXP.cost", "median.ssopt.EXPrate", "sim"),
names_to = "Method",
values_to = "NR")
Bp.plot$Method <- factor(Bp.plot$Method)
level_order <- c(grep("NR.min", levels(Bp.plot$Method)),
grep("NR.fs", levels(Bp.plot$Method)),
grep("NR.ssopt", levels(Bp.plot$Method)),
grep("NR.ffopt", levels(Bp.plot$Method)),
grep("NR.act", levels(Bp.plot$Method)),
grep("NR.opp", levels(Bp.plot$Method)))
Bp.plot$Method <- factor(Bp.plot$Method,
levels(Bp.plot$Method)[level_order])
mgmt_labels <- c("Min. Rate",
"FS",
"SS.Opt",
"FF.Opt",
"Actual",
"Alt. Price")
yMIN <- DescTools::RoundTo(min(Bp.plot$NR,na.rm=T),5,floor)
yMAX <- DescTools::RoundTo(max(Bp.plot$NR,na.rm=T),5,ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <-
ggplot2::ggplot(Bp.plot) +
ggplot2::geom_boxplot(ggplot2::aes(x = Method, y = NR),
fill = "green3",
notch = FALSE) +
ggplot2::scale_y_continuous(name = paste0("Average Net Return (", nr_lab, ")"),
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP)) +
ggplot2::scale_x_discrete(name = "Management Strategy", labels = mgmt_labels) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
if (any(Bp.plot$NR < 0, na.rm = TRUE)) {
p <- p +
ggplot2::geom_hline(yintercept = 0, color = "red", linetype = 2)
}
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/NR/NRboxPlots/",
fieldname, "_avgNR_box_",
fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Plot the variation in the difference in net returns over different price years.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/NR/NRboxplots' if selected.
plotDiffNRbox = function(Bp.var,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
nr_lab <- "$/ha"
} else {
nr_lab <- "$/ac"
}
Bp.plot <- as.data.frame(Bp.var)
Bp.plot$`SS.Opt - Min. Rate` <- Bp.plot$NR.ssopt - Bp.plot$NR.min
Bp.plot$`SS.Opt - FS` <- Bp.plot$NR.ssopt - Bp.plot$NR.fs
Bp.plot$`SS.Opt - FF.Opt` <- Bp.plot$NR.ssopt - Bp.plot$NR.ffopt
Bp.plot$`SS.Opt - Actual` <- Bp.plot$NR.ssopt - Bp.plot$NR.act
Bp.plot$`SS.Opt - Alt. Price` <- Bp.plot$NR.ssopt - Bp.plot$NR.opp
Bp.plot <- Bp.plot[, c("SS.Opt - Min. Rate", "SS.Opt - FS", "SS.Opt - FF.Opt", "SS.Opt - Actual", "SS.Opt - Alt. Price")]
Bp.plot <- tidyr::pivot_longer(Bp.plot,
c("SS.Opt - Min. Rate", "SS.Opt - FS", "SS.Opt - FF.Opt", "SS.Opt - Actual", "SS.Opt - Alt. Price"),
names_to = "Comparison",
values_to = "Diff.NR")
Bp.plot$Comparison <- factor(Bp.plot$Comparison,
levels = c("SS.Opt - Min. Rate", "SS.Opt - FS", "SS.Opt - FF.Opt", "SS.Opt - Actual", "SS.Opt - Alt. Price"))
yMIN <- DescTools::RoundTo(min(Bp.plot$Diff.NR,na.rm=T),5,floor)
yMAX <- DescTools::RoundTo(max(Bp.plot$Diff.NR,na.rm=T),5,ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <-
ggplot2::ggplot(Bp.plot) +
ggplot2::geom_boxplot(ggplot2::aes(x = Comparison, y = Diff.NR),
fill = "green3",
notch = FALSE) +
ggplot2::scale_y_continuous(name = paste0("Difference in Net Return (", nr_lab, ")"),
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP)) +
ggplot2::scale_x_discrete(name = "Strategy Comparison", labels = levels(Bp.plot$Comparison)) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
if (any(Bp.plot$Diff.NR < 0, na.rm = TRUE)) {
p <- p +
ggplot2::geom_hline(yintercept = 0, color = "red", linetype = 2)
}
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/NR/NRboxPlots/",
fieldname, "_diffAvgNR_box_",
fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png", scale = 1,
width = 10, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Calculate the probability that a given management method is more profitable
#' than the other management strategies. For example, compare the probability that
#' the site-specific optimization strategy yielded a higher net-retrun than the
#' other management strategies. This probability is generated from across all
#' iterations of the Monte Carlo simulation, so reflects the variation in economic
#' conditions.
#' @param mgmt Select the management strategy to compare, select from; 'NR.ssopt',
#' 'NR.min', 'NR.fs', 'NR.ffopt', 'NR.act', 'NR.opp'. See 'SimClass' documentation
#' or this class' description for acronym definitions.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param rt Number of times to sample the data with replacement to generate the
#' probability that the net-return of the selected strategy was higher than each
#' other strategy.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A data.frame and data saved in 'Outputs/NR/NRProbabilities' if selected.
mgmtNRprobTable = function(mgmt,
Bp.var,
rt,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
Bp.var <- as.data.frame(Bp.var)
mgmt_opts <- c("NR.min", "NR.fs", "NR.ssopt", "NR.ffopt", "NR.act", "NR.opp")
mgmt_opts <- mgmt_opts[-grep(mgmt, mgmt_opts)]
mgmt_prob <- rep(list(0), length(mgmt_opts)) %>%
`names<-`(mgmt_opts)
mgmt_col <- grep(mgmt, names(Bp.var))
for (nn in 1:rt) {
rp <- as.integer(runif(1, 1, nrow(Bp.var)))
for (j in 1:length(mgmt_prob)) {
alt_col <- grep(names(mgmt_prob)[j], names(Bp.var))
mgmt_prob[[j]] <- ifelse(Bp.var[rp, mgmt_col] > Bp.var[rp, alt_col],
1 + mgmt_prob[[j]],
0 + mgmt_prob[[j]])
}
}
mgmt_prob <- lapply(mgmt_prob, function(x) x / rt)
Table.p <- do.call(rbind, mgmt_prob) %>%
t() %>%
as.data.frame()
colnames(Table.p) <- paste0(mgmt, ">", names(Table.p))
row.names(Table.p) <- c("Probability")
if (SAVE) {
data.table::fwrite(Table.p,
paste0(out_path, "/Outputs/NR/NRProbabilities/",
fieldname, "_mgmtProbs_", fxn,
"_", sim_year, "_", opt, ".txt"),
sep="\t")
}
return(Table.p)
},
#' @description
#' Plot the mean net return for each management strategy as a
#' bar plot and a quarter standard deviation around each mean.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/NR/NRbarplots' if selected.
plotNRbar = function(Bp.var,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
nr_lab <- "$/ha"
} else {
nr_lab <- "$/ac"
}
Bp.var <- as.data.frame(Bp.var)
#Average net return over the field for the different N application methods without organic
TF3 <- data.frame(NR.ssopt = mean(Bp.var[, 'NR.ssopt'], na.rm = TRUE),
NR.min = mean(Bp.var[, 'NR.min'], na.rm = TRUE),
NR.fs = mean(Bp.var[, 'NR.fs'], na.rm = TRUE),
NR.ffopt = mean(Bp.var[, 'NR.ffopt'], na.rm = TRUE),
NR.act = mean(Bp.var[, 'NR.act'], na.rm = TRUE),
NR.opp = mean(Bp.var[, 'NR.opp'], na.rm = TRUE))
#If plotting metric values don't forget to relabel y-axis in barplot
#Put std. error bars on
TF3 <- tidyr::pivot_longer(TF3,
dplyr::everything(),
names_to = "Method",
values_to = "MeanNR")
TF3$Method <- factor(TF3$Method)
level_order <- c(grep("NR.min", levels(TF3$Method)),
grep("NR.fs", levels(TF3$Method)),
grep("NR.ssopt", levels(TF3$Method)),
grep("NR.ffopt", levels(TF3$Method)),
grep("NR.act", levels(TF3$Method)),
grep("NR.opp", levels(TF3$Method)))
TF3$Method <- factor(TF3$Method,
levels(TF3$Method)[level_order])
TF3$SDNR <- c(sd(Bp.var[, 'NR.ssopt'], na.rm = TRUE),
sd(Bp.var[, 'NR.min'], na.rm = TRUE),
sd(Bp.var[, 'NR.fs'], na.rm = TRUE),
sd(Bp.var[, 'NR.ffopt'], na.rm = TRUE),
sd(Bp.var[, 'NR.act'], na.rm = TRUE),
sd(Bp.var[, 'NR.opp'], na.rm = TRUE))
TF3$SDNR <- ifelse(is.na(TF3$SDNR), 0, TF3$SDNR)
TF3$max <- TF3$MeanNR + (1.96 * TF3$SDNR / 10)
TF3$min <- TF3$MeanNR - (1.96 * TF3$SDNR / 10)
if (any(any(na.omit(TF3$min) < 0) | any(na.omit(TF3$max) > 0))) { # if means above and below zero
MIN <- DescTools::RoundTo((min(TF3$MeanNR, na.rm = T) -
(1.96 * min(TF3$SDNR, na.rm = T) / 10)),
5, floor)
MAX <- DescTools::RoundTo((max(TF3$MeanNR, na.rm = T) +
(1.96 * max(TF3$SDNR, na.rm = T) / 10)),
5, ceiling)
STEP <- (MAX - MIN) / 10
}
if (all(na.omit(TF3$min) < 0)) { # if means below zero
MIN <- DescTools::RoundTo((min(TF3$MeanNR, na.rm = T) -
(1.96 * min(TF3$SDNR, na.rm = T) / 10)),
5, floor)
MAX <- 0
STEP <- (MAX - MIN) / 10
}
if (all(na.omit(TF3$min > 0))) { # if means above zero
MIN <- 0
MAX <- DescTools::RoundTo((max(TF3$MeanNR, na.rm = T) +
(1.96 * max(TF3$SDNR, na.rm = T) / 10)),
5, ceiling)
STEP <- (MAX - MIN) / 10
}
mgmt_labels <- c("Min. Rate",
"FS",
"SS.Opt",
"FF.Opt",
"Actual",
"Alt. Price")
p <-
ggplot2::ggplot(TF3, ggplot2::aes(x = Method, y = MeanNR)) +
ggplot2::geom_bar(stat = "identity", fill = "green3") +
ggplot2::geom_errorbar(ggplot2::aes(ymin = MeanNR - (1.96 * SDNR / 10),
ymax = MeanNR + (1.96 * SDNR / 10)),
width = .2) + #
ggplot2::scale_y_continuous(name = paste0("Average net return (", nr_lab, ")"),
limits = c(MIN, MAX),
breaks = seq(MIN, MAX, STEP)) +
ggplot2::scale_x_discrete(name = "Management Strategy",
labels = mgmt_labels) +
ggplot2::geom_text(data = TF3,
ggplot2::aes(x = Method, y = MeanNR, label = round(MeanNR, 2)),
hjust = 1.15,
vjust = 1.25) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/NR/NRbarPlots/",
fieldname, "_avgNR_bar_",
fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Create a bar plot of the total experimental input applied with
#' each strategy. The amount of the site-specific input applied
#' is averaged across the simulation for each point and summed to
#' get a total amount of experimental input applied with the site-
#' specific method. The total amount of experimental input applied
#' with all other strategies is calculated the same way.
#' @param TF4 Data.frame with columns for the amount of experimental rates applied
#' for each strategy.
#' @param expvar Experimental variable optimized, select/input
#' 'As-Applied Nitrogen' or 'As-Applied Seed Rate'. This is the type of
#' input that was experimentally varied across the field as part of the
#' on-farm experimentation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/EXP/EXPapplied' if selected.
plotTotExpAppl = function(TF4,
expvar = self$expvar,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
exp_lab <- "kg"
} else {
exp_lab <- "lbs"
}
if (TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] == 0) {
txt <- "FF.Opt Rate = 0"
} else {
if (TF4[which(TF4$Method == "EXP.ssopt"),"EXP"] >
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"]) {
txt <- paste0(as.integer(
(TF4[which(TF4$Method == "EXP.ssopt"), "EXP"] -
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"]) /
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] * 100),
"% more ",
ifelse(expvar == "aa_n", "N", "Seed"), " used in SS.Opt than FF.Opt")
}else{
txt <- paste0(as.integer(
(TF4[which(TF4$Method == "EXP.ssopt"), "EXP"] -
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"]) /
TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] * -100),
"% less ",
ifelse(expvar == "aa_n", "N", "Seed"), " used in SS.Opt than FF.Opt")
}
}
y_lab <- paste0(ifelse(expvar == "aa_n", "Nitrogen", "Seed"),
" used on field (", exp_lab, "/field)")
yMIN <- 0
yMAX <- DescTools::RoundTo(max(TF4$EXP, na.rm = TRUE), 5, ceiling)
ySTEP <- (yMAX - yMIN) / 10
gg_title <- paste0(ifelse(expvar == "aa_n", "Nitrogen", "Seed"),
" used with each application strategy")
mgmt_labels <- c("Min. Rate", "FS", "SS.Opt", "FF.Opt")
p <-
ggplot2::ggplot(TF4, ggplot2::aes(x = Method, y = EXP)) +
ggplot2::geom_bar(stat = "identity", fill = "blue") +
ggplot2::scale_y_continuous(name = y_lab,
limits=c(yMIN, yMAX),
breaks=seq(yMIN, yMAX, ySTEP)) +
ggplot2::scale_x_discrete(name = "Management Strategy",
labels = mgmt_labels) +
ggplot2::geom_text(data = TF4,
ggplot2::aes(x = Method, y = EXP, label = round(EXP, 2)),
vjust = -.25) +
ggplot2::theme_bw() +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14)) +
ggplot2::ggtitle(gg_title, subtitle = txt)
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/EXP/EXPapplied/",
fieldname, "_tot",
ifelse(expvar == "aa_n", "N", "SR"),
"applied_bar_", fxn, "_", sim_year, "_",
opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Generate a histogram of site specific experimental inptut rates within
#' a field. Uses the average site-specific rate for each point across the
#' simulation, indicating the mean optimum rate at each point across all
#' economic conditions.
#' @param NRopt Data frame with the net-returns and experimental optimums
#' for every point for every simulation iteration.
#' @param TF4 Data.frame with columns for the amount of experimental rates applied
#' for each strategy.
#' @param expvar Experimental variable optimized, select/input
#' 'As-Applied Nitrogen' or 'As-Applied Seed Rate'. This is the type of
#' input that was experimentally varied across the field as part of the
#' on-farm experimentation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/EXP/ssoptEXP' if selected.
plotSSOPThist = function(NRopt,
TF4,
expvar = self$expvar,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
exp_lab <- "kg/ha"
} else {
exp_lab <- "lbs/acre"
}
NRopt <- as.data.frame(NRopt)
gg_title <- paste0("Site-specific profit maximizing ",
ifelse(expvar=="aa_n", "nitrogen", "seed"), " rates")
sub_title <- paste0("Mean SS.Opt ",
ifelse(expvar == "aa_n", "N", "Seed"),
" Rate = ",
round(mean(NRopt$EXP.rate.ssopt, na.rm = TRUE), 2),
" ", exp_lab, ", SD = ",
round(sd(NRopt$EXP.rate.ssopt, na.rm = TRUE), 2),
" ", exp_lab)
if (TF4[which(TF4$Method == "EXP.ssopt"), "EXP"] == 0) {
p <-
ggplot2::ggplot(NRopt) +
ggplot2::geom_histogram(ggplot2::aes(x = EXP.rate.ssopt + 1),
bins = 1,
col = "white",
fill = "blue",
na.rm = TRUE) +
ggplot2::labs(x = exp_lab) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle = sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
} else {
bin_num <- DescTools::RoundTo(max(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 5, ceiling) / 5
xMIN <- DescTools::RoundTo(min(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 5, floor)
xMAX <- DescTools::RoundTo(max(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 5, ceiling)
xSTEP <- (xMAX - xMIN) / 10
p <-
ggplot2::ggplot(NRopt) +
ggplot2::geom_histogram(ggplot2::aes(x = EXP.rate.ssopt + 1),
bins = bin_num,
col = "white",
fill = "blue",
na.rm = TRUE) +
ggplot2::scale_x_continuous(name = exp_lab,
limits = c(xMIN, xMAX),
breaks = seq(xMIN, xMAX, xSTEP),
oob = scales::squish) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle = sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
}
yMIN <- 0
yMAX <- DescTools::RoundTo(max(ggplot2::ggplot_build(p)$data[[1]]$count), 5, ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <- p +
ggplot2::geom_vline(ggplot2::aes(xintercept = round(mean(NRopt$EXP.rate.ssopt,
na.rm = TRUE), 0)),
col = "red",
linetype = 2) +
ggplot2::scale_y_continuous(name = "Frequency",
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP))
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/EXP/ssoptEXP/",
fieldname, "_ssopt",
ifelse(expvar == "aa_n", "N", "SR"),
"_hist_", fxn, "_", sim_year, "_",
opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' Generate a histogram of full-field experimental input rates identified
#' in the simulation. Shows the distribution of full-field optimum rates
#' across all economic scenarios in the simulation.
#' @param Bp.var Data.table containing the average net-returns from the field for each
#' management type for every iteration of the simulation.
#' @param TF4 Data.frame with columns for the amount of experimental rates applied
#' for each strategy.
#' @param expvar Experimental variable optimized, select/input
#' 'As-Applied Nitrogen' or 'As-Applied Seed Rate'. This is the type of
#' input that was experimentally varied across the field as part of the
#' on-farm experimentation.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param fxn The functional form of the models used for analysis.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param opt The optimization method used in the simulation.
#' @param SAVE Logical, whether to save figure.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @return A 'ggplot' object and data saved in 'Outputs/EXP/ffoptEXP' if selected.
plotFFOPThist = function(Bp.var,
TF4,
expvar = self$expvar,
fieldname = self$unique_fieldname,
fxn = self$unique_fxn,
sim_year,
opt = self$opt,
SAVE = self$SAVE,
out_path = self$out_path) {
if (self$SI) {
exp_lab <- "kg/ha"
} else {
exp_lab <- "lbs/acre"
}
Bp.var <- as.data.frame(Bp.var)
x_lab <- paste0("Profit Maximizing Top-Dress ",
ifelse(expvar == "aa_n", "N", "Seed"),
" Rate (", exp_lab, ")")
gg_title <- paste0("Variation in best full-field uniform ",
ifelse(expvar == "aa_n", "N", "seed"),
" rate over years")
sub_title <- paste0("Mean FF.Opt ",
ifelse(expvar == "aa_n", "N", "Seed"),
" Rate = ",
round(mean(Bp.var[, "ffopt.EXPrate"], na.rm = TRUE), 2),
" ", exp_lab, ", SD = ",
round(sd(Bp.var[, "ffopt.EXPrate"], na.rm = TRUE), 2),
" ", exp_lab)
if (TF4[which(TF4$Method == "EXP.ffopt"), "EXP"] == 0) {
p <-
ggplot2::ggplot(Bp.var) +
ggplot2::geom_histogram(ggplot2::aes(x = ffopt.EXPrate + 1),
bins=1,
col="white",
fill="blue",
na.rm=TRUE) +
ggplot2::labs(x=x_lab) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle=sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
} else {
bin_num <- DescTools::RoundTo(max(Bp.var$ffopt.EXPrate,
na.rm = TRUE), 5, ceiling) / 5
xMIN <- DescTools::RoundTo(min(Bp.var$ffopt.EXPrate,
na.rm = TRUE), 5, floor)
xMAX <- DescTools::RoundTo(max(Bp.var$ffopt.EXPrate,
na.rm = TRUE), 5, ceiling)
xSTEP <- (xMAX - xMIN) / 10
p <-
ggplot2::ggplot(Bp.var) +
ggplot2::geom_histogram(ggplot2::aes(x = ffopt.EXPrate),
bins = bin_num,
col = "white",
fill = "blue",
na.rm = TRUE) +
ggplot2::scale_x_continuous(name = x_lab,
limits = c(xMIN, xMAX),
breaks = seq(xMIN, xMAX, xSTEP),
oob = scales::squish) +
ggplot2::theme_bw() +
ggplot2::ggtitle(gg_title, subtitle = sub_title) +
ggplot2::theme(axis.text = ggplot2::element_text(size = 12),
axis.title = ggplot2::element_text(size = 14))
}
yMIN <- DescTools::RoundTo(min(ggplot2::ggplot_build(p)$data[[1]]$count), 5, floor)
yMAX <- DescTools::RoundTo(max(ggplot2::ggplot_build(p)$data[[1]]$count), 5, ceiling)
ySTEP <- (yMAX - yMIN) / 10
p <- p +
ggplot2::geom_vline(ggplot2::aes(xintercept = round(mean(Bp.var[, "ffopt.EXPrate"],
na.rm = TRUE), 0)),
col = "red",
linetype = 2) +
ggplot2::scale_y_continuous(name = "Frequency",
limits = c(yMIN, yMAX),
breaks = seq(yMIN, yMAX, ySTEP))
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(paste0(out_path, "/Outputs/EXP/ffoptEXP/",
fieldname, "_ffopt",
ifelse(expvar == "aa_n", "N", "SR"),
"_hist_", fxn, "_", sim_year, "_",
opt, ".png"),
plot = p, device = "png", scale = 1,
width = 7.5, height = 7.5, units = "in")
}
return(p)
},
#' @description
#' This method is for plotting maps of simulation outcomes. These
#' include the net-returns from the management strategies, predicted
#' responses, and optimized rates. Is a wrapper to OFPE::plotMaps().
#' @param dat Data frame with the net-returns and experimental optimums
#' for every point for every simulation iteration.
#' @param var The label of the variable to map. Used in figure name.
#' @param var_col_name The name of the column of the variable in the
#' supplied data ('dat').
#' @param var_label The label to be applied to the legend of the map
#' corresponding to the variable mapped.
#' @param var_main_label The main label to apply to the map.
#' @param fieldname Unique field name corresponding to all fields used in the simulation.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used.
#' @param fxn The functional form of the models used for analysis.
#' @param opt The optimization method used in the simulation. For labeling.
#' @param SAVE Logical, whether to save figure.
#' @param farmername The name of the farmer that manages the field.
#' @param out_path The path to the folder in which to store and
#' save outputs from the simulation.
#' @param db Connection to the OFPE database to identify UTM zone.
#' @param utm_fieldname Name of the field for identifying the UTM zone.
#' @return A 'ggmap' object and maps saved in 'Outputs/Maps' folder if selected.
plotSimMaps = function(dat,
var,
var_col_name,
var_label,
var_main_label,
fieldname = self$unique_fieldname,
sim_year,
fxn = self$unique_fxn,
opt = self$opt,
SAVE = self$SAVE,
farmername = self$farmername,
out_path = self$out_path,
db = self$db,
utm_fieldname = self$utm_fieldname) {
utm_zone <- OFPE::findUTMzone(db,
fieldname = utm_fieldname)
p <- OFPE::plotMaps(dat,
var_col_name,
var_label,
var_main_label,
fieldname,
farmername,
utm_zone) %>%
suppressMessages() %>%
suppressWarnings()
if (SAVE) {
try({dev.off()}, silent = TRUE)
ggplot2::ggsave(
file = paste0(out_path, "/Outputs/Maps/", sim_year, "/",
fieldname, "_", tolower(var),
"_map_", fxn, "_", sim_year, "_", opt, ".png"),
plot = p, device = "png",
width = 7.5, height = 7.5, units = "in"
)
}
return(p)
},
#' @description
#' Gets the Bp.var output table from the 'SimClass' simulation. The Bp.var table contains means
#' across each field for each iteration of the simulation. This records base
#' price received, the cost of the input, the mean net-return per acre for
#' each management strategy, and the full-field optimum rate. Arguments are
#' used to fill in the file name in the same format as how the file was
#' saved from the simulation.
#' @param dat_path The path to the folder where the data is saved.
#' @param unique_fieldname Unique field name corresponding to the fieldname
#' in the file name.
#' @param unique_fxn The functional form of the models used for analysis
#' found in the file name.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used. Found in the file name.
#' @param opt The optimization method used in the simulation. Found in file
#' name.
#' @return Data frame with 10 columns and rows equal to the number of
#' iterations in the simulation.
loadBpVar = function(dat_path = self$dat_path,
unique_fieldname = self$unique_fieldname,
unique_fxn = self$unique_fxn,
sim_year,
opt = self$opt) {
Bp.var <- data.table::fread(
paste0(dat_path, sim_year, "/",
unique_fieldname, "_BpVar_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv")
)
return(Bp.var)
},
#' @description
#' Gets the NRffmax table that is an output from the 'SimClass' simulation with rows for
#' each iteration and the full-field optimum experimental input rate and the
#' total net-return from the full-field. Arguments are used to fill in the file
#' name in the same format as how the file was saved from the simulation.
#' @param dat_path The path to the folder where the data is saved.
#' @param unique_fieldname Unique field name corresponding to the fieldname
#' in the file name.
#' @param unique_fxn The functional form of the models used for analysis
#' found in the file name.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used. Found in the file name.
#' @param opt The optimization method used in the simulation. Found in file
#' name.
#' @return Data frame with 3 columns and rows equal to the number of
#' iterations in the simulation.
loadNRffmax = function(dat_path = self$dat_path,
unique_fieldname = self$unique_fieldname,
unique_fxn = self$unique_fxn,
sim_year,
opt = self$opt) {
NRffmax <- data.table::fread(
paste0(dat_path, sim_year, "/",
unique_fieldname, "_NRffMaxData_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv")
)
return(NRffmax)
},
#' @description
#' Create the table TF4, which is used for the output functions 'plotFFOPThist',
#' 'plotTotExpAppl', and 'plotSSOPThist'. This contains mean values for each
#' of the strategies to create bar charts of experimental inputs used.
#' @param NRopt 'SimClass' output table that contains every point with
#' the net-return and experimenal rates for each management strategy, and
#' the estimated crop responses for each location for the simulation year that
#' are the means across all simulations.
#' @param NRffmax 'SimClass' output table that contains the optimum full-field
#' experimental rate, the total net-return for each iteration of the simulation.
#' @param fieldsize Size of the field in acres.
#' @param fs The farmer selected experimental rate used in the simulation.
#' @return Data frame with SSOPT and FFOPT experimental input results for
#' each management strategy.
makeTF4 = function(NRopt,
NRffmax,
fieldsize = self$fieldsize,
fs = self$fs) {
TF4 <- data.frame(EXP.ssopt = (sum(NRopt[,'EXP.rate.ssopt'], na.rm = T) *
fieldsize) / nrow(NRopt),
EXP.min = self$AAmin * fieldsize,
EXP.fs = fs * fieldsize,
EXP.ffopt = mean(NRffmax$EXP.rate, na.rm = T) * fieldsize)
TF4 <- tidyr::pivot_longer(TF4,
dplyr::everything(),
names_to = "Method",
values_to = "EXP")
TF4$Method <- factor(TF4$Method)
mgmt_order <- c(grep("EXP.min", levels(TF4$Method)),
grep("EXP.fs", levels(TF4$Method)),
grep("EXP.ssopt", levels(TF4$Method)),
grep("EXP.ffopt", levels(TF4$Method)))
TF4$Method <- factor(TF4$Method, levels(TF4$Method)[mgmt_order])
return(TF4)
},
#' @description
#' Gets the NRopt table that is an output from the 'SimClass' simulation with rows for
#' each in the field. Columns contain the means across all iterations of the
#' simulation. Columns in this file are the price received, cost of nitrogen,
#' SSOPT and FFOPT experimental rates and the estimated net-return, yield, and
#' protein for each management strategy at each location in the field.
#' Arguments are used to fill in the file name in the same format as how the file was saved from the simulation.
#' @param dat_path The path to the folder where the data is saved.
#' @param unique_fieldname Unique field name corresponding to the fieldname
#' in the file name.
#' @param unique_fxn The functional form of the models used for analysis
#' found in the file name.
#' @param sim_year Year that the simulation was performed for. Indicates the
#' weather conditions used. Found in the file name.
#' @param opt The optimization method used in the simulation. Found in file
#' name.
#' @return Data frame with 22 columns and rows equal to the number of
#' observations in the dataset.
loadNRopt = function(dat_path = self$dat_path,
unique_fieldname = self$unique_fieldname,
unique_fxn = self$unique_fxn,
sim_year,
opt = self$opt) {
on.exit(closeAllConnections())
## get col names
NRopt_names <- read.csv(file = paste0(dat_path, sim_year, "/",
unique_fieldname, "_NRopt_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv"),
nrows = 1,
header = TRUE) %>%
names()
## average all cols by x and y
NRopt_names <- NRopt_names[!1:length(NRopt_names) %in%
grep("^x$|^y$|^field$", NRopt_names)]
avg_query <- paste0("AVG([", NRopt_names, "])") %>%
paste(collapse = ", ")
NRopt <- sqldf::read.csv.sql(
paste0(dat_path, sim_year, "/",
unique_fieldname, "_NRopt_",
unique_fxn, "_",
sim_year, "_",
opt, ".csv"),
sql = paste0("SELECT x, y, field, ", avg_query,
"FROM file GROUP BY x, y")
) %>%
`names<-`(c("x", "y", "field", NRopt_names))
NRcols <- grep("NR", names(NRopt))
for (i in 1:length(NRcols)) {
NRopt <- NRopt[NRopt[, NRcols[i]] > -50000 &
NRopt[, NRcols[i]] < 50000, ]
}
return(NRopt)
}
),
private = list(
.setupOP = function() {
cwd <- paste0(self$out_path, "/Outputs") # outputs working directory
if (!file.exists(cwd)) {
dir.create(cwd)
#dir.create(paste0(cwd, "/", "Exploratory"))
dir.create(paste0(cwd, "/", "Maps"))
for (i in 1:length(self$sim_years)) {
dir.create(paste0(cwd, "/", "Maps/", self$sim_years[i]))
}
dir.create(paste0(cwd, "/", "EXP"))
dir.create(paste0(cwd, "/", "EXP/ffoptEXP"))
dir.create(paste0(cwd, "/", "EXP/ssoptEXP"))
dir.create(paste0(cwd, "/", "EXP/EXPapplied"))
dir.create(paste0(cwd, "/", "NR"))
dir.create(paste0(cwd, "/", "NR/NRbarPlots"))
dir.create(paste0(cwd, "/", "NR/NRboxPlots"))
dir.create(paste0(cwd, "/", "NR/NRprobabilities"))
} else {
# if (!file.exists(paste0(cwd, "/", "Exploratory"))) {
# dir.create(paste0(cwd, "/", "Exploratory"))
# }
if (!file.exists(paste0(cwd, "/", "Maps"))) {
dir.create(paste0(cwd, "/", "Maps"))
for (i in 1:length(self$sim_years)) {
dir.create(paste0(cwd, "/", "Maps/", self$sim_years[i]))
}
} else {
for (i in 1:length(self$sim_years)) {
if (!file.exists(paste0(cwd, "/", "Maps/", self$sim_years[i]))) {
dir.create(paste0(cwd, "/", "Maps/", self$sim_years[i]))
}
}
}
if (!file.exists(paste0(cwd, "/", "EXP"))) {
dir.create(paste0(cwd, "/", "EXP"))
}
if (!file.exists(paste0(cwd, "/", "EXP/ffoptEXP"))) {
dir.create(paste0(cwd, "/", "EXP/ffoptEXP"))
}
if (!file.exists(paste0(cwd, "/", "EXP/ssoptEXP"))) {
dir.create(paste0(cwd, "/", "EXP/ssoptEXP"))
}
if (!file.exists(paste0(cwd, "/", "EXP/EXPapplied"))) {
dir.create(paste0(cwd, "/", "EXP/EXPapplied"))
}
if (!file.exists(paste0(cwd, "/", "NR"))) {
dir.create(paste0(cwd, "/", "NR"))
}
if (!file.exists(paste0(cwd, "/", "NR/NRbarPlots"))) {
dir.create(paste0(cwd, "/", "NR/NRbarPlots"))
}
if (!file.exists(paste0(cwd, "/", "NR/NRboxPlots"))) {
dir.create(paste0(cwd, "/", "NR/NRboxPlots"))
}
if (!file.exists(paste0(cwd, "/", "NR/NRprobabilities"))) {
dir.create(paste0(cwd, "/", "NR/NRprobabilities"))
}
}
}
)
)
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.