SARIMAXで株価予測してみた

こんにちは。
株価って一概にグラフ理論だけじゃ未来予測できないだろうな、、と思いつつ敢えて今回は挑戦してみました。

SARIMAXは、パラメータがorder =(p,d,q)、seasonal_order = (P,D,Q,s)があり、pは自己回帰係数、qは移動平均係数(厳密にはどちらも次数)とdは差分の距離です。
数2の数列を思い出しますね、等差数列ってd使ってましたよね。

そして、時系列データは自然要因による季節に特有の変化と、社会的慣習や制度による季節 変化を含んでいます。夏にはエアコンが売れたり、ボーナス時期には売り上げが増えたりすることはよく知られていますね。これらの季節要因を考慮して調整を入れるものを季節調整と言います。P,D,Qは季節調整を入れたp,d,qの値で、sは季節調整の周期を決定する要素です。
p,d,qは0から2、P,D,Qは0か1がいいようです。
だいたい1年ごとに周期があるのでsは12が通常なようです。

自己回帰モデルは時点 t におけるモデル出力が時点 t 以前のモデル出力に依存する確率過程で、移動平均モデルは将来の予測値は過去の予測値と実績値との誤差により決まるとしたときの分析対象予測過程です。これらを組み合わせたのがARMA、季節調整を入れたものがSARIMAXとなります。

数学的な説明は省きましたが、さらに数式から理解したい方はこのサイトを参考にするとよいでしょう。
to-kei.net

今回はstatsmodelsというライブラリを主に用いるのでpip installを済ませておく必要があります。

import pandas as pd
import numpy as np
import io
from statsmodels import api as sm
from datetime import datetime
from statsmodels.tsa.statespace.sarimax import SARIMAX
from statsmodels.tsa.arima_model import ARIMA
import warnings
import itertools

次に、株価のcsvファイルをダウンロードしてきます。
こちらのサイトから引っ張ってきました。
欠損値処理をdropnaでおこないます。
kabuoji3.com

kabu_data = pd.read_csv("csv")
kabu_ = kabu_data.drop(columns = ["date","start","high","low","close","dekidaka","end"])
kabu_.dropna()

次に、ADF検定をします。
ADF検定とは、単位根過程(yt が非定常過程、また差分系列yt−yt−1=Δytが定常過程である時、ytは単位根過程である)に従うxt,yt を回帰分析を行うと、まったく関係のないxt,ytの間に有意な相関を見出してしまう「見せかけの回帰」が行われます。「見せかけの回帰」は、例えば株価と二酸化炭素濃度という関係のない２つの間に相関性を導き出してしまいます。

そのため時系列データに対して回帰分析を行う時は「見せかけの回帰」に気を付ける必要があります。

adf_result = sm.tsa.stattools.adfuller(kabu_data["ave"],autolag='AIC')
adf = pd.Series(adf_result[0:4], index=['Test Statistic','p-value','#Lags Used','Number of Observations Used'])
print(adf)

結果はこちら

注目すべきはp値とAIC値で、単位根過程の帰無仮説を棄却できる範囲にp値があるのでok、一方でAIC値はちょっと大きめ。一桁台なら理想だけど。

あとは学習させてしばらく先を予測して見ます。

100個学習して100個先を予測、ですね。
forecast( )関数の引数に100を入れます。

N = 100

test = kabu_data.loc[:,'ave']
kabu = kabu_data.loc[:N,'ave']
print(kabu)

SARIMA_kabu = sm.tsa.statespace.SARIMAX(kabu.astype(float),order=(1,0,1),seasonal_order = (1,1,1,12), enforce_stationarity = False, enforce_invertibility = False,trend = "n").fit(trend='nc',disp=False)
print(SARIMA_kabu.summary())
pred = SARIMA_kabu.predict()
pred2 = SARIMA_kabu.forecast()


kabu_diff = kabu_.diff()
kabu_diff.index = kabu_data.index
kabu_diff = kabu_diff.dropna()

plt.plot(kabu_,color = "b")
plt.plot(pred,color = "r")
plt.plot(pred2,color="y")
plt.ylim([2000,3000])
plt.show()

結果はこちら。

青が元データ、赤が学習過程、黄色が予測です。

ぜんっぜんグラフの形似てませんね

次数をいじったりなんらかのパラメータを変えていく必要はありそうですね、、、とりあえず予測モデルとして理論上可能ということはわかりました。
もう少し外乱の少ないネタならSARIMAXは予測機能としては遜色ないのでしょう。

スクレイピングで手っ取り早く欲しい画像を集める

こんにちは、寝正月で元旦を消費してしまいました。
2日の今日はブログを書こうと心に決めました、、、

今回は、スクレイピング技術試してみようかなという記事です。
スクレイピングとは、ウェブサイトから情報を抽出するコンピュータソフトウェア技術のことです(Wikipediaより)

プログラミングを使った自動化技術のうちの一つですね、作業効率化にたいへん寄与してくれるはずです。

pythonを使ったコードですが、その前にいくつかライブラリをpipコマンドでインストールしておきましょう。

必要なものとしては、

bs4：BeautifulSoupと呼ばれるHTML や XML から狙ったデータを抽出するためのライブラリです。

requests : HTMLの取得に使います。

urlib : URL を扱う幾つかのモジュールを集めたパッケージです

argparse : 自動的にヘルプと使用方法メッセージを生成し、ユーザーが不正な引数をプログラムに指定したときにエラーを発生させます。

このくらいが新規にpipインストールすべきであろうライブラリで、他はデフォルトであるのではないかと思います。

さっそく必要なライブラリのimportです。

import bs4
import requests
import re
import urllib.request, urllib.error
import os
import argparse
import sys
import json

次に、BeautifulSoupを用いてURLを引っ張り出してきます。
そしてterminal上で -s:検索ワード(英字表記のみ) -n:数字で何枚取るか -o:自分のローカルのどこに保存するかのパス
を記載します(後述)

def get_soup(url,header):
    return bs4.BeautifulSoup(urllib.request.urlopen(urllib.request.Request(url,headers=header)),'html.parser')

def main(args):
    parser = argparse.ArgumentParser(description='Options for scraping Google images')
    parser.add_argument('-s', '--search', default='kensakuword', type=str, help='search term')
    parser.add_argument('-n', '--num_images', default=10, type=int, help='num of images to scrape')
    parser.add_argument('-o', '--directory', default=' 保存したい場所のパス', type=str, help='output directory')
    args = parser.parse_args()

    query = args.search.split()
    query = '+'.join(query)
    max_images = args.num_images

続いて、画像をフォルダでグルーピングします。

save_directory = args.directory + '/' + query
    if not os.path.exists(save_directory):
        os.makedirs(save_directory) 

そして最後にお決まりごとと画像を取得している様子がわかるようなコードを追加していきます。基本的にjpg形式で取得してます。

url="https://www.google.co.jp/search?q="+query+"&source=lnms&tbm=isch"
    header={'User-Agent':"Mozilla/5.0 (Windows NT 6.1; WOW64) AppleWebKit/537.36 (KHTML, like Gecko) Chrome/43.0.2357.134 Safari/537.36"}
    soup = get_soup(url,header)
    ActualImages=[]

    for a in soup.find_all("div",{"class":"rg_meta"}):
        link , Type =json.loads(a.text)["ou"]  ,json.loads(a.text)["ity"]
        ActualImages.append((link,Type))
    for i , (img , Type) in enumerate( ActualImages[0:max_images]):
        try:
            Type = Type if len(Type) > 0 else 'jpg'
            print("Downloading image {} ({}), type is {}".format(i, img, Type))
            raw_img = urllib.request.urlopen(img).read()
            f = open(os.path.join(save_directory , "img_"+str(i)+"."+Type), 'wb')
            f.write(raw_img)
            f.close()
        except Exception as e:
            print ("could not load : "+img)
            print (e)

if __name__ == '__main__':
    from sys import argv
    try:
        main(argv)
    except KeyboardInterrupt:
        pass
    sys.exit()

上記のコードにて、私はこのように検索しました、ちょうど壁紙欲しかったので、、、

取得画面はこんな具合

そして、自分の保存したいローカルにはこのように自動的に画像が保存されていきました。

以上、スクレイピングを一通り経験してみたと言った感じです。
単純に楽なのと、ちょっと深層学習いろいろテストしたいなという野望があるので学習データ画像取ってくるにはもってこいの技術ですよね。
論文検索でキーワードからURLのみ引っ張ってくる、などにも応用できそうです。

家でラーメン二郎つくった-どん二郎-

こんにちは。

今回は地方でなかなか話題の本場のラーメン二郎に有り付けない人必見の自作ラーメン二郎の記事です！

あの中毒性のあるラーメン、体が欲しても店舗まで遠いし家でできないかなーと思い調べて実際に作ってみました！

まずは材料

世界一美味しい「どん二郎」の作り方 誰も思いつかなかった激ウマ! B級フードレシピ

どん兵衛(うどん)1つ
キャベツ　お好みで
もやし　お好みで
チャーシュー　市販のひとつ
牛脂　１つ

せっかく自分で二郎作るんだから、お好みでええねん！

特に牛脂が二郎っぽくするのには必須らしい、ほう。
牛脂は業務用スーパーなどにいけば無料でもらえるし、お肉コーナーに涼しい顔して鎮座してたりするので手に入りやすいですよ

もやしとキャベツを茹でて、その熱湯でをどん兵衛に注ぎ待つこと3分。
牛脂も同時に入れておきます。

所要時間10分から15分ほどで完成！

見た目なんとなく二郎だ

さっそく実食！

一口めは二郎っぽさがちゃんとでてました、牛脂パワーか

麺もわりかし二郎に近い気がする、合格点

しかし問題がひとつ

味うっっっっす

これ醤油とかいれるべきでしたね、二郎でいうからめってやつですね

しかし材料費は200~300円で済みますし、待ち時間もさほどでもなし、利点はたくさんあります。
もっと本格化したいなら背脂を手に入れて野菜マシアブラマシが実現できるような準備をするといいんじゃないかな

以上、お手軽二郎でした！

ラーメン二郎にまなぶ経営学 ―大行列をつくる26(ジロー)の秘訣

主成分分析をゲキ簡単な手法でやってみた

こんにちは。

今回は統計学にカテゴライズされる記事で、主成分分析というものをやってみます。

主成分分析は詳しい説明は省きますが、簡単に言うと、ごちゃごちゃしたいろんな要素で構成される結果を、少ない要素で表すことはできないかという試みのこと。例えば、身長が伸びれば体重はだいたい大きくなって正の相関があるわけだから、この2つの要素はぶっちゃけどっちか一つで良くね？という要領で要素を減らしていきます。

御託はこのくらいにして、今回はContrastive PCA(
GitHub - abidlabs/contrastive: Contrastive PCA)
を用いて,
こちらのデータ
https://raw.githubusercontent.com/maromato/RNA_data_analysis_1/master/TNBC10vNormal10_Counts_2.csv

を解析しました。このデータはハーバード大学医学部バイオインフォマティックスコースで利用されている乳がん細胞と正常細胞の遺伝子発現データです。

from sklearn.cluster import AffinityPropagation, KMeans, DBSCAN, SpectralClustering
from sklearn.manifold import MDS, TSNE, Isomap
from sklearn.metrics import silhouette_score

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from scipy.linalg import logm, expm
from contrastive import CPCA

これらのライブラリを用意します。
特にcontrastiveモジュールは特有のモジュールなのでgithubから必ずcloneしてそのフォルダがある位置でプログラム実行してくださいね。

それでまあ学習用データ、バリデーションデータと分けて解析したわけですが、matplotlib使わなくてもcontrastiveモジュールにグラフ表示機能が内包されているのか勝手に結果が表示されました。

黒が正常細胞、赤が乳がん細胞のようです。しっかりと分類されてるのがみて取れます。2軸で今回は表されてるわけですね。ほんとは多次元なところを平面2次元グラフに表すところがミソなようです。

それにしてもCPCAが便利すぎるでしょ(笑) もう少し原理に基づいてプログラムを書くためにscikit-learnのPCAで次はなにかしらやってみようと思います。

スラスラプログラム書きたい場合、プログラミングは動画教材使ったり、なんならオンラインスクールのぞいてみるのもいいですね、成長は早いです。

参考にしたサイト
https://note.mu/tokitky/n/nf2cadefc4633

pythonでEMアルゴリズムによる混合ガウス分布推定やってみた

こんにちは。

今回は、pythonを使ってEMアルゴリズムによる混合ガウス分布推定やってみました。

まず、EMアルゴリズムってなんすか？って方はこちらの記事を参考にしてみてください。

qiita.com

ちなみに数学的に理解するには尤度関数とか偏微分とか果てはイェンセンの不等式(私もこれは知りませんでした)とかが身についてないと理解できないです、それが苦痛な人は数学的理解はひとまず置いといて、概念的になにやってるのかを自分の口で説明できればそれでいいんじゃない？と思います。

自分的に説明するならば、正規分布がとある割合で混ざってて、それの比率を取り出すためにパラメータ推定するもの、ってイメージですね。

ひとまず結果見せた方が早いと思うのでやってみますよ

まずは正規分布まぜこぜなものを作成

import numpy as np
import matplotlib.pyplot as plt
from sklearn import mixture

dist1 = np.random.normal(0,1,75000)
dist2 = np.random.normal(3,0.5,25000)
dist = np.r_[dist1, dist2]

plt.hist(dist,bins = 100)
plt.title('正規分布 N(0,1)とN(3,0.5)を7500対2500で混合した分布')

こんなかんじっすね

ではGaussianMixtureっていうモジュールを使って解析してみましょう

dist1 = np.reshape(dist,(np.shape(dist)[0],1))
max_n_components = 8
lowest_bic = np.infty
lowest_bic_ix = 0

for i in range(1,max_n_components):
    clf = mixture.GaussianMixture(n_components = i, covariance_type = 'full')
    clf.fit(dist1)
    bic = clf.bic(dist1)
    print('bic',i,bic)
    if bic < lowest_bic:
        lowest_bic = bic
        lowest_bic_ix = i
print ('lowest  bic case no', lowest_bic_ix)

clf = mixture.GaussianMixture(n_components = lowest_bic_ix, covariance_type = 'full')
clf.fit(dist1)
print('weights:', clf.weights_)
print('means:', clf.means_)
print('std dev :', np.sqrt(clf.covariances_))

これで得た結果がこちら

bicっていうのはベイズ情報量規準ってやつでこれが一番小さい要素分布数を選べばいいわけです。すると、書いてある通り2が最適解となります。
つまりふたつの正規分布が混ざってるよってことを言ってるわけですね。

かなり恣意的なやりかたでしたが、ゆくゆくは実際の医療データとかでやってみたいなとは思います。練習ってことで！

機械学習で家の値段の決定要因を調べてみた

こんにちは。

今回は、機械学習プログラムを用いて重回帰分析を実際にやってみました。
こういうのって一度それなりのプログラム書けばあとはおんなじようなプログラムでマイナーチェンジで済むので財産としてローカルに残しておく次第です。

さて、今回はKaggleにあるDatasetsのうちのこちらを用いました

www.kaggle.com

CSVファイルをダウンロードして、いざスタート。

その都度プログラムを確認したかったので今回はJupyter Notebookを使って一行一行確認しながらやっていきました。

ひとまずデータの中身を確認しましょう

import pandas as pd
house = pd.read_csv('train.csv'
                   )
house.head()

すると、

こんな具合でどんなデータが揃ってるか見れます。

次に、こんな便利な関数があることも記しておきます。

house.describe()

なんとこの一言だけで平均、最大最小、標準偏差などなど一気に出してくれます。

とりあえず必要なモジュールをインポート

import pandas as pd 
import matplotlib.pyplot as plt 
from sklearn.cross_validation import train_test_split 
from sklearn import linear_model 
from sklearn.metrics import r2_score 
import warnings 
warnings.filterwarnings('ignore') 

%matplotlib inline

データを訓練用と検証用にわける

path = 'train.csv' 
house= pd.read_csv(path)

# 目的変数と説明変数を作成
y = house['SalePrice']


X = house[['OverallQual', 'YearBuilt', 'TotalBsmtSF', 'GrLivArea']]

# 学習用、検証用データに分割
(X_train, X_test, y_train, y_test) = train_test_split(X, y, test_size = 0.2, random_state = 5)

今回のデータでは、全部の項目を上げると煩わしいので
目的変数はSalePriceとして、
説明変数は
'OverallQual':Overall material and finish quality
'YearBuilt' :Original construction date
'GarageCars':Size of garage in car capacity
'GrLivArea':Above grade (ground) living area square feet

の4要素を説明変数としてみて、どれが一番家の値段を引き上げる要因なのかを検討してみました。

#回帰モデル
model = linear_model.LinearRegression()

#モデルの学習
model.fit(X_train, y_train)

fit関数でサクッと学習です

ではどの要素が相関が強いかみてみましょう

print(pd.DataFrame({"Name":['OverallQual', 'YearBuilt', 'GarageCars', 'GrLivArea'],"Coefficients":model.coef_}).sort_values(by = 'Coefficients').round(4))

この表示結果がこちら

どうやらOverallQual:Overall material and finish qualityが主要因で、今回選択した中だとGrLivArea:Above grade (ground) living area square feetが家の値段引き上げ要因にはもっともならないという結論に至りました

ちなみにバリデーションデータとの比較を最後にきちっとやって

#予測値の算出
y_train_pred = model.predict(X_train) 
y_test_pred = model.predict(X_test) 

#決定係数を表示
print('R2(Train) : % 0.3f' % r2_score(y_train,  y_train_pred))
print('R2(Test) : % 0.3f' % r2_score(y_test, y_test_pred))

#グラフを描画する
plt.scatter(y_train, y_train_pred)
plt.scatter(y_test, y_test_pred)
plt.show()

トレーニングデータと検証用データでそんなに相違はないようです

こんな感じで重回帰分析をぶん回して今後も遊んでいこうと思います！

旅する大学生ならEduroamを使おう！

こんにちはっす。

今回は旅をしたりよく各地に移動する大学生におすすめの、お金のかからないノマド生活推奨記事です。

私も実家は通ってる大学からは遠い場所にあり、移動先でwifiが使いまくれる場所ないかなと常々思っておりました。

そこで、Eduroamを使います。

Eduroamが使用可能な機関はこちらでチェックです

www.eduroam.jp

ただし、注意点としては2点あります。

1. 認証時に必要な中間CA証明書がスマートフォンや端末にインストールされている必要がある．

2. IDには，各々の大学IDの後ろに「@"各々の大学名アルファベット".ac.jp」を付加する．

これらの条件さえ満たせば移動先にある大学、研究所で優雅にwifiを用いることができます。

簡単に申し込みして利用できる大学図書館も多々ありますし、学生という立場を存分に利用しましょう！
スタバなどのカフェも良いのですがちょっと高いですもんね学生にとっては。