diff options
Diffstat (limited to 'python/analytics')
| -rw-r--r-- | python/analytics/__init__.py | 1 | ||||
| -rw-r--r-- | python/analytics/black.py | 15 | ||||
| -rw-r--r-- | python/analytics/index.py | 313 | ||||
| -rw-r--r-- | python/analytics/option.py | 120 | ||||
| -rw-r--r-- | python/analytics/utils.py | 7 |
5 files changed, 456 insertions, 0 deletions
diff --git a/python/analytics/__init__.py b/python/analytics/__init__.py new file mode 100644 index 00000000..90d1b1fe --- /dev/null +++ b/python/analytics/__init__.py @@ -0,0 +1 @@ +from .index import Index diff --git a/python/analytics/black.py b/python/analytics/black.py new file mode 100644 index 00000000..647c92b2 --- /dev/null +++ b/python/analytics/black.py @@ -0,0 +1,15 @@ +import numpy as np +from scipy.stats import norm + +def d1(F, K, sigma, T): + return (np.log(F / K) + sigma**2 * T / 2) / (sigma * np.sqrt(T)) + +def d2(F, K, sigma, T): + return d1(F, K, sigma, T) - sigma * np.sqrt(T) + +def black(F, K, T, sigma, option_type = "payer"): + chi = 1 if option_type == "payer" else -1 + if option_type == "payer": + return F * norm.cdf(d1(F, K, sigma, T)) - K * norm.cdf(d2(F, K, sigma, T)) + else: + return K * norm.cdf(- d2(F, K, sigma, T)) - F * norm.cdf(- d1(F, K, sigma, T)) diff --git a/python/analytics/index.py b/python/analytics/index.py new file mode 100644 index 00000000..85ba3474 --- /dev/null +++ b/python/analytics/index.py @@ -0,0 +1,313 @@ +import array +import datetime +import math +import pandas as pd + +from pyisda.legs import ContingentLeg, FeeLeg +from quantlib.settings import Settings +from quantlib.time.api import ( + Date, Schedule, WeekendsOnly, CDS, Following, + Unadjusted, Period, pydate_from_qldate ) +from termcolor import colored +from pandas.tseries.offsets import BDay +from dates import prev_immdate +from db import dbconn +from psycopg2 import DataError +from pyisda.curve import SpreadCurve +from yieldcurve import YC, ql_to_jp, roll_yc, rate_helpers +from quantlib.time.api import Actual365Fixed + +serenitasdb = dbconn('serenitasdb') + +class Index(): + """ minimal class to represent a credit index """ + def __init__(self, start_date, end_date, recovery, fixed_rate, + notional = 10e6): + """ + start_date : :class:`datetime.date` + index start_date (Could be issue date, or last imm date) + end_date : :class:`datetime.date` + index last date + recovery : + recovery rate (between 0 and 1) + fixed_rate : + fixed coupon (in bps) + """ + self.fixed_rate = fixed_rate + self.notional = notional + self._sched = Schedule(Date.from_datetime(start_date), + Date.from_datetime(end_date), + Period("3M"), + WeekendsOnly(), + Following, + Unadjusted, + CDS) + self._start_date = start_date + self._end_date = end_date + self.recovery = recovery + + self._fee_leg = FeeLeg(self._start_date, end_date, True, 1, 1) + self._default_leg = ContingentLeg(self._start_date, end_date, 1) + self._trade_date = None + self._yc = None + self._risky_annuity = None + self._spread = None + self.name = None + + @property + def start_date(self): + return self._start_date + + @property + def end_date(self): + return self._end_date + + @start_date.setter + def start_date(self, d): + self._fee_leg = FeeLeg(d, self.end_date, True, 1, 1) + self._default_leg = ContingentLeg(d, self.end_date, 1) + self._start_date = d + self._sched = Schedule(Date.from_datetime(d), + Date.from_datetime(self.end_date), + Period("3M"), + WeekendsOnly(), + Following, + Unadjusted, + CDS) + + @end_date.setter + def end_date(self, d): + self._fee_leg = FeeLeg(self.start_date, d, True, 1, 1) + self._default_leg = ContingentLeg(self.start_date, d, 1) + self._end_date = d + self._sched = Schedule(self.start_date, + d, + Period("3M"), + WeekendsOnly(), + Following, + Unadjusted, + CDS) + + def survival_probability(self, d): + if d > self.trade_date: + return 1 + else: + return math.exp( - self.flat_hazard * (d - self.trade_date)/365) + + def forward_pv(self, exercise_date): + """This is default adjusted forward price at time exercise_date""" + step_in_date = exercise_date + datetime.timedelta(days=1) + a = self._fee_leg.pv(self.trade_date, step_in_date, self._value_date, + self._yc, self._sc, False) + Delta = self._fee_leg.accrued(step_in_date) + value_date = (pd.Timestamp(exercise_date) + 3* BDay()).date() + df = self._yc.discount_factor(value_date) + q = self.survival_probability(exercise_date) + clean_forward_annuity = a - Delta * df * q + dl_pv = self._default_leg.pv( + self.trade_date, step_in_date, self._value_date, + self._yc, self._sc, self.recovery) + forward_price = self.notional * (dl_pv - clean_forward_annuity * self.fixed_rate*1e-4) + fep = self.notional * (1 - self.recovery) * (1 - q) + return forward_price * self._yc.discount_factor(self._value_date) / df + fep + + @property + def spread(self): + return self._spread * 1e4 + + def _update(self): + self._sc = SpreadCurve(self.trade_date, self._yc, self.start_date, + self._step_in_date, self._value_date, + [self.end_date], array.array('d', [self._spread]), + self.recovery) + self._risky_annuity = self._fee_leg.pv(self.trade_date, self._step_in_date, + self._value_date, self._yc, + self._sc, False) + self._dl_pv = self._default_leg.pv( + self.trade_date, self._step_in_date, self._value_date, + self._yc, self._sc, self.recovery) + self._pv = self._dl_pv - self._risky_annuity * self.fixed_rate * 1e-4 + self._clean_pv = self._pv + self._accrued * self.fixed_rate * 1e-4 + self._price = 100 * (1 - self._clean_pv) + + @spread.setter + def spread(self, s: float): + """ s: spread in bps """ + self._spread = s * 1e-4 + self._update() + + @property + def flat_hazard(self): + sc_data = self._sc.inspect()['data'] + ## conversion to continuous compounding + return math.log(1 + sc_data[0][1]) + + @property + def pv(self): + return self.notional * self._pv + + @property + def accrued(self): + return - self.notional * self._accrued * self.fixed_rate * 1e-4 + + @property + def days_accrued(self): + return int(self._accrued * 360) + + @property + def clean_pv(self): + return self.notional * self._clean_pv + + @property + def price(self): + return self._price + + @price.setter + def price(self, val): + pass + + @property + def DV01(self): + old_pv = self.pv + self.spread += 1 + dv01 = self.pv - old_pv + self.spread -= 1 + return dv01 + + @property + def IRDV01(self): + old_pv = self.pv + old_yc = self._yc + for rh in self._helpers: + rh.quote += 1e-4 + self._yc = ql_to_jp(self._ql_yc) + self._update() ## to force recomputation + new_pv = self.pv + for r in self._helpers: + r.quote -= 1e-4 + self._yc = old_yc + self._update() + return new_pv - old_pv + + @property + def rec_risk(self): + old_pv = self.pv + old_recovery = self.recovery + self.recovery = old_recovery - 0.01 + self._update() + pv_minus = self.pv + self.recovery = old_recovery + 0.01 + self._update() + pv_plus = self.pv + self.recovery = old_recovery + self._update() + return (pv_plus - pv_minus)/2 + + @property + def jump_to_default(self): + return self.notional * (1 - self.recovery) - self.clean_pv + + @property + def risky_annuity(self): + return self._risky_annuity - self._accrued + + @property + def trade_date(self): + if self._trade_date is None: + raise AttributeError('Please set trade_date first') + else: + return self._trade_date + + @trade_date.setter + def trade_date(self, d): + settings = Settings() + settings.evaluation_date = Date.from_datetime(d) + self.start_date = pydate_from_qldate( + self._sched.previous_date(settings.evaluation_date)) + self._helpers = rate_helpers(self.currency) + self._ql_yc = YC(self._helpers) + self._yc = ql_to_jp(self._ql_yc) + self._trade_date = d + self._step_in_date = self.trade_date + datetime.timedelta(days=1) + self._accrued = self._fee_leg.accrued(self._step_in_date) + self._value_date = (pd.Timestamp(self._trade_date) + 3* BDay()).date() + if self._spread is not None: + self._update() + + @classmethod + def from_name(cls, index, series, tenor, trade_date = datetime.date.today(), + notional = 10e6): + try: + with serenitasdb.cursor() as c: + c.execute("SELECT maturity, coupon FROM index_maturity " \ + "WHERE index=%s AND series=%s AND tenor = %s", + (index.upper(), series, tenor)) + maturity, coupon = next(c) + except DataError as e: + raise + else: + recovery = 0.4 if index.lower() == "ig" else 0.3 + instance = cls(trade_date, maturity, recovery, coupon) + instance.name = "MARKIT CDX.NA.{}.{} {:%m/%y} ".format( + index.upper(), + series, + maturity) + if index.upper() in ["IG", "HY"]: + instance.currency = "USD" + else: + instance.currency = "EUR" + instance.notional = notional + instance.trade_date = trade_date + return instance + + def __repr__(self): + if self.days_accrued > 1: + accrued_str = "Accrued ({} Days)".format(self.days_accrued) + else: + accrued_str = "Accrued ({} Day)".format(self.days_accrued) + s = ["{:<20}\t{:>15}".format("CDS Index", colored(self.name, attrs = ['bold'])), + "", + "{:<20}\t{:>15}".format("Trade Date", ('{:%m/%d/%y}'. + format(self.trade_date))), + "{:<20}\t{:>15.2f}\t\t{:<20}\t{:>10,.2f}".format("Trd Sprd (bp)", + self.spread, + "Coupon (bp)", + self.fixed_rate), + "{:<20}\t{:>15.2f}\t\t{:<20}\t{:>10}".format("1st Accr Start", + self.spread, + "Payment Freq", + "Quarterly"), + "{:<20}\t{:>15}\t\t{:<20}\t{:>10.2f}".format("Maturity Date", + ('{:%m/%d/%y}'. + format(self.end_date)), + "Rec Rate", + self.recovery), + "{:<20}\t{:>15}\t\t{:<20}\t{:>10}".format("Bus Day Adj", + "Following", + "Day Count", + "ACT/360"), + "", + colored("Calculator", attrs = ['bold']), + "{:<20}\t{:>15}".format("Valuation Date", ('{:%m/%d/%y}'. + format(self.trade_date))), + "{:<20}\t{:>15}".format("Cash Settled On", ('{:%m/%d/%y}'. + format(self._value_date))), + "", + "{:<20}\t{:>15.8f}\t\t{:<20}\t{:>10,.2f}".format("Price", + self.price, + "Spread DV01", + self.DV01), + "{:<20}\t{:>15,.0f}\t\t{:<20}\t{:>10,.2f}".format("Principal", + self.clean_pv, + "IR DV01", + self.IRDV01), + "{:<20}\t{:>15,.0f}\t\t{:<20}\t{:>10,.2f}".format(accrued_str, + self.accrued, + "Rec Risk (1%)", + self.rec_risk), + "{:<20}\t{:>15,.0f}\t\t{:<20}\t{:>10,.0f}".format("Cash Amount", + self.pv, + "Def Exposure", + self.jump_to_default) + ] + return "\n".join(s) diff --git a/python/analytics/option.py b/python/analytics/option.py new file mode 100644 index 00000000..32f4f947 --- /dev/null +++ b/python/analytics/option.py @@ -0,0 +1,120 @@ +from .black import black +from .utils import GHquad +from yieldcurve import roll_yc +from pandas.tseries.offsets import BDay + +class Option: + def __init__(self, index, exercise_date, strike, option_type="payer"): + self.index = index + self._exercise_date = exercise_date + self._forward_yc = roll_yc(self.index._yc, self.exercise_date) + self.exercise_date_settle = (pd.Timestamp(self.exercise_date) + 3* BDay()).date() + self._T = None + self.strike = strike + self.option_type = option_type.lower() + self._Z, self._w = GHquad(50) + self.notional = 1 + + @property + def exercise_date(self): + return self._exercise_date + + @exercise_date.setter + def exercise_date(self, d : datetime.date): + self._exercise_date = d + self.exercise_date_settle = (pd.Timestamp(d) + 3* BDay()).date() + self._forward_yc = roll_yc(self.index._yc, self.exercise_date) + + @property + def pv(self): + fp = self.index.forward_pv(self.exercise_date) / self.index.notional + T = self.T + tilt = np.exp(-self.sigma**2/2 * T + self.sigma * self._Z * math.sqrt(T)) + rolled_curve = roll_yc(self.index._yc, self.exercise_date) + args = (fp, self.exercise_date, self.exercise_date_settle, + self.index, self._forward_yc, tilt, self._w) + eta = 1.1 + a = self.index.spread + b = self.index.spread * eta + while True: + if calib(*((b,) + args)) > 0: + break + b *= eta + + S0 = brentq(calib, a, b, args) + + G = g(self.index, self.strike, self.exercise_date) + if T == 0: + pv = self.notional * (g(self.index, self.index.spread, self.exercise_date) - G) + if self.option_type == "payer": + return pv if self.index.spread > self.strike else 0 + else: + return - pv if self.index.spread < self.strike else 0 + + Zstar = (math.log(self.strike/S0) + self.sigma**2/2 * T) / \ + (self.sigma * math.sqrt(T)) + + if self.option_type == "payer": + Z = Zstar + np.logspace(0, 1.5, 300) - 1 + elif self.option_type == "receiver": + Z = Zstar - np.logspace(0, 1.5, 300) + 1 + else: + raise ValueError("option_type needs to be either 'payer' or 'receiver'") + S = S0 * np.exp(-self.sigma**2/2 * T + self.sigma * Z * math.sqrt(T)) + a, b = strike_vec(S * 1e-4, rolled_curve, self.exercise_date, + self.exercise_date_settle, + self.index.start_date, self.index.end_date, self.index.recovery) + val = ((a - b * self.index.fixed_rate*1e-4) - G) * 1/math.sqrt(2*math.pi) * np.exp(-Z**2/2) + df_scale = self.index._yc.discount_factor(self.exercise_date_settle) + return self.notional * simps(val, Z) * df_scale + + @property + def pv2(self): + G = g(self.index, self.strike, self.exercise_date) + fp = self.index.forward_pv(self.exercise_date) / self.index.notional + forward_annuity = self.index.forward_annuity(self.exercise_date) + DA_forward_spread = fp / forward_annuity + self.index.fixed_rate * 1e-4 + strike_tilde = self.index.fixed_rate * 1e-4 + G / forward_annuity + return forward_annuity * black(DA_forward_spread, + strike_tilde, + self.T, + self.sigma, + self.option_type) * self.notional + + @property + def delta(self): + old_index_pv = self.index.pv + old_pv = self.pv + self.index.spread += 0.1 + notional_ratio = self.index.notional/self.notional + delta = (self.pv - old_pv)/(self.index.pv - old_index_pv) * notional_ratio + self.index.spread -= 0.1 + return delta + + + @property + def T(self): + if self._T: + return self._T + else: + return year_frac(self.index.trade_date, self.exercise_date) + 1/365 + + @property + def gamma(self): + pass + + @property + def theta(self): + old_pv = self.pv + self._T = self.T - 1/365 + theta = self.pv - old_pv + self._T = None + return theta + + @property + def vega(self): + old_pv = self.pv + self.sigma += 0.01 + vega = self.pv - old_pv + self.sigma -= 0.01 + return vega diff --git a/python/analytics/utils.py b/python/analytics/utils.py new file mode 100644 index 00000000..a67d6f32 --- /dev/null +++ b/python/analytics/utils.py @@ -0,0 +1,7 @@ +import numpy as np +from scipy.special import h_roots + +def GHquad(n : int): + """Gauss-Hermite quadrature weights""" + Z, w = h_roots(n) + return Z*np.sqrt(2), w/np.sqrt(np.pi) |