diff --git a/pyfeng/__init__.py b/pyfeng/__init__.py index 83b5f11..a1944ac 100644 --- a/pyfeng/__init__.py +++ b/pyfeng/__init__.py @@ -12,6 +12,7 @@ SabrChoiWu2021H, SabrChoiWu2021P, ) +from .garch import GarchMcTimeStep, GarchUncorrBaroneAdesi2004 from .sabr_int import SabrUncorrChoiWu2021 from .sabr_mc import SabrMcCond from .nsvh import Nsvh1, NsvhMc diff --git a/pyfeng/data/heston_benchmark.xlsx b/pyfeng/data/heston_benchmark.xlsx index 0e5c59d..382af59 100644 Binary files a/pyfeng/data/heston_benchmark.xlsx and b/pyfeng/data/heston_benchmark.xlsx differ diff --git a/pyfeng/ex.py b/pyfeng/ex.py index a3497f1..d279ce5 100644 --- a/pyfeng/ex.py +++ b/pyfeng/ex.py @@ -4,7 +4,6 @@ from .heston_mc import HestonMcAndersen2008, HestonMcGlassermanKim2011, HestonMcTseWan2013, HestonMcChoiKwok2023 from .sv32_mc import Sv32McCondQE, Sv32McAe2 from .sv32_mc2 import Sv32McTimeStep, Sv32McExactBaldeaux2012, Sv32McExactChoiKwok2023 -from .garch import GarchMcTimeStep, GarchUncorrBaroneAdesi2004 from .subord_bm import VarGammaQuad, ExpNigQuad # SABR / OUSV models for research diff --git a/pyfeng/garch.py b/pyfeng/garch.py index 767310e..d2921f8 100644 --- a/pyfeng/garch.py +++ b/pyfeng/garch.py @@ -81,7 +81,7 @@ class GarchMcTimeStep(sv.SvABC, sv.CondMcBsmABC): var_process = True scheme = 1 # - def set_mc_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, scheme=1): + def set_num_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, scheme=1): """ Set MC parameters @@ -95,7 +95,7 @@ def set_mc_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, sc References: - Andersen L (2008) Simple and efficient simulation of the Heston stochastic volatility model. Journal of Computational Finance 11:1–42. https://doi.org/10.21314/JCF.2008.189 """ - super().set_mc_params(n_path, dt, rn_seed, antithetic) + super().set_num_params(n_path, dt, rn_seed, antithetic) self.scheme = scheme def var_step_euler(self, var_0, dt, milstein=True): diff --git a/pyfeng/heston_mc.py b/pyfeng/heston_mc.py index 4273e57..de0db85 100644 --- a/pyfeng/heston_mc.py +++ b/pyfeng/heston_mc.py @@ -162,7 +162,7 @@ class HestonMcAndersen2008(HestonMcABC): >>> spot = 100 >>> sigma, vov, mr, rho, texp = 0.04, 1, 0.5, -0.9, 10 >>> m = pfex.HestonMcAndersen2008(sigma, vov=vov, mr=mr, rho=rho) - >>> m.set_mc_params(n_path=1e5, dt=1/8, rn_seed=123456) + >>> m.set_num_params(n_path=1e5, dt=1/8, rn_seed=123456) >>> m.price(strike, spot, texp) >>> # true price: 44.330, 13.085, 0.296 array([44.31943535, 13.09371251, 0.29580431]) @@ -170,7 +170,7 @@ class HestonMcAndersen2008(HestonMcABC): psi_c = 1.5 # parameter used by the Andersen QE scheme scheme = 4 - def set_mc_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, scheme=4): + def set_num_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, scheme=4): """ Set MC parameters @@ -184,7 +184,7 @@ def set_mc_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, sc References: - Andersen L (2008) Simple and efficient simulation of the Heston stochastic volatility model. Journal of Computational Finance 11:1–42. https://doi.org/10.21314/JCF.2008.189 """ - super().set_mc_params(n_path, dt, rn_seed, antithetic) + super().set_num_params(n_path, dt, rn_seed, antithetic) self.scheme = scheme def var_step_qe(self, var_0, dt): @@ -305,7 +305,7 @@ class HestonMcGlassermanKim2011(HestonMcABC): kk = 1 # K for series truncation. tabulate_x2_z = False - def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, kk=1): + def set_num_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, kk=1): """ Set MC parameters @@ -317,7 +317,7 @@ def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, kk=1): kk: truncation index """ - super().set_mc_params(n_path, dt, rn_seed, antithetic=False) + super().set_num_params(n_path, dt, rn_seed, antithetic=False) self.scheme = scheme self.kk = kk @@ -828,14 +828,14 @@ class HestonMcTseWan2013(HestonMcGlassermanKim2011): >>> spot = 100 >>> sigma, vov, mr, rho, texp = 0.04, 1, 0.5, -0.9, 10 >>> m = pfex.HestonMcTseWan2013(sigma, vov=vov, mr=mr, rho=rho) - >>> m.set_mc_params(n_path=1e4, rn_seed=123456) + >>> m.set_num_params(n_path=1e4, rn_seed=123456) >>> m.price(strike, spot, texp) >>> # true price: 44.330, 13.085, 0.296 array([12.08981758, 0.33379748, 42.28798189]) # not close so far """ dist = 'ig' - def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, dist=None): + def set_num_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, dist=None): """ Set MC parameters @@ -847,7 +847,7 @@ def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, dist=None dist: distribution to use for approximation. 'ig' for inverse Gaussian (default), 'ga' for Gamma, 'ln' for LN """ - super().set_mc_params(n_path, dt, rn_seed, scheme=scheme) + super().set_num_params(n_path, dt, rn_seed, scheme=scheme) if dist is not None: self.dist = dist @@ -891,7 +891,7 @@ class HestonMcChoiKwok2023(HestonMcGlassermanKim2011): dist = 'ig' - def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, kk=0, dist=None): + def set_num_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, kk=0, dist=None): """ Set MC parameters @@ -903,7 +903,7 @@ def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, scheme=3, kk=0, dis dist: distribution to use for approximation. 'ig' for inverse Gaussian (default), 'ga' for Gamma, 'ln' for LN """ - super().set_mc_params(n_path, dt, rn_seed, scheme=scheme, kk=kk) + super().set_num_params(n_path, dt, rn_seed, scheme=scheme, kk=kk) if dist is not None: self.dist = dist diff --git a/pyfeng/ousv.py b/pyfeng/ousv.py index 086f2ed..61d3733 100644 --- a/pyfeng/ousv.py +++ b/pyfeng/ousv.py @@ -156,6 +156,26 @@ def cond_spot_sigma(self, vol_0, texp): sigma_cond = np.sqrt((1 - self.rho**2) * var_mean) / vol_0 return spot_cond, sigma_cond + def avgvar_mv(self, var_0, texp): + """ + Mean and variance of the variance V(t+dt) given V(0) = var_0 + (variance is not implemented yet) + + Args: + var_0: initial variance + texp: time step + + Returns: + mean, variance + """ + + mr_t = self.mr * texp + e_mr = np.exp(-mr_t) + x0 = var_0 - self.theta + vv = self.vov**2/2/self.mr + self.theta**2 + \ + ((x0**2 - self.vov**2/2/self.mr)*(1 + e_mr) + 4*self.theta * x0)*(1 - e_mr)/(2*self.mr*texp) + return vv + class OusvMcTimeStep(OusvMcABC): """ @@ -211,7 +231,7 @@ class OusvMcChoi2023(OusvMcABC): n_sin = 2 - def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, antithetic=True, n_sin=2): + def set_num_params(self, n_path=10000, dt=None, rn_seed=None, antithetic=True, n_sin=2): """ Set MC parameters @@ -224,7 +244,7 @@ def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, antithetic=True, n_ assert n_sin % 2 == 0 self.n_sin = n_sin - super().set_mc_params(n_path, dt, rn_seed, antithetic) + super().set_num_params(n_path, dt, rn_seed, antithetic) @classmethod def _a2sum(cls, mr_t, ns=0, odd=None): @@ -444,7 +464,7 @@ def cond_states_step(self, vol_0, dt, zn=None): cosh = np.cosh(mr_t) vovn = self.vov * np.sqrt(dt) # normalized vov - x_0 = self.sigma - self.theta + x_0 = vol_0 - self.theta if zn is None: x_t = self.vol_step(vol_0, dt) - self.theta else: @@ -512,6 +532,11 @@ def cond_states_step(self, vol_0, dt, zn=None): return x_t, vv_t, uu_t + def unexplained_var_ratio(self, mr_t, ns): + ns = ns or self.n_sin + rv = self._a4sum(mr_t, ns=ns) / self._a4sum(mr_t) + return rv + def vol_path_sin(self, tobs, zn=None): """ vol path composed of sin terms @@ -548,3 +573,21 @@ def vol_path_sin(self, tobs, zn=None): + self.vov * np.sqrt(dt) * (an*sin) @ zn[1:,:] return sigma_path + + def price_var_option(self, strike, texp, cp=1): + """ + Price of variance option + + Args: + strike: + texp: + cp: + + Returns: + + """ + df = np.exp(-self.intr * texp) + vol_t, vv_t, uu_t = self.cond_states(self.sigma, texp) + # vv_t is the average variance + price = df * np.fmax(np.sign(cp)*(vv_t[:, None] - strike), 0).mean(axis=0) + return price \ No newline at end of file diff --git a/pyfeng/sabr_mc.py b/pyfeng/sabr_mc.py index 137f290..6f099e9 100644 --- a/pyfeng/sabr_mc.py +++ b/pyfeng/sabr_mc.py @@ -117,7 +117,7 @@ class SabrMcExactCai2017(sabr.SabrABC, sv.CondMcBsmABC): comb_coef = None nn = None - def set_mc_params(self, n_path=10000, m_inv=20, m_euler=20, n_euler=35, rn_seed=None, antithetic=True): + def set_num_params(self, n_path=10000, m_inv=20, m_euler=20, n_euler=35, rn_seed=None, antithetic=True): """ Set MC parameters diff --git a/pyfeng/sv32_mc2.py b/pyfeng/sv32_mc2.py index 21335de..e909d2d 100644 --- a/pyfeng/sv32_mc2.py +++ b/pyfeng/sv32_mc2.py @@ -57,7 +57,7 @@ class Sv32McTimeStep(Sv32McABC): """ scheme = 1 # Milstein - def set_mc_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, scheme=1): + def set_num_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, scheme=1): """ Set MC parameters @@ -71,13 +71,13 @@ def set_mc_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True, sc References: - Andersen L (2008) Simple and efficient simulation of the Heston stochastic volatility model. Journal of Computational Finance 11:1–42. https://doi.org/10.21314/JCF.2008.189 """ - super().set_mc_params(n_path, dt, rn_seed, antithetic) + super().set_num_params(n_path, dt, rn_seed, antithetic) self.scheme = scheme mr = self.mr * self.theta theta = (self.mr + self.vov**2)/mr self._m_heston = heston_mc.HestonMcAndersen2008(1/self.sigma, self.vov, self.rho, mr, theta) - self._m_heston.set_mc_params(n_path, dt, rn_seed, antithetic, scheme=scheme) + self._m_heston.set_num_params(n_path, dt, rn_seed, antithetic, scheme=scheme) def var_step_euler(self, var_0, dt, milstein=True): """ @@ -155,7 +155,7 @@ class Sv32McExactBaldeaux2012(Sv32McABC): is_fwd: Bool, true if asset price is forward """ - def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, antithetic=True, scheme=1): + def set_num_params(self, n_path=10000, dt=None, rn_seed=None, antithetic=True, scheme=1): """ Set MC parameters @@ -169,13 +169,13 @@ def set_mc_params(self, n_path=10000, dt=None, rn_seed=None, antithetic=True, sc References: - Andersen L (2008) Simple and efficient simulation of the Heston stochastic volatility model. Journal of Computational Finance 11:1–42. https://doi.org/10.21314/JCF.2008.189 """ - super().set_mc_params(n_path, dt, rn_seed, antithetic) + super().set_num_params(n_path, dt, rn_seed, antithetic) self.scheme = scheme mr = self.mr * self.theta theta = (self.mr + self.vov**2)/mr self._m_heston = heston_mc.HestonMcAndersen2008(1/self.sigma, self.vov, self.rho, mr, theta) - self._m_heston.set_mc_params(n_path, dt, rn_seed, antithetic, scheme=scheme) + self._m_heston.set_num_params(n_path, dt, rn_seed, antithetic, scheme=scheme) def laplace(self, bb, var_0, var_t, dt): phi, _ = self._m_heston.phi_exp(dt) diff --git a/pyfeng/sv_abc.py b/pyfeng/sv_abc.py index 8e685be..6dc6eff 100644 --- a/pyfeng/sv_abc.py +++ b/pyfeng/sv_abc.py @@ -116,7 +116,7 @@ class CondMcBsmABC(smile.OptSmileABC, abc.ABC): correct_fwd = True result = {} - def set_mc_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True): + def set_num_params(self, n_path=10000, dt=0.05, rn_seed=None, antithetic=True): """ Set MC parameters @@ -152,7 +152,7 @@ def tobs(self, texp): if self.dt is None: return np.array([texp]) else: - n_dt = np.ceil(texp / (2 * self.dt)) * 2 + n_dt = np.ceil(texp / self.dt) tobs = np.arange(1, n_dt + 1) / n_dt * texp return tobs diff --git a/tests/test_heston.py b/tests/test_heston.py index 70326f9..f15e8af 100644 --- a/tests/test_heston.py +++ b/tests/test_heston.py @@ -55,7 +55,7 @@ def test_price_mc(self): """ for no in [1, 2, 3]: m, p, rv = pfex.HestonMcAndersen2008.init_benchmark(no) - m.set_mc_params(n_path=1e5, dt=1/8, rn_seed=123456) + m.set_num_params(n_path=1e5, dt=1/8, rn_seed=123456) m.correct_fwd = False vol0 = pf.Bsm(None, intr=m.intr, divr=m.divr).impvol(rv['val'], **rv['args_pricing']) @@ -65,14 +65,14 @@ def test_price_mc(self): np.testing.assert_allclose(m.result['spot error'], 0, atol=2e-3) m, *_ = pfex.HestonMcGlassermanKim2011.init_benchmark(no) - m.set_mc_params(n_path=1e5, rn_seed=123456, kk=10) + m.set_num_params(n_path=1e5, rn_seed=123456, kk=10) m.correct_fwd = False vol1 = m.vol_smile(**rv['args_pricing']) np.testing.assert_allclose(vol0, vol1, atol=5e-3) np.testing.assert_allclose(m.result['spot error'], 0, atol=2e-3) m, *_ = pfex.HestonMcTseWan2013.init_benchmark(no) - m.set_mc_params(n_path=1e5, rn_seed=123456, dt=1) + m.set_num_params(n_path=1e5, rn_seed=123456, dt=1) m.correct_fwd = False vol1 = m.vol_smile(**rv['args_pricing']) np.testing.assert_allclose(vol0, vol1, atol=5e-3) @@ -80,12 +80,12 @@ def test_price_mc(self): m, *_ = pfex.HestonMcChoiKwok2023.init_benchmark(no) m.correct_fwd = False - m.set_mc_params(n_path=1e5, rn_seed=123456, kk=10, dt=None) + m.set_num_params(n_path=1e5, rn_seed=123456, kk=10, dt=None) vol1 = m.vol_smile(**rv['args_pricing']) np.testing.assert_allclose(vol0, vol1, atol=5e-3) np.testing.assert_allclose(m.result['spot error'], 0, atol=2e-3) - m.set_mc_params(n_path=1e5, rn_seed=123456, kk=1, dt=1/4) + m.set_num_params(n_path=1e5, rn_seed=123456, kk=1, dt=1/4) vol1 = m.vol_smile(**rv['args_pricing']) np.testing.assert_allclose(vol0, vol1, atol=5e-3) np.testing.assert_allclose(m.result['spot error'], 0, atol=2e-3) diff --git a/tests/test_sabr.py b/tests/test_sabr.py index 160acb5..7e15bd8 100644 --- a/tests/test_sabr.py +++ b/tests/test_sabr.py @@ -64,7 +64,7 @@ def test_UnCorrChoiWu2021(self): def test_CondMc(self): for k in [19, 20]: # can test 22 (Korn&Tang) also, but difficult to pass m, df, rv = pf.SabrMcCond.init_benchmark(k) - m.set_mc_params(n_path=5e4, dt=0.05, rn_seed=1234) + m.set_num_params(n_path=5e4, dt=0.05, rn_seed=1234) p = m.price(**rv["args_pricing"]) np.testing.assert_almost_equal(p, rv["val"], decimal=4)