Started writing tests for calibration of historical volatility.

quantdsl/priceprocess/blackscholes.py

@@ -2,10 +2,11 @@
 
 import datetime
 import math
+from collections import defaultdict
 
 import scipy
 import scipy.linalg
-from collections import defaultdict
+from pandas_datareader.data import DataReader
 from scipy.linalg import LinAlgError
 
 from quantdsl.exceptions import DslError
@@ -150,11 +151,47 @@ def get_correlation_from_calibration(self, market_calibration, name_i, name_j):
         return correlation
 
 
-def calc_sigma(curve):
+def calc_sigma_from_curve(curve):
     dates = scipy.array([i[0] for i in curve])
     prices = scipy.array([i[1] for i in curve])
+    return historical_volatility_std_over_mean(prices, dates)
+
+
+def historical_volatility_std_over_mean(prices, dates):
     volatility = prices.std() / prices.mean()
-    duration = max(dates) - min(dates)
-    years = duration.total_seconds() / datetime.timedelta(365).total_seconds()
+    max_date = max(dates)
+    min_date = min(dates)
+    if isinstance(max_date, np.datetime64):
+        max_date = max_date.astype('O') / 1e9
+        min_date = min_date.astype('O') / 1e9
+    duration = max_date - min_date
+    if isinstance(duration, datetime.timedelta):
+        duration_seconds = duration.total_seconds()
+    else:
+        duration_seconds = duration
+    years = duration_seconds / datetime.timedelta(365).total_seconds()
     assert years > 0, "Can't calculate volatility for price series with zero days duration"
     return float(volatility) / math.sqrt(years)
+
+
+import numpy as np
+
+
+def historical_volatility(sym, days):
+    quotes = get_yahoo_data(sym, days)
+    "Return the annualized stddev of daily log returns of `sym`."
+    # return historical_volatility_log_returns(quotes)
+    return historical_volatility_std_over_mean(quotes.values, list(quotes.index.values))
+
+
+def historical_volatility_log_returns(quotes):
+    logreturns = np.log(quotes / quotes.shift(1))
+    return np.sqrt(252 * logreturns.var())
+
+
+def get_yahoo_data(sym, days):
+    try:
+        return DataReader(sym, 'yahoo')['Close'][-days:]
+    except Exception as e:
+        msg = "Error getting data for symbol '{}': {}".format(sym, e)
+        raise Exception(msg)

quantdsl/tests/test_calc_and_plot.py

@@ -2,8 +2,9 @@
 
 from eventsourcing.domain.model.events import assert_event_handlers_empty
 
-from quantdsl.exceptions import TimeoutError, CallLimitError
-from quantdsl.interfaces.calcandplot import calc_print
+from quantdsl.application.base import Results
+from quantdsl.exceptions import CallLimitError
+from quantdsl.interfaces.calcandplot import calc_print, calc
 
 
 class TestCalcPrint(TestCase):
@@ -209,3 +210,51 @@ def test_dependency_graph_size_limit(self):
                 },
                 max_dependency_graph_size=1,
             )
+
+    def test_results_dataframes(self):
+        source_code = """
+from quantdsl.lib.storage2 import GasStorage        
+GasStorage(Date('2011-1-1'), Date('2011-4-1'), 'GAS', 0, 0, 50000, TimeDelta('1m'), 1)
+        """
+
+        results = calc(
+            source_code=source_code,
+            observation_date='2011-1-1',
+            interest_rate=2.5,
+            periodisation='monthly',
+            price_process={
+                'name': 'quantdsl.priceprocess.blackscholes.BlackScholesPriceProcess',
+                'market': ['GAS'],
+                'sigma': [0.5],
+                'curve': {
+                    'GAS': (
+                        ('2011-1-1', 13.5),
+                        ('2011-2-1', 16.5),
+                        ('2011-3-1', 19.5),
+                        ('2011-4-1', 17.5),
+                    )
+                }
+            },
+        )
+        assert isinstance(results, Results)
+        # self.assertIsInstance(results.cash_mean, DataFrame)
+
+
+        import matplotlib.pyplot as plt
+        #
+        plt.ioff()
+
+        results.prices_mean.plot(title='Prices')
+
+        plt.pause(1)
+
+        results.hedges_mean.plot(title='Hedges')
+
+        plt.pause(1)
+
+        results.cash.plot(title='Cash')
+
+        plt.show()
+
+        # import matplotlib; import matplotlib.pyplot; print(matplotlib.backends.backend)
+

quantdsl/tests/test_price_processes.py

@@ -3,7 +3,8 @@
 
 import scipy
 
-from quantdsl.priceprocess.blackscholes import BlackScholesPriceProcess, calc_sigma
+from quantdsl.priceprocess.blackscholes import BlackScholesPriceProcess, calc_sigma_from_curve, historical_volatility, \
+    get_yahoo_data, historical_volatility_std_over_mean, historical_volatility_log_returns
 
 
 class TestBlackScholesPriceProcess(unittest.TestCase):
@@ -198,10 +199,20 @@ def test_simulate_future_prices_from_longer_curve(self):
 
 class TestCalcSigma(unittest.TestCase):
     def test(self):
-        sigma = calc_sigma([
+        sigma = calc_sigma_from_curve([
             (datetime.datetime(2011, 1, 1), 10),
             (datetime.datetime(2011, 2, 1), 11),
             (datetime.datetime(2011, 3, 1), 9),
             (datetime.datetime(2011, 4, 1), 10),
         ])
         self.assertAlmostEqual(sigma, 0.14, places=2)
+
+
+class TestHistoricalVolatility(unittest.TestCase):
+    def test(self):
+        quotes = get_yahoo_data('GOOG', 30)
+        "Return the annualized stddev of daily log returns of `sym`."
+        vol_log_returns = historical_volatility_log_returns(quotes)
+        vol_std_over_mean = historical_volatility_std_over_mean(quotes.values, list(quotes.index.values))
+
+        raise Exception((vol_log_returns, vol_std_over_mean))