僕がデータ分析者として覚醒するまで

しがない会社員がデータ分析者として覚醒するまでのブログ

不偏推定量

ブログを書き始めたころに、覚えたはずの不偏推定量。 今日統計検定2級の問題を解いていて、不偏推定量の意味が分からず、悶絶。 ブログに記録したことは覚えてたのに記憶に定着していなかった。。。 下記のURLのような意味なんだけど、結局のところ自分で書き換えられないと覚えたことにはならないよなぁ。。。 bellcurve.jp

そもそも、確率分布におけるE[X],V[X]とかの意味がすっと出てこない時点でだめだなぁ。。。

層化抽出法

層化抽出法とは

  • 母集団をあらかじめいくつかの層(グループ)に分けておき、各層の中から必要な数の調査対象を無作為に抽出する方法

bellcurve.jp

統計検定2級の試験問題に出てきたんだけど、母集団をいくつかのグループに分ける必要がある。その分けるグループは当然ながら偏った性質を持っていてはいけない。

各層はランダムかつ母集団とよく似た性質にわけるようにする。 実際に標本集団を作るときは偏った性質の層を作らないようにする必要がある。

はやく統計検定2級受けたい。。。

自己相関関数と偏自己相関関数

コレログラム

コレログラムは、時系列データにおいて、周期性や自己相関がどの時点にあるかを 判断できるグラフである。
時系列データを分析する上で、どこに相関があるかは重要で、ARモデルなどで係数を推定するために必要なものである。

自己相関

自己相関関数は、原系列に対して、遅れのあるデータを作成し、それらの1次遅れ・2次遅れ・・・n次遅れのデータを作成し、それぞれ相関をとることで周期性や何次遅れのデータが原系列にとって相関があるかを知ることができる。
bellcurve.jp

偏自己相関

偏自己相関は、自己相関が時系列データのn次遅れのデータそれぞれに行うに対して、偏自己相関は自己相関係数から、時間によって受ける影響を除去した自己相関となる。

r_{xy\cdot z} = \frac{r_{xy} - r_{xz}r_{yz}}
{\sqrt{1 - r_{xz}^2} \sqrt{1 - r_{yz}^2}} 

算術平均と幾何平均

算術平均と幾何平均

算術平均

  • 学生時代に学ぶ平均値のこと
    [1,2,3,4]の数値があった場合、(1+2+3+4)/4 = 2.5 となる。 これは、特に違和感なく理解できる。

幾何平均

oto-suu.seesaa.net

  • 幾何平均違和感ありまくり。
  • 何に使えるかは成長”率”とか寄与”率”とかの割合の平均を出したいときに使える。
  • 預金を例にとると、この幾何平均を使うと複利を含めた計算で平均を出すことができる。

あぁ、そういうことか。 割合の平均を求めるには、単に足して個数分で割っただけでは、複利みたいな増加分の計算ができないってことか。

ちょっと納得。これで明日から使える。

VARモデルについてのメモ

VARモデルについてのメモ

概要

  • VARとは、VectorAutoRegression(ベクトル自己回帰)モデルである。
  • VARはARモデル(単変量)を拡張したものである。
  • VARは、多変量のモデルで、目的変数が自己のみではなく、他の変量にも依存する変量を解析するために用いる。

建築設備のエネルギーへの適用

  • 室内負荷の予測を行うために、目的変数:室内負荷熱量 、説明変数:外気温として、回帰分析を行う。
  • 空調機電力量の予測を行うために、目的変数:空調機電力量、説明変数:室内負荷熱量 or 外気温度として、回帰分析を行う。

Pythonでの実装するには

Statsmodels(VARモデルの使い方)
※英語

説明変数の取得方法

外気温

  • 各建物で蓄積しているデータを用いる。
  • Webページ上からスクレイピングでデータを取得する。
  • 有料のサービスを使用する。 n-kishou.com

室内負荷熱量

  • VARで予測した結果を用いる。
  • 翌日の設定温度からおおよその室内負荷を演算子使用する。

まとめ

室内負荷は、外気温度と相関が大きくあるため、単純なARモデルで予測するよりも、 VARモデルなど、説明変数を用いた予測の方が精度は高くなると考えられる。 データの分析はしていない。 室内負荷の予測方法なにかいいものはないものか。。。

建築設備のエネルギーデータをpandasで処理するときの注意

Pandasで処理するときの注意

建築設備のエネルギーデータ

一般に、建築設備におけるエネルギーデータは、カラム名が日本語でついます。 f:id:snuow:20180124115522j:plain

pandasで読み込む際も、下記のようにencoding = 'Shift-JIS'と指定してあげないとエラーを吐きます。

# header:0行目,index_col:0行目,パースする日時:0行目
df= pd.read_csv(r'Data\\' + file + '.csv',encoding="Shift-JIS",header=0,index_col=0,parse_dates=[0])

もちろん、時系列データとしてタイムスタンプもついていますので、それぞれ指定してあげないといけません。

詳しい処理は、省きますが、元のヘッダ情報を使用して、四則演算論理演算を行ったうえで、ファイルを作成するという処理を行いました。

問題

何も考えずに、下記のコードで出したら、エラーも吐かずに出力できました。が、よくよく中身を見ると文字化けパレードだったので、別システムで読み込みができずにエラーログを大量に出していました。

df.to_csv(r'add_items_Result\\' + file + '.csv')

解決策

結果として、下記のように読込の時と同じように出力時もencoding = 'Shift-JIS'を指定してあげないとダメみたいです。

df.to_csv(r'add_items_Result\\' + file + '.csv',encoding='Shift-JIS')

まとめ

文字コードって難しい。

時系列データ分析SARIMAモデル 「学習データのinputについて」【statsmodels】

stasmodelsのSARIMAX

statsmodelsのSARIMAXに学習データを読み込む際の注意点についてのメモです。

SARIMAXのモデル作成

学習データ

statsmodelsではSARIMAモデルを扱うのに下記のクラスを使用します。 学習データは、endogに入力します。(arrayやDataframeに対応)←今回の問題箇所

 class statsmodels.tsa.statespace.sarimax.SARIMAX(endog, exog=None,
 order=(1, 0, 0), seasonal_order=(0, 0, 0, 0), trend=None, 
measurement_error=False, time_varying_regression=False, 
mle_regression=True, simple_differencing=False, 
enforce_stationarity=True, enforce_invertibility=True, 
hamilton_representation=False, **kwargs)

私の場合は下記のようにpandasのDataframeでフォルダ内のcsvを全部読み込んで結合してendogに入力していました。 もっとスマートに作りたかったのは内緒です。

import pandas as pd
import glob,os
def CombineFiletoDF(folderpath):
    csvs = glob.glob(folderpath + '\*.csv')
    for csv in csvs:
        if csv == csvs[0]:
            df = pd.read_csv(csvs[0],dtype=np.float64,encoding='SHIFT-JIS',index_col=[0],parse_dates=[0])
        else:
            df = df.append(pd.read_csv(csv,encoding='SHIFT-JIS',index_col=[0],parse_dates=[0],skiprows=0)) 
    
    df = df.resample('H').asfreq()
    return df

で、下記のようにtrainDFつくって入力してフィッティングして、実データと比較してました。 ソフトウェアエンジニアではないので、私の場合はプロトタイプ作成を行った形です。 出力結果としては、タイムスタンプ付きのDataframe(Series?)として出力されます。

mod = sm.tsa.statespace.SARIMAX(endog = trainDF,#Dataframeで入力
                                       enforce_stationarity=False,
                                       enforce_invertibility=False,
                                       order=SARIMAorder,
                                       seasonal_order=SARIMAORDER,
                                       freq='H'
                                      )

prediction_result = mod.fit(disp=False)

prediction_result.predict(start=49,end=72)

# Output
2017-01-24 01:00:00       3.074818
2017-01-24 02:00:00       3.000000
2017-01-24 03:00:00       3.925182
2017-01-24 04:00:00       3.074818
2017-01-24 05:00:00       3.925182
2017-01-24 06:00:00       3.074818
2017-01-24 07:00:00       3.925182
2017-01-24 08:00:00       3.074818
2017-01-24 09:00:00    1426.854329
2017-01-24 10:00:00    1479.087658
2017-01-24 11:00:00    1384.257394
2017-01-24 12:00:00    1287.257394
2017-01-24 13:00:00    1279.074818
2017-01-24 14:00:00    1165.202675
2017-01-24 15:00:00    1131.945281
2017-01-24 16:00:00    1111.496370
2017-01-24 17:00:00    1027.658846
2017-01-24 18:00:00     802.657166
2017-01-24 19:00:00     312.232253
2017-01-24 20:00:00     313.782805
2017-01-24 21:00:00     300.346733
2017-01-24 22:00:00      26.996639
2017-01-24 23:00:00       2.224455
2017-01-25 00:00:00       2.925182
Freq: H, dtype: float64

他方、ソフトウェア開発を行っている協力会社の方は、numpyのarrayでendogに入力していました。numpyのarrayにはindexがついていません。要するにどのデータを実際に使ったのかをindexで確認できない状況でした。結果としては、同じ値が出力されています。

mod = sm.tsa.statespace.SARIMAX(endog = np.array(trainDF),#numpyのarrayで入力
                                       enforce_stationarity=False,
                                       enforce_invertibility=False,
                                       order=SARIMAorder,
                                       seasonal_order=SARIMAORDER,
                                       freq='H'
                                      )

prediction_result = mod.fit(disp=False)

prediction_result.predict(start=49,end=72)

# Output
array([    3.0748185 ,     3.        ,     3.9251815 ,     3.0748185 ,
           3.9251815 ,     3.0748185 ,     3.9251815 ,     3.0748185 ,
        1426.85432885,  1479.08765798,  1384.25739448,  1287.25739448,
        1279.0748185 ,  1165.20267547,  1131.94528099,  1111.49637   ,
        1027.65884648,   802.65716595,   312.23225339,   313.78280451,
         300.346733  ,    26.99663893,     2.2244555 ,     2.9251815 ])

start,endの固定値は、協力会社さんの作ったプログラムから出力した結果を決め打ちで入れています。

何が問題って、本来予測したかったのは、2017/1/25であったこと。

Dataframe入力しておけば、すぐ気づけた間違いも、np.arrayにしたことで見えなくなってしまってました。 データの操作が色々あったので、配列のほうが都合がよかったのかと思いますが、僕はDataframeで入力することをお勧めします。

だって、時系列データ分析してるのに、タイムスタンプとったらわけわからんと思うんだ。