Nothing
#' Fit beta distribution for one expert
#'
#' @description
#' Fit beta distribution to weights elicited from one expert via the roulette method.
#'
#' @param df A dataframe generated by \code{get_model_input_1exp}.
#'
#' @details This function is based on \code{SHELF::fitdist} and yields identical results.
#'
#' @return Parameters of beta fit.
#'
#' @examples
#' chips <- c(0, 2, 3, 2, 1, 1, 1, 0, 0, 0)
#' x <- get_cum_probs_1exp(chips)
#' y <- get_model_input_1exp(x)
#' fit_beta_1exp(df = y)$par
#'
fit_beta_1exp <- function (df)
{
w <- df$w
cum_probs <- df$cum_probs
### exclude cum_probs that are 0 or 1
inc <- (cum_probs>0) & (cum_probs<1)
min_cum_probs <- min(cum_probs[inc])
max_cum_probs <- max(cum_probs[inc])
min_weight <- min(w[inc])
max_weight <- max(w[inc])
###
min_q <- qnorm(min_cum_probs)
max_q <- qnorm(max_cum_probs)
m <- (min_weight * max_q - max_weight * min_q) / (max_q - min_q)
v <- ((max_weight - min_weight)/(max_q - min_q))^2
### set starting values
alpha <- abs(m^3/v * (1/m - 1) - m)
beta <- abs(alpha/m - alpha)
if(identical(cum_probs[inc], w[inc])){
alpha <- beta <- 1
}
### function to be optimized
beta_error <- function(params, weight, probs) {
sum((pbeta(q = weight, shape1 = exp(params[1]), shape2 = exp(params[2])) -
probs)^2)
}
### optimization
beta_fit <- optim(par = c(log(alpha), log(beta)), fn = beta_error,
weight = w[inc], probs = cum_probs[inc])
if (beta_fit$convergence != 0) {
warning("Algorithm did not converge.")
}
beta_fit$par <- exp(beta_fit$par)
names(beta_fit$par) <- c("alpha", "beta")
return(beta_fit)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.