[R-sig-finance] fBonds svensson fitting
Diethelm Wuertz
wuertz at itp.phys.ethz.ch
Thu Feb 2 16:23:18 CET 2006
# Here comes a program to compute the Nelson-Siegel approach for the
# Yield Curve. It can be modified for forward rates (straightforward
# but a little more lengthy), and to other approaches including the Svenson
# Approach and/or several Spline Interpolations.
# The data set used here are an extension of the McCulloch U.S. Treasury
# term structure data (appearing in the Handbook of Monetary Economics)
# to February 1991 as used by S-Plus.
# The results of my program are compared with those obtained from
# SPlus/Finmetrics. They agree to more than 6 digets!
NelsonSiegel =
function(Yield, Maturity)
{ # A function written by Diethelm Wuertz
# Description:
# Fit the Yield Curve by the Nelson-Siegel Method
#
# Details:
# This function finds a global solution. The start values for the
# betas are solved exactly as a function of tau using OLS.
#
# Copyright:
# Diethelm Wuertz, (c) 2004 fBonds
#
# Source:
# Partial copy from 'fBonds' from 'Rmetrics' (unpublished).
# FUNCTION:
# Find Optimal Start Solution by OLS of beta's vs. Yields:
n = length(Maturity)
gmin = 1.0e99
for (i in 1:n) {
tau = Maturity[i]
x = Maturity/tau
a = matrix(rep(NA, times = 9), nrow = 3)
a[1,1] = 1
a[1,2] = a[2,1] = mean((1-exp(-x))/x)
a[1,3] = a[3,1] = mean((1-exp(-x))/x - exp(-x))
a[2,2] = mean( ((1-exp(-x))/x)^2 )
a[2,3] = a[3,2] = mean(((1-exp(-x))/x)*((1-exp(-x))/x-exp(-x)))
a[3,3] = mean(((1-exp(-x))/x - exp(-x))^2)
b = c(
mean ( Yield ),
mean ( Yield * ((1-exp(-x))/x)),
mean ( Yield * (((1-exp(-x))/x - exp(-x)))))
beta = solve(a, b)
yfit = beta[1] + beta[2]*exp(-x) + beta[3]*x*exp(-x)
fmin = sum( (Yield-yfit)^2 )
if (fmin < gmin) {
gmin = fmin
gvec = c(beta, tau)
}
}
# Function to be optimized:
fx <- function(Maturity, x) {
x[1] + x[2] * (1-exp(-Maturity/x[4]))/(Maturity/x[4]) +
x[3] *
((1-exp(-Maturity/x[4]))/(Maturity/x[4])-exp(-Maturity/x[4]))
}
func <- function(x) { sum( (Yield - fx(Maturity, x))^2 ) }
# Optimize:
fit = nlminb(objective = func, start = gvec)
fit$start = gvec
names(fit$par) = c("beta1", "beta2", "beta3", "tau")
# Plot Curve:
yfit = fx(Maturity, gvec)
plot(Maturity, Yield, ylim = c(min(c(Yield, yfit)), max(c(Yield,
yfit))),
pch = 19, cex = 0.5, main = "Nelson-Siegel" )
lines(Maturity, yfit, col = "steelblue")
# Return Value:
fit
}
Yield = c(
0.04984, 0.05283, 0.05549, 0.05777, 0.05961, 0.06102, 0.06216, 0.06314,
0.06403,
0.06488, 0.06568, 0.06644, 0.06717, 0.06786, 0.06852, 0.06913, 0.06969,
0.07020,
0.07134, 0.07205, 0.07339, 0.07500, 0.07710, 0.07860, 0.08011, 0.08114,
0.08194,
0.08274, 0.08355, 0.08434, 0.08512, 0.08588, 0.08662, 0.08731, 0.08794,
0.08851,
0.08900, 0.08939, 0.08967, 0.08980, 0.08976, 0.08954, 0.08910, 0.08843,
0.08748,
0.08626, 0.08474, 0.08291)
Maturity = c(
0.083, 0.167, 0.250, 0.333, 0.417, 0.500, 0.583, 0.667,
0.750, 0.833,
0.917, 1.000, 1.083, 1.167, 1.250, 1.333, 1.417, 1.500,
1.750, 2.000,
2.500, 3.000, 4.000, 5.000, 6.000, 7.000, 8.000, 9.000, 10.000,
11.000,
12.000, 13.000, 14.000, 15.000, 16.000, 17.000, 18.000, 19.000, 20.000,
21.000,
22.000, 23.000, 24.000, 25.000, 26.000, 27.000, 28.000, 29.000)
NelsonSiegel(Yield, Maturity)
# $par
# beta1 beta2 beta3 tau
# 8.885514e-02 -3.671098e-02 1.398381e-11 1.025116e+00
# Compare with S-Plus - (only runs under SPlus):
# ans = term.struct(Yield, Maturity, na.rm = T, method = "ns", input =
"spot")
# c(ans$coefficients, tau = ans$tau)
# b0 b1 b2 tau
# 0.08885513 -0.03671097 5.254415e-008 1.02511
RMSE =
function(Yield, Maturity, coef)
{
sum(Yield - coef[1] +
coef[2]*(1-exp(-Maturity/coef[4]))/(Maturity/coef[4]) +
coef[3] * ( (1-exp(-Maturity/coef[4]))/(Maturity/coef[4]) -
exp(-Maturity/coef[4])) )^2
}
coef = c(8.885514e-02, -3.671098e-02, 1.398381e-11, 1.025116e+00)
RMSE(Yield, Maturity, coef)
# 1.492437
coef = c(0.08885513, -0.03671097, 5.254415e-008, 1.02511)
RMSE(Yield, Maturity, coef)
# 1.492431
DW
Thomas Steiner wrote:
>Dear Kris,
>excellent, thank you.
>
>
>
>>r = nlm(f=nelson,p=c(1,10,10,10),y=fwdcrv,steptol=1e-10,iterlim=500)
>>
>>
>
>Where is the difference to nls? I tried something like
>
>nls(formula=yc$yield~nelsonsiegel(yc$time,init), data = yc, start=init)
>
>where yc is my dataframe with time and forward rates.
>Anyway, even good starting estimates seem always to produce a singular
>gradient matrix.
>
>
>
>>Note svensson is double humped so you have an additional parameter (
>>x[5] )
>>
>>
>
>I think Svensson is 6-dimensional.
>Some interpretations of the parameters of eg Nelson-Siegel can be used
>to get better initial estimates. But this is as well mentioned in the
>paper (pdf) you linked.
>
>
>
>>Hope this helps,
>>
>>
>
>Yes, a lot. BTW: How do you calculate forward rates out of SWAP (Libor) rates?
>
>Thomas
>
>_______________________________________________
>R-sig-finance at stat.math.ethz.ch mailing list
>https://stat.ethz.ch/mailman/listinfo/r-sig-finance
>
>
>
More information about the R-sig-finance
mailing list