VGAM (version 1.1-2)

levy: Levy Distribution Family Function

Description

Estimates the scale parameter of the Levy distribution by maximum likelihood estimation.

Usage

levy(location = 0, lscale = "loglink", iscale = NULL)

Arguments

location

Location parameter. Must have a known value. Called \(a\) below.

lscale

Parameter link function for the (positive) scale parameter \(b\). See Links for more choices.

iscale

Initial value for the \(b\) parameter. By default, an initial value is chosen internally.

Value

An object of class "vglmff" (see vglmff-class). The object is used by modelling functions such as vglm, and vgam.

Details

The Levy distribution is one of three stable distributions whose density function has a tractable form. The formula for the density is $$f(y;b) = \sqrt{\frac{b}{2\pi}} \exp \left( \frac{-b}{2(y - a)} \right) / (y - a)^{3/2} $$ where \(a<y<\infty\) and \(b>0\). Note that if \(a\) is very close to min(y) (where y is the response), then numerical problem will occur. The mean does not exist. The median is returned as the fitted values.

References

Nolan, J. P. (2005) Stable Distributions: Models for Heavy Tailed Data.

See Also

The Nolan article was at http://academic2.american.edu/~jpnolan/stable/chap1.pdf.

Examples

# NOT RUN {
nn <- 1000; loc1 <- 0; loc2 <- 10
myscale <- 1  # log link ==> 0 is the answer
ldata <- data.frame(y1 = loc1 + myscale/rnorm(nn)^2,  # Levy(myscale, a)
                    y2 = rlevy(nn, loc = loc2, scale = exp(+2)))

# Cf. Table 1.1 of Nolan for Levy(1,0)
with(ldata, sum(y1 > 1) / length(y1))  # Should be 0.6827
with(ldata, sum(y1 > 2) / length(y1))  # Should be 0.5205

fit1 <- vglm(y1 ~ 1, levy(location = loc1), data = ldata, trace = TRUE)
coef(fit1, matrix = TRUE)
Coef(fit1)
summary(fit1)
head(weights(fit1, type = "work"))

fit2 <- vglm(y2 ~ 1, levy(location = loc2), data = ldata, trace = TRUE)
coef(fit2, matrix = TRUE)
Coef(fit2)
c(median = with(ldata, median(y2)), fitted.median = head(fitted(fit2), 1))
# }