複数の株価推移の傾向がわかる!?動的因子分析とは?

このエントリーをはてなブックマークに追加

こんにちは、データデザイン部で職場体験中の新入社員、衞藤と申します。

大学院時代には神経科学の研究室で統計や機械学習についても学んでいました(神経科学とデータサイエンスとの関りは非常に深いです)。

学生時代から気分も一新。導入研修を終えて部署での職場体験に臨んだ新入社員衞藤には、トレーナーHさんからタスクが与えられます。

Hさん
「今案件で来客数予測をしてるんだけど、もうちょっと精度を上げなきゃいけないんだよね。」

衞藤
「そうなんですね。新しいモデルを検討したりするんですか?」

Hさん
「それもあるんだけど、今のモデルの予測誤差について店舗間で共通した傾向が生じているのか調べてほしいんだよね。」

衞藤
「なるほど?」

Hさん
「そのために動的因子分析で予測誤差を説明する因子を見つけて欲しいんだわ。誤差を説明する因子を考察すれば、予測精度を上げるために何をすればいいかが分かるかもしれないしね。」

衞藤
「(動的因子分析、初めて聞くなあ…。)」

初めてのタスクとして与えられたのは動的因子分析を用いた来客数予測誤差の解析でした!新入社員衞藤は果たしてどうなってしまうのか(どうにかした)。

こちらの記事では、動的因子分析で予測誤差を解析するタスクで得た知見をもとに、因子分析とは何なのかから実データ(株価)を用いた動的因子分析の実施までをpythonのコードと共にまとめています。本記事では、それぞれの分析手法のアルゴリズムや数式にはあまり着目せず、それぞれの手法がどのような特徴を持っているのかをコードと共に説明しています。

【実行環境】
 OS:Ubuntu 18.04.1 LTS
 言語:python 3.7.3 on Jupyter notebook

 
 
 
【図:本記事に出てくる手法と疑問点】

目次

動的因子分析を理解する前に

そもそも因子分析とは

衞藤
「(因子分析、聞いたことあるけどどんな手法だったっけ…)」

動的因子分析を理解する前に、そもそも因子分析とはなんなのでしょうか。因子分析(Factor Analysis)とは複数のデータ(観測データ)を説明する隠れた共通変数(因子)を見つけるための統計的手法です。例として挙げるならば、教室の生徒40人の5教科(数理英国社)のテストの点数を因子分析することで、記憶力や論理的思考といったような生徒たちの点数を左右する因子の抽出が期待できます。

因子分析の解析上の特性

さて、因子分析の概要をおさらいした所で次にエンジニアとして気になるのは手法の特徴です。どのような場面で因子分析を用いればいいのかを判断するためには、手法の特徴を正しく理解している必要があります。他の手法と比べて因子分析を選択する理由は何が挙げられるのでしょう。

因子分析と類似した手法として挙げられるのは主成分分析(principal components analysis, PCA)や独立成分分析(independent components analysis, ICA)です。読者の方の中には、これらの手法を使ったことがあるといった方がいらっしゃると思います。私も、因子分析を調べ始めた際に、

???「それってPCAじゃ駄目なんですか?」

と自問していました。PCAやICAは特に次元圧縮の文脈で使われることの多い手法です。多次元のデータを2次元や3次元のデータに縮約することで、データの分布を可視化することが可能になります(2018年にはUMAPと呼ばれる次元圧縮手法も提案され話題になりました)(https://arxiv.org/abs/1802.03426)。

因子分析の特徴と、他の次元圧縮の手法との違いを確認するために2次元空間上に相関を持つ2つの基底ベクトルからなる点群データをプロットしました。解析者としては、データを説明する因子(or 成分)を各手法が抽出することを期待します。

 


# 各種ライブラリのインポート
import matplotlib.pyplot as plt
from sklearn.decomposition import PCA
from sklearn.decomposition import FastICA
from sklearn.decomposition import FactorAnalysis
import numpy as np
import seaborn as sns

# Jupyter notebook上でのグラフ描画
%matplotlib inline

# 2次元平面上に直交しない2つのデータの分布を作成
# 自由度1.5のt分布からランダムにを20000×2点サンプリング
rng = np.random.RandomState(42)
S = rng.standard_t(1.5, size=(20000, 2))
S[:, 0] *= 2.

# 2つの基底ベクトルをもとにデータを作成
A = np.array([[1, 1], [0, 2]])
X = np.dot(S, A.T)
Xp = X / X.std()

# 点群を散布図にプロット
plt.figure(figsize=(5, 5))
sns.scatterplot(Xp[:, 0], Xp[:, 1])
plt.xlim([-3, 3])
plt.ylim([-3, 3])
plt.xlabel("x")
plt.ylabel("y")


【相関を持つ2つの基底ベクトルによって生成された点群データ。横軸にxの値、縦軸にyの値を示している。青丸で示しているのがそれぞれのデータ点】

こちらのデータ対して、PCA、ICA、因子分析を適用した結果が以下の図になります。今回は、ライブラリとしてscikit-learnを用いました。


# PCAの実行
pca = PCA().fit(X)
PCA(copy=True, whiten=False)
components_pca = pca.components_ / pca.components_.std()
print(components_pca)

# ICAの実行
ica = FastICA(random_state=rng).fit(X)
components_ica = ica.mixing_ / ica.mixing_.std()
print(components_ica)

# 因子分析の実行
fa = FactorAnalysis(random_state=rng).fit(X)
components_fa = fa.components_ / fa.components_.std()
print(components_fa)

# 散布図上にそれぞれの成分分解の手法で抽出された成分(基底ベクトル)を描画
plt.figure(figsize=(15, 5))
plt.subplot(1, 3, 1)
sns.scatterplot(Xp[:, 0], Xp[:, 1])
plt.xlabel("x")
plt.ylabel("y")
plt.xlim([-3, 3])
plt.ylim([-3, 3])
plt.quiver(0, 0, components_pca[0, 0], components_pca[0, 1],
color="red", width=0.01, scale=6)
plt.quiver(0, 0, components_pca[1, 0], components_pca[1, 1],
color="blue", width=0.01, scale=6)
plt.title("PCA")

plt.subplot(1, 3, 2)
sns.scatterplot(Xp[:, 0], Xp[:, 1])
plt.xlabel("x")
plt.ylabel("y")
plt.xlim([-3, 3])
plt.ylim([-3, 3])
plt.quiver(0, 0, components_ica[0, 0], components_ica[1, 0],
color="red", width=0.01, scale=6)
plt.quiver(0, 0, components_ica[0, 1], components_ica[1, 1],
color="blue", width=0.01, scale=6)
plt.title("ICA")

plt.subplot(1, 3, 3)
sns.scatterplot(Xp[:, 0], Xp[:, 1])
plt.xlabel("x")
plt.ylabel("y")
plt.xlim([-3, 3])
plt.ylim([-3, 3])
plt.quiver(0, 0, components_fa[0, 0], components_fa[0, 1],
color="red", width=0.01, scale=6)
plt.quiver(0, 0, components_fa[1, 0], components_fa[1, 1],
color="blue", width=0.01, scale=6)
plt.title("FA")


【PCA,ICA、因子分析が抽出した成分(基底ベクトル)。横軸にxの値、縦軸にyの値を示しており、赤矢印で示すのが1つ目の因子(成分)、青矢印で示すのが2つ目の因子(成分)。左から、主成分分析、独立成分分析、因子分析の結果となっている。】

結果を確認すると、PCAの結果抽出された新しい次元の軸と、因子分析で抽出された因子の軸がほぼ同じであることがわかります。これは、どちらの手法もデータの分散が最も大きくなる基底を因子(成分)として選択するアルゴリズムであるためです。また、1つ目の次元と2つ目の次元が直交する性質があることもこちらの図から理解することができるでしょう。また、ICAはこのようなデータの分布に対して、最も適切にデータを説明する基底ベクトルが選択できている印象を受けます。これは、ICAが直交性の仮定を持っておらず、最もデータの分布が正規分布から離れるように新しい基底を求めるためです。

この結果を踏まえると、読者の方からは因子分析と主成分分析は何も変わらないじゃないかという声が挙がってくることが容易に想像ができます。このような議論は長らく続いているようですが、明確な回答の1つとして、因果の関係がこれらの手法は逆転している点が挙げられます。主成分分析は、主成分の構成要素として観測変数を定義していますが、因子分析では観測変数の構成要素として因子を定義しています。そのため、主成分分析は主としてデータの要約という目的、因子分析はデータの構成要素の抽出という目的で用いられます。

衞藤「(観測変数の構成要素を抽出する際に因子分析を用いればいいんだな。)」

動的因子分析

動的因子分析とは

衞藤「(結局動的因子分析ってなんなんだ。)」

さて、因子分析の概要をつかんだ上で動的因子分析とは結局なんなのでしょうか。結論から言うと、動的因子分析は因子分析の時系列データへの応用版です。これまではテストの点数などある時点でのデータに着目しましたが、経済データなどは多くの場合時系列の変化を伴います。複数の時系列データの変動を説明する因子を見つける手法が動的因子分析というわけです。

例えば、複数の店舗の来客数の変動に対して動的因子分析を実施することで、主婦客の来客数変動とサラリーマン客の来客数変動と特徴づけられるような因子を抽出できるかもしれません。このことから、それぞれの店舗の来客数は主婦とサラリーマンの客層の因子の重ね合わせで表現されていると考えることができます。


【動的因子分析の活用例。複数店舗の来客数変動から2つの因子を抽出することで、主婦層とサラリーマン層によって来客数が変動していると解釈したケース。】

実データへの適用(株価推移)

それでは実際に株価推移のデータを使って、動的因子分析を行ってみましょう。今回利用する株価データは[株価データ]( https://kabuoji3.com/stock/)さんからダウンロードしました。選択した企業(A, B, C, D)の2018年の株価データを読み込んでプロットします。


# 各種ライブラリのインポート
import pandas as pd
import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
import scipy

# Jupyter notebook上でのグラフ描画
%matplotlib inline

# 企業毎の各日付の株価データを取得
df_stock_price = pd.DataFrame()
companies = ["A", "B", "C", "D"]

for idx, company in enumerate(companies):
df = pd.read_csv(company + "_2018.csv", encoding="SHIFT-JIS", header=1)
df.rename(columns={"終値調整値": company}, inplace=True)
if idx == 0:
df_stock_price = df[["日付", company]]
else:
df_stock_price = pd.merge(df_stock_price, df[["日付", company]], on="日付")
df_stock_price.set_index("日付", inplace=True)

# 成型後のデータを確認
df_stock_price.head()

# 日付の抽出
dates = df_stock_price.index

# 4社の株価をプロット
num_companies = len(companies)
TICK_DAYS = 30
len_dates = len(dates)
range_dates = range(0, len_dates, TICK_DAYS)
fig, ax = plt.subplots(num_companies, 1, figsize=(13, 10))
fig.text(0.04, 0.5, "Stock price [JPY]", va="center", rotation="vertical", fontsize=18)
idx = 0
for company, price in df_stock_price.iteritems():
ax[idx].plot(dates, price, linewidth=3.0)
if idx == 3:
ax[idx].set_xticks(range_dates)
ax[idx].set_xticklabels(dates[range_dates], rotation=90)
else:
ax[idx].set_xticks(range_dates)
ax[idx].set_xticklabels([], rotation=90)
ax[idx].set_xlim(0, len_dates)
ax[idx].set_title("Company " + company, fontsize=18)
ax[idx].spines["right"].set_color("none")
ax[idx].spines["top"].set_color("none")
idx += 1
plt.xlabel("Dates", fontsize=18)


【4企業の2018年における株価変動。縦軸にある日の株価、横軸に日付を示している。上から、企業A、B、C、Dそれぞれの株価変動を示している。】

それぞれの企業の株価データを確認すると上昇傾向や下降傾向があることがわかります。このようなトレンドを抽出したい場合もありますが、日毎の変化を比較したい場合もあります。その場合は時系列データの対数差分を取ることで日毎の変化率を捉えることが可能になります。こちらの変化率をz変換することで、動的因子分析の入力に適した値に変換します。その結果もプロットしてみましょう。


# 株価を企業ごとに対数差分データにしたのちz値に変換
df_logdiff = df_stock_price.copy()
df_logdiff = df_logdiff.applymap(np.log)
df_logdiff = df_logdiff.diff()
df_logdiff = df_logdiff.iloc[1:, :]
df_logdiff = df_logdiff.apply(scipy.stats.zscore, axis=0)

# 対数差分のz値に変換した株価データを確認
df_logdiff.head(10)

# 4社の株価の対数差分のz値をプロット
dates_logdiff = df_logdiff.index
len_logdiff = len(dates_logdiff)
range_logdiff = range(0, len_logdiff, TICK_DAYS)
idx = 0
fig, ax = plt.subplots(num_companies, 1, figsize=(13, 10))
fig.text(0.04, 0.5, "Stock price fluctuations", va="center", rotation="vertical", fontsize=18)
for company, logdiff in df_logdiff.iteritems():
ax[idx].plot(dates_logdiff, logdiff, linewidth=3.0)
if idx == 3:
ax[idx].set_xticks(range_logdiff)
ax[idx].set_xticklabels(dates_logdiff[range_logdiff], rotation=90)
else:
ax[idx].set_xticks(range_logdiff)
ax[idx].set_xticklabels([], rotation=90)
ax[idx].set_xlim(0, len_logdiff)
ax[idx].set_title("Company " + company, fontsize=18)
ax[idx].spines["right"].set_color("none")
ax[idx].spines["top"].set_color("none")
idx += 1
plt.xlabel("Dates", fontsize=18)


【4企業の2018年における株価の対数差分のz値。縦軸にある日の株価変動率のz値、横軸に日付を示している。上から、企業A、B、C、Dそれぞれの株価の対数差分のz値を示している。】

こちらを実施したことで、先ほど確認されたトレンド成分が除去されたことがわかります。また、各企業の株価変化のグラフを確認すると、同じタイミングで株価が変動していたりすることも分かります。このような観測変数から動的因子を抽出することはできるのでしょうか?

今回は、pythonでの動的因子分析の実施にあたりライブラリとしてstatsmodelsを用いました。Statsmodelで動的因子モデルを実施する際にはいくつかパラメータ(因子数:2、AR次元:1、誤差のAR次元:1)を設定しました。因子数は抽出する因子の数(観測変数の数以下)、AR次元は観測変数の値が何個前のデータから影響を受けるか、誤差のAR次元は誤差が何個前の誤差から影響を受けるかという値になります。これらのパラメータは、基本的に解析対象のデータの特性や、入力次元数に依存して適宜変更する必要があります。赤池情報量基準(AIC)などを用いて適切なモデルを選択するプロセスを挟むことも多いですが、今回は省略します。

それでは4企業の株価変化のデータを用いて動的因子分析を行っていきます。


# 動的因子分析のモデル宣言
# 因子数:2 因子が影響を受けるデータ範囲: 1 誤差が影響を受けるデータ範囲: 1
mod = sm.tsa.DynamicFactor(df_logdiff,
k_factors=2,
factor_order=1,
price_order=1)
# モデルの学習、methodは最適化手法を設定
initial_res = mod.fit(method="powell", disp=False)
res = mod.fit(initial_res.params)

# 推定結果
print(res.summary(separate_params=False))

# 因子の時間変化による変動
num_factor = len(res.factors.filtered)
fig, ax = plt.subplots(num_factor, 1, figsize=(10, 6))
fig.text(0.04, 0.5, "Fluctuations", va="center", rotation="vertical", fontsize=18)
for idx, factor in enumerate(res.factors.filtered):
ax[idx].plot(dates_logdiff, factor, linewidth=3.0)
ax[idx].set_xticks(range_logdiff)
ax[idx].set_xlim(0, len_logdiff)
ax[idx].set_title("Factor {}".format(idx + 1), fontsize=18)
ax[idx].spines["right"].set_color("none")
ax[idx].spines["top"].set_color("none")
if idx == 0:
ax[idx].set_xticklabels([])
else:
ax[idx].set_xticklabels(dates_logdiff[range_logdiff], rotation=90)
ax[idx].set_xlabel("Dates", fontsize=18)

 

動的因子分析によって抽出された因子が以下です。


【動的因子分析で抽出された2つの因子。縦軸にある日の変動、横軸に日付を示している。上から、因子1、因子2それぞれの因子の変動を示している。】

無事に因子を抽出することができました!これら2つの因子を確認すると、10月の終わりにのみ因子1ではピークが立ち上がっているが因子2では期間を通して大きく変動するなど違いがあることが分かります。2018年の10月末には携帯通信料金の値下げの流れが生じ、携帯事業の収益が悪化する可能性が見込まれたため携帯キャリア関連株が暴落しました。このことから、今回抽出された因子1は「携帯事業に大きく影響を受ける因子」、因子2は「携帯事業に影響を受けず期間を通して大きく値が変動する因子」などの解釈が可能です(株価に関する知識をお持ちの方であればより詳細な解釈が可能かもしれません)。

また、これらの因子がどの程度それぞれの観測変数を説明できるかは、各因子の各観測変数に対する寄与率を計算することで判断ができます。


# 各店舗に対する動的因子の寄与率
res.plot_coefficients_of_determination(figsize=(8, 10));


【因子ごとの観測変数に対する寄与率。縦軸に寄与率、横軸にそれぞれの観測変数を示している。上から因子1の結果、因子2の結果を示している。】

結果を確認すると、因子1は主に企業B、Dの分散を説明する因子であり、因子2が主に企業A、Cの分散をある程度説明する因子であったことが分かります。また、因子1、2だけで説明できる分散は十分ではなく、これらの因子の他に様々な要因が株価の変動に影響していることがわかるでしょう。

4企業の株価変化という少数の観測変数から動的因子分析を実施したところ、結果として抽出された因子はある程度解釈はできたものの、(筆者が株について詳しくないこともあり)細かな解釈が難しい内容でした。因子分析においても同様ですが、抽出された因子の解釈については、分析者の恣意性がある程度求められます。抽出された因子に〇〇因子と解釈を行う場合は、何故その解釈を行ったのかの根拠を明確にする必要性があると考えられます。

まとめ

  • 因子分析は複数の観測変数を構成する成分を抽出する手法
  • 動的因子分析は因子分析の時系列応用版

本記事では、新入社員衛藤が動的因子分析を用いた分析タスクを任されてから実際に解析を行うまでに得た知見を共有させていただきました。複数のデータ(観測変数)をいくつかの構成要素に分解するニーズは、データ解析業務においては往々にしてあると思います。その際の解析手法として、今回ご紹介した動的因子分析やその他の手法がご参考になれば幸いです。

 

参考文献

経済・ファイナンスデータの軽量時系列分析、沖本竜義、2010、朝倉書店

https://kabuoji3.com/stock/
https://qiita.com/supersaiakujin/items/138c0d8e6511735f1f45
https://scikit-learn.org/stable/auto_examples/decomposition/plot_ica_vs_pca.html#sphx-glr-auto-examples-decomposition-plot-ica-vs-pca-py
http://www.chadfulton.com/topics/dfm_coincident.html
https://limo.media/articles/-/8275

サラリーマンのためのデータサイエンス基礎講座

データサイエンスの基礎からビジネスに使えるフレームワーク、データを直接操作してAI開発を体験できるハンズオンまで網羅した、半日集中講座を毎月開催しております。

データサイエンス基礎講座はこちら

資料の無料ダウンロード

貴社の事業課題に向けたAI活用の各種資料をご案内いたします。

資料請求はこちら

WRITER
Shota Eto

Shota Eto

衛藤 祥太Shota Eto

2019年度富士通クラウドテクノロジーズ株式会社新入社員。修士時代はヒトの脳活動データを機械学習を用いて解析。 取りたい資格は狩猟免許(一種・わな) 。

最新記事

ページTOPへ