diff options
Diffstat (limited to 'python/analytics')
| -rw-r--r-- | python/analytics/black.pxd | 2 | ||||
| -rw-r--r-- | python/analytics/black.pyx | 2 | ||||
| -rw-r--r-- | python/analytics/cms_spread.py | 96 | ||||
| -rw-r--r-- | python/analytics/cms_spread_utils.pyx | 68 |
4 files changed, 75 insertions, 93 deletions
diff --git a/python/analytics/black.pxd b/python/analytics/black.pxd new file mode 100644 index 00000000..e3b9fafa --- /dev/null +++ b/python/analytics/black.pxd @@ -0,0 +1,2 @@ +#cython: language_level=3 +cdef double cnd_erf(double d) nogil diff --git a/python/analytics/black.pyx b/python/analytics/black.pyx index 80d13a1a..9f733282 100644 --- a/python/analytics/black.pyx +++ b/python/analytics/black.pyx @@ -3,7 +3,7 @@ from libc.math cimport log, sqrt, erf from scipy.stats import norm import cython -cpdef double cnd_erf(double d): +cdef double cnd_erf(double d) nogil: """ 2 * Phi where Phi is the cdf of a Normal """ cdef double RSQRT2 = 0.7071067811865475 return 1 + erf(RSQRT2 * d) diff --git a/python/analytics/cms_spread.py b/python/analytics/cms_spread.py index 24c45e73..fa61c7cf 100644 --- a/python/analytics/cms_spread.py +++ b/python/analytics/cms_spread.py @@ -1,12 +1,12 @@ -import atexit +from . import cms_spread_utils +from .cms_spread_utils import h_call, h_put import datetime import matplotlib.pyplot as plt import numpy as np import pandas as pd import re from math import exp, sqrt, log, pi -from .black import bachelier, cnd_erf -from numba import cfunc, types, float64, vectorize +from .black import bachelier from quantlib.time.api import ( Date, Period, @@ -45,95 +45,7 @@ from utils.db import dawn_engine, serenitas_pool __all__ = ["CmsSpread"] -@vectorize( - [ - float64( - float64, - float64, - float64, - float64, - float64, - float64, - float64, - float64, - float64, - ) - ], - cache=True, - nopython=True, -) -def h_call(z, K, S1, S2, mu_x, mu_y, sigma_x, sigma_y, rho): - # conditionned on S2, integral wrt S1 - # z = (y - mu_y) / sigma_y - u1 = mu_x + rho * sigma_x * z - Ktilde = K + S2 * exp(mu_y + sigma_y * z) - u2 = log(Ktilde / S1) - - v = sigma_x * sqrt(1 - rho * rho) - v2 = sigma_x * sigma_x * (1 - rho * rho) - x = (u1 - u2) / v - return 0.5 * (S1 * exp(u1 + 0.5 * v2) * cnd_erf(x + v) - Ktilde * cnd_erf(x)) - - -@vectorize( - [ - float64( - float64, - float64, - float64, - float64, - float64, - float64, - float64, - float64, - float64, - ) - ], - cache=True, - nopython=True, -) -def h_put(z, K, S1, S2, mu_x, mu_y, sigma_x, sigma_y, rho): - # z = (y - mu_y) / sigma_y - u1 = mu_x + rho * sigma_x * z - Ktilde = K + S2 * exp(mu_y + sigma_y * z) - u2 = log(Ktilde / S1) - - v = sigma_x * sqrt(1 - rho * rho) - v2 = sigma_x * sigma_x * (1 - rho * rho) - x = (u2 - u1) / v - return 0.5 * (Ktilde * cnd_erf(x) - S1 * exp(u1 + 0.5 * v2) * cnd_erf(x - v)) - - -_sig = types.double(types.intc, types.CPointer(types.double)) - - -@cfunc(_sig, cache=True, nopython=True) -def _h1(n, args): - # z = (y - mu_y) / sigma_y - z = args[0] - K = args[1] - S1 = args[2] - S2 = args[3] - mu_x = args[4] - mu_y = args[5] - sigma_x = args[6] - sigma_y = args[7] - rho = args[8] - u1 = mu_x + rho * sigma_x * z - Ktilde = K + S2 * exp(mu_y + sigma_y * z) - u2 = log(Ktilde / S1) - - v = sigma_x * sqrt(1 - rho * rho) - v2 = sigma_x * sigma_x * (1 - rho * rho) - x = (u1 - u2) / v - return ( - 0.5 - * (S1 * exp(u1 + 0.5 * v2) * cnd_erf(x + v) - Ktilde * cnd_erf(x)) - * exp(-0.5 * z * z) - ) - - -_call_integrand = LowLevelCallable(_h1.ctypes) +_call_integrand = LowLevelCallable.from_cython(cms_spread_utils, "h1") def get_fixings(conn, tenor1, tenor2, fixing_date=None): diff --git a/python/analytics/cms_spread_utils.pyx b/python/analytics/cms_spread_utils.pyx new file mode 100644 index 00000000..17396054 --- /dev/null +++ b/python/analytics/cms_spread_utils.pyx @@ -0,0 +1,68 @@ +# cython: language_level=3, cdivision=True +from libc.math cimport log, exp, sqrt +from .black cimport cnd_erf +cimport cython +import numpy as np + +cdef api double h1(int n, double* data) nogil: + # z = (y - mu_y) / sigma_y + cdef: + double z = data[0] + double K = data[1] + double S1 = data[2] + double S2 = data[3] + double mu_x = data[4] + double mu_y = data[5] + double sigma_x = data[6] + double sigma_y = data[7] + double rho = data[8] + double u1, u2, Ktilde, v, v2, x + + u1 = mu_x + rho * sigma_x * z + Ktilde = K + S2 * exp(mu_y + sigma_y * z) + u2 = log(Ktilde / S1) + + v = sigma_x * sqrt(1 - rho * rho) + v2 = sigma_x * sigma_x * (1 - rho * rho) + x = (u1 - u2) / v + return ( + 0.5 + * (S1 * exp(u1 + 0.5 * v2) * cnd_erf(x + v) - Ktilde * cnd_erf(x)) + * exp(-0.5 * z * z) + ) + + +cpdef double[::1] h_call(double[::1] z, double K, double S1, double S2, double mu_x, double mu_y, double sigma_x, double sigma_y, double rho): + # conditionned on S2, integral wrt S1 + # z = (y - mu_y) / sigma_y + cdef double[::1] r = np.empty_like(z) + cdef: + double u1, u2, Ktilde, x + double v = sigma_x * sqrt(1 - rho * rho) + double v2 = sigma_x * sigma_x * (1 - rho * rho) + int i + + for i in range(z.shape[0]): + u1 = mu_x + rho * sigma_x * z[i] + Ktilde = K + S2 * exp(mu_y + sigma_y * z[i]) + u2 = log(Ktilde / S1) + x = (u1 - u2) / v + r[i] = 0.5 * (S1 * exp(u1 + 0.5 * v2) * cnd_erf(x + v) - Ktilde * cnd_erf(x)) + return r + +cpdef double[::1] h_put(double[::1] z, double K, double S1, double S2, double mu_x, double mu_y, double sigma_x, double sigma_y, double rho): + # z = (y - mu_y) / sigma_y + cdef: + double[::1] r = np.empty_like(z) + double u1, u2, Ktilde, x + double v = sigma_x * sqrt(1 - rho * rho) + double v2 = sigma_x * sigma_x * (1 - rho * rho) + int i + + for i in range(z.shape[0]): + u1 = mu_x + rho * sigma_x * z[i] + Ktilde = K + S2 * exp(mu_y + sigma_y * z[i]) + u2 = log(Ktilde / S1) + x = (u2 - u1) / v + r[i] = 0.5 * (Ktilde * cnd_erf(x) - S1 * exp(u1 + 0.5 * v2) * cnd_erf(x - v)) + return r |