Fix FX variance swap fair strike and the realised variance normalisation

financepy/products/equity/equity_variance_swap.py

@@ -256,7 +256,13 @@ def realised_variance(self, close_prices: float, use_logs=True):
                 x = (close_prices[i] - close_prices[i - 1]) / close_prices[i - 1]
                 cum_x2 += x * x
 
-        var = cum_x2 * 252.0 / num_observations
+        # N prices give N - 1 returns: the market-standard realised variance
+        # is 252 / (N - 1) times the sum of the squared returns.
+        num_returns = num_observations - 1
+        if num_returns < 1:
+            raise FinError("At least two prices are needed for a realised variance")
+
+        var = cum_x2 * 252.0 / num_returns
         return var
 
     ###########################################################################

financepy/products/fx/fx_variance_swap.py

@@ -11,8 +11,7 @@
 from ...utils.global_vars import ONE_MILLION
 from ...utils.global_vars import G_DAYS_IN_YEAR
 from ...utils.global_types import OptionTypes
-from .fx_vanilla_option import FXVanillaOption
-from ...models.black_scholes import BlackScholes
+from ...models.black_scholes_analytic import european_value
 from ...market.curves.discount_curve import DiscountCurve
 
 from ...utils.helpers import check_argument_types
@@ -165,7 +164,9 @@ def fair_strike(
         if min_strike < strike_spacing:
             k = sstar
             klist = [sstar]
-            while k >= strike_spacing:
+            # keep the put strikes strictly positive: the log-contract weight is
+            # undefined at a zero strike
+            while k - strike_spacing > 0.0:
                 k -= strike_spacing
                 klist.append(k)
             put_k = np.array(klist)
@@ -180,10 +181,17 @@ def fair_strike(
 
         self.call_strikes = call_k
 
-        option_total = 2.0 * (r * t_mat - (s0 * g / sstar - 1.0) - np.log(sstar / s0)) / t_mat
+        # Demeterfi, Derman, Kamal and Zhou (1999), eq. (29), with the foreign
+        # rate q in the role of a dividend yield: the drift of the log contract is
+        # (r - q) and the option portfolio is compounded at the domestic rate r.
+        option_total = (
+            2.0 * ((r - q) * t_mat - (s0 * g / sstar - 1.0) - np.log(sstar / s0)) / t_mat
+        )
 
-        self.call_wts = np.zeros(num_call_options)
-        self.put_wts = np.zeros(num_put_options)
+        # The number of puts may have been reduced above when the requested
+        # strikes would have gone below zero, so use the stored counts throughout.
+        self.call_wts = np.zeros(self.num_call_options)
+        self.put_wts = np.zeros(self.num_put_options)
 
         def f(x):
             return (2.0 / t_mat) * ((x - sstar) / sstar - np.log(x / sstar))
@@ -203,25 +211,21 @@ def f(x):
             sum_wts += self.call_wts[n]
 
         pi_put = 0.0
-        for n in range(0, num_put_options):
+        for n in range(0, self.num_put_options):
             k = put_k[n]
             vol = volatility_curve.volatility(k)
-            opt = FXVanillaOption(self.maturity_dt, k, put_type)
-            model = BlackScholes(vol)
-            v = opt.value(value_dt, s0, discount_curve, dividend_curve, model)
+            v = european_value(s0, t_mat, k, r, q, vol, put_type.value)
             pi_put += v * self.put_wts[n]
 
         pi_call = 0.0
-        for n in range(0, num_call_options):
+        for n in range(0, self.num_call_options):
             k = call_k[n]
             vol = volatility_curve.volatility(k)
-            opt = FXVanillaOption(self.maturity_dt, k, call_type)
-            model = BlackScholes(vol)
-            v = opt.value(value_dt, s0, discount_curve, dividend_curve, model)
+            v = european_value(s0, t_mat, k, r, q, vol, call_type.value)
             pi_call += v * self.call_wts[n]
 
         pi = pi_call + pi_put
-        option_total += g * pi
+        option_total += np.exp(r * t_mat) * pi
         var = option_total
 
         return var
@@ -249,7 +253,13 @@ def realised_variance(self, close_prices: float, use_logs=True):
                 x = (close_prices[i] - close_prices[i - 1]) / close_prices[i - 1]
                 cum_x2 += x * x
 
-        var = cum_x2 * 252.0 / num_observations
+        # N prices give N - 1 returns: the market-standard realised variance
+        # is 252 / (N - 1) times the sum of the squared returns.
+        num_returns = num_observations - 1
+        if num_returns < 1:
+            raise FinError("At least two prices are needed for a realised variance")
+
+        var = cum_x2 * 252.0 / num_returns
         return var
 
     ###########################################################################

tests/unit/test_FinEquityVarianceSwap.py

@@ -136,3 +136,22 @@ def test_fair_strike_with_dividend_yield_and_zero_strike_grid():
         )
         assert np.isfinite(fair_var), (r, q)
         assert abs(fair_var - sigma**2) / sigma**2 < 0.01, (r, q, fair_var)
+
+
+def test_realised_variance_of_constant_log_returns():
+    """N prices give N - 1 returns. A series growing by the same log return every
+    day has a realised variance of 252 times that squared return."""
+    import numpy as np
+    import pytest
+
+    from financepy.products.equity.equity_variance_swap import EquityVarianceSwap
+    from financepy.utils.date import Date
+    from financepy.utils.error import FinError
+
+    swap = EquityVarianceSwap(Date(1, 1, 2026), Date(1, 1, 2027), 0.09)
+    daily = 0.01
+    closes = 100.0 * np.exp(daily * np.arange(11))
+    assert np.isclose(swap.realised_variance(closes), 252.0 * daily**2)
+    assert np.isclose(swap.realised_variance(closes[:2]), 252.0 * daily**2)
+    with pytest.raises(FinError):
+        swap.realised_variance(closes[:1])

tests/unit/test_FinFXVarianceSwap.py

@@ -0,0 +1,69 @@
+########################################################################################
+# Copyright (C) 2018, 2019, 2020 Dominic O'Kane
+########################################################################################
+
+import numpy as np
+
+from financepy.market.curves.flat_discount_curve import FlatDiscountCurve
+from financepy.market.volatility.equity_vol_curve import EquityVolCurve
+from financepy.products.fx.fx_variance_swap import FinFXVarianceSwap
+from financepy.utils.date import Date
+
+value_dt = Date(1, 1, 2026)
+maturity_dt = Date(1, 1, 2027)
+t_exp = (maturity_dt - value_dt) / 365.0
+sigma = 0.12
+spot_fx_rate = 1.10
+strikes = np.linspace(0.3, 2.5, 221)
+
+########################################################################################
+
+
+def _fair_strike(r_dom, r_for, num_call_options=100, num_put_options=100, spacing=0.005):
+    domestic_curve = FlatDiscountCurve(value_dt, r_dom)
+    foreign_curve = FlatDiscountCurve(value_dt, r_for)
+    vol_curve = EquityVolCurve(
+        strikes, np.full(len(strikes), sigma), spot_fx_rate, t_exp, r_dom, r_for
+    )
+    swap = FinFXVarianceSwap(value_dt, maturity_dt, sigma**2)
+    fair_var = swap.fair_strike(
+        value_dt,
+        spot_fx_rate,
+        foreign_curve,
+        vol_curve,
+        num_call_options,
+        num_put_options,
+        spacing,
+        domestic_curve,
+    )
+    return swap, fair_var
+
+
+def test_fair_strike_equals_flat_variance_for_any_rates():
+    """Under a flat volatility the fair variance strike is the squared volatility
+    whatever the domestic and foreign rates. The foreign rate plays the role of a
+    dividend yield: the log-contract drift is r_d - r_f and the option portfolio
+    compounds at r_d."""
+    for r_dom, r_for in [(0.03, 0.0), (0.03, 0.05), (0.01, 0.04), (0.04, 0.04), (0.0, 0.0)]:
+        _, fair_var = _fair_strike(r_dom, r_for)
+        assert np.isfinite(fair_var), (r_dom, r_for)
+        assert abs(fair_var - sigma**2) / sigma**2 < 0.01, (r_dom, r_for, fair_var)
+
+
+def test_fair_strike_with_more_puts_than_positive_strikes():
+    """A put grid that would reach below zero is truncated to positive strikes and
+    the weights and pricing loop use the reduced count."""
+    swap, fair_var = _fair_strike(0.03, 0.05, num_put_options=300)
+    assert swap.num_put_options < 300
+    assert len(swap.put_wts) == swap.num_put_options
+    assert np.all(np.asarray(swap.put_strikes) > 0.0)
+    assert abs(fair_var - sigma**2) / sigma**2 < 0.01
+
+
+def test_realised_variance_of_constant_log_returns():
+    """A series growing by the same log return every day has a realised variance of
+    252 times that squared return."""
+    swap = FinFXVarianceSwap(value_dt, maturity_dt, sigma**2)
+    daily = 0.01
+    closes = spot_fx_rate * np.exp(daily * np.arange(11))
+    assert np.isclose(swap.realised_variance(closes), 252.0 * daily**2)