日本株の最新情報からPythonで株価分析！ | とりあえずやってみよう！C++, Python

こんな方におすすめ

Pythonで株価を分析したい
Pythonで日本の株価の最新情報を分析したい
モンテカルロシミュレーションからValue at Riskを算出したい

今回は、Pythonで株価を分析していきたいと思います。日本株の最新情報を無料サイトから取得する方法も理解できると思います。データ取得の方法だけ利用して、自分なりの株価分析を行うのもいいと思います。

この記事では、過去データからモンテカルロシミュレーションを用いてVaR(Value of Risk)を算出まで行います。

<br />

プログラムの内容

必要なモジュールをインポートします。

import pandas as pd
from pandas import Series,DataFrame
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_style('whitegrid')
%matplotlib inline

# webからデータを読み込むためのモジュール
from pandas_datareader import data, wb

# 日付と時刻を扱うためのモジュール
from datetime import datetime

from pandas.plotting import *

import pandas as pd

from pandas import Series,DataFrame

import numpy as np

import matplotlib.pyplot as plt

import seaborn as sns

sns.set_style('whitegrid')

%matplotlib inline

# webからデータを読み込むためのモジュール

from pandas_datareader import data, wb

# 日付と時刻を扱うためのモジュール

from datetime import datetime

from pandas.plotting import *

stooqから日本株の最新情報を読み込みます。今回は例としてオリエンタルランドの株を読み込み、直近１年間のデータで解析していきます。

OLC_df=data.DataReader('4661.JP','stooq')

#直近１年間のデータ
end = datetime.now()
start = datetime(end.year - 1,end.month,end.day)
#取得したデータを昇順にならびかえ
OLC_df_s=test_df.sort_index(ascending=True)
#１年間の区間を設ける
OLC_df_s_one_year=OLC_df_s[start:end]
OLC_df_s_one_year

OLC_df=data.DataReader('4661.JP','stooq')

#直近１年間のデータ

end = datetime.now()

start = datetime(end.year - 1,end.month,end.day)

#取得したデータを昇順にならびかえ

OLC_df_s=test_df.sort_index(ascending=True)

#１年間の区間を設ける

OLC_df_s_one_year=OLC_df_s[start:end]

OLC_df_s_one_year

次のような直近１年間のデータが得られます。（PCのみ途中まで表示）

Open	High	Low	Close	Volume
Date
2019-06-24	13231.40	13445.77	13191.51	13430.81	477990
2019-06-25	13370.99	13525.54	13361.02	13430.81	413403
2019-06-26	13356.03	13430.81	13281.25	13356.03	353127
2019-06-27	13321.13	13351.05	13181.54	13306.18	440581
2019-06-28	13266.29	13346.06	13211.46	13301.20	486716
2019-07-01	13420.85	13470.70	13311.16	13470.70	433762
2019-07-02	13450.76	13495.63	13400.90	13490.64	379705

このデータをデータフレームのdescribe()関数を用いて、概観を掴んでみます。（PCのみ表示）

OLC_df_s_one_year.describe()

1	OLC_df_s_one_year.describe()

	Open	High	Low	Close	Volume
count	243.000000	243.000000	243.000000	243.000000	2.430000e+02
mean	14782.910370	14951.739918	14628.067531	14802.666420	9.841943e+05
std	1060.187515	1002.291284	1071.264795	1028.682657	7.518405e+05
min	11481.930000	12420.450000	11232.320000	12051.030000	2.325760e+05
25%	14080.000000	14230.000000	13893.805000	14084.460000	5.237235e+05
50%	14876.590000	15046.320000	14826.660000	14931.440000	7.520000e+05
75%	15475.640000	15597.205000	15342.925000	15472.500000	1.155901e+06
max	16913.380000	16953.320000	16688.730000	16953.320000	4.874360e+06

終値でグラフを作成してみます。

plt.plot(OLC_df_s_one_year['Close']

1	plt.plot(OLC_df_s_one_year['Close']

次に出来高でグラフを作成してみます。

plt.plot(OLC_df_s_one_year['Volume']

1	plt.plot(OLC_df_s_one_year['Volume']

終値から移動平均線を描いてみます。

# 10日、20日、50日の移動平均を描いてみます。
ma_day = [10,20,50]

for ma in ma_day:
    column_name = "MA {}".format(str(ma))
    OLC_df_s_one_year[column_name]=pd.rolling_mean(OLC_df_s_one_year['Close'],ma)

plt.plot(OLC_df_s_one_year[['Close','MA 10','MA 20','MA 50']])

# 10日、20日、50日の移動平均を描いてみます。

ma_day = [10,20,50]

for ma in ma_day:

column_name = "MA {}".format(str(ma))

OLC_df_s_one_year[column_name]=pd.rolling_mean(OLC_df_s_one_year['Close'],ma)

plt.plot(OLC_df_s_one_year[['Close','MA 10','MA 20','MA 50']])

日々の変化率を算出し、グラフにします。

# pct_changeを使います
OLC_df_s_one_year['Daily Return'] = OLC_df_s_one_year['Close'].pct_change()
# 変化率をプロット
plt.plot(OLC_df_s_one_year['Daily Return'])

# pct_changeを使います

OLC_df_s_one_year['Daily Return'] = OLC_df_s_one_year['Close'].pct_change()

# 変化率をプロット

plt.plot(OLC_df_s_one_year['Daily Return'])

変化率のヒストグラムを作成します。dropna()はNaNを取り除く関数です。

sns.distplot(OLC_df_s_one_year['Daily Return'].dropna(),bins=100,color='skyblue')

1	sns.distplot(OLC_df_s_one_year['Daily Return'].dropna(),bins=100,color='skyblue')

変動幅の平均値と偏差からリスクを算出します。

rets = OLC_df_s_one_year['Daily Return'].dropna()

area = np.pi*20

plt.scatter(rets.mean(), rets.std(),alpha = 0.5,s =area)

# Set the x and y limits of the plot
plt.ylim([0.01,0.030])
plt.xlim([-0.005,0.01])

#Set the plot axis titles
plt.xlabel('Expected returns')
plt.ylabel('Risk')
plt.annotate(
    label, 
    xy = (x, y), xytext = (0, 100),
    textcoords = 'offset points', ha = 'right')
    #rrowprops = dict(arrowstyle='-', connectionstyle= 'arc3'))

rets = OLC_df_s_one_year['Daily Return'].dropna()

area = np.pi*20

plt.scatter(rets.mean(), rets.std(),alpha = 0.5,s =area)

# Set the x and y limits of the plot

plt.ylim([0.01,0.030])

plt.xlim([-0.005,0.01])

#Set the plot axis titles

plt.xlabel('Expected returns')

plt.ylabel('Risk')

plt.annotate(

label,

xy = (x, y), xytext = (0, 100),

textcoords = 'offset points', ha = 'right')

#rrowprops = dict(arrowstyle='-', connectionstyle= 'arc3'))

5パーセンタイルの位置にある変動率を算出する。この意味は95%の確率でこれ以上日々の変動率はないことを示す。下の結果だと-3.1%なので、100万円の株価であったら、96.9万円（-3.1万円）は下回らないということが言えます。これがVaRの考え方です。

rets.quantile(0.05)

1	rets.quantile(0.05)

-0.03196286211312713

1	-0.03196286211312713

１年間の株価の情報から、日々の変動をモンテカルロ法を用いてシミュレーションを行い、VaRを算出してみる。これにより、１年間保持した時に、どれくらいの損失に収まるかが予測できます。

＃基準を１年にする。
days = 365
dt = 1/days

# 日々の変動の平均
mu = rets.mean()

# 変動の標準偏差
sigma = rets.std()

#モンテカルロシミュレーション関数の作成
def stock_monte_carlo(start_price,days,mu,sigma):
    
    # 戻り値となる価格のリストを返す
    price = np.zeros(days)
    price[0] = start_price
    # driftは変動の平均と時間のかけ算と shockはランダム関数でどちらかに動くのを表現する
    shock = np.zeros(days)
    drift = np.zeros(days)
    
    # 指定された日数で、計算
    for x in range(1,days):
        #  shockを算出
        shock[x] = np.random.normal(loc=mu * dt, scale=sigma * np.sqrt(dt))
        # Driftを算出
        drift[x] = mu * dt
        # これらを使って価格を計算します。
        price[x] = price[x-1] + (price[x-1] * (drift[x] + shock[x]))
        
    return price

＃基準を１年にする。

days = 365

dt = 1/days

# 日々の変動の平均

mu = rets.mean()

# 変動の標準偏差

sigma = rets.std()

#モンテカルロシミュレーション関数の作成

def stock_monte_carlo(start_price,days,mu,sigma):

# 戻り値となる価格のリストを返す

price = np.zeros(days)

price[0] = start_price

# driftは変動の平均と時間のかけ算と shockはランダム関数でどちらかに動くのを表現する

shock = np.zeros(days)

drift = np.zeros(days)

# 指定された日数で、計算

for x in range(1,days):

# shockを算出

shock[x] = np.random.normal(loc=mu * dt, scale=sigma * np.sqrt(dt))

# Driftを算出

drift[x] = mu * dt

# これらを使って価格を計算します。

price[x] = price[x-1] + (price[x-1] * (drift[x] + shock[x]))

return price

ランダム回数を5回でモンテカルロの動きについて確認してみる。最初の終値から始めてみます。

start_price = rohm_df.iloc[0,3]

for run in range(5):
    plt.plot(stock_monte_carlo(start_price,days,mu,sigma))
plt.xlabel("Days")
plt.ylabel("Price")  
plt.title('Monte Carlo Analysis for OLC')

start_price = rohm_df.iloc[0,3]

for run in range(5):

plt.plot(stock_monte_carlo(start_price,days,mu,sigma))

plt.xlabel("Days")

plt.ylabel("Price")

plt.title('Monte Carlo Analysis for OLC')

シミュレーションの回数を増やして、VaRを算出します。10000回でやってみます。

# 回数
runs = 10000

# 結果保持用
simulations = np.zeros(runs)

# 表示のオプション
np.set_printoptions(threshold=5)

for run in range(runs):    
    # シミュレーション
    simulations[run] = stock_monte_carlo(start_price,days,mu,sigma)[days-1];

# 回数

runs = 10000

# 結果保持用

simulations = np.zeros(runs)

# 表示のオプション

np.set_printoptions(threshold=5)

for run in range(runs):

# シミュレーション

simulations[run] = stock_monte_carlo(start_price,days,mu,sigma)[days-1];

# モンテカルロで得られた結果からヒストグラムを表示
plt.hist(simulations,bins=200)

# 1パーセンタイルの位置を設定
q = np.percentile(simulations, 1)

# 最初の株価
plt.figtext(0.6, 0.8, s="Start price: {:0.2f}".format(start_price))
# 最終的な株価の平均値
plt.figtext(0.6, 0.7, "Mean final price: {:0.2f}".format(simulations.mean()))
# VaR (信頼区間99%で計算）
plt.figtext(0.6, 0.6, "VaR(0.99): {:0.2f}".format(start_price - q))

# 1パーセンタイル
plt.figtext(0.15, 0.6, "q(0.99): {:0.2f}".format(q))

# 1% クォンタイルに線を入れます
plt.axvline(x=q, linewidth=4, color='r')

# タイトル
plt.title("Final price distribution for OLC Stock after {} days".format(days), weight='bold');

# モンテカルロで得られた結果からヒストグラムを表示

plt.hist(simulations,bins=200)

# 1パーセンタイルの位置を設定

q = np.percentile(simulations, 1)

# 最初の株価

plt.figtext(0.6, 0.8, s="Start price: {:0.2f}".format(start_price))

# 最終的な株価の平均値

plt.figtext(0.6, 0.7, "Mean final price: {:0.2f}".format(simulations.mean()))

# VaR (信頼区間99%で計算）

plt.figtext(0.6, 0.6, "VaR(0.99): {:0.2f}".format(start_price - q))

# 1パーセンタイル

plt.figtext(0.15, 0.6, "q(0.99): {:0.2f}".format(q))

# 1% クォンタイルに線を入れます

plt.axvline(x=q, linewidth=4, color='r')

# タイトル

plt.title("Final price distribution for OLC Stock after {} days".format(days), weight='bold');

この結果から、オリエンタルランドの株を１年間保有した場合、99%の確率で543円の損失以内に収まるということが言えます。シミュレーションを行う期間や回数によって結果は変わりますが、リスク解析の目安として使用できると思います。

まとめ

日本の株価の最新情報の入手法と、Pythonを使ったリスク解析法の一例ついて説明してみました。オリエンタルランドの株価を用いて試しましたが、自分が保有している株や自分の会社でトライしてみると楽しいと思います。