import numpy as np
[docs]def Galfa(M_Obs: np.ndarray, r_Obs: np.ndarray, ufr, alpha, Tau):
"""Calculates the gap at the convergence point
Calculates the gap at the convergence point between the alowable tolerance Tau and the curve extrapolated using the Smith-Wilson algorithm.
interpolation and extrapolation of rates.
Args:
M_Obs: n x 1 ndarray of maturities of bonds, that have rates provided in input (r). Ex. u=[[1], [3]]
r_Obs: n x 1 ndarray of rates, for which you wish to calibrate the algorithm. Each rate belongs to an observable zero coupon bond with a known maturity. Ex. r = [[0.0024], [0.0034]]
ufr: 1 x 1 floating number, representing the ultimate forward rate. Ex. ufr = 0.042
alpha: 1 x 1 floating number representing the convergence speed parameter alpha. Ex. alpha = 0.05
Tau: 1 x 1 floating number representing the allowed difference between ufr and actual curve. Ex. Tau = 0.00001
Returns:
1 x 1 floating number representing the distance between ufr input and the maximum allowed discrepancy Tau
Example:
>>> import numpy as np
>>> M_Obs = np.transpose(np.array([1, 2, 4, 5, 6, 7]))
>>> r_Obs = np.transpose(np.array([0.01, 0.02, 0.03, 0.032, 0.035, 0.04]))
>>> alfa = 0.15
>>> ufr = 0.04
>>> Precision = 0.0000000001
>>> Tau = 0.0001
>>> Galfa(M_Obs, r_Obs, ufr, alfa, Tau)
-8.544212205612438e-05
For more information see https://www.eiopa.europa.eu/sites/default/files/risk_free_interest_rate/12092019-technical_documentation.pdf
Implemented by Gregor Fabjan from Qnity Consultants on 17/12/2021.
"""
from smith_wilson_funcs import SWCalibrate
U = max(M_Obs) # Find maximum liquid maturity from input
T = max(U + 40, 60) # Define the convergence point as defined in paragraph 120 and again in 157
C = np.identity(M_Obs.size) # Construct cash flow matrix described in paragraph 137 assuming the input is ZCB bonds with notional value of 1
d = np.exp(-np.log(1 + ufr) * M_Obs) # Calculate vector d described in paragraph 138
Q = np.diag(d) @ C # Matrix Q described in paragraph 139
b = SWCalibrate(r_Obs, M_Obs, ufr, alpha) # Calculate the calibration vector b using the equation from paragraph 149
K = (1+alpha * M_Obs @ Q@ b) / (np.sinh(alpha * M_Obs.transpose())@ Q@ b) # Calculate kappa as defined in the paragraph 155
return( alpha/np.abs(1 - K*np.exp(alpha*T))-Tau) # Size of the gap at the convergence point between the alowable tolerance Tau and the actual curve. Defined in paragraph 158
[docs]def BisectionAlpha(xStart, xEnd, M_Obs, r_Obs, ufr, Tau, Precision, maxIter):
"""Bisection root finding algorithm for finding the root of a function.
The function here is the allowed difference between the ultimate forward rate and the extrapolated curve using Smith & Wilson.
Args:
xStart: 1 x 1 floating number representing the minimum allowed value of the convergence speed parameter alpha. Ex. alpha = 0.05
xEnd: 1 x 1 floating number representing the maximum allowed value of the convergence speed parameter alpha. Ex. alpha = 0.8
M_Obs: n x 1 ndarray of maturities of bonds, that have rates provided in input (r). Ex. u=[[1], [3]]
r_Obs: n x 1 ndarray of rates, for which you wish to calibrate the algorithm. Each rate belongs to an observable zero coupon bond with a known maturity. Ex. r = [[0.0024], [0.0034]]
ufr: 1 x 1 floating number, representing the ultimate forward rate. Ex. ufr = 0.042
Tau: 1 x 1 floating number representing the allowed difference between ufr and actual curve. Ex. Tau = 0.00001
Precision: 1 x 1 floating number representing the precision of the calculation. Higher the precision, more aqurate the estimation of the root
maxIter: 1 x 1 positive integer representing the maximum number of iterations allowed. This is to prevent an infinite loop in case the method does not converge to a solution
Returns:
1 x 1 floating number representing the optimal value of the parameter alpha
Example:
>>> import numpy as np
>>> M_Obs = np.transpose(np.array([1, 2, 4, 5, 6, 7]))
>>> r_Obs = np.transpose(np.array([0.01, 0.02, 0.03, 0.032, 0.035, 0.04]))
>>> xStart = 0.05
>>> xEnd = 0.5
>>> maxIter = 1000
>>> alfa = 0.15
>>> ufr = 0.042
>>> Precision = 0.0000000001
>>> Tau = 0.0001
>>> BisectionAlpha(xStart, xEnd, M_Obs, r_Obs, ufr, Tau, Precision, maxIter)
0.11549789285636511
For more information see https://www.eiopa.europa.eu/sites/default/files/risk_free_interest_rate/12092019-technical_documentation.pdf and https://en.wikipedia.org/wiki/Bisection_method
Implemented by Gregor Fabjan from Qnity Consultants on 17/12/2021.
"""
yStart = Galfa(M_Obs, r_Obs, ufr, xStart, Tau) # Check if the initial point is a solution
yEnd = Galfa(M_Obs, r_Obs, ufr, xEnd, Tau) # Check if the final point is a solution
if np.abs(yStart) < Precision:
return xStart # If initial point already satisifes the conditions return start point
if np.abs(yEnd) < Precision:
return xEnd # If final point already satisifes the conditions return end point
iIter = 0
while iIter <= maxIter:
xMid = (xEnd+xStart)/2 # calculate mid point
yMid = Galfa(M_Obs, r_Obs, ufr, xMid, Tau) # What is the solution at midpoint
if (yMid == 0 or (xEnd-xStart)/2 < Precision): # Solution found
return xMid
else: # Solution not fund
iIter += 1
if np.sign(yMid) == np.sign(yStart): # If the start point and the middle point have the same sign, then the root must be in the second half of the interval
xStart = xMid
else: # If the start point and the middle point have a different sign than by mean value theorem the interval must contain at least one root
xEnd = xMid
print("Method failed to converge")