https://github.com/cran/fGarch
Revision e48fa5936fa9f105a043cce59083c44659f1d3b4 authored by Rmetrics Core Team on 09 November 2009, 00:00:00 UTC, committed by Gabor Csardi on 09 November 2009, 00:00:00 UTC
1 parent 9bb4ec8
Tip revision: e48fa5936fa9f105a043cce59083c44659f1d3b4 authored by Rmetrics Core Team on 09 November 2009, 00:00:00 UTC
version 2110.80.1
version 2110.80.1
Tip revision: e48fa59
dist-sged.R
# This library is free software; you can redistribute it and/or
# modify it under the terms of the GNU Library General Public
# License as published by the Free Software Foundation; either
# version 2 of the License, or (at your option) any later version.
#
# This library is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Library General Public License for more details.
#
# You should have received a copy of the GNU Library General
# Public License along with this library; if not, write to the
# Free Foundation, Inc., 59 Temple Place, Suite 330, Boston,
# MA 02111-1307 USA
# Copyrights (C)
# for this R-port:
# 1999 - 2008, Diethelm Wuertz, Rmetrics Foundation, GPL
# Diethelm Wuertz <wuertz@itp.phys.ethz.ch>
# info@rmetrics.org
# www.rmetrics.org
# for the code accessed (or partly included) from other R-ports:
# see R's copyright and license files
# for the code accessed (or partly included) from contributed R-ports
# and other sources
# see Rmetrics's copyright file
################################################################################
# FUNCTION: GED DISTRIBUTION:
# dged Density for the Generalized Error Distribution
# pged Probability function for the GED
# qged Quantile function for the GED
# rged Random Number Generator for the GED
# FUNCTION: SKEW GED DISTRIBUTION:
# dsged Density for the skewed GED
# psged Probability function for the skewed GED
# qsged Quantile function for the skewed GED
# rsged Random Number Generator for the skewed GED
# FUNCTION: PARAMETER ESTIMATION:
# gedFit Fit the parameters for a GED distribution
# sgedFit Fit the parameters for a skew GED distribution
# FUNCTION: SLIDER:
# sgedSlider Displays Generalized Error Distribution and RVS
################################################################################
dged <-
function(x, mean = 0, sd = 1, nu = 2)
{
# A function imlemented by Diethelm Wuertz
# Description:
# Compute the density for the
# generalized error distribution.
# FUNCTION:
# Compute Density:
z = (x - mean ) / sd
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
g = nu / ( lambda * (2^(1+1/nu)) * gamma(1/nu) )
result = g * exp (-0.5*(abs(z/lambda))^nu) / sd
# Return Value
result
}
# ------------------------------------------------------------------------------
pged <-
function(q, mean = 0, sd = 1, nu = 2)
{
# A function implemented by Diethelm Wuertz
# Description:
# Compute the probability for the
# generalized error distribution.
# FUNCTION:
# Compute Probability:
q = (q - mean ) / sd
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
g = nu / ( lambda * (2^(1+1/nu)) * gamma(1/nu) )
h = 2^(1/nu) * lambda * g * gamma(1/nu) / nu
s = 0.5 * ( abs(q) / lambda )^nu
result = 0.5 + sign(q) * h * pgamma(s, 1/nu)
# Return Value:
result
}
# ------------------------------------------------------------------------------
qged <-
function(p, mean = 0, sd = 1, nu = 2)
{
# A function implemented by Diethelm Wuertz
# Description:
# Compute the quantiles for the
# generalized error distribution.
# FUNCTION:
# Compute Quantiles:
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
q = lambda * (2*qgamma((abs(2*p-1)), 1/nu))^(1/nu)
result = q*sign(2*p-1) * sd + mean
# Return Value:
result
}
# ------------------------------------------------------------------------------
rged <-
function(n, mean = 0, sd = 1, nu = 2)
{
# A function implemented by Diethelm Wuertz
# Description:
# Generate GED random deviates. The function uses the
# method based on the transformation of a Gamma random
# variable.
# FUNCTION:
# Generate Random Deviates:
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
# print(lambda)
r = rgamma(n, 1/nu)
z = lambda * (2*r)^(1/nu) * sign(runif(n)-1/2)
result = z * sd + mean
# Return Value:
result
}
# ------------------------------------------------------------------------------
.dsged <-
function(x, nu, xi)
{
# A function implemented by Diethelm Wuertz
# Description:
# Internal Function
# FUNCTION:
# Standardize:
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
g = nu / ( lambda * (2^(1+1/nu)) * gamma(1/nu) )
m1 = 2^(1/nu) * lambda * gamma(2/nu) / gamma(1/nu)
mu = m1*(xi-1/xi)
sigma = sqrt((1-m1^2)*(xi^2+1/xi^2) + 2*m1^2 - 1)
z = x*sigma + mu
# Compute:
Xi = xi^sign(z)
g = 2 / (xi + 1/xi)
Density = g * dged(x = z/Xi, nu=nu)
# Return Value:
Density * sigma
}
dsged <-
function(x, mean = 0, sd = 1, nu = 2, xi = 1.5)
{
# A function implemented by Diethelm Wuertz
# Description:
# Compute the density function of the
# skewed generalized error distribution
# FUNCTION:
# Shift and Scale:
result = .dsged(x = (x-mean)/sd, nu = nu, xi = xi) / sd
# Return Value:
result
}
# ------------------------------------------------------------------------------
.psged <-
function(q, nu, xi)
{
# A function implemented by Diethelm Wuertz
# Description:
# Internal Function
# FUNCTION:
# Standardize:
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
g = nu / ( lambda * (2^(1+1/nu)) * gamma(1/nu) )
m1 = 2^(1/nu) * lambda * gamma(2/nu) / gamma(1/nu)
mu = m1*(xi-1/xi)
sigma = sqrt((1-m1^2)*(xi^2+1/xi^2) + 2*m1^2 - 1)
z = q*sigma + mu
# Compute:
Xi = xi^sign(z)
g = 2 / (xi + 1/xi)
Probability = Heaviside(z) - sign(z) * g * Xi * pged(q = -abs(z)/Xi, nu=nu)
# Return Value:
Probability
}
psged <-
function(q, mean = 0, sd = 1, nu = 2, xi = 1.5)
{
# A function implemented by Diethelm Wuertz
# Description:
# Compute the distribution function of the
# skewed generalized error distribution
# FUNCTION:
# Shift and Scale:
result = .psged(q = (q-mean)/sd, nu = nu, xi = xi)
# Return Value:
result
}
# ------------------------------------------------------------------------------
.qsged <-
function(p, nu, xi)
{
# A function implemented by Diethelm Wuertz
# Description:
# Internal Function
# FUNCTION:
# Standardize:
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
g = nu / ( lambda * (2^(1+1/nu)) * gamma(1/nu) )
m1 = 2^(1/nu) * lambda * gamma(2/nu) / gamma(1/nu)
mu = m1*(xi-1/xi)
sigma = sqrt((1-m1^2)*(xi^2+1/xi^2) + 2*m1^2 - 1)
# Compute:
g = 2 / (xi + 1/xi)
sig = sign(p-1/2)
Xi = xi^sig
p = (Heaviside(p-1/2)-sig*p) / (g*Xi)
Quantile = (-sig*qged(p=p, sd=Xi, nu=nu) - mu ) / sigma
# Return Value:
Quantile
}
qsged <-
function(p, mean = 0, sd = 1, nu = 2, xi = 1.5)
{
# A function implemented by Diethelm Wuertz
# Description:
# Compute the quantile function of the
# skewed generalized error distribution
# FUNCTION:
# Shift and Scale:
result = .qsged(p = p, nu = nu, xi = xi) * sd + mean
# Return Value:
result
}
# ------------------------------------------------------------------------------
.rsged <-
function(n, nu, xi)
{
# A function implemented by Diethelm Wuertz
# Description:
# Internal Function
# FUNCTION:
# Generate Random Deviates:
weight = xi / (xi + 1/xi)
z = runif(n, -weight, 1-weight)
Xi = xi^sign(z)
Random = -abs(rged(n, nu=nu))/Xi * sign(z)
# Scale:
lambda = sqrt ( 2^(-2/nu) * gamma(1/nu) / gamma(3/nu) )
g = nu / ( lambda * (2^(1+1/nu)) * gamma(1/nu) )
m1 = 2^(1/nu) * lambda * gamma(2/nu) / gamma(1/nu)
mu = m1*(xi-1/xi)
sigma = sqrt((1-m1^2)*(xi^2+1/xi^2) + 2*m1^2 - 1)
Random = (Random - mu ) / sigma
# Return value:
Random
}
# ------------------------------------------------------------------------------
rsged <-
function(n, mean = 0, sd = 1, nu = 2, xi = 1.5)
{
# A function implemented by Diethelm Wuertz
# Description:
# Generate random deviates from the
# skewed generalized error distribution
# FUNCTION:
# Shift and Scale:
result = .rsged(n = n, nu = nu, xi = xi) * sd + mean
# Return Value:
result
}
################################################################################
gedFit <-
function(x, ...)
{
# A function implemented by Diethelm Wuertz
# Description:
# Fit the parameters for a skew Normal distribution
# FUNCTION:
# Start Value:
start = c(mean = mean(x), sd = sqrt(var(x)), nu = 2)
# Log-likelihood Function:
loglik = function(x, y = x){
f = -sum(log(dged(y, x[1], x[2], x[3])))
f }
# Minimization:
fit = nlminb(start = start, objective = loglik,
lower = c(-Inf, 0, 0), upper = c(Inf, Inf, Inf), y = x, ...)
# Add Names to $par
names(fit$par) = c("mean", "sd", "nu")
# Return Value:
fit
}
# ------------------------------------------------------------------------------
sgedFit <-
function(x, ...)
{
# A function implemented by Diethelm Wuertz
# Description:
# Fit the parameters for a skew Normal distribution
# FUNCTION:
# Start Value:
start = c(mean = mean(x), sd = sqrt(var(x)), nu = 2, xi = 1)
# Log-likelihood Function:
loglik = function(x, y = x){
f = -sum(log(dsged(y, x[1], x[2], x[3], x[4])))
f }
# Minimization:
fit = nlminb(start = start, objective = loglik,
lower = c(-Inf, 0, 0, 0), upper = c(Inf, Inf, Inf, Inf), y = x, ...)
# Add Names to $par
names(fit$par) = c("mean", "sd", "nu", "xi")
# Return Value:
fit
}
# ------------------------------------------------------------------------------
sgedSlider <-
function(type = c("dist", "rand"))
{
# A function implemented by Diethelm Wuertz
# Description:
# Displays interactively skew GED distribution
# Note:
# dsged(x, mean = 0, sd = 1, nu = 5, xi = 1.5)
# FUNCTION:
# Internal Function:
refresh.code = function(...)
{
# Sliders:
N = .sliderMenu(no = 1)
mean = .sliderMenu(no = 2)
sd = .sliderMenu(no = 3)
nu = .sliderMenu(no = 4)
xi = .sliderMenu(no = 5)
invert = .sliderMenu(no = 6)
# Compute Data:
if (invert == 1) xi = round(1/xi, digits = 4)
xmin = round(qsged(0.01, mean, sd, nu, xi), digits = 2)
xmax = round(qsged(0.99, mean, sd, nu, xi), digits = 2)
s = seq(xmin, xmax, length = N)
y1 = dsged(s, mean, sd, nu, xi)
y2 = psged(s, mean, sd, nu, xi)
main1 = paste("Skew GED Density\n",
"mean = ", as.character(mean), " | ",
"sd = ", as.character(sd), " | ",
"nu = ", as.character(nu), " | ",
"xi = ", as.character(xi) )
main2 = paste("Skew GED Probability\n",
"xmin [0.01] = ", as.character(xmin), " | ",
"xmax [0.99] = ", as.character(xmax) )
# Random Numbers:
if (type[1] == "rand") {
x = rsged(N, mean, sd, nu, xi)
}
# Frame:
par(mfrow = c(2, 1), cex = 0.7)
# Density:
if (type[1] == "rand") {
hist(x, probability = TRUE, col = "steelblue", border = "white",
breaks = "FD",
xlim = c(xmin, xmax), ylim = c(0, 1.1*max(y1)), main = main1 )
lines(s, y1, col = "orange")
} else {
plot(s, y1, type = "l", xlim = c(xmin, xmax), col = "steelblue")
abline (h = 0, lty = 3)
title(main = main1)
grid()
}
# Probability:
plot(s, y2, type = "l", xlim = c(xmin, xmax), ylim = c(0, 1),
col = "steelblue" )
abline (h = 0, lty = 3)
title(main = main2)
grid()
# Frame:
par(mfrow = c(1, 1), cex = 0.7)
}
# Open Slider Menu:
.sliderMenu(refresh.code,
names = c( "N", "mean", "sd", "nu", "xi", "xi.inv"),
minima = c( 10, -5.0, 0.1, 2.1, 1.0, 0 ),
maxima = c( 1000, +5.0, 5.0, 10.0, 10.0, 1 ),
resolutions = c( 10, 0.1, 0.1, 0.1, 0.1, 1 ),
starts = c( 100, 0.0, 1.0, 5.0, 1.0, 0 )
)
}
################################################################################
![swh spinner](/static/img/swh-spinner.gif)
Computing file changes ...