こんにちは。システムトレーダーの卵ことKenKenです。本日も「トレンドのモデル化」について、「Python3ではじめるシステムトレード ──環境構築と売買戦略」で勉強したことやそれに関連することを調べたのでまとめていきたいと思います。
前回の記事では、1949年から現在までの期間で日経平均株価のトレンドのモデル化を試みました。結果としては、決定係数の高い回帰係数は作ることができたものの、残差項の分布に着目してみると、短い期間でトレンドが発生している様子が伺えました。したがって今回は、期間を複数に分けて、期間ごとの日経平均株価のトレンドをモデル化してみようと思います。期間の分け方は、以前も紹介した，内閣府が公表している景気の分類をベースに以下のように分けたものを使用します。
| 循環期 | 景気 | 期間－始点 | 終点 |
| 1, 2 | 戦後復興期(recover) | 1949/5/16 | 1954/11/30 |
| 3, 4, 5, 6 | 高度経済成長気(growth) | 1954/12/1 | 1971/12/31 |
| 7, 8, 8, 10 | 安定期(stable) | 1972/1/1 | 1986/11/30 |
| 11 | バブル経済期(bubble) | 1986/12/1 | 1993/10/31 |
| 12, 13, 14, 15,16 | 経済変革期(reform) | 1993/11/1 | 2019/11/30 |
| 独自設定 | 新型コロナウィルス期(covid-19) | 2019/12/1 | 現在(2020/11/25) |
戦後復興期
import pandas_datareader.data as pdr
import statsmodels.api as sm
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
# データ取得
N225 = pdr.DataReader('NIKKEI225', 'fred', '1949/5/16').dropna()
ln_N225 = np.log(N225)
# 戦後復興期(～1954年)
y = ln_N225[:'1954/11/30'] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
決定係数は0.762と全期間の時の値(0.756)より改善している。また、p－値から回帰係数と切片についてゼロであるという帰無仮説を棄却するため、この2つの推定値については意味がある。グラフを見てみると、山と谷が一致していない為、さらに期間を分割する必要がありそうだ。
高度経済成長期
# 高度経済成長期(1954年～1971年)
y = ln_N225['1954/12/1':'1971/12/31'] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
証券不況の間、日経平均株価は低迷していたことがわかる。上図をさらに、証券不況まで、証券不況下、その後の回復期間に分けることで時間トレンドを日経平均株価の動きで説明できそうである。
安定期
# 安定期(1972年～1986年)
y = ln_N225['1972/1/1':'1986/11/30'] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
決定係数は、0.91とかなり良好である。チャートを見てみると、安定した上昇トレンドがわかる。
バブル成長期
# バブル成長期(1986年～1993年)
y = ln_N225['1986/12/1':'1993/10/31'] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
決定係数は、0.215とこれまでで一番低い。バブルによる上昇相場とバブル終戦による暴落相場を含んでいるため、低くなっているようである。チャートを見てみると、バブルのピーク前後でトレンド形成の仕方が異なるのがわかる。したがって、これらの期間を分割する必要があるようだ。そこで二つを分けてみる。
バブル成長期(日経平均株価ピークまで)
回帰係数は、これまでで一番高く強い上昇トレンドが在るように見える。グラフからもその様子が伺える。さらに、残差についてみてみる。
# 残差のヒストグラム
plt.figure()
ax = plt.subplot(1, 1, 1)
ax.hist(results.resid, bins=10, color='lightgray')
plt.ylabel('frequency')
残差は一つの分布に従い、平均・分散は一定しているように見える。
バブル暴落時(日経平均株価のピーク以降)
上記と同様の手順で、バブル暴落時の様子を見てみる。
# バブル成長期(ピークまで)
y = ln_N225['1990/1/1':'1992/8/31'] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
# 残差のヒストグラム
plt.figure()
ax = plt.subplot(1, 1, 1)
ax.hist(results.resid, bins=10, color='lightgray')
plt.ylabel('frequency')
決定係数は、0.816と悪くない。グラフを見ても、回帰モデルとしては良好に見える。しかし、残差についてみてみると、一つの分布に従っているとは思えない。
経済変革期
# 経済変革期(1993年～2019年11/30)
y = ln_N225['1993/11/1':'2019/11/30'] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
# 残差のヒストグラム
plt.figure()
ax = plt.subplot(1, 1, 1)
ax.hist(results.resid, bins=10, color='lightgray')
plt.ylabel('frequency')
結果はこれまでで一番ひどいものとなった。決定係数は、0.0と作成したモデルは一切説明できていないことを表している。チャートを見ても期間の取り方が適切でないことがわかる。残差のヒストグラムでもいくつかの分布から構成されていると考えられる。したがって、期間をさらに短くすることで、時間トレンドで日経平均株価を説明できる可能性はある。
新型コロナウイルス期
# 新型コロナウイルス期(2019年12/1～現在)
y = ln_N225['2019/12/1':] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
# 残差のヒストグラム
plt.figure()
ax = plt.subplot(1, 1, 1)
ax.hist(results.resid, bins=10, color='lightgray')
plt.ylabel('frequency')
この期間についても、先ほど同様、良い結果は得られなかった。チャートを見てみると、2020年3月19日までの大暴落とそれ以降で分割する必要がありそうである。
新型コロナウイルス期(暴落まで)
# 新型コロナウイルス期(2019年12/1～2020/3/19)
y = ln_N225['2019/12/1':'2020/3/19'] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
# 残差のヒストグラム
plt.figure()
ax = plt.subplot(1, 1, 1)
ax.hist(results.resid, bins=10, color='lightgray')
plt.ylabel('frequency')
決定係数は、0.48といまいち。チャートを見ても、暴落前と暴落時でさらに分割が必要そう。残差のヒストグラムも分布は一つではなく、平均・分散も安定して異なことがわかる。
新型コロナウイルス期(暴落後～現在)
# 新型コロナウイルス期(2020/3/19～現在)
y = ln_N225['2020/3/19':] # 被説明変数
x = range(len(y)) # 説明変数
x = sm.add_constant(x) # 切片
model = sm.OLS(y, x) # 最小二乗法
results = model.fit()
print(results.summary())
# グラフを図示
plt.figure()
ax = plt.subplot(1, 1, 1)
plt.plot(y, label='ln_N225', color='darkgray')
plt.plot(y.index, results.fittedvalues, label='predict')
plt.ylabel('log(N225 index)')
plt.legend()
# 残差のヒストグラム
plt.figure()
ax = plt.subplot(1, 1, 1)
ax.hist(results.resid, bins=10, color='lightgray')
plt.ylabel('frequency')
結果は良好であるように見える。決定係数は、0.815である。回帰係数、切片については、ゼロである帰無仮説を棄却し、意味のある推定値であることがわかる。また、回帰係数は9.87と高く、バブルのころに匹敵するほどの上昇トレンドが見られる。残差のヒストグラムもひとつの分布から形成されており、残差の平均・分散は一定している。したがって、安定した上昇トレンドが確認された。
まとめ
今回は、日経平均株価を景気循環における期間に分け、それぞれの期間におけるトレンドをモデル化してみた。期間幅を適切に設定することで、トレンドをモデル化することができることがよくわかった。また、決定係数が高くても、残差の分布についても確認する必要があった。
シストレ開発において、トレンドを戦略に入れる場合は、この辺の知識は役に立つと思われる。数日間、終値を更新した場合を上昇トレンドのような単純な戦略より、統計的なデータを組み入れていったほうが良いシステムが構築できるかもしれない。近々、統計的モデルに基づくトレンドを組み入れた戦略でトレーディングシステムの開発してみようかな。。。作成した場合は、また記事にまとめたいと思います。
以上
