見出し画像

データを準備する

数理科学エージェント

この章の目的

手元に便利なデータとツールをおいて何か思いつたときすぐに試せるようしておきたい、ということで入手が容易なデータとツールを準備します．  


1. データを入手する

身近な題材として株価関連のデータを金融マーケット情報を提供するWebサイトからダウンロードする方法，そして公開APIを使ったプログラムでダウンロードする方法をそれぞれ簡単に紹介しましょう．

1.1. Yahoo finance USから株価データをダウンロードする

Yahoo financeのサイトhttps://finance.yahoo.com/にアクセスします．ページ上部の検索バーに証券コード（例えば「極洋」なら1301）を入力すると以下のように候補がプルダウンリストに表示されるので選択します．

yahoo finance US

Yahoo financeでは日本株を指定するために証券コードに'.T'を付与しています．極洋なら'1301.T'です．そのほかマーケット指標，為替には下表のように固有のシンボルが定義されています．以下サイトから検索することもできます．yahoo finance us/ ticker lookup


$$
\begin{array}{cc}
\hline
マーケット指標，為替 & シンボル名 \\
\hline
日経225 & \hat{} N225 \\
\hline
SP 500 & \hat{} GSPC \\
\hline
USD/JPY（米ドル円レート） & JPY=X \\
\hline
USD/EUR（米ドルユーロレート）& EURUSD=X \\
\hline
\end{array}

$$

遷移した先の個別株価情報のページで，'Historical Data'のメニュータブを選びます．そして表示したい期間などを設定し'Download'をクリックすればデータがcsvファイルとしてダウンロードされます．

株価データダウンロード画面

1.2. pandas_datareaderを使ってデータを取得する

pandas_datareaderはPythonプログラムに直接データを取り込むためのAPIを提供するポピュラーなツールです．python環境と，2次元データを操作するための定番ツールであるpandasが必要となります．事前にpython環境を準備したうえで，以下のコマンドをコマンドプロンプトやターミナルで実行してツールをインストールします．

pip install pandas pandas_datareader

以下のpythonコードは株価データのダウンロードとCSVファイルへの保存を行います．

import pandas_datareader.data as web
# ダウンロードするデータの情報を設定
start_date = '2023-01-01'  # データの開始日を指定
end_date = '2023-05-30'  # データの終了日を指定
stock_symbol = '1301.JP'  # 株価データをダウンロードする企業のシンボルを指定
# データプロバイダとしてを'stook'を指定して株価データをダウンロード
stock_data = web.DataReader(stock_symbol, data_source='stooq', start=start_date, end=end_date)

# CSVファイルとして保存する．シリアル番号のインデクスを付与しない．
stock_data.to_csv('stock_data.csv',index=False)

stooqのほかにも多数のデータプロバイダが定義されており多様なデータを取得できることがdata_readerの特徴です．ツールのサイトhttps://pandas-datareader.readthedocs.io/en/latest/を参照すると定義済のプロバイダを確認できる．ただし設定するとエラーが返るプロバイダ名もあるようです（例：'yahoo'）．

1.3. yahoo_ finance_api2を使ってデータを取得する

もうひとつツールを紹介しておきましょう．yahoo_finance_api2はpandas_datareaderとくらべてより柔軟なダウンロード設定ができるところが良い．コマンドプロンプトやターミナルで次のコマンドを実行してツールをインストールします．

pip install yahoo_finance_api2

以下のpythonコードは日経平均株価データのダウンロードとCSVファイルへの保存を行います．

  • ダウンロードされるデータは次のような辞書型で返されます．扱いやすくするためにpandasデータフレーム形式へ変換します．

{
  'timestamp': [...], 'open': [...],'high': [...],'low': [...],'close': [...],'volume': [...]
}

  • yahoo_finance_api2ではダウンロードするデータの期間は現在から数えた期間を指定します．過去の特定期間を指定したダウンロードはできません．

  • ダウンロードされたtimestamp列は文字型なのでto_datetime()を使ってdatetime64型へ変換します．

  • ダウンロードされるデータの日付は協定世界時（UTC）のため日本ローカルタイムに変換します.

  • カラム名を'Date', 'Open', 'High', 'Low', 'Close', 'Volume'に変え，さに'Date'列をindexに設定します．

# pre_yfapi2.py
# 
import pandas as pd
from yahoo_finance_api2 import share

my_share = share.Share("^N225") # 個別企業の場合には証券コードに.Tを付加する
period_type = share.PERIOD_TYPE_YEAR # 取得する期間を年単位とする
period = 1 # 取得する期間を１年とする
frequency_type = share.FREQUENCY_TYPE_DAY # 日次データを取得する
frequency = 1 # 日次データを取得する間隔を毎日とする．
symbol_data = my_share.get_historical(period_type,period,frequency_type,frequency)
df = pd.DataFrame(symbol_data.values(), index=symbol_data.keys()).T
# timestamp列をdatetime64型へ変換
df.timestamp = pd.to_datetime(df.timestamp, unit='ms')
# timestamp列をローカルタイム（日本標準時）へ変換
df.timestamp =df.timestamp.dt.tz_localize('UTC').dt.tz_convert('Asia/Tokyo')
# カラム名を変更
df.columns = ['Date', 'Open', 'High', 'Low', 'Close', 'Volume']
# Date列のフォーマットを日付のみに変更
df['Date'] = df['Date'].dt.strftime('%Y-%m-%d')
# Date列をindexに変更
df.set_index('Date', inplace=True)
# csv出力
df.to_csv("^N225.csv")

作成した株価csvファイルはpythonのチャート作図用ライブラリmatplotlibを使うと簡単に表示することができます．
残念なことにmatplotlibは既定値では日本語表示に対応していません．日本語を表示する設定はネット情報でもいくつか探すことができるのですが，かなり環境依存性が強いので注意が必要です．ここで使用したjapanize_matplotlib は，設定ファイルを修正したりプログラムでの設定も特に必要なく日本語表示ができます．

# pre_nikkei.py
# 
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.dates as mdates
import japanize_matplotlib

# CSVを読み込みデータフレームを返す
data = pd.read_csv('^N225.csv')

# 'Date'列をdatetime64型に変換する
data['Date'] = pd.to_datetime(data['Date'])

# 'Date'列をindexに設定する
data.set_index('Date', inplace=True)

# 'Close'列を折れ線チャートで描画する
plt.plot(data.index, data['Close'])

# チャートのタイトルと軸ラベルを設定する
plt.title('Nikkei 225 Index')
plt.xlabel('Date')
plt.ylabel('Closing Price')

# x軸ラベルの日付フォーマットを'年-月'に設定する
date_format = mdates.DateFormatter('%y-%m')
plt.gca().xaxis.set_major_formatter(date_format)

#見やすくするためにx軸目盛りラベルを左に45度回転する
plt.xticks(rotation=45)

# グリッド線を表示する
plt.grid(True)

# 指定した価格の高さで水平線を描画する
specified_price1 = 28814  # Close price at 2021-06-01
specified_price2 = 33485  # Close price at 2023-06-15
plt.axhline(specified_price1, color='blue')
plt.axhline(specified_price2, color='red')
# 水平線の近くにアノテーションをつける
memo1 = "28,814yen"
plt.annotate(memo1, xy=(data.index[0], specified_price1), xytext=(data.index[0], specified_price1 + 2000),
             arrowprops=dict(arrowstyle="->"), color='red')
memo1 = "33,485yen"
plt.annotate(memo1, xy=(data.index[-1], specified_price2), xytext=(data.index[-1], specified_price2 - 2000),
             arrowprops=dict(arrowstyle="->"), color='red')

# チャートを表示する
plt.show()


日経平均株価の推移


少しの余談　トレード戦略の評価について

私たちが個人として株式の売買を行うときは，様々な情報を収集し複数の案を比較したうえで取引を実行に移していくでしょう．これは売買を能動的（アクティブ）に行うという意味で「アクティブ運用」という方法をとっているといえます．
これに対し金融サービス事業者は，金融マーケットのよく知られた指標（例えば日経平均株価）の上がり下がりと同期したパフォーマンスが得られる商品を販売しています．これはパッシブ運用型商品とかインデックス運用型商品と呼ばれています．

日付ｔでの日経平均株価を $${price(t)}$$ とし，日経平均株価に追従するように設計されたパッシブ運用型商品Ａをa円で購入したとします．すると１年後$${(t+1year)}$$にＡを売却したときの利益  $${ profit(t+1year)}$$は売却手数料を$${fee}$$として，以下で計算できます．

$$
profit(t+1year) =  a * (price(t+1year) / p(t) - 1) - fee
$$

日経平均が１年後に首尾よく上昇してくれれば何の苦労もなく利益を得ることができるわけです．別の見方をすれば私たちがアクティブ運用にトライする場合，上記のprofitを大きく上回る利益を獲得することを目標とすることは自然と思います．
本稿で説明する予測技術や分類技術を使った際のパフォーマンスの評価はパッシブ運用で得られるパフォーマンスとの比較で良し悪しを判断する必要があるでしょう．
このチャートを見ると2年間で１６％ほど値上がりしているのでこの期間に日経平均株価と連動するパッシブ運用型商品を買った人にはうれしい状況ではないでしょうか．タイミングよく売却していたならばですが．

2. 「statsmodels」の短いまとめ

本稿で統計分析用ソフトの「statsmodels」，数値解析用ソフトの「SciPy」，機械学習用のソフトの「scikit-learn」をしばしば活用します．この節では「statsmodels」の使い方とサンプルデータについて簡単にまとめておきます．SciPy,scikit-learnについては特にまとめていませんがサンプルコードのなかで適宜使い方についてコメントを入れるようにします．

2.1. オープンソースソフトウェアを使うことについて

statsmodelsは修正 BSD License（3 条項 BSD ライセンス）に従うオープンソースソフトウェアです（statsmodelsのライセンス条項）．
オープンソースの利用は私的・商用問わずメリットだけでなくデメリットも理解したうえで利用することが必要でしょう．以下に参考としてstatsmodelsが従う修正 BSD Licenseの特徴と注意点を整理しておきます．

特徴

再配布の自由： ソフトウェアを再配布する自由を保障する．

  • 商用利用可： 商用利用は許可される．つまりこのライセンスの下で開発したソフトウェアを商業的なプロダクトに組み込むことができる．

  • 修正可： ソフトウェアを修正し，派生作品を作成することが許可される．派生作品も同じライセンスの下で再配布できる．

  • 著作権表示： ソフトウェアのコピーには，元の著作権表示とライセンス条件が含まれている必要がある．

  • 無保証： ソフトウェアは無保証である．ソフトウェアの使用に関連するリスクは利用者が自己負担する必要がある．

注意点

再配布に制約がないが故にソフトウェアが不正確であったり，他のソフトウェアと競合する可能性がある場合でも，再配布ができる．このため，派生作品がソフトウェアの品質を損なう可能性がある．

  • 特許関連の制約なし： 特許関連の保護や制約が含まれていないため，特許侵害に対する保護が提供されない．これは特に特許関連の法的問題に対処する必要がある場合に問題となる可能性がある．

  • 変更の追跡不要： 派生作品を元のプロジェクトに戻す必要がないため，変更の追跡が不要である．これにより，派生プロジェクトがソフトウェアを変更しても，元のプロジェクトに変更を報告する必要がない．

2.2. インストールとカリフォルニア州住宅データの入手

サポート対象のPythonのバージョンはWebサイトを参照してください
Anaconda環境でのインストール

conda install -c conda-forge statsmodels

pipでのインストール

python -m pip install statsmodels


次に本稿でしばしば使われる有用なデータについて説明します．

カリフォルニア住宅データ（California Housing Dataset）は，機械学習やデータ分析の教育や実践に広く使用されているデータセットです．これはカリフォルニア州の住宅価格に関連する様々な統計情報から構成されたデータであり、以下のような使い方に好適です．

  1. 回帰分析の実験: 住宅価格の予測に使用できる典型的な回帰分析のデータセットである．これを用いて，住宅価格を他の特徴（部屋数，人口，所得など）に基づいて予測するモデルを構築し，回帰分析の手法を試すことができる．

  2. 特徴量エンジニアリング: データセットには住宅価格以外にもさまざまな特徴が含まれている．これらの特徴を活用して特徴量エンジニアリングを試すことができる．例えば，新しい特徴を作成し，カテゴリカルデータをエンコードする方法をトライするなど．

  3. データ可視化: データセットには地理情報が含まれており，地図上で住宅価格の分布を可視化することができる．データ可視化の手法を試すことができる．

  4. 機械学習アルゴリズムの比較: カリフォルニア住宅データを使用して，異なる機械学習アルゴリズム（線形回帰，決定木，ランダムフォレストなど）の性能を比較できる．ハイパーパラメータの調整やモデルの評価にも役立つ．

  5. 実世界データの扱い: 次世界のデータセットであるため，欠損値の処理，外れ値の検出，データクリーニングなど，実際のデータに対処する手法を試すことができる．

データの取得

StatLibのwebサイトにアクセスして「Houses.zip」ダウンロードします．

これを展開して出力されるcadata.txtをテキストエディタを使って開き，冒頭数行の説明書きをすべて削除し数値データのみとし保存します．

csvファイルの準備

以下のコードを実行してヘッダ付きのcsv形式のファイルをつくります．

# ref_housetocsv.py
# convert original text data to csv
# ダウンロードしたzipファイルを展開してできるcadata.txtで
# 数値データより上のコメントを全て削除した後で実行する

import pandas as pd

# テキストファイルを読み取り、データをDataFrameに変換
data = []
# delimiterがひとつ以上の空白のため正規表現（'\s+'）を使う
df = pd.read_csv('cadata.txt',delimiter=r'\s+',header=None,
    names=['median_house_value', 'median_income', 'housing_median_age',
           'total_rooms', 'total_bedrooms', 'population', 'households',
           'latitude', 'longitude'])
# データ型を数値に変換
df = df.apply(pd.to_numeric, errors='coerce')
print(df)
# CSVファイルに保存
df.to_csv('cadata.csv', index=False)


回帰予測を実行する

入手したCalifornia Housingデータセットを用いて回帰予測によるパラメータ推定を行ってみましょう．なお回帰や帰無仮説など統計学の技法については別の章で扱うためここでは説明は省略します．

以下のコードではサンプルデータの「median_house_value」列を目的変数，そのほかの列を説明変数として重回帰予測を行いレポートを出力するとともにデータの散布図を表示します．

# resf_sm_ols.py
# 最小二乗法を用いた重回帰推定
#
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm
import japanize_matplotlib
import pandas as pd

# サンプルデータ生成
df = pd.read_csv("cadata.csv")
df = df.dropna()
column_header = ["median_income", "housing_median_age", "total_rooms",
                 "total_bedrooms", "population", "households", "latitude",
                 "longitude"]
# データ数が大きいので最初の500件のみ使用する
df = df.iloc[0:500, ]
df_X = df.loc[:, column_header].astype(float)
X = df_X.values
df_y = df.loc[:, "median_house_value"].astype(float)
y = df_y.values

X = sm.add_constant(X)  # 切片を含むために説明変数に定数項を追加

model = sm.OLS(y, X).fit()  # 最小二乗法による線形回帰モデル

# パラメータ推定結果の表示
print(model.summary())
print("Parameters: ", model.params)

# プロット
fig = plt.figure(figsize=(10, 9))
column_header_jp = ["世帯収入", "住宅築年数", "全部屋数",
                    "寝室数", "ブロック内人口", "ブロック内世帯数", "ブロック中心緯度",
                    "ブロック中心経度"]
title_jp = "住宅価格 vs "
clrs = ['b', 'g', 'r', 'c', 'm', 'y', 'b', 'g', 'r', 'c', 'm', 'y']
# 4行2列のサブプロットに散布図scatter()で描画
# marker: "o","+","x",etc
# s:マーカーのサイズ
for i in range(8):
    pos = 421 + i
    ax = fig.add_subplot(pos)
    x = X[:, i + 1]

    ax.scatter(x, y, color=clrs[i], marker='o', s=10)
    ax.set_title(title_jp + column_header_jp[i])

plt.tight_layout()
plt.show()

実行すると以下のように各目的変数と各説明変数との相関を示す散布図が表示されます．

2.3. 処理結果レポートの見方

上記のサンプルコードを実行したときの処理結果レポート

results.summary()で出力される情報は，statsmodelsを使用して推定した統計モデルの詳細な統計情報を提供しています．この情報には，回帰モデルの適合度，係数の推定値，仮説検定の結果，モデルの性能評価などが含まれます．以下に，表示される結果の各セクションを説明します．

  1. Dep. Variable（従属変数）：
    モデルで設定した従属変数（目的変数）の名前．

  2. Model（モデル）：
    使用された統計モデルの名前や仕様．

  3. Method（方法）：
    使用された推定方法（最尤法，最小二乗法など）．

  4. Date（日付）：
    モデルを実行した日付．

  5. Time（時間）：
    モデルの実行にかかった時間．

  6. No. Observations（観測数）：
    モデルに使用された観測データの総数．

  7. Df Residuals（残差の自由度）：
    残差の自由度（残差の数から推定されるパラメータの数を差し引いたもの）．

  8. Df Model（モデルの自由度）：
    モデルの自由度（推定されたパラメータの数）．

  9. Covariance Type（共分散の種類）：
    共分散行列のタイプ（Homoscedastic（等分散），Heteroscedastic（非等分散）．

  10. R-squared（決定係数）：
    決定係数（R-squared）は，モデルがデータをどれだけよく説明しているかを示す指標で，0から1の値を取る．1に近いほどモデルの適合度が高いことを示す．

  11. Adj. R-squared（調整決定係数）：
    調整決定係数（Adj. R-squared）は，説明変数の数を考慮した決定係数で，過剰適合（オーバーフィッティング）を考慮した指標．

  12. F-statistic（F統計量）：
    F統計量（F-statistic）は，回帰分析において，回帰モデル全体の統計的有意性を評価するための統計量である．統計的有意性は帰無仮説の検証によって判定される．
    帰無仮説  H0（null hypothesis） ： 全ての説明変数（独立変数）を含まない，つまり単純な平均値モデルであり，帰無仮説モデルの下では，説明変数が目的変数に対して影響を持たないと仮定する．
    対立仮説 H1（Alternative hypothesis）： 対立仮説モデルの下では，少なくとも一つ以上の説明変数が目的変数に対して影響を持つと仮定する．
    F統計量が大きい場合，帰無仮説モデルと対立仮説モデルとの間に統計的に有意な違いがあることを示し，説明変数が目的変数に対して有意な影響を持っている可能性が高いことを示唆する．一般的に，F統計量の値が一定の有意水準を超える場合，帰無仮説を棄却し，回帰モデル全体が有意であると結論できる．

  13. Prob (F-statistic)（F統計量の確率値）：
    F統計量の確率値は，F分布の下でのモデルの有意性を示します．小さい値ほどモデルが有用であることを示す．

  14. AIC（赤池情報量基準）および BIC（ベイズ情報量基準）：
    AICとBICは，モデルの適合度とモデルの複雑さをトレードオフで考慮する情報量基準．モデルの選択に役立つ．

  15. coef（係数）：
    各説明変数の係数（回帰係数）の推定値．
    以下15.から19.までの予測式のパラメータ推定値と評価値が１行ごとに表示される．

  16. std err（標準誤差）：
    係数の推定値の標準誤差．

  17. t（t値）：
    各係数のt値が表示される．通常，t値が絶対値が2以上で，p値（t分布の下での帰無仮説を棄却する確率）が0.05より小さい場合，係数は統計的に有意であるとみなされる．この数値の絶対値が大きいほど目的変数の予測への貢献が大きいと考えられる．

  18. P>|t|（p値）：
    この値が有意水準以下であれば検定で設定された帰無仮説は棄却される．

  19. [0.025 0.975]（信頼区間）：
    各係数の信頼区間が表示されます．信頼区間は，推定値の信頼性を示します．

  20. Omnibus（オムニバス検定）および Prob(Omnibus)（オムニバス検定のｐ値）：
    オムニバス検定は，モデル全体の適合度を評価する統計検定である．

  21. Durbin-Watson（ダービン・ワトソン統計量）：
    ダービン・ワトソン統計量は，残差の自己相関を評価する指標で，時系列データの場合に特に重要

  22. Jarque-Bera (JB)：
    ジャーク・ベラ統計量（正規性の検定）．

  23. Prob(JB)：
    JB統計量のｐ値．データが正規分布に従うかどうかを評価する．

  24. Skew：
    歪度（スキュー）はデータの分布が左右対称からどれだけずれているかを示す．

  25. Kurtosis：
    尖度（カートーシス）はデータの尖度を評価し，正規分布と比較する．

下部には残差の統計情報や回帰診断統計量が表示されることがあります．これらの情報は，モデルの適合度や残差の性質を評価するために役立つでしょう．

前書き　数理的手法を使っていろいろ予測したい  

次の章　仮説検定を活用する


この記事が気に入ったらサポートをしてみませんか？