確率的トレンドとは何か?

一般にトレンドには何か発生の原因があるのではないかとかんえられがちである。しかし、乱数の和からもでもトレンドは発生する。このような動きを確率的トレンドといい、ランダムウォークの特徴です。

1900年にフランスの数学者、ルイ・バシュリエは、最初にランダムウォーク理論を提唱 し、株価が酔っ払いの千鳥足のように動く様子を示した。 株式、債券、通貨などの取引が活発な金融商品の価格は、ランダム

ウォークとしてモデル化される。 価格の動きがランダムである理由は、それが本質的にランダムであったり、外生的なショックがランダムに発生 し、それが価格に影響を与えているとか、または金融市場の効率性を反映し

た結果であるとか、さまざまな考え方がある。ランダムウォークはこのランダムな価格の動きが累積した結果です。

トレンドの発生に理由はいらない!  

確率過程  

ランダムに発生した変数を時間の順に並べた集まりを確率過程とよびます。 この確率過程は、株式、債券、または通貨などの価格の動きをモデル化するため利用される。

時系列データを扱う際には、得られたデータを実現値といいあす。なぜ標本とよばずに実現値とよぶかというと、それは特定の確率過程から抽出され、そして実現した特別な一連の値だからです。

この実現値がある区間内の任意の実数であれば、それを連続型の確率変数といいます。とびとびの値を取る場合には離散型の確率変数といいます。実際の金融関連の価格は、とびとびの値を取るので離散型の確率変数です。

また金融市場で取引可能な価格を提示したり、実際の取引を行う行為は連続ではなく、これもとびとびに起こるのでこのような確率過程を離散時間の確率過程といいます。一方、価格の理論を構築するような場合には、連続

時間の確率過程が便利です。ブラウン運動とかウィーナー過程は連続時間の確率過程です。

これらの確率過程を定常な確率過程と非定常な確率過程に分けることができる。 定常確率過程とは、確率変数の値が初期値から大きく離れることなく、その周辺をうろうろとする状態を継続することができる確率過程です。

その際に、確率変数の平均、分散、共分散は時間に対して一定である。一方、非定常確率過程とは時間の経過と共に、確率変数の平均、分散、共分散などが変化する確率過程である。

確定的トレンドは非定常確率過程の例です。

確率変数の値が時間の経過と共に初期値から徐々に離れて行きトレンドを形成しながら、その周辺をうろうろするからである。そのために確率変数の平均値は時間の経過と共に変化するが、分散は変化しない。

この確率過程から確定的トレンドを取り除くと、残った確率過程は定常過程となる。従って確定的トレンドはトレンド定常性ともよば れる。

また、トレンドを形成し初期値から離れた後に、再び原点に復帰する再帰性をともなう確率過程もある。ランダムウォークはこのような確率過程の1つです。

ランダムウォーク  

資産価格の分散が時間とともに変化するのであれば、その確率過程はランダムウォークとしてモデル化できる。 ランダムウォークには幾つかのバリエーションがあり、ここでは2つ紹介します。  

ドリフト無しランダムウォーク  

ドリフト無しランダムウォークは確率変数の和として定義できます。

t時の金融資産の価格を$W_{t}$とすると  

$$ W_{t}=W_{0}+\sum_{i=1}^{t} w_{i} $$

となる。

ここで$w_{i}$は平均ゼロ、分散$\sigma ^{2}$の正規分布に従う確率変数であり、かく乱項とよばれる。 このような確率過程においては1つ1つの w の影響は減衰することなく永久に残り続ける(permanent effect)。

そしてその効果が確率的トレンドを発生させる。

このモデルは  

$$ W_{t}=W_{t-1}w_{i} $$

と書き直すこともできる。

ここで w の期待値はゼロ$E(w_{t})= 0)$である。

従って、$W_{t}$(価格)の分散は  

$$var(w_{t}) = E(w^{2}) - (E(w))^{2} = (E(w^{2})) = \sigma ^{2} $$

となります。

よって、$w_{t}$の平均、分散は  

$$E(W_{t}) = W_{0}$$

 

$$var(W_{t}) = tvar(w_{t}) = t\sigma ^{2} $$

となります。期待値は定数で、分散は時間の関数に成り立っています。  

ドリフト付きモデル  

ドリフト無しランダムウォークのモデルに切片を加えたモデルがドリフト付きモデルです。   $$W_{t} = a + W_{t-1} + w_{t}$$

ここで a は定数でありドリフト率とよばれる。このモデルの特性は

$$W_{t} = at + W_{0} + \sum_{i=1}^{t} w_{i}$$

であらわされる。ドリフト項がトレンドの役割を担っている。

$W_{t}$の期待値と分散は、

$$E(W_{t}) = W_{0} + at$$

$$var(W_{t}) = t var(W_{t}) = t\sigma ^{2} $$

です。

どちらも時間の関数である。また、Wの差分の期待値は  

$$ E(W_{t}-W_{t-1}) = E(\bigtriangleup W_{t}) = E(a) + E(W_{t}) = a $$

となります。

確率的トレンドとドリフト

確率的トレンドは$W_{1}$,…,+$W_{i}$=$ΣW_{i}$の結果として生じるトレンドです。

$W_{i}$は平均がゼロ、分散が$ \sigma ^{2}$の正規分布に従う確率変数として定義される。

従って、$W_{i}$と$W_{i+j}$ , j =1,… との相関はゼロである。

価格が上昇の後に続けて上昇したとしても、その後に下落をして、価格が最初の値に戻る必要はない。

現在の価格の動きは過去の価格の動きとは独立です。

従って、上昇と下落を繰り返しながら上昇トレンドと、下落トレンドの発生が可能になる。 このトレンドはランダムな現象の結果であるので、この価格の動

きから売買のタイミングを図ることで利益を上げることはできない。  

しかし、この確率的トレンドを最小二乗法で線形回帰を行うと、確定的トレンドと判断されてしまうことがある。 これを見せかけの確定的トレンドという。  

ドリフト付きランダムウォークの場合はどうであろうか

$$W_{t} = W_{0} + at + \sum_{i=1}^{t} x_{w}^{2}$$

上式では、残差項$W_{i}$が平均ゼロ、分散$aw^{2}$と定められているが、それはドリフト率が存在するからです。

ドリフト率をゼロとしてそれと同じ効果をもつ平均値を$w_{i}$がもてば同じ確率過程を生成できます。

このドリフトから発生するトレンドは確定的トレンドか? 

それとも確率的トレンドか? 

ドリフト付きランダムウォークは時間の関数となっている。しかし、ここで注意が必要です。このドリフト率の役割は平均値をもつ残差項の平均を表現し

ているのです。従って、必ずしも単位時間当たり、aが加算されるのではなく、平均として単位時間当たりaが加算されていくので、確率過程の一部であ

る。従ってドリフト付きランダムウォークが作るトレンドは確率的トレンドです。  

ランダムウォークの判定  

与えらえた時系列のデータがランダムウォークであるかどうかの判定問題は、長い歴史をもち、経済時系列解析の中でも最も注目され続けている課題の内

の1つです。それを単位根問題という。定常過程においては、母集団から抽出された標本は、その標本の数が増えれば標本平均は真の平均に近づき、こ

の差が従う分布は標本数の増加に従い正規分布に従うという中心極限定理が成り立つ。しかし、単位根をもつ非定常な過程ではこの中心極限定理が成り立

たない。それが問題となるのです。そのために膨大な研究結果が報告されている。そのような問題には一切触れることなく、ADF検定の方法のみを説明します。

ADFとは拡張ディッキー・フラー検定(Augmented Dicky-Fuller Test)の略称です。

単位根  

金融資産の価格(W)を線形回帰式  

$$W_{t} = KW_{t-1} + W_{t}$$

で表現し、K=1であるならば、このモデルは単位根をもつという。

単位根をもつとその時系列は非定常確率過程になり、分散は時間の経過と共に発散する。このような非定常な確率過程に対しては最小二乗法が使えない。し

たがって、もし1次の階差$( \bigtriangleup W_{t} = W_{t} - W{t-1})$を取って時系列$W_{t}$が定常になれば$W_{t}$はランダムウォークであり、最小二乗法がを使えます。

そこで  

$$( \bigtriangleup W_{t} = W_{t} - W{t-1})$$

というモデルを立て、この回帰係数を最小二乗法で推定し、γ= 0であれば、

上式は$W_{t} = γW_{t-1} - W_{t}$となり、$W_{t}$がランダムウォークである可能性ある。

従って帰無仮説をγ= 0として仮説検定を行えばよいのです。この帰無仮説が棄却できなければ時系列はランダムウォークである可能性があります。このよう

な検定の仕方をディッキー・フラー(DF)検定といいます。  

拡張ディッキー・フラー検定  

ADF検定はDF検定を次のように拡張したものです。

$$ \bigtriangleup W{t} = a+ \beta t+ \gamma W{t-1} + \sum_{i=1}^{t} \sigma _{i}W_{t-1} + w_{t} $$

I = 0とすれば、DF検定になります。

statsmodelでは

$$1. \bigtriangleup W{t} = \gamma W{t-1} + \sum_{i=1}^{t} \sigma _{i}\bigtriangleup W_{t-1} + w_{t} $$

$$2. \bigtriangleup W{t} = a+ \gamma W{t-1} + \sum_{i=1}^{t} \sigma _{i}\bigtriangleup W_{t-1} + w_{t} $$

$$3. \bigtriangleup W{t} = a+ \beta t + \gamma W{t-1} + \sum_{i=1}^{t} \sigma _{i}\bigtriangleup W_{t-1} + w_{t} $$

とさらに4.ドリフト付き、時間トレンド付き、かつ時間の2乗に比例するトレンド付きモデルについて検定を行うことができる。

帰無仮説を棄却できなければデータはランダムウォークと判定される。

In [86]:
%matplotlib inline
import matplotlib.pyplot as plt
import statsmodels.api as sm
import pandas as pd
import pandas_datareader.data as web
import numpy as np
end='2017/7/31'
lnn225 = np.log(pdr.DataReader("NIKKEI225", 'fred',"1949/5/16",end)).dropna()
In [42]:
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='nc')[0])#検定統計量
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='nc')[1])#p-値
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='nc')[2])#ラグの数
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='nc')[3])#データの数
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='nc')[4])#臨界値
2.17369686757
0.994169731135
38
16730
{'1%': -2.5658736531258919, '5%': -1.9410160670068219, '10%': -1.6168041340187047}

結果の解釈

日経平均株価の単位根検定の結果は単純ランダムウォークではバブルのピークから谷を除くとすべての期間で帰無仮設を棄却できない。

従って、単純ランダムウォークという結果になった。しかし、最小二乗法でkを推定するとバブルのピークから谷と経済変革期をの除くとプラスであるから、

単純ランダムウォークの可能性のあるのは経済変革期だけである。経済変化期の日経平均株価の動きには確定トレンドは無く、トレンドに見えるすべての部分

が確率的現象である可能性があることをこの結果はあらわしています。

プログラムコードの解説: ADF検定

statsmodels.tsa. Stattools.adfuller(x,maxlag=None,regression='c',auto='AIC',store=False, regreults=False)

x:時系列データ maxlag:Iの数を指定

regression='c'がデフォルトであり、cは上述の式2に相当し、ドリフト付きランダムウォークです。

ncはドリフト無しランダムウォーク(式1)である。ctはドリフト付ランダムウォーク(式3)です。

auto=AIC,BIC,t-stat,無しの4つの選択が可能

store:結果の保存

regression:デフォルトは False(偽)、Trueはすべての結果を返す。

プログラムコードの解説

z=lnn225

対数価格のデータの作成

y=z.diff().dropna()

対数収益率のデータの作成

x=z.shift(1).dropna()

1日前の対数価格からなる時系列データの作成

model=sm.OLS(y,x)

Statsmodelの線形回帰のモデルを使用

results=model.fit()

ritを用いてモデルを最適化

print("without drift ",results.params[0])

回帰係数の出力

x=sm.add_constant(x)

モデルにドリフト項(切片)を加えるために準備

model=sm.OLS(y,x)

results=model.fit()

print("with drift ",results.params[0],results.params[1])

[0]はドリフト率、[1]は回帰係数

x["t"]=range(len(y))

時間経過の時系列を作成

model=sm.OLS(y,x)

results=model.fit()

print("with drift + time trend ",results.params[0],results.params[1],results.params[2])

[0]はドリフト率の係数、[1]は回帰係数、[2]は時間トレンドの係数

In [43]:
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='ct')[1])
0.845163622052
In [44]:
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='c')[1])
0.0888254028253
In [45]:
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='c')[3])
16730
In [46]:
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='c')[4])
{'1%': -3.4307409326737046, '5%': -2.861712776642042, '10%': -2.5668619646086057}
In [47]:
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='c')[:2])
(-2.6205426848019742, 0.088825402825320532)
In [48]:
print(sm.tsa.adfuller(lnn225.NIKKEI225,regression='nc')[1])
0.994169731135
In [84]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix[:'1954/11/30'],regression='ct')[1])
0.504338667965
In [50]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix[:'1954/11/30'],regression='c')[1])
0.903086924527
In [51]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix[:'1954/11/30'],regression='nc')[1])
0.895822001571
In [52]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1954/12/1':'1971/12/31'],regression='ct')[1])
0.615588208256
In [53]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1954/12/1':'1971/12/31'],regression='c')[1])
0.480043863331
In [54]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1954/12/1':'1971/12/31'],regression='nc')[1])
0.999689754599
In [55]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1972/1/1':'1986/11/30'],regression='ct')[1])
0.840240937298
In [56]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1972/1/1':'1986/11/30'],regression='c')[1])
0.96175537736
In [57]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1972/1/1':'1986/11/30'],regression='nc')[1])
0.999996991887
In [58]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1986/12/1':'1989/12/31'],regression='ct')[1])
0.312378263486
In [59]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1986/12/1':'1989/12/31'],regression='c')[1])
0.768377939279
In [60]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1986/12/1':'1989/12/31'],regression='nc')[1])
0.999525130045
In [61]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1986/12/1':'1993/10/30'],regression='ct')[1])
0.408168347822
In [62]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1986/12/1':'1993/10/30'],regression='ct')[1])
0.408168347822
In [63]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1986/12/1':'1993/10/30'],regression='c')[1])
0.614194333311
In [64]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1986/12/1':'1993/10/30'],regression='nc')[1])
0.71409077222
In [65]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix[:'1993/10/30'],regression='ct')[1])
0.620086823254
In [66]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix[:'1993/10/30'],regression='c')[1])
0.652008364051
In [67]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix[:'1993/10/30'],regression='nc')[1])
0.999821116797
In [68]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1993/10/30':],regression='ct')[1])
0.735623826062
In [69]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1993/10/30':],regression='c')[1])
0.326699911361
In [70]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1990/1/1':'1992/8/31'],regression='nc')[1])
0.0739555163473
In [71]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1990/1/1':'1992/8/31'],regression='c')[1])
0.442806176355
In [72]:
print(sm.tsa.adfuller(lnn225.NIKKEI225.ix['1990/1/1':'1992/8/31'],regression='ct')[1])
0.323414633113
景気 期間始点   終点 ADF p-値 y 判定
全期間              0.99 0.0000 X
戦後復興期(recover) 1949/5/16 1954/11/30 0.90 0.0001 X
高度経済成長期(growth) 1954/12/1 1971/12/31 1.00 0.0001 X
安定期(stable) 1972/1/1 1986/11/30 1.00 0.0001 X
バブル経済期(bubble) 1986/12/0 1993/10/31 0.71 0.0000 X
バブルのピークまで   1986/12/1   1989/12/31   1.00 0.0001 X
バブルのピークから谷 1990/1/1 1992/8/31   0.07 -0.0000 X
経済変革期(reform) 1993/11/1 現在   0.63 -0.0001 OK

ADF 検定の後にモデルの回帰係数を求め、回帰係数がマイナスであるかどうかをチェックする。回帰係数がプラスであれば、κが1より大きくなり、モデルは発散してしまう

ドリフト付きモデル  次にドリフト付きモデルをテストしてみよう。

regression=cに設定を変更するだけでADF検定ができる

また、pd.olsではinterept=Trueにすればよい。

全期間を除いたすべての期間でドリフト付きモデルでランダムウォークの可能性がある。

ドリフトは時間トレンドを説明しているわけではなく、与えられた期間の中の確率項の平均値を表していて確率トレンドが存在することを示唆している。ドリフトの存在は確率的トレンドが存在しないことを意味しないことに注意してほしいです。

ドリフト+時間トレンド付きモデル  

すべての期間においてドリフト+時間トレンド付きモデルでランダムウォークである可能性があります。

In [73]:
z=lnn225
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("without drift  ",results.params[0])
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("with drift  ",results.params[0],results.params[1])
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("with drift + time trend  ",results.params[0],results.params[1],results.params[2])
without drift   2.64035301275e-05
with drift   0.00159089854801 -0.000156562803399
with drift + time trend   0.00155046216659 -0.000149176839762 -2.60661370915e-09
In [82]:
z=lnn225.ix[:"1954/11/30"]
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("nc  ",results.params)
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("c  ",results.params[0],results.params[1])
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("ct  ",results.params[0],results.params[1],results.params[2])
nc   NIKKEI225    0.000082
dtype: float64
c   0.00101585631908 -0.00010764780611
ct   0.015455834379 -0.00348332776496 5.00736004787e-06
In [75]:
z=lnn225.ix["1954/11/30":'1971/12/31']
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("nc  ",results.params)
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("c  ",results.params[0],results.params[1])
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("ct  ",results.params[0],results.params[1],results.params[2])
nc   NIKKEI225    0.000068
dtype: float64
c   0.00354522118342 -0.000434452045716
ct   0.00776380402623 -0.00113940523217 3.41862728775e-07
In [76]:
z=lnn225.ix['1971/12/31':'1986/11/30']
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("nc  ",results.params)
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("c  ",results.params[0],results.params[1])
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("ct  ",results.params[0],results.params[1],results.params[2])
nc   NIKKEI225    0.000057
dtype: float64
c   0.000594490997009 -1.03330061863e-05
ct   0.0126779199481 -0.00151484490346 6.08743200288e-07
In [77]:
z=lnn225.ix['1986/11/30':'1993/10/31']
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("nc  ",results.params)
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("c  ",results.params[0],results.params[1])
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("ct  ",results.params[0],results.params[1],results.params[2])
nc   NIKKEI225    0.000003
dtype: float64
c   0.022443821649 -0.00221648372153
ct   0.0437321302035 -0.00415890963348 -1.94402463757e-06
In [78]:
z=lnn225.ix['1986/11/30':'1989/12/31']
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("nc  ",results.params)
print("nc  ",results.aic)
print("nc  ",results.bic)
print("nc  ",results.rsquared)
print("nc  ",results.rsquared_adj)
print("nc  ",results.mse_resid)
print("nc  ",results.eigenvals)
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("c  ",results.params[0],results.params[1])
print("c  ",results.aic)
print("c  ",results.bic)
print("c  ",results.rsquared)
print("c  ",results.rsquared_adj)
print("c  ",results.mse_resid)
print("c  ",results.eigenvals)
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("ct  ",results.params[0],results.params[1],results.params[2])
print("ct  ",results.aic)
print("ct  ",results.bic)
print("ct  ",results.rsquared)
print("ct  ",results.rsquared_adj)
print("ct  ",results.mse_resid)
print("ct  ",results.eigenvals)
nc   NIKKEI225    0.000096
dtype: float64
nc   -4690.71922738
nc   -4686.08065959
nc   0.00754526995726
nc   0.0062445429192
nc   0.000126045202902
nc   [ 79769.67617582]
c   0.0273577963315 -0.00258119955419
c   -4690.19492165
c   -4680.91778607
c   0.00179377601062
c   0.000483794089372
c   0.000125967071361
c   [  8.05334306e+04   2.45603153e-01]
ct   0.227868128943 -0.0228692184611 1.77248395794e-05
ct   -4696.34982225
ct   -4682.43411889
ct   0.0123919203254
ct   0.00979636689658
ct   0.000124793426035
ct   [  1.48417739e+08   1.89889943e+04   2.27656103e-02]
In [79]:
z=lnn225.ix['1993/10/31':]
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("nc  ",results.params)
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("c  ",results.params[0],results.params[1])
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("ct  ",results.params[0],results.params[1],results.params[2])
nc   NIKKEI225   -0.000004
dtype: float64
c   0.01438084477 -0.00151219089951
ct   0.014378777831 -0.00151200061398 9.0106282207e-11
In [80]:
z=lnn225.ix['1990/1/1':'1992/8/31']
y=z.diff().dropna()
x=z.shift(1).dropna()
model=sm.OLS(y,x)
results=model.fit()
print("nc  ",results.params)
x=sm.add_constant(x)
model=sm.OLS(y,x)
results=model.fit()
print("c  ",results.params[0],results.params[1])
x["t"]=range(len(y))
model=sm.OLS(y,x)
results=model.fit()
print("ct  ",results.params[0],results.params[1],results.params[2])
nc   NIKKEI225   -0.000118
dtype: float64
c   0.0559706709658 -0.00566027578103
ct   0.207015979612 -0.0200296525262 -1.83524258751e-05
景気 ADF p値   y    切片 時間 判定
全期間        0.84 -0.00012 0.0016 -0.0 OK
戦後復興期(recover) 0.50 -0.0035 0.00155 0.0 OK
高度経済成長期(growth) 0.62 -0.0011 0.0078 0.0 OK
安定期(stable) 0.84 -0.0015 0.0128 0.0 OK
バブル経済期(bubble) 0.41 -0.0042 0.0437 -0.0 OK
バブルのピークまで   0.31 -0.0229 0.2278 0.0 OK
バブルのピークから谷 0.32 -0.0015 0.0142 -0.0 OK
経済変革期(reform) 0.76 -0.0200 0.2070 -0.0 OK

結果の解説

$$1. \bigtriangleup W{t} = \gamma W_{t-1} + w_{t} $$

$$2. \bigtriangleup W{t} = a+ \gamma W_{t-1} + w_{t} $$

$$3. \bigtriangleup W{t} = a+ \beta t + \gamma W_{t-1} + w_{t} $$ モデル2はモデル1の残差にまだ別の要因により説明できる部分があるのではないかとモデル1にドリフト項を加えたと考えることができます。

同様に、モデル3はモデル2の残差にまだ別の理由で説明できる要素があるのではないかと時間トレンドの項をモデル2に加えたととらえることできます。

さて、このなかでどれがこの場合良いモデルなのだろうか? 

どのようにモデルを選択したらよいのだろうか? 

これは難しい問題で今回の範囲を超えるので専門書を読んでいただきたいが、簡単に説明してみると、真のモデルの推定が目的の場合には、

BIC(Bayesian Information Criteria)を用い、また予測が目的の場合には、AIC(Akaike Information Critria)を用いる。

これはstatsmodelから得ることができます。

AIC: results.aic

BIC: results.bic です。

また、同様にMSE(Mean Squared Error)やR2(R-Squared)は

MSE: results.mse

R-Squared: results.rsquared で得ることができます。

ランダムウォークについて  

金融市場における価格の動きがランダムウォークである理由は大きく3つ考えられる。

1つ目は価格の動きが本質的に確率的であり、それは理想気体の中の分子の動きと同じように、どうすることもできな確率的な動きであるとする考え方です。  

2つ目は、人々が自らの暮らしを豊かにするために多くの活動を行うが、経済が順調に成長すれば、そこから新たなリスクが生じ、停滞すればそこから抜け出すためには革新的な行動が必要となる。 いつでもそのような経済の成長の一歩一歩は予測不可能だからとする考え方です。  

3つ目は、富の分配が公平になるように金融システムを構築したり、市場が効率的になるように裁定取引に多くの人々が参加したり、できるだけの多くの情報を集めて理論価格を計算したりといった人々の活動が効率的で確率的な価格の変化をもたらすとする考え方です。

確定的トレンドと確率的トレンド  

確定的トレンドは時間トレンドともよばれ、時間の経過にともない価格が継続して上昇または下落する現象です。

Material 5ではこのような現象には何らかの原因があるはずであるという説明をしました。

一方、Material 6では確率的トレンドを説明し、これは確率変数の和が作り上げるトレンドであるので、継続性のおいては予測が不可能であるという話をしました。

日本の場合にはバブル崩壊前には継続的な上昇トレンドがあるようにも思えるが、バブル崩壊後には価格は大きな波のような動きを続け、長期の上昇トレンドをみることはできない

そこで米国の株式市場をみてみよう。FREDから最も長い期間にわたり手に入る株式指数はWilshine5000 であるので、DLして確定トレンドの有無を確認します。

In [87]:
%matplotlib inline
import matplotlib.pyplot as plt
import statsmodels.api as sm
import pandas_datareader.data as pdr
import numpy as np
end='2016/9/30'
lnw5000 = np.log(pdr.DataReader("WILL5000INDFC", 'fred',"1949/5/16",end)).dropna()
lnw5000.columns=['Close']
plt.plot(lnw5000.Close,color='hotpink')
lnw5000["t"]=range(len(lnw5000))
model=sm.OLS(lnw5000.Close,lnw5000.t)
results=model.fit()
results.fittedvalues.plot(label='prediction',style='--')
plt.ylabel('log(n225 index)')
Out[87]:
<matplotlib.text.Text at 0x7f6a9a9e8c18>

また、残差を目視で確かめてみる。

In [88]:
results.resid.hist(bins=50,color='lightgreen')
plt.ylabel('frequency')
plt.xlabel('residual')
Out[88]:
<matplotlib.text.Text at 0x7f6ad0655a58>
In [89]:
print(results.summary())
                            OLS Regression Results                            
==============================================================================
Dep. Variable:                  Close   R-squared:                       0.967
Model:                            OLS   Adj. R-squared:                  0.967
Method:                 Least Squares   F-statistic:                 2.770e+05
Date:                Wed, 02 Aug 2017   Prob (F-statistic):               0.00
Time:                        20:12:46   Log-Likelihood:                -7929.6
No. Observations:                9397   AIC:                         1.586e+04
Df Residuals:                    9396   BIC:                         1.587e+04
Df Model:                           1                                         
Covariance Type:            nonrobust                                         
==============================================================================
                 coef    std err          t      P>|t|      [0.025      0.975]
------------------------------------------------------------------------------
t              0.0006   1.07e-06    526.352      0.000       0.001       0.001
==============================================================================
Omnibus:                     2185.621   Durbin-Watson:                   0.000
Prob(Omnibus):                  0.000   Jarque-Bera (JB):             1204.052
Skew:                          -0.744   Prob(JB):                    3.50e-262
Kurtosis:                       2.071   Cond. No.                         1.00
==============================================================================

Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
In [ ]: