#' Fit a single phenology model
#'
#' General purpose model fit function, called by higher order
#' pr_fit_comaprison() and pr_cross_validate() functions. In itself a wrapper
#' around the pr_fit_parameter() function.
#'
#' This function together with pr_predict() is in many ways
#' the core of the package.
#'
#' @param model the model name to be used in optimizing the model
#' @param data dataset generated by the format_phenocam() or
#' format_modis() routines, or adhering to the general model optimization
#' input format.
#' @param method optimization method to use (default = GenSA)
#' - GenSA : Generalized Simulated Annealing algorithm
#' - genoud : GENetic Optimization Using Derivatives
#' - BayesianTools: various bayesian based optimization tools
#' @param par_ranges a vector of starting parameter values (function specific)
#' defaults to the parameter ranges as provided with the current models
#' and set forth by Basler (2016)
#' @param control list of control parameters to be passed to the optimizer
#' @param plot TRUE / FALSE, plot model fit
#' @param ... additional control parameters to be passed
#' @keywords phenology, model, validation
#' @export
#' @examples
#'
#' \dontrun{
#' pr_fit(model,
#' par_ranges = "parameter_ranges.csv")
#' }
#This is a test
pr_fit <- function(
model = "TT",
data = phenor::phenocam_DB,
method = "GenSA",
control = list(max.call = 2000),
par_ranges = system.file(
"extdata",
"parameter_ranges.csv",
package = "phenor",
mustWork = TRUE
),
plot = TRUE,
...
){
# convert to a flat format for speed
data <- pr_flatten(data)
# read parameter ranges
d <- pr_parameters(model = model,
par_ranges = par_ranges)
# optimize paramters
optim_par = pr_fit_parameters(
par = NULL,
data = data,
model = model,
method = method,
lower = as.numeric(d[1,]),
upper = as.numeric(d[2,]),
control = control,
...
)
# estimate model output using optimized
# parameters
out <- pr_predict(
data = data,
model = model,
par = optim_par$par
)
# basic statistics
RMSE_NULL <- sqrt(mean((data$transition_dates - null(data)) ^ 2, na.rm = T))
RMSE <- rmse(par = optim_par$par, data = data, model = model)
Ac <- AICc(measured = data$transition_dates,
predicted = out,
k = length(optim_par$par))
output <- list(
"model" = model,
"par" = optim_par$par,
"measured" = data$transition_dates,
"predicted" = out,
"rmse" = RMSE,
"rmse_null" = RMSE_NULL,
"aic" = Ac,
"opt_out" = optim_par$opt_out
)
# assign class
class(output) <- "pr_fit"
# return optimized parameters and stats
return(output)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.