【第1章】pythonで「経済・ファイナンスデータの計量時系列分析」の章末問題を解く

自分の勉強用に、[沖本竜義「経済・ファイナンスデータの計量経済分析」](https://www.asakura.co.jp/G_12.php?isbn=ISBN978-4-254-12792-8)の章末問題をpythonで解いてみました。
jupyter notebook上で記述したものをほぼそのまま載せてあります。
ソースコードの細かい説明は省いてありますが、今後余裕があれば説明を加えようと思います。


こちらの記事は第1章の章末問題を解いたものになります。
[第2章](https://qiita.com/mckeeeen/items/a0126a20116dd27ecba9)、[第4章](https://qiita.com/mckeeeen/items/4afa1167008c1f315d0b)、[第5章](https://qiita.com/mckeeeen/items/d293f8005fee5452ecbb)、[第6章](https://qiita.com/mckeeeen/items/55cb0130e5a7fa035e1e)、[第7章](https://qiita.com/mckeeeen/items/ce000b00b0c725503072)も公開しています。

各種設定・モジュールのインポート

```python
%matplotlib inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from statsmodels.tsa.stattools import acf
from statsmodels.stats.diagnostic import acorr_ljungbox


import matplotlib as mpl
font = {"family":"IPAexGothic"}
mpl.rc('font', **font)
plt.rcParams["font.size"] = 12
```


# 1.3

(1.8)式のモデルからデータを生成する

```python
gwn_0_1 = pd.Series(np.random.randn(100),name=u'\u03bc=0,\u03c3=1')
gwn_2_1 = pd.Series(2+np.random.randn(100),name=u'\u03bc=2,\u03c3=1')
gwn__2_1 = pd.Series(-2+np.random.randn(100),name=u'\u03bc=-2,\u03c3=1')
gwn_0_2 = pd.Series(2*np.random.randn(100),name=u'\u03bc=0,\u03c3=2')
gwn_0_3 = pd.Series(3*np.random.randn(100),name=u'\u03bc=0,\u03c3=3')
gwn_2_2 = pd.Series(2+2*np.random.randn(100),name=u'\u03bc=2,\u03c3=2')
```

プロットする

```python
fig, ax = plt.subplots(3,2,figsize=(20,20))
ax[0,0].plot(gwn_0_1)
ax[0,0].set_title(u'\u03bc=0,\u03c3=1')
ax[0,1].plot(gwn_2_1)
ax[0,1].set_title(u'\u03bc=2,\u03c3=1')
ax[1,0].plot(gwn__2_1)
ax[1,0].set_title(u'\u03bc=-2,\u03c3=1')
ax[1,1].plot(gwn_0_2)
ax[1,1].set_title(u'\u03bc=0,\u03c3=2')
ax[2,0].plot(gwn_0_3)
ax[2,0].set_title(u'\u03bc=0,\u03c3=3')
ax[2,1].plot(gwn_2_2)
ax[2,1].set_title(u'\u03bc=2,\u03c3=2')

plt.show()
```


![output_3_0.png](https://qiita-image-store.s3.amazonaws.com/0/235450/6f220a2a-29a6-f0cc-16ea-204971351ca9.png)


# 1.5
##(1)

ファイルをインポートする

```python
eco = pd.read_csv('./input/economicdata.csv')
print(eco.shape)
eco.head()
```

    (364, 7)





<div>

<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>date</th>
      <th>topix</th>
      <th>exrate</th>
      <th>indprod</th>
      <th>cpi</th>
      <th>saunemp</th>
      <th>intrate</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>Jan-75</td>
      <td>276.09</td>
      <td>29.13</td>
      <td>47.33</td>
      <td>52.625</td>
      <td>1.7</td>
      <td>12.67</td>
    </tr>
    <tr>
      <th>1</th>
      <td>Feb-75</td>
      <td>299.81</td>
      <td>29.70</td>
      <td>46.86</td>
      <td>52.723</td>
      <td>1.8</td>
      <td>13.00</td>
    </tr>
    <tr>
      <th>2</th>
      <td>Mar-75</td>
      <td>313.50</td>
      <td>29.98</td>
      <td>46.24</td>
      <td>53.114</td>
      <td>1.8</td>
      <td>12.92</td>
    </tr>
    <tr>
      <th>3</th>
      <td>Apr-75</td>
      <td>320.57</td>
      <td>29.80</td>
      <td>47.33</td>
      <td>54.092</td>
      <td>1.8</td>
      <td>12.02</td>
    </tr>
    <tr>
      <th>4</th>
      <td>May-75</td>
      <td>329.65</td>
      <td>29.79</td>
      <td>47.33</td>
      <td>54.385</td>
      <td>1.8</td>
      <td>11.06</td>
    </tr>
  </tbody>
</table>
</div>



date列をdatetime型に変換してecoのDatetimeIndexとする

```python
eco.index=pd.to_datetime(eco.date.values,format='%b-%y')
eco.drop('date',axis=1,inplace=True)
eco.head()
```




<div>

<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>topix</th>
      <th>exrate</th>
      <th>indprod</th>
      <th>cpi</th>
      <th>saunemp</th>
      <th>intrate</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>1975-01-01</th>
      <td>276.09</td>
      <td>29.13</td>
      <td>47.33</td>
      <td>52.625</td>
      <td>1.7</td>
      <td>12.67</td>
    </tr>
    <tr>
      <th>1975-02-01</th>
      <td>299.81</td>
      <td>29.70</td>
      <td>46.86</td>
      <td>52.723</td>
      <td>1.8</td>
      <td>13.00</td>
    </tr>
    <tr>
      <th>1975-03-01</th>
      <td>313.50</td>
      <td>29.98</td>
      <td>46.24</td>
      <td>53.114</td>
      <td>1.8</td>
      <td>12.92</td>
    </tr>
    <tr>
      <th>1975-04-01</th>
      <td>320.57</td>
      <td>29.80</td>
      <td>47.33</td>
      <td>54.092</td>
      <td>1.8</td>
      <td>12.02</td>
    </tr>
    <tr>
      <th>1975-05-01</th>
      <td>329.65</td>
      <td>29.79</td>
      <td>47.33</td>
      <td>54.385</td>
      <td>1.8</td>
      <td>11.06</td>
    </tr>
  </tbody>
</table>
</div>



プロットする

```python
title=['TOPIX','実効為替レート','鉱工業生産指数','CPI','失業率','コールレート']
ylim=[3500,120,120,120,6,14]

fig,ax = plt.subplots(nrows=3,ncols=2,figsize=[10,15])
for i in range(3):
    for j in range(2):
        ax[i,j].plot(eco.iloc[:,2*i+j])
        ax[i,j].set_title(title[2*i+j])
        ax[i,j].set_ylim(0,ylim[2*i+j])
        ax[i,j].set_xlim(eco.index[0],eco.index[-1])
        ax[i,j].hlines(ax[i,j].get_yticks(),ax[i,j].get_xlim()[0],ax[i,j].get_xlim()[1],linewidth=0.5)
        ax[i,j].spines['right'].set_visible(False)
        ax[i,j].set_xticks(eco.index[eco.index.is_year_start][::3])
        ax[i,j].set_xticklabels(eco.index[eco.index.is_year_start][::3].strftime('%y'))
        ax[i,j].set_xlabel('年')
plt.subplots_adjust(wspace=0.2,hspace=0.3)
plt.show()
```


![output_7_0.png](https://qiita-image-store.s3.amazonaws.com/0/235450/8079cb4d-dfd5-9ab1-1fc1-2f1ce088cbba.png)


## (2)

対数差分系列を計算する

```python
eco['topix_dlog'] = np.log(eco.topix).diff()
eco['exrate_dlog'] = np.log(eco.exrate).diff()
eco['indprod_dlog'] = np.log(eco.indprod).diff()
eco.head()
```




<div>

<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>topix</th>
      <th>exrate</th>
      <th>indprod</th>
      <th>cpi</th>
      <th>saunemp</th>
      <th>intrate</th>
      <th>topix_dlog</th>
      <th>exrate_dlog</th>
      <th>indprod_dlog</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>1975-01-01</th>
      <td>276.09</td>
      <td>29.13</td>
      <td>47.33</td>
      <td>52.625</td>
      <td>1.7</td>
      <td>12.67</td>
      <td>NaN</td>
      <td>NaN</td>
      <td>NaN</td>
    </tr>
    <tr>
      <th>1975-02-01</th>
      <td>299.81</td>
      <td>29.70</td>
      <td>46.86</td>
      <td>52.723</td>
      <td>1.8</td>
      <td>13.00</td>
      <td>0.082422</td>
      <td>0.019378</td>
      <td>-0.009980</td>
    </tr>
    <tr>
      <th>1975-03-01</th>
      <td>313.50</td>
      <td>29.98</td>
      <td>46.24</td>
      <td>53.114</td>
      <td>1.8</td>
      <td>12.92</td>
      <td>0.044650</td>
      <td>0.009383</td>
      <td>-0.013319</td>
    </tr>
    <tr>
      <th>1975-04-01</th>
      <td>320.57</td>
      <td>29.80</td>
      <td>47.33</td>
      <td>54.092</td>
      <td>1.8</td>
      <td>12.02</td>
      <td>0.022301</td>
      <td>-0.006022</td>
      <td>0.023299</td>
    </tr>
    <tr>
      <th>1975-05-01</th>
      <td>329.65</td>
      <td>29.79</td>
      <td>47.33</td>
      <td>54.385</td>
      <td>1.8</td>
      <td>11.06</td>
      <td>0.027931</td>
      <td>-0.000336</td>
      <td>0.000000</td>
    </tr>
  </tbody>
</table>
</div>



## (3)

プロットする

```python
title=['TOPIX','実効為替レート','鉱工業生産指数']
ylim_u=[15,12,5]
ylim_l=[-20,-8,-5]
yticks_int=[5,2,1]

fig,ax = plt.subplots(nrows=3,ncols=1,figsize=[10,15])
for i in range(3):
    ax[i].plot(eco.iloc[:,6+i]*100)
    ax[i].set_title(title[i])
    ax[i].set_ylim(ylim_l[i],ylim_u[i])
    ax[i].set_xlim(eco.index[0],eco.index[-1])
    ax[i].spines['right'].set_visible(False)
    ax[i].set_xticks(eco.index[eco.index.is_year_start][::3])
    ax[i].set_xticklabels(eco.index[eco.index.is_year_start][::3].strftime('%y'))
    ax[i].set_yticks(range(ylim_l[i],ylim_u[i]+1)[::yticks_int[i]])
    ax[i].hlines(ax[i].get_yticks(),ax[i].get_xlim()[0],ax[i].get_xlim()[1],linewidth=0.5)
    ax[i].set_xlabel('年')
plt.subplots_adjust(wspace=0.2,hspace=0.3)
plt.show()
```


![output_11_0.png](https://qiita-image-store.s3.amazonaws.com/0/235450/f067b897-a45d-70a2-49ae-4889e0a59ef8.png)


## (4)

indprodの変化率の標本自己相関とQ(m)とP値を求める

```python
autocorr_indprod=acf(eco.indprod_dlog[1:])[1:]
q_m_indprod=acorr_ljungbox(eco.indprod_dlog[1:])[0]
p_value_indprod=acorr_ljungbox(eco.indprod_dlog[1:])[1]
lag=range(1,41)
test_indprod=pd.DataFrame({'lag':lag,'autocorr':autocorr_indprod,'Q(m)':q_m_indprod,'P_value':p_value_indprod})
test_indprod[['lag','autocorr','Q(m)','P_value']].iloc[0:10,:].T
```




<div>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
      <th>3</th>
      <th>4</th>
      <th>5</th>
      <th>6</th>
      <th>7</th>
      <th>8</th>
      <th>9</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>lag</th>
      <td>1.000000e+00</td>
      <td>2.000000e+00</td>
      <td>3.000000e+00</td>
      <td>4.000000e+00</td>
      <td>5.000000e+00</td>
      <td>6.000000e+00</td>
      <td>7.000000e+00</td>
      <td>8.000000e+00</td>
      <td>9.000000e+00</td>
      <td>1.000000e+01</td>
    </tr>
    <tr>
      <th>autocorr</th>
      <td>-3.056815e-01</td>
      <td>1.231154e-01</td>
      <td>1.699805e-01</td>
      <td>6.173203e-02</td>
      <td>1.677640e-02</td>
      <td>8.291975e-02</td>
      <td>-1.082355e-03</td>
      <td>8.959420e-02</td>
      <td>3.584451e-02</td>
      <td>-4.913541e-02</td>
    </tr>
    <tr>
      <th>Q(m)</th>
      <td>3.420025e+01</td>
      <td>3.976335e+01</td>
      <td>5.039731e+01</td>
      <td>5.180377e+01</td>
      <td>5.190793e+01</td>
      <td>5.445974e+01</td>
      <td>5.446017e+01</td>
      <td>5.745610e+01</td>
      <td>5.793698e+01</td>
      <td>5.884316e+01</td>
    </tr>
    <tr>
      <th>P_value</th>
      <td>4.972263e-09</td>
      <td>2.320056e-09</td>
      <td>6.574768e-11</td>
      <td>1.516132e-10</td>
      <td>5.634492e-10</td>
      <td>5.958688e-10</td>
      <td>1.907042e-09</td>
      <td>1.466964e-09</td>
      <td>3.342017e-09</td>
      <td>5.995413e-09</td>
    </tr>
  </tbody>
</table>
</div>



標本自己相関をプロットする

```python
fig,ax=plt.subplots(1,1,figsize=[8,5])
ax.bar(test_indprod.lag[:20],test_indprod.autocorr[:20])
ax.set_ylim(-0.5,0.5)
ax.set_xticks([1,5,10,15,20])
ax.set_yticks(np.asarray(range(-5,6))*0.1)
ax.hlines([1.96/np.sqrt(eco.shape[0]-1),-1.96/np.sqrt(eco.shape[0]-1)],0,21,linestyles='dashed',linewidth=1)
ax.set_xlabel('ラグ')
ax.set_ylabel('自己相関')
ax.set_title('鉱工業生産指数の成長率')
plt.show()
```


![output_14_0.png](https://qiita-image-store.s3.amazonaws.com/0/235450/79b24c8d-564c-d902-9889-f25c4bab9171.png)


## (5)

TOPIXについて自己相関検定を行う

```python
autocorr_topix=acf(eco.topix_dlog[1:])[1:]
q_m_topix=acorr_ljungbox(eco.topix_dlog[1:])[0]
p_value_topix=acorr_ljungbox(eco.topix_dlog[1:])[1]
lag=range(1,41)
test_topix=pd.DataFrame({'lag':lag,'autocorr':autocorr_topix,'Q(m)':q_m_topix,'P_value':p_value_topix})
test_topix[['lag','autocorr','Q(m)','P_value']].iloc[0:10,:].T
```




<div>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
      <th>3</th>
      <th>4</th>
      <th>5</th>
      <th>6</th>
      <th>7</th>
      <th>8</th>
      <th>9</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>lag</th>
      <td>1.000000e+00</td>
      <td>2.000000e+00</td>
      <td>3.000000e+00</td>
      <td>4.000000e+00</td>
      <td>5.000000e+00</td>
      <td>6.000000e+00</td>
      <td>7.000000</td>
      <td>8.000000</td>
      <td>9.000000</td>
      <td>10.000000</td>
    </tr>
    <tr>
      <th>autocorr</th>
      <td>3.077999e-01</td>
      <td>6.824984e-02</td>
      <td>3.529433e-02</td>
      <td>6.801588e-02</td>
      <td>3.130508e-02</td>
      <td>-5.079917e-02</td>
      <td>0.015837</td>
      <td>0.025024</td>
      <td>0.106257</td>
      <td>0.064149</td>
    </tr>
    <tr>
      <th>Q(m)</th>
      <td>3.467592e+01</td>
      <td>3.638552e+01</td>
      <td>3.684399e+01</td>
      <td>3.855135e+01</td>
      <td>3.891405e+01</td>
      <td>3.987178e+01</td>
      <td>39.965130</td>
      <td>40.198841</td>
      <td>44.424647</td>
      <td>45.969211</td>
    </tr>
    <tr>
      <th>P_value</th>
      <td>3.894192e-09</td>
      <td>1.255984e-08</td>
      <td>4.964713e-08</td>
      <td>8.622952e-08</td>
      <td>2.471338e-07</td>
      <td>4.827198e-07</td>
      <td>0.000001</td>
      <td>0.000003</td>
      <td>0.000001</td>
      <td>0.000001</td>
    </tr>
  </tbody>
</table>
</div>



標本自己相関をプロットする

```python
fig,ax=plt.subplots(1,1,figsize=[8,5])
ax.bar(test_topix.lag[:20],test_topix.autocorr[:20])
ax.set_ylim(-0.5,0.5)
ax.set_xticks([1,5,10,15,20])
ax.set_yticks(np.asarray(range(-5,6))*0.1)
ax.hlines([1.96/np.sqrt(eco.shape[0]-1),-1.96/np.sqrt(eco.shape[0]-1)],0,21,linestyles='dashed',linewidth=1)
ax.set_xlabel('ラグ')
ax.set_ylabel('自己相関')
ax.set_title('TOPIXの成長率')
plt.show()
```


![output_17_0.png](https://qiita-image-store.s3.amazonaws.com/0/235450/454e327e-7ddb-872d-927c-18a3ef610257.png)


実効為替レートについて自己相関検定を行う

```python

autocorr_exrate=acf(eco.exrate_dlog[1:])[1:]
q_m_exrate=acorr_ljungbox(eco.exrate_dlog[1:])[0]
p_value_exrate=acorr_ljungbox(eco.exrate_dlog[1:])[1]
lag=range(1,41)
test_exrate=pd.DataFrame({'lag':lag,'autocorr':autocorr_exrate,'Q(m)':q_m_exrate,'P_value':p_value_exrate})
test_exrate[['lag','autocorr','Q(m)','P_value']].iloc[0:10,:].T
```




<div>

<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
      <th>3</th>
      <th>4</th>
      <th>5</th>
      <th>6</th>
      <th>7</th>
      <th>8</th>
      <th>9</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>lag</th>
      <td>1.000000e+00</td>
      <td>2.000000e+00</td>
      <td>3.000000e+00</td>
      <td>4.000000e+00</td>
      <td>5.000000e+00</td>
      <td>6.000000e+00</td>
      <td>7.000000e+00</td>
      <td>8.000000e+00</td>
      <td>9.000000e+00</td>
      <td>1.000000e+01</td>
    </tr>
    <tr>
      <th>autocorr</th>
      <td>3.217763e-01</td>
      <td>4.860159e-02</td>
      <td>9.068437e-02</td>
      <td>3.667455e-02</td>
      <td>-3.474773e-02</td>
      <td>-8.257040e-02</td>
      <td>-3.047480e-02</td>
      <td>1.070266e-01</td>
      <td>1.015502e-01</td>
      <td>-1.826474e-02</td>
    </tr>
    <tr>
      <th>Q(m)</th>
      <td>3.789649e+01</td>
      <td>3.876344e+01</td>
      <td>4.179009e+01</td>
      <td>4.228649e+01</td>
      <td>4.273335e+01</td>
      <td>4.526370e+01</td>
      <td>4.560934e+01</td>
      <td>4.988452e+01</td>
      <td>5.374427e+01</td>
      <td>5.386948e+01</td>
    </tr>
    <tr>
      <th>P_value</th>
      <td>7.459915e-10</td>
      <td>3.824946e-09</td>
      <td>4.445457e-09</td>
      <td>1.454942e-08</td>
      <td>4.184633e-08</td>
      <td>4.148156e-08</td>
      <td>1.041334e-07</td>
      <td>4.300973e-08</td>
      <td>2.111351e-08</td>
      <td>5.131560e-08</td>
    </tr>
  </tbody>
</table>
</div>


標本自己相関をプロットする

```python
fig,ax=plt.subplots(1,1,figsize=[8,5])
ax.bar(test_exrate.lag[:20],test_exrate.autocorr[:20])
ax.set_ylim(-0.5,0.5)
ax.set_xticks([1,5,10,15,20])
ax.set_yticks(np.asarray(range(-5,6))*0.1)
ax.hlines([1.96/np.sqrt(eco.shape[0]-1),-1.96/np.sqrt(eco.shape[0]-1)],0,21,linestyles='dashed',linewidth=1)
ax.set_xlabel('ラグ')
ax.set_ylabel('自己相関')
ax.set_title('実効為替レートの成長率')
plt.show()
```


![output_19_0.png](https://qiita-image-store.s3.amazonaws.com/0/235450/2b72e781-e994-a913-c7b5-f248db818099.png)