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
garch-Hessian.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
################################################################################
# FUNCTION: PARAMETER ESTIMATION:
# .garchRoptimhess Uses R internal optimhess function
# .garchRCDAHessian Computes R coded CDA Hessian matrix
# .garchRTSHessian Computes R coded Two Sided Hessian matrix
# .hessian2sided Function called from .garchRTSHessian
################################################################################
.garchRoptimhess <-
function(par, .params, .series, eps = 1.0e-4)
{
# A function implemeted by Diethelm Wuertz
# Description:
# Arguments:
# par -
# .params -
# .series -
# eps -
# FUNCTION:
# Take Time:
.StartHessian <- Sys.time()
# Control Parameters:
con <- list(trace = 0, fnscale = 1, parscale = rep.int(1,
length(par)), ndeps = rep.int(0.001, length(par)), maxit = 100L,
abstol = -Inf, reltol = sqrt(.Machine$double.eps), alpha = 1,
beta = 0.5, gamma = 2, REPORT = 10, type = 1, lmm = 5,
factr = 1e+07, pgtol = 0, tmax = 10, temp = 10)
# Compute Hessian:
H <- optimHess(par, .garchLLH, NULL, con)
H <- 0.5 * (H + t(H))
nm <- names(par)
dimnames(H) <- list(nm, nm)
# Elapsed Time:
time = Sys.time() - .StartHessian
attr(H, "time") = time
# Return Value:
H
}
# ------------------------------------------------------------------------------
.garchRCDAHessian <-
function(par, .params, .series, eps = 1.0e-4)
{
# A function implemented by Diethelm Wuertz
# Description:
# Compute CDA (central difference approximated) Hessian
# Arguments:
# par -
# .params -
# .series -
# eps -
# Reference:
# http://status.sph.umich.edu/computing/manuals/sas8/stat/chap46/sect26.htm
# FUNCTION:
# Take start time
.StartHessian <- Sys.time()
# Algorithm:
algorithm = .params$control$algorithm[1]
.trace = FALSE
# Compute Hessian:
eps = eps * par
n = length(par)
H = matrix(0, ncol = n, nrow = n)
for (i in 1:n) {
for (j in 1:n) {
x1 = x2 = x3 = x4 = par
x1[i] = x1[i] + eps[i]
x1[j] = x1[j] + eps[j]
x2[i] = x2[i] + eps[i]
x2[j] = x2[j] - eps[j]
x3[i] = x3[i] - eps[i]
x3[j] = x3[j] + eps[j]
x4[i] = x4[i] - eps[i]
x4[j] = x4[j] - eps[j]
H[i, j] = (
.garchLLH(x1, .trace) -
.garchLLH(x2, .trace) -
.garchLLH(x3, .trace) +
.garchLLH(x4, .trace) ) / (4*eps[i]*eps[j])
}
}
colnames(H) = rownames(H) = names(par)
# Attribute execution time:
time = Sys.time() - .StartHessian
attr(H, "time") = time
# Return Value:
H
}
# ------------------------------------------------------------------------------
.garchTSHessian <-
function(par, .params, .series, eps = NA)
{
# A function implemented by Diethelm Wuertz
# Description:
# Compute two sided (TS) approximated Hessian
# Arguments:
# par -
# .params -
# .series -
# eps - not used
# FUNCTION:
# Take start time
.StartHessian <- Sys.time()
# Algorithm:
algorithm = .params$control$algorithm[1]
# Compute Hessian:
H <- .hessian2sided(f = .garchLLH, x = par, trace = FALSE, fGarchEnv = FALSE)
colnames(H) = rownames(H) = names(par)
# Attribute execution time:
time = Sys.time() - .StartHessian
attr(H, "time") = time
# Return Value:
H
}
# ------------------------------------------------------------------------------
.hessian2sided <-
function(f, x, ...)
{
# A function adapted from Kevin Sheppard's Matlab garch toolbox
# ... implemented by Alexios Ghalanos in his rgarch package
# ... R port for Rmetrics' fGarch by Diethelm Wuertz
# Description:
# Computes two sided (TS) approximated Hessian
# Arguments:
# f -
# x -
# Notes:
# requires package Matrix (added as suggestion)
# FUNCTION:
# Settings:
n = length(x)
fx <- f(x, ...)
eps = .Machine$double.eps
# Compute the stepsize (h)
h = eps^(1/3) *
apply( as.data.frame(x), 1, FUN = function(z) max(abs(z), 1.0e-2))
xh = x + h
h = xh - x
ee = Matrix(diag(h), sparse = TRUE)
# Compute forward and backward steps:
gp = vector(mode = "numeric", length = n)
for(i in 1:n) gp[i] <- f(x + ee[, i], ...)
gm = vector(mode = "numeric", length = n)
for(i in 1:n) gm[i] <- f(x - ee[, i], ...)
H = h %*% t(h)
Hm = H
Hp = H
# Compute "double" forward and backward steps:
for(i in 1:n){
for(j in i:n){
Hp[i, j] <- f(x + ee[, i] + ee[, j], ...)
Hp[j, i] <- Hp[i, j]
Hm[i, j] <- f(x - ee[, i] - ee[, j], ...)
Hm[j, i] <- Hm[i, j]
}
}
# Compute the hessian:
for(i in 1:n){
for(j in i:n){
H[i, j] = ((Hp[i, j] - gp[i] - gp[j] + fx + fx - gm[i] - gm[j] +
Hm[i, j]) / H[i, j]) / 2
H[j, i] = H[i, j]
}
}
# Return Value:
return(H)
}
################################################################################
![swh spinner](/static/img/swh-spinner.gif)
Computing file changes ...