from pyisda.legs import ContingentLeg, FeeLeg from pyisda.flat_hazard import strike_vec from pyisda.curve import YieldCurve, BadDay, SpreadCurve import math from scipy.optimize import brentq from scipy.integrate import simps import datetime from tranche_functions import GHquad class Index(): """ minimal class to represent a credit index """ def __init__(self, start_date, end_date, recovery, fixed_rate): self.fixed_rate = fixed_rate self.notional = 1 self._start_date = start_date self._end_date = end_date self.recovery = recovery self._fee_leg = FeeLeg(start_date, end_date, True, 1, 1) self._default_leg = ContingentLeg(start_date, end_date, 1) @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 @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 def pv(self, trade_date, exercise_date, yc, sc, h = None): step_in_date = exercise_date + datetime.timedelta(days=1) a = self._fee_leg.pv(trade_date, step_in_date, trade_date, yc, sc, False) Delta = self._fee_leg.accrued(step_in_date) if exercise_date > trade_date: Delta *= math.exp(-h * year_frac(trade_date, exercise_date)) clean_forward_annuity = a - Delta dl_pv = self._default_leg.pv( trade_date, step_in_date, trade_date, yc, sc, self.recovery) return dl_pv - clean_forward_annuity * self.fixed_rate, clean_forward_annuity def year_frac(d1, d2, day_count_conv = "Actual/365"): """ compute the year fraction between two dates """ if day_count_conv.lower() in ["actual/365", "act/365"]: return (d2-d1).days/365 elif day_count_conv.lower() in ["actual/360", "act/360"]: return (d2-d1).days/360 def DAforward_price(spread, trade_date, exercise_date, yield_curve, index): cash_settle_date = (pd.Timestamp(trade_date) + 3* BDay()).date() h, sc = flat_hazard(spread, yc, trade_date, cash_settle_date, start_date = trade_date, end_date = index.end_date, recovery_rate = index.recovery) exercise_date_settle = (pd.Timestamp(exercise_date) + 3* BDay()).date() forward_price = index.pv(trade_date, exercise_date, yield_curve, sc, h)[0] fep = (1 - index.recovery) * (1-math.exp(-h*year_frac(trade_date, exercise_date))) df = yc.discount_factor(exercise_date_settle) price = fep + 1/df * forward_price return price def flat_hazard(spread, yc, trade_date=datetime.date.today(), cash_settle_date = None, start_date = datetime.date.today(), end_date = datetime.date(2021, 6, 20), recovery_rate = 0.4): step_in_date = start_date + datetime.timedelta(days=1) if cash_settle_date is None: cash_settle_date = (pd.Timestamp(trade_date) + 3* BDay()).date() sc = SpreadCurve(trade_date, yc, trade_date, step_in_date, cash_settle_date, [end_date], array.array('d', [spread]), recovery_rate) sc_data = sc.inspect()['data'] ## conversion to continuous compounding hazard_rate = math.log(1 + sc_data[0][1]) return (hazard_rate, SpreadCurve.from_flat_hazard(trade_date, hazard_rate)) def calib(S0, fp, forward_yield_curve, exercise_date_settle, index, tilt, w): S = S0 * tilt a, b = strike_vec(S, forward_yield_curve, exercise_date_settle, index.start_date, index.end_date, index.recovery) vec = a - index.fixed_rate * b df = forward_yield_curve.discount_factor(exercise_date_settle) return 1/df*np.inner(vec - fp, w) def g(spread, forward_yield_curve, index): """ computes the strike price using the expected forward yield curve """ exercise_date = forward_yield_curve.base_date step_in_date = exercise_date + datetime.timedelta(days=1) exercise_date_settle = (pd.Timestamp(exercise_date) + 3* BDay()).date() sc = SpreadCurve(exercise_date, forward_yield_curve, exercise_date, step_in_date, exercise_date_settle, [index.end_date], array.array('d', [spread]), index.recovery) a = index._fee_leg.pv(exercise_date, step_in_date, exercise_date, forward_yield_curve, sc, True) dl_pv = index._default_leg.pv( exercise_date, step_in_date, exercise_date, forward_yield_curve, sc, index.recovery) df = forward_yield_curve.discount_factor(exercise_date_settle) return 1/df * (dl_pv - a * index.fixed_rate) def ATMstrike(spread, trade_date, exercise_date, yield_curve, index): fp = DAforward_price(spread, trade_date, exercise_date, yc, index) yc_forward = yc.expected_forward_curve(exercise_date) closure = lambda S: g(S, yc_forward, index) - fp eta = 1.1 a = spread b = spread * eta while True: if closure(b) > 0: break return brentq(closure, a, b) def option(index, ref, trade_date, exercise_date, yield_curve, sigma, K, option_type="payer"): """ computes the pv of an option using Pedersen's model """ fp = DAforward_price(ref, trade_date, exercise_date, yield_curve, index) forward_yc = yield_curve.expected_forward_curve(exercise_date) #expiry is end of day (not sure if this is right) T = year_frac(trade_date, exercise_date) Z, w = GHquad(50) tilt = np.exp(-sigma**2/2 * T + sigma * Z * math.sqrt(T)) exercise_date_settle = (pd.Timestamp(exercise_date) + 3* BDay()).date() args = (fp, forward_yc, exercise_date_settle, index, tilt, w) ## atm forward is greater than spread eta = 1.1 a = ref b = ref * eta while True: if calib(*((b,) + args)) > 0: break b *= eta S0 = brentq(calib, a, b, args) S = S0 * tilt G = g(K, forward_yc, index) handle = lambda Z: g(S0 * math.exp(-sigma**2/2 * T + sigma * Z * math.sqrt(T)), forward_yc, index) - G Zstar = brentq(handle, -3, 3) if option_type.lower() == "payer": Z = Zstar + np.logspace(0, 1.1, 300) - 1 elif option_type.lower() == "receiver": Z = Zstar - np.logspace(0, 1.1, 300) + 1 else: raise ValueError("option_type needs to be either 'payer' or 'receiver'") S = S0 * np.exp(-sigma**2/2 * T + sigma * Z * math.sqrt(T)) df = forward_yc.discount_factor(exercise_date_settle) a, b = strike_vec(S, forward_yc, exercise_date_settle, index.start_date, index.end_date, index.recovery) val = ((a - b * index.fixed_rate)/df - G) * 1/math.sqrt(2*math.pi) * np.exp(-Z**2/2) return simps(val, Z) * yield_curve.discount_factor(exercise_date_settle)