-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy path[2]自我回歸模型.py
120 lines (87 loc) · 2.19 KB
/
[2]自我回歸模型.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
# -*- coding: utf-8 -*-
"""[2]自我回歸模型.ipynb
Automatically generated by Colaboratory.
Original file is located at
https://colab.research.google.com/drive/1zni84oomRNu0tTKKy9lDKTVPuOWhMifN
"""
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
# Commented out IPython magic to ensure Python compatibility.
# %matplotlib inline
import ffn
import statsmodels.api as sm
"""# 自我回歸模型
- 模擬不同beta下的`AR(1)`序列
- 利用匯率資料估計`AR(1)`模型參數
- 購買力困惑 (PPP puzzle)
## 模擬不同beta下的AR(1)序列
- $ \beta_1 = 0.2$
- $\beta_1 = 0.95$
- $\beta_0 = 0$
- 模型設定
$$
y_t = \beta_0 + \beta_1 y_{t-1} + \epsilon_t
$$
- 設定初始值 $y_0 = 0$
- 生成隨機數 $\epsilon_1\sim N(0,1)$
- 得到序列值 $y_1$
- 依此類推可得所有序列值
"""
len(y)
y = np.zeros(1001)
range(1,1000)
for i in range(1,10):
print(i)
for i in range(1,10):
print(i)
# set beta
b_0 = 0
b_11 = 0.2
b_12 = 0.95
sigma = 1
x = np.zeros(1000)
y = np.zeros(1000)
for i in range(1,1000):
x[i] = b_0 + b_11*x[i-1] + np.random.randn()
y[i] = b_0 + b_12*y[i-1] + np.random.randn()
plt.plot(x,'--')
plt.plot(y,'-')
plt.title('beta=0.2 & 0.95')
plt.xlabel('Time')
plt.ylabel('Series')
data = ffn.get('EURUSD=X',start='2010-01-01',end='2020-12-31')
data.info()
# AR(1): eurusd
y = data.eurusdx
x = sm.add_constant(data.eurusdx.shift(1))
res = sm.OLS(y,x,missing='drop').fit()
print(res.summary())
res.resid.plot()
#res.fittedvalues.plot()
# generate log return
data['ret'] = np.log(data.eurusdx).diff()
data.head()
data.ret
# AR(1): ret_eurusd
x = sm.add_constant(data.ret.shift(1))
y = data.ret
res = sm.OLS(y,x,missing='drop').fit()
print(res.summary())
# test PPP puzzle
test = pd.read_csv('re_us.csv',index_col='DATE')
test.info()
# set date index
test.index = pd.to_datetime(test.index)
x = sm.add_constant(test.RBUSBIS.shift())
y = test.RBUSBIS
res = sm.OLS(y,x,missing='drop').fit()
print(res.summary())
# calculate the half-life (month frequency)
np.log(0.5)/np.log(0.9888)
# import package
from statsmodels.tsa.ar_model import ar_select_order
# data: y
# ic='aic'/'bic'
mod = ar_select_order(y, ic='aic', maxlag=10)
mod.ar_lags