|
| 1 | + |
| 2 | + |
| 3 | +############### |
| 4 | +# Authored by Weisheng Jiang |
| 5 | +# Book 6 | From Basic Arithmetic to Machine Learning |
| 6 | +# Published and copyrighted by Tsinghua University Press |
| 7 | +# Beijing, China, 2022 |
| 8 | +############### |
| 9 | + |
| 10 | +import pandas as pd |
| 11 | +import numpy as np |
| 12 | +import matplotlib.pyplot as plt |
| 13 | +import pandas_datareader |
| 14 | +import scipy.stats as stats |
| 15 | +import pylab |
| 16 | + |
| 17 | +df_price = pandas_datareader.data.DataReader(['sp500'], |
| 18 | + data_source='fred', |
| 19 | + start='08-01-2018', end='08-01-2021') |
| 20 | +df_price = df_price.dropna() |
| 21 | + |
| 22 | +#%% Rolling max, min |
| 23 | + |
| 24 | +df_max_100 = df_price.rolling(100).max() |
| 25 | +df_min_100 = df_price.rolling(100).min() |
| 26 | + |
| 27 | +fig, ax = plt.subplots() |
| 28 | +# sp500 price |
| 29 | +ax.plot(df_price['sp500'], label = 'Price') |
| 30 | +ax.plot(df_max_100, label = 'Max, 100') |
| 31 | +ax.plot(df_min_100, label = 'Min, 100') |
| 32 | + |
| 33 | +plt.ylabel('Price level') |
| 34 | +plt.legend(loc='upper left') |
| 35 | +plt.xticks(rotation = 45) |
| 36 | +# Rotates X-Axis Ticks by 45-degrees |
| 37 | +fig.tight_layout() |
| 38 | + |
| 39 | +#%% Rolling mean |
| 40 | + |
| 41 | +df_mean_50 = df_price.rolling(50).mean() |
| 42 | +df_mean_100 = df_price.rolling(100).mean() |
| 43 | +df_mean_250 = df_price.rolling(250).mean() |
| 44 | + |
| 45 | +fig, ax = plt.subplots() |
| 46 | +# sp500 price |
| 47 | +ax.plot(df_price['sp500'], label = 'Price') |
| 48 | +ax.plot(df_mean_50, label = 'Mean, 50') |
| 49 | +ax.plot(df_mean_100,label = 'Mean, 100') |
| 50 | +ax.plot(df_mean_250,label = 'Mean, 150') |
| 51 | + |
| 52 | +plt.ylabel('Price level') |
| 53 | +plt.legend(loc='upper left') |
| 54 | +plt.xticks(rotation = 45) |
| 55 | +# Rotates X-Axis Ticks by 45-degrees |
| 56 | +fig.tight_layout() |
| 57 | + |
| 58 | + |
| 59 | + |
| 60 | +#%% daily log return |
| 61 | + |
| 62 | +daily_log_r = df_price.apply(lambda x: np.log(x) - np.log(x.shift(1))) |
| 63 | + |
| 64 | +daily_log_r = daily_log_r.dropna() |
| 65 | + |
| 66 | +#%% Rolling mean, skew, kurtosis on returns |
| 67 | + |
| 68 | +daily_log_r_mean = daily_log_r.rolling(50).mean() |
| 69 | +daily_log_r_std = daily_log_r.rolling(50).std() |
| 70 | +daily_log_r_skew = daily_log_r.rolling(50).skew() |
| 71 | +daily_log_r_kurt = daily_log_r.rolling(50).kurt() |
| 72 | + |
| 73 | +fig, axes = plt.subplots(5,1) |
| 74 | + |
| 75 | +# plot daily log returns |
| 76 | +axes[0].plot(daily_log_r[daily_log_r_mean.first_valid_index():daily_log_r_mean.index[-1]], |
| 77 | + marker='.', markersize = 3, |
| 78 | + color = 'b',linestyle='None', |
| 79 | + label = 'Daily log r') |
| 80 | +axes[0].axhline(y = 0, color='r', linestyle='-') |
| 81 | +axes[0].set_ylabel ('Log return') |
| 82 | +axes[0].set_xticks([]) |
| 83 | + |
| 84 | +# plot first moment, mean |
| 85 | +axes[1].plot(daily_log_r_mean,label = 'Mean', color = 'b') |
| 86 | +axes[1].axhline(y = 0, color='r', linestyle='-') |
| 87 | +axes[1].set_ylabel ('Mean') |
| 88 | +axes[1].set_xticks([]) |
| 89 | + |
| 90 | +# plot second moment, std |
| 91 | +axes[2].plot(daily_log_r_std,label = 'Std', color = 'b') |
| 92 | +axes[2].axhline(y = 0, color='r', linestyle='-') |
| 93 | +axes[2].set_ylabel ('Volatility') |
| 94 | +axes[2].set_xticks([]) |
| 95 | + |
| 96 | +# plot third moment, skew |
| 97 | +axes[3].plot(daily_log_r_skew,label = 'Skew', color = 'b') |
| 98 | +axes[3].axhline(y = 0, color='r', linestyle='-') |
| 99 | +axes[3].set_ylabel ('Skew') |
| 100 | +axes[3].set_xticks([]) |
| 101 | + |
| 102 | +# plot fourth moment, kurtosis |
| 103 | +axes[4].plot(daily_log_r_kurt,label = 'Kurtosis', color = 'b') |
| 104 | +axes[4].set_ylabel ('Kurt') |
| 105 | +axes[4].axhline(y = 3, color='r', linestyle='-') |
| 106 | + |
| 107 | +plt.xticks(rotation = 45) |
| 108 | +# Rotates X-Axis Ticks by 45-degrees |
| 109 | +fig.tight_layout() |
| 110 | + |
| 111 | + |
| 112 | +#%% Moving quantile |
| 113 | + |
| 114 | +alpha_95 = 0.95 |
| 115 | + |
| 116 | +df_95_percent = daily_log_r.rolling(100).quantile(alpha_95) |
| 117 | +df_05_percent = daily_log_r.rolling(100).quantile(1 - alpha_95) |
| 118 | + |
| 119 | +fig, ax = plt.subplots() |
| 120 | +ax.plot(daily_log_r,marker='.', color = 'b',linestyle='None') |
| 121 | +ax.plot(df_95_percent,label = '95% percentile') |
| 122 | +ax.plot(df_05_percent,label = '5% percentile') |
| 123 | + |
| 124 | +plt.axhline(y=0, color='r', linestyle='-') |
| 125 | + |
| 126 | +plt.xlabel('Date') |
| 127 | +plt.ylabel('Daily log return') |
| 128 | +plt.xticks(rotation = 45) # Rotates X-Axis Ticks by 45-degrees |
| 129 | +fig.tight_layout() |
| 130 | +plt.legend(loc='upper left') |
| 131 | + |
| 132 | +#%% Moving average volatility |
| 133 | + |
| 134 | +MA_vol_50 = daily_log_r.rolling(50).std() |
| 135 | +MA_vol_100 = daily_log_r.rolling(100).std() |
| 136 | +MA_vol_250 = daily_log_r.rolling(250).std() |
| 137 | + |
| 138 | +# plot daily log returns |
| 139 | +fig, axes = plt.subplots(2,1) |
| 140 | +# sp500 daily log returns |
| 141 | +axes[0].plot((daily_log_r**2)[MA_vol_250.first_valid_index():MA_vol_250.index[-1]], |
| 142 | + color = 'b') |
| 143 | + |
| 144 | +axes[0].set_xticks([]) |
| 145 | +axes[0].axhline(y=0, color='r', linestyle='-') |
| 146 | +axes[0].set_ylabel('Daily log return squared') |
| 147 | + |
| 148 | +# Moving average volatility |
| 149 | +axes[1].plot(MA_vol_50[MA_vol_250.first_valid_index():MA_vol_250.index[-1]], label = 'Window = 50') |
| 150 | +axes[1].plot(MA_vol_100[MA_vol_250.first_valid_index():MA_vol_250.index[-1]],label = 'Window = 100') |
| 151 | +axes[1].plot(MA_vol_250[MA_vol_250.first_valid_index():MA_vol_250.index[-1]],label = 'Window = 250') |
| 152 | +axes[1].set_xlabel("Date") |
| 153 | +plt.xticks(rotation = 45) # Rotates X-Axis Ticks by 45-degrees |
| 154 | +axes[1].set_ylabel("MA vol") |
| 155 | +fig.tight_layout() |
| 156 | +plt.legend(loc='upper left') |
| 157 | + |
| 158 | +list_df = [MA_vol_50, MA_vol_100, MA_vol_250] |
| 159 | + |
| 160 | +for data in zip(list_df): |
| 161 | + |
| 162 | + data = data[0] |
| 163 | + |
| 164 | + fig, ax = plt.subplots() |
| 165 | + |
| 166 | + # daily return of selected date range |
| 167 | + plt.plot(daily_log_r[data.first_valid_index():data.index[-1]], |
| 168 | + marker='.', markersize = 3, |
| 169 | + color = 'b',linestyle='None') |
| 170 | + |
| 171 | + upper_bound = 2*data[data.first_valid_index():data.index[-1]]; |
| 172 | + lower_bound = -upper_bound; |
| 173 | + |
| 174 | + ax.plot(upper_bound, color = 'r') |
| 175 | + ax.plot(lower_bound, color = 'r') |
| 176 | + ax.fill_between(upper_bound.index, upper_bound['sp500'], |
| 177 | + lower_bound['sp500'], color = '#DBEEF3') |
| 178 | + ax.axhline(y=0, color='k', linestyle='-') |
| 179 | + |
| 180 | + ax.set_ylabel('Daily log return') |
| 181 | + plt.xticks(rotation = 45) # Rotates X-Axis Ticks by 45-degrees |
| 182 | + plt.tight_layout() |
| 183 | + |
| 184 | +#%% squared root of time |
| 185 | + |
| 186 | +# plot daily log returns |
| 187 | +fig, ax = plt.subplots() |
| 188 | + |
| 189 | +# Moving average volatility |
| 190 | +ax.plot(MA_vol_50*np.sqrt(250)*100, label = 'Window = 50') |
| 191 | +ax.plot(MA_vol_100*np.sqrt(250)*100, label = 'Window = 100') |
| 192 | +ax.plot(MA_vol_250*np.sqrt(250)*100, label = 'Window = 250') |
| 193 | +ax.set_xlabel("Date") |
| 194 | +plt.xticks(rotation = 45) # Rotates X-Axis Ticks by 45-degrees |
| 195 | +ax.set_ylabel("MA annualized volatility, %") |
| 196 | +fig.tight_layout() |
| 197 | +plt.legend(loc='upper left') |
| 198 | + |
| 199 | +#%% EWMA mean |
| 200 | + |
| 201 | +df_mean_99 = df_price.ewm(alpha=0.01, adjust=False).mean() |
| 202 | +df_mean_975 = df_price.ewm(alpha=0.025, adjust=False).mean() |
| 203 | +df_mean_94 = df_price.ewm(alpha=0.06, adjust=False).mean() |
| 204 | + |
| 205 | +fig, ax = plt.subplots() |
| 206 | +# sp500 price |
| 207 | +ax.plot(df_price['sp500']) |
| 208 | +ax.plot(df_mean_99, label = 'Lambda = 0.99') |
| 209 | +ax.plot(df_mean_975,label = 'Lambda = 0.975') |
| 210 | +ax.plot(df_mean_94, label = 'Lambda = 0.94') |
| 211 | + |
| 212 | +plt.ylabel('Price level') |
| 213 | +plt.legend(loc='upper left') |
| 214 | +plt.xticks(rotation = 45) |
| 215 | +# Rotates X-Axis Ticks by 45-degrees |
| 216 | +fig.tight_layout() |
| 217 | +plt.legend(loc='upper left') |
| 218 | + |
| 219 | + |
| 220 | +#%% EWMA rolling volatility |
| 221 | + |
| 222 | +EWMA_vol_99 = daily_log_r.ewm(alpha=0.01, adjust=False).std() |
| 223 | +EWMA_vol_975 = daily_log_r.ewm(alpha=0.025, adjust=False).std() |
| 224 | +EWMA_vol_94 = daily_log_r.ewm(alpha=0.06, adjust=False).std() |
| 225 | + |
| 226 | +# plot daily log returns |
| 227 | +fig, axes = plt.subplots(2,1) |
| 228 | +# sp500 daily log returns |
| 229 | +axes[0].plot((daily_log_r**2)[MA_vol_250.first_valid_index():MA_vol_250.index[-1]], |
| 230 | + color = 'b') |
| 231 | + |
| 232 | +axes[0].set_xticks([]) |
| 233 | +axes[0].axhline(y=0, color='r', linestyle='-') |
| 234 | +axes[0].set_ylabel('Daily log return') |
| 235 | + |
| 236 | +# Moving average volatility |
| 237 | + |
| 238 | +axes[1].plot(EWMA_vol_99[MA_vol_250.first_valid_index():MA_vol_250.index[-1]], label = 'Lambda = 0.99') |
| 239 | +axes[1].plot(EWMA_vol_975[MA_vol_250.first_valid_index():MA_vol_250.index[-1]],label = 'Lambda = 0.975') |
| 240 | +axes[1].plot(EWMA_vol_94[MA_vol_250.first_valid_index():MA_vol_250.index[-1]], label = 'Lambda = 0.94') |
| 241 | +axes[1].set_xlabel("Date") |
| 242 | +plt.xticks(rotation = 45) # Rotates X-Axis Ticks by 45-degrees |
| 243 | +axes[1].set_ylabel("EWMA volatility") |
| 244 | +fig.tight_layout() |
| 245 | +plt.legend(loc='upper left') |
| 246 | + |
| 247 | +list_df = [EWMA_vol_99, EWMA_vol_975, EWMA_vol_94] |
| 248 | + |
| 249 | +for data in zip(list_df): |
| 250 | + |
| 251 | + data = data[0] |
| 252 | + |
| 253 | + fig, ax = plt.subplots() |
| 254 | + |
| 255 | + # daily return of selected date range |
| 256 | + plt.plot(daily_log_r[data.first_valid_index():data.index[-1]], |
| 257 | + marker='.', markersize = 3, |
| 258 | + color = 'b',linestyle='None') |
| 259 | + |
| 260 | + upper_bound = 2*data[data.first_valid_index():data.index[-1]]; |
| 261 | + lower_bound = -upper_bound; |
| 262 | + |
| 263 | + ax.plot(upper_bound, color = 'r') |
| 264 | + ax.plot(lower_bound, color = 'r') |
| 265 | + ax.fill_between(upper_bound.index, upper_bound['sp500'], |
| 266 | + lower_bound['sp500'], color = '#DBEEF3') |
| 267 | + ax.axhline(y=0, color='k', linestyle='-') |
| 268 | + |
| 269 | + ax.set_ylabel('Daily log return') |
| 270 | + plt.xticks(rotation = 45) # Rotates X-Axis Ticks by 45-degrees |
| 271 | + plt.tight_layout() |
| 272 | + |
| 273 | +#%% Compare volatilities |
| 274 | + |
| 275 | +fig, ax = plt.subplots() |
| 276 | +plt.plot(MA_vol_50*np.sqrt(250)*100, label = 'MA, 50', linestyle = '--') |
| 277 | +plt.plot(MA_vol_100*np.sqrt(250)*100, label = 'MA, 100', linestyle = '--') |
| 278 | +plt.plot(MA_vol_250*np.sqrt(250)*100, label = 'MA, 250', linestyle = '--') |
| 279 | +plt.plot(EWMA_vol_99*np.sqrt(250)*100, label = 'EWMA, $\lambda$ = 0.99') |
| 280 | +plt.plot(EWMA_vol_975*np.sqrt(250)*100,label = 'EWMA, $\lambda$ = 0.975') |
| 281 | +plt.plot(EWMA_vol_94*np.sqrt(250)*100, label = 'EWMA, $\lambda$ = 0.94') |
| 282 | + |
| 283 | +ax.set_xlabel("Date") |
| 284 | +plt.xticks(rotation = 45) |
| 285 | +# Rotates X-Axis Ticks by 45-degrees |
| 286 | +ax.set_ylabel("Annualized volatility, %") |
| 287 | +fig.tight_layout() |
| 288 | +plt.legend(loc='upper left') |
0 commit comments