Modified volatility predicts DJI returns

1from AlgorithmImports import *
2from pandas.core.frame import DataFrame
3import pandas as pd
4import statsmodels.api as sm
5# endregion
6
7class ModifiedVolatilityPredictsDJIReturns(QCAlgorithm):
8
9def Initialize(self):
10self.SetStartDate(2010, 1, 1)
11self.SetCash(100000)
12
13self.leverage:int = 3
14self.equity:Symbol = self.AddEquity("DIA", Resolution.Minute, leverage=self.leverage).Symbol
15self.cash:Symbol = self.AddEquity("BIL", Resolution.Minute, leverage=self.leverage).Symbol
16
17self.gamma:float = 3.
18self.estimation_period:int = 60
19self.std_period:int = 24
20self.months_in_year:float = 12.
21self.allocation_cap:List[float] = [-0.5, 1.5]
22
23self.recent_month:int = -1
24
25def OnData(self, data: Slice) -> None:
26if not(data.ContainsKey(self.equity) and data.ContainsKey(self.cash) and data[self.equity] and data[self.cash]):
27    return
28
29# monthly rebalance
30if self.Time.month == self.recent_month:
31    return
32self.recent_month = self.Time.month
33
34# get monthly closes
35period:int = (self.estimation_period + self.std_period) * 2 + 1
36history:DataFrame = self.History(self.equity, period * 31, Resolution.Daily)['close'].unstack(level=0)
37history = history.groupby(pd.Grouper(freq='M')).last()[self.equity]
38
39if len(history) >= period:
40    history = history.iloc[-period:]
41    returns:DataFrame = history.pct_change().iloc[1:]
42
43    # long-term volatility 
44    lv:DataFrame = returns.rolling(self.std_period).std() * np.sqrt(self.months_in_year)
45    lv = lv.dropna()
46    
47    # short-term volatility 
48    # SVAR is stock variance as in Welch and Goyal (2008)
49    svar:DataFrame = (returns ** 2).rolling(self.std_period).sum()
50    svar = svar.dropna()
51    sv:DataFrame = np.sqrt(svar / self.months_in_year)
52
53    lv_mean:DataFrame = lv.rolling(self.estimation_period).mean()
54    lv_mean = lv_mean.dropna()
55    sv_mean:DataFrame = sv.rolling(self.estimation_period).mean()
56    sv_mean = sv_mean.dropna()
57
58    # beta estimation
59    beta_factor:np.ndarray = np.array([((lv.iloc[i-self.estimation_period:i] - lv_mean.iloc[i]) * (sv.iloc[i-self.estimation_period:i] - sv_mean.iloc[i])).sum() for i in range(-1, -(len(lv_mean)+1), -1)])
60    beta_devisor:np.ndarray = np.array([((sv.iloc[i-self.estimation_period:i] - sv_mean.iloc[i]) ** 2).sum() for i in range(-1, -(len(lv_mean)+1), -1)])
61    beta:np.ndarray = beta_factor / beta_devisor
62    adj_lv:np.ndarray = (lv.iloc[-len(beta):].values - (beta * sv.iloc[-len(beta):].values))[-self.estimation_period:]
63
64    # regression to predict return
65    model = self.multiple_linear_regression(adj_lv[:-1], returns.iloc[-len(adj_lv):].values[1:])
66    ret_pred:float = model.predict(adj_lv[-1])
67    variance_pred:float = svar.iloc[-1] #lv.iloc[-1] ** 2
68
69    # allocation
70    equity_allocation:float = (ret_pred / (self.gamma * variance_pred))[-1]
71    equity_allocation = min(max(equity_allocation, min(self.allocation_cap)), max(self.allocation_cap))
72    cash_allocation:float = 1. - equity_allocation
73
74    # trade execution
75    self.SetHoldings(self.equity, equity_allocation)
76    self.SetHoldings(self.cash, cash_allocation)
77else:
78    if self.Portfolio.Invested:
79        self.Liquidate()
80
81def multiple_linear_regression(self, x:np.ndarray, y:np.ndarray):
82x:np.ndarray = np.array(x).T
83# x = sm.add_constant(x)
84result = sm.OLS(endog=y, exog=x).fit()
85return result