LambertW (version 0.6.4)

loglik-LambertW-utils: Log-Likelihood for Lambert W\(\times\) F RVs

Description

Evaluates the log-likelihood for \(\theta\) given observations y.

loglik_LambertW computes the log-likelihood of \(\theta\) for a Lambert W \(\times\) F distribution given observations y.

loglik_input computes the log-likelihood of various distributions for the parameter \(\boldsymbol \beta\) given the data x. This can be used independently of the Lambert W x F framework to compute the log-likelihood of parameters for common distributions.

loglik_penalty computes the penalty for transforming the data back to the input (see Goerg 2016). This penalty is independent of the distribution specified by distname, but only depends on \(\tau\). If type = "s" then the penalty term exists if the distribution is non-negative (see get_distname_family) and gamma >= 0; otherwise, it returns NA.

Usage

loglik_LambertW(theta, y, distname, type, return.negative = FALSE,
  flattened.theta.names = names(theta), use.mean.variance = TRUE)
loglik_input(beta, x, distname, dX = NULL, log.dX = function(x, beta)
  log(dX(x, beta)))
loglik_penalty(tau, y, type = c("h", "hh", "s"), is.non.negative = FALSE)

Arguments

theta

list; a (possibly incomplete) list of parameters alpha, beta, gamma, delta. complete_theta fills in default values for missing entries.

y

a numeric vector of real values (the observed data).

distname

character; name of input distribution; see get_distnames.

type

type of Lambert W \(\times\) F distribution: skewed "s"; heavy-tail "h"; or skewed heavy-tail "hh".

return.negative

logical; if TRUE it returns the negative log-likelihood as a scalar (which is useful for numerical minimization algorithms for maximum likelihood estimation); otherwise it returns a list of input log-likelihood, penalty, and their sum = full likelihood. Default: FALSE.

flattened.theta.names

vector of strings with names of flattened theta; this is necessary for optimization functions since they drop the names of a vector, but all functions in this package use names to select elements of (the flattened) theta.

use.mean.variance

logical; if TRUE it uses mean and variance implied by \(\boldsymbol \beta\) to do the transformation (Goerg 2011). If FALSE, it uses the alternative definition from Goerg (2016) with location and scale parameter.

beta

numeric vector (deprecated); parameter \(\boldsymbol \beta\) of the input distribution. See check_beta on how to specify beta for each distribution.

x

a numeric vector of real values (the input data).

dX

optional; density function of x. Common distributions are already built-in (see distname). If you want to supply your own density, you must supply a function of (x, beta) and set distname = "user".

log.dX

optional; a function that returns the logarithm of the density function of x. Often -- in particular for exponential families -- the \(\log\) of \(f_X(x)\) has a simpler form (and is thus faster to evaluate).

tau

named vector \(\tau\) which defines the variable transformation. Must have at least 'mu_x' and 'sigma_x' element; see complete_tau for details.

is.non.negative

logical; by default it is set to TRUE if the distribution is not a location but a scale family.

Value

loglik_input and loglik_penalty return a scalar; loglik_LambertW returns a list with 3 values:

loglik.input

loglikelihood of beta given the transformed data,

loglik.penalty

penalty for transforming the data,

loglik.LambertW

total log-likelihood of theta given the observed data; if the former two values exist this is simply their sum.

Details

For heavy-tail Lambert W\(\times\) F distributions (type = "h" or type = "hh") the log-likelihood decomposes into an input log-likelihood plus a penalty term for transforming the data.

For skewed Lambert W \(\times\) F distributions this decomposition only exists for non-negative input RVs (e.g., "exp"onential, "gamma", "f", …). If negative values are possible ("normal", "t", "unif", "cauchy", …) then loglik_input and loglik_penalty return NA, but the value of the output log-likelihood will still be returned correctly as loglik.LambertW.

See Goerg (2016) for details on the decomposition of the log-likelihood into a log-likelihood on the input parameters plus a penalty term for transforming the data.

Examples

# NOT RUN {
set.seed(1)
yy <- rLambertW(n = 1000, distname = "normal", 
                theta = list(beta = c(0, 1), delta = 0.2))
loglik_penalty(tau = theta2tau(list(beta = c(1, 1), delta = c(0.2, 0.2)),
                               distname = "normal"), 
               y = yy, type = "hh")
# For a type = 's' Lambert W x F distribution with location family input
# such a decomposition doesn't exist; thus NA.
loglik_penalty(tau = theta2tau(list(beta = c(1, 1), gamma = 0.03), 
                               distname = "normal"),
               is.non.negative = FALSE,
               y = yy, type = "s") 
# For scale-family input it does exist
loglik_penalty(tau = theta2tau(list(beta = 1, gamma = 0.01), 
                               distname = "exp"),
               is.non.negative = TRUE,
               y = yy, type = "s") 
               
# evaluating the Gaussian log-likelihood
loglik_input(beta = c(0, 1), x = yy, distname = "normal") # built-in version
# or pass your own log pdf function
loglik_input(beta = c(0, 1), x = yy, distname = "user", 
             log.dX = function(xx, beta = beta) { 
                dnorm(xx, mean = beta[1], sd = beta[2], log = TRUE)
             })
# }
# NOT RUN {
# you must specify distname = 'user'; otherwise it does not work
loglik_input(beta = c(0, 1), x = yy, distname = "mydist", 
             log.dX = function(xx, beta = beta) { 
                dnorm(xx, mean = beta[1], sd = beta[2], log = TRUE)
                })
# }
# NOT RUN {
### loglik_LambertW returns all three values
loglik_LambertW(theta = list(beta = c(1, 1), delta = c(0.09, 0.07)), 
                y = yy, type = "hh", distname ="normal")

# can also take a flattend vector; must provide names though for delta
loglik_LambertW(theta = flatten_theta(list(beta = c(1, 1), 
                                          delta = c(delta_l = 0.09, 
                                                    delta_r = 0.07))), 
                y = yy, type = "hh", distname ="normal")

# }