Nothing
# ============================== texmex::evmt =============================== #
# Methods for class texmex_evmOpt
#' @export
logLikVec.texmex_evmOpt <- function(object, pars = NULL, ...) {
if (!missing(...)) {
warning("extra arguments discarded")
}
# If the parameter estimates have not been provided in pars then extract
# them from the fitted object
if (is.null(pars)) {
pars <- coef(object)
}
n_pars <- length(pars)
response_data <- object$data$y
# Determine the model and whether or not there are covariates
if (object$family$name == "GEV") {
if (n_pars == 3) {
mu <- pars["mu: (Intercept)"]
phi <- pars["phi: (Intercept)"]
xi <- pars["xi: (Intercept)"]
} else {
mu_mat <- object$data$D$mu
phi_mat <- object$data$D$phi
xi_mat <- object$data$D$xi
n_mu <- ncol(mu_mat)
n_phi <- ncol(phi_mat)
n_xi <- ncol(xi_mat)
mu_pars <- pars[1:n_mu]
phi_pars <- pars[(n_mu + 1):(n_mu + n_phi)]
xi_pars <- pars[(n_mu + n_phi + 1):n_pars]
mu <- as.vector(mu_mat %*% mu_pars)
phi <- as.vector(phi_mat %*% phi_pars)
xi <- as.vector(xi_mat %*% xi_pars)
}
# Calculate the loglikelihood contributions
sigma <- exp(phi)
if (any(sigma <= 0)) {
val <- -Inf
} else {
val <- revdbayes::dgev(response_data, loc = mu, scale = sigma,
shape = xi, log = TRUE)
}
} else if (object$family$name == "GPD") {
if (n_pars == 2) {
phi <- pars["phi: "]
xi <- pars["xi: "]
} else {
phi_mat <- object$data$D$phi
xi_mat <- object$data$D$xi
n_phi <- ncol(phi_mat)
n_xi <- ncol(xi_mat)
phi_pars <- pars[1:n_phi]
xi_pars <- pars[(n_phi + 1):n_pars]
phi <- as.vector(phi_mat %*% phi_pars)
xi <- as.vector(xi_mat %*% xi_pars)
}
sigma <- exp(phi)
# Calculate the loglikelihood contributions
if (any(sigma <= 0)) {
val <- -Inf
} else {
val <- revdbayes::dgp(response_data, loc = object$threshold,
scale = sigma, shape = xi, log = TRUE)
}
} else if (object$family$name == "EGP3") {
if (n_pars == 2) {
lambda <- pars["lambda: "]
phi <- pars["phi: "]
xi <- pars["xi: "]
} else {
lambda_mat <- object$data$D$lambda
phi_mat <- object$data$D$phi
xi_mat <- object$data$D$xi
n_lambda <- ncol(lambda_mat)
n_phi <- ncol(phi_mat)
n_xi <- ncol(xi_mat)
lambda_pars <- pars[1:n_lambda]
phi_pars <- pars[(n_lambda + 1):(n_lambda + n_phi)]
xi_pars <- pars[(n_lambda + n_phi + 1):n_pars]
lambda <- as.vector(lambda_mat %*% lambda_pars)
phi <- as.vector(phi_mat %*% phi_pars)
xi <- as.vector(xi_mat %*% xi_pars)
}
sigma <- exp(phi)
kappa <- exp(lambda)
# Calculate the loglikelihood contributions
if (any(sigma <= 0)) {
val <- -Inf
} else {
val <- texmex::degp3(response_data, kappa = kappa, sigma = sigma,
xi = xi, u = object$threshold, log.d = TRUE)
}
}
# Return the usual attributes for a "logLik" object
attr(val, "nobs") <- nobs(object)
attr(val, "df") <- length(pars)
class(val) <- "logLikVec"
return(val)
}
#' @export
nobs.texmex_evmOpt <- function(object, ...) {
return(length(object$data$y))
}
#' @export
vcov.texmex_evmOpt <- function(object, ...) {
vc <- object$cov
par_names <- names(coef(object))
dimnames(vc) <- list(par_names, par_names)
return(vc)
}
#' @export
coef.texmex_evmOpt <- function(object, ...) {
return(object$coefficients)
}
#' @export
logLik.texmex_evmOpt <- function(object, ...) {
return(logLik(logLikVec(object)))
}
# See texmex_methods.R for nobs and vcov methods for class "evmOpt"
# (texmex already has coef and logLik methods)
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.