pandas+matplotlibで世界時価総額トップ銘柄のアルファ値とベータ値をアニメ化する

#1. はじめに

　機械学習と金融工学の両方に興味があったので，一石二鳥を狙ってMachine learning for trading @ Udacityを受講しています．株価の評価指標である$\alpha$と$\beta$について学んだので，pandas+matplotlibで時変化をアニメ化して遊んでみました．対象は，2017年2月時点の世界時価総額ランキングトップ10銘柄です．Githubはこちらです．
　ご覧の通りド素人ですので，誤りや不備があればご指摘いただけると幸いです．

#2. 環境

##2-1. 株価データの入手

　日足の株価データを使いますので，事前に入手する必要があります．例えば，YAHOO! FINANCEからは，次のような手順で入手できます．なお，スクレイピングは規約違反ですのでご注意ください．

1. 目的の銘柄のサマリーページにアクセス（例：google（GOOG.csv））
2. ページ上部のHistrical Dataをクリック．
3. Time Period:が所望の期間であること，またFrequency:がDailyであることを確認し，Download Dataをクリック．
4. 必要に応じて，ファイル名を変更（例：table.csv→GOOG.csv）．

なお，本記事では，GOOG.csv等の株価データがdata/ディレクトリに格納されていることを前提に話を進めます．

##2-2. python

　データ成形にpandas，アニメ化にmatplotlibを使いますので，事前にインストールする必要があります．以下の環境で動作を確認しました．

• Ubuntu 14.04.5 LTS
• python 2.7.11
• NumPy 1.10.4
• matplotlib 1.5.1
• pandas 0.18.0

#3. 評価指標

3-1. daily-return

　daily-return $r$とは，簡単に言うと株価の前日比です．厳密には下式で定義されます．ここで，$x[t]$は日付$t$における調整後終値を表します．

r = \frac{x[t]}{x[t-1]} - 1 

pandasを使うと，簡単にdaily-returnを計算することができます．詳細はMachine learning for trading @ Udacityをご参照ください．ご参考までに，2006年12月1日から2016年12月1日のGOOGのdaily-returnを以下に示します．

3-2. アルファ値とベータ値

　市場平均とdaily-returnを比較することで，その銘柄の特徴を評価することができます．ここでは，SPYを市場平均と仮定し，SPYとGOOGのdaily-returnの散布図を書き，回帰直線を引いてみます．期間は，2006年12月1日から2016年12月1日です．

　この回帰直線$y=\beta x + \alpha$の切片をアルファ値，傾きをベータ値と呼びます．銘柄が市場平均と一致する時，回帰直線は直線$y=x$と一致する，つまり$\alpha=0$かつ$\beta=1$となるはずです．$\alpha$が大きいほど，市場平均に対する超過リターンが大きく，$\beta$が大きいほど，市場平均に対する連動性が大きいことを示します．ちなみに上図の場合は，$\alpha=0.000309329688614 $，$\beta=0.917720842929$です．
　$\alpha$および$\beta$は，個別銘柄のアクティブリターンとパッシブリターンの評価に利用される指標だそうです．

#4. 評価

　matplotlibを使って，アルファ値とベータ値の推移をアニメ化します．対象は，以下の2017年2月時点の世界時価総額ランキングトップ10銘柄です．業種はWikipediaを参考にしました．
　知識が足りないので，特に考察はありません．結果を見ていろいろ妄想するだけです．

#	symbol	企業名	業種
1	AAPL	Apple	電気機器
2	GOOG	Alphabet	コングロマリット
3	MSFT	Microsoft	情報・通信業
4	BRK-A	Berkshire Hathaway	保険業
5	AMZN	Amazon.com	小売業
6	FB	Facebook	インターネット
7	XON	Exxon Mobil	石油・石炭製品
8	JNJ	Johnson & Johnson	サービス業
9	JPM	JPMorgan Chase	その他金融業
10	WFC	Wells Fargo	金融業

##4-1. 前処理

　csvファイルから株価情報を読みだして，daily-returnを計算します．本記記事では，Udacityで紹介された以下の関数¹を使います．

• get_data(symbol)：symbolのCSVファイルから調整後終値を読みだして，pandas.DataFrame形式で出力する関数です．
• fill_missing_values(stock_data)：stock_data（pandas.DataFrame）の欠損を穴埋めする関数です．
• compute_daily_returns(stock_data)：stock_data（pandas.DataFrame）から，daily-returnを計算し，pandas.DataFrameで出力する関数です．

##4-2. アニメ化

　$\alpha$と$\beta$をアニメ化する関数を作りました．デフォルト設定では，2006年12月1日（start_date）から2016年12月1日（end_date）の上記10銘柄（symbols）について，直近2年間（period）のデータから値を計算します．10日（interval）ずつ計算範囲をずらしながら，animate_polyfit()で1フレームずつ描画し，ab.gifを出力します．凡例の大きさはSPYと当該銘柄の相関係数の大きさを表しています．グラフ右上の日付は，計算範囲の最終日を表します．つまり，右上の日付から過去period日分のデータをもとに計算した値を表示しています．
　なお，各銘柄のcsvファイルをdataディレクトリに置く必要があります．また，4.1節の関数を定義したudacity.pyを同じディレクトリに置く必要がありますので，ご注意ください．

import os
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.animation as ani
import matplotlib.cm as cm
import udacity  # Functions defined in "Machine learning for trading"


def animate_a_b(
        symbols=["AAPL", "GOOG", "MSFT", "BRK-A", "AMZN",
                 "FB", "XOM", "JNJ", "JPM", "WFC"],
        start_date="2006-12-01", end_date="2016-12-01",
        period=252 * 2):

    """ --- Preprocess: You have to take Udacity! --- """
    # Read data
    dates = pd.date_range(start_date, end_date)  # date range as index
    stock_data = udacity.get_data(symbols, dates)  # get data for each symbol

    # Fill missing values
    udacity.fill_missing_values(stock_data)
    # Daily returns
    daily_returns = udacity.compute_daily_returns(stock_data)

    """ --- Make animation --- """
    interval = 10  # Nframe interval
    frames = (len(stock_data) - period) / interval  # Num of frames
    markers = ["o", "^", "s"]

    def animate_polyfit(nframe):

        plt.clf()
        daily_returns_p = daily_returns[
            nframe * interval: nframe * interval + period]
        corr = daily_returns_p.corr(method="pearson")

        xmin, xmax = -0.003, 0.003
        ymin, ymax = 0.0, 2.0

        plt.plot([0, 0], [ymin, ymax], '--', color='black')
        plt.plot([xmin, xmax], [1, 1], '--', color='black')
        for n, symbol in enumerate(symbols[1:]):
            beta, alpha = np.polyfit(daily_returns_p["SPY"],
                                     daily_returns_p[symbol], 1)
            plt.plot(alpha, beta, markers[n % len(markers)], alpha=0.7,
                     label=symbol, color=cm.jet(n * 1. / len(symbols)),
                     ms=np.absolute(corr.ix[0, n + 1]) * 25)
        plt.xlim([xmin, xmax])
        plt.ylim([ymin, ymax])
        plt.xlabel("Alpha")
        plt.ylabel("Beta")
        plt.text(xmax, ymax, str(daily_returns_p.index[-1]),
                 ha="right", va="bottom")
        plt.legend(loc="upper left")

    fig = plt.figure(figsize=(8, 8))
    anim = ani.FuncAnimation(fig, animate_polyfit, frames=frames)
    anim.save("ab.gif", writer="imagemagick", fps=18)

　2008年ごろのリーマン・ショック周辺で，金融業界の２銘柄（JPM，WFC）がエキサイトしていますね…．ちなみに，FBは比較的若い企業なので，アニメの途中から登場します．

#5. おわりに
　$\alpha$および$\beta$の推移をpandas+matplotlibでアニメ化してみました．使い道はさておき，二次元以上の指標の時変化を直感的に表現する方法として，アニメ化は有りかもしれません²．今後も受講を続け，記事にできるようなことがあれば，ご紹介したいと思います．
　最後まで読んで下さり，ありがとうございました！

#Appendix. Trendalyzer
　コメント欄でご指摘の通り，本記事の元ネタはHans RoslingさんのTrendalyzerです．Trendalyzerは，Googleスプレッドシートのモーショングラフとして利用できますので，使ってみました．こちらでスプレッドシートを公開しておりますので，自由に使ってください．シート1にデータ本体，シート2にモーショングラフを置いてあります．

A-1. 準備

　googleスプレッドシートで開けるよう，計算結果をエクセル形式で保存しました．モーショングラフ用のフォーマット³に合わせて，pandas.DataFrameを作成します．

save_a_b()

def save_a_b(
        symbols=["AAPL", "GOOG", "MSFT", "BRK-A", "AMZN",
                 "FB", "XOM", "JNJ", "JPM", "WFC"],
        start_date="2006-12-01", end_date="2016-12-01",
        period=252 * 2):

    """ --- Preprocess: You have to take Udacity! --- """
    # Read data
    dates = pd.date_range(start_date, end_date)  # date range as index
    stock_data = udacity.get_data(symbols, dates)  # get data for each symbol

    # Fill missing values
    udacity.fill_missing_values(stock_data)
    # Daily returns
    daily_returns = udacity.compute_daily_returns(stock_data)

    """ --- Calculate spreadsheet of alpha and beta ---"""
    sheet = pd.DataFrame(columns=["Symbol", "Date", "Alpha", "Beta", "Color",
                                  "Size"])
    interval = 10  # Nframe interval
    frames = (len(stock_data) - period) / interval  # Num of frames

    for nframe in range(frames):

        daily_returns_p = daily_returns[
            nframe * interval: nframe * interval + period]
        corr = daily_returns_p.corr(method="pearson")

        for n, symbol in enumerate(symbols[1:]):
            beta, alpha = np.polyfit(daily_returns_p["SPY"],
                                     daily_returns_p[symbol], 1)

            new_row = pd.DataFrame(
                [[symbol, daily_returns_p.index[-1].strftime("%Y/%m/%d"),
                  alpha, beta, n, np.absolute(corr.ix[0, n + 1]) * 25]],
                columns=sheet.columns)
            sheet = sheet.append(new_row, ignore_index=True)

    sheet.to_excel("ab.xlsx")

A-2. 遊んでみた

　スタートボタンを押すと，アニメーションを開始します．アニメーションの速度を変更することも可能です．

　タイムバーを使って，早送りや巻き戻しも思いのままです．特定の銘柄を選択すると，軌跡を表示することができます．すごい！

#参考

　本記事の作成にあたっては，以下を参考にさせて頂きました．ありがとうございました！

• Machine learning for trading @ Udacity：株価データの入手，加工について
• 投資信託のアルファとベータってなに？：アルファ値とベータ値について
• Pythonを使った回帰分析の概念の解説 番外編１：アニメーションの作成について
• Hans Rosling:
  New insights on poverty：本記事の元ネタ．

1. 権利上の問題を回避するため，ソースコードは掲載しません．無料講座ですので，ぜひ受講をば！ ↩

2. 実は昔TEDで，これの100倍格好いいアニメを使ったプレゼンを見たのを思い出したので，真似してみました．詳しくは覚えてませんが，GDPと何かを軸にとって，主要国の動きを解説してました． ↩

3. %m/%d/%Y形式で出力したのですが，google先生がモーショングラフに変換してくださらず，ハマりました．結局，%Y/%m/%d形式で出力し，googleスプレッドシートの表示形式>数字で日付を指定すると，解決しました．よくわかってません… ↩