1. アメダス気象データ分析チャレンジ！（Python版）　診断的分析

ここでは，過去の気象データと過去の電力消費量データを分析することで気温と電力消費量との関係性を分析します．このような分析を診断的分析（Diagnostic Analytics）と言います．診断的分析は，過去に蓄積されたデータに見られる特徴的な変化の要因（因果関係）を見つけ出すことによって，そのような変化がなぜ起きたのかを明らかにします．どのように読み解きどのような行動に移すのかを判断するのは分析者自身となります．

• 1.1 気象データ分析とは？
• 1.2 使用するデータ
  • 1.2.1 気象庁AMeDASの気温データ
  • 1.2.2 東京電力の電力消費量データ
• 1.3 Pythonプログラムの構成と実行
  • 1.3.1 ライブラリのインポート
  • 1.3.2 気象庁AMeDASの気温データ（csv形式）を読む
  • 1.3.3 東京電力の電力消費量データ（csv形式）を読む
  • 1.3.4 2つのデータを結合する
  • 1.3.5 データの整理（NaNを取り除く）
  • 1.3.6 時系列図を作成する
  • 1.3.7 散布図と回帰直線を作成する
  • 1.3.8 データの整理と分析（仮説１：夏の期間のデータの抽出）
  • 1.3.9 データの整理と分析（仮説２：12～15時のデータの抽出）
  • 1.3.10 データの整理と分析（仮説３：休日のデータの削除）
  • 1.3.11 データの整理と分析（仮説４：お盆・祝日のデータの削除）
• 1.4 まとめ
• 1.5 著作権について

1.1 気象データ分析とは？

ビジネスで活用したい気象データは、過去に観測された気象データではなく、未来の予測された気象データかもしれません．しかし、予測データを適切に活用するためには、過去データから学ばなくてはなりません． そのため，この講義では，まず未来の気象予測データではなく，過去の気象観測データを扱います．

図: データ分析の目的

データ分析を行う上で重要なことは、「そこにデータがあるから」データ分析をするのではなく，何らかの活動に役立てる目的があり，その目的を達成するためにデータ分析をすることになります．そこでまずは，データ分析の目的を言葉で表現してみましょう．分析者はこれまでの知識や経験から「なにかありそうだ」と感じたこと，すなわち直感から仮説を立てます．

図: データ分析は仮説と検証の繰り返し

仮説が立てば，それに応じて，必要なデータと分析の手法が定まってきます．そして，必要なデータを入手して，データを分析できる形へと整理します．得られたデータをグラフや統計値などから分析することで，仮説が正しいかどうかを検証していきます．

ほとんどの場合，1回の分析でよい結果が得られるわけではありません．何度も仮説を立てなおして，データを整理し直し，分析し直すことで，仮説と検証を繰り返す必要があるのです．これを繰り返すことによって，データの中に潜んでいる本質を理解することができるようになります．

このようなデータ分析を行うことで，ビジネスの現場で，データに基づいた適切な意思決定ができるようになります．また，分析者の知識と経験が増えることによって，仮説を立てる際の「なにかありそうだ」という直感の精度があがっていくことになります．

図: データ分析の流れ

データ分析の大まかな流れについて見てみましょう．気象データもビジネスデータもビッグデータともいうべき大量なデータが存在しています．大量のデータをそのまま分析することはできません．まず，はじめに必要な情報だけを並び替えたり，集計したりすることによってデータを要約する必要があります．これにより，データ全体の傾向や特徴が分かるようになります．次に，データの傾向や特徴を視覚的に分かりやすくするためにグラフで表します．数字の羅列だけではよく分からないことも，グラフを作成することによって新たな特徴を発見できるようになります．そして，データやグラフから読み取った定性的な理解に対して，相関係数などの統計値を算出してデータの傾向や特徴を数値で示します．これにより，直感的な理解に根拠を与えることができ，説得力が増すことになります．

図：アメダス気象データ分析チャレンジ（Excel版）https://www.wxbc.jp/weather-challange/

このようなデータ分析の基礎を学習したい場合には，まずはExcelから入ることをおすすめします．気象ビジネス推進コンソーシアムではExcelによる気象データ分析の動画を配信していますので，こちらで学習してみて下さい．

1.2 使用するデータ

1.2.1 気象庁AMeDASの気温データ

早速，Pythonにより気象データとビジネスデータを掛け合わせたデータ分析をしてみましょう．繰り返しになりますが，この講義では，過去の気象データと過去の電力消費量データを分析することで気温と電力消費量との関係性を分析します．気温と電力消費量との間には関係性があるのではないか？ との仮説を立てたわけです．

このプログラムでは，気象データとして，気象庁AMeDASの2017年1月～12月の1時間毎の1年間分の東京の気温データを使用しています．ファイル名は，amedas2017.csvとします．過去の気象データとビジネスデータとの関係性を調べることでデータに潜む本質にせまり，最終的には，予測へとつなげることが目標です．データのダウンロードの方法は，以下の通りとなります （注意：以下の例では仙台のデータをダウンロードしています） ．

気象庁のホームページのリンク： http://www.data.jma.go.jp/gmd/risk/obsdl/

図: まず気象庁ホームページから必要とする地点を選びます

図: 必要とする地点を指定します

図: 次に気象要素を選びます

図: ここでは時間間隔と気象要素を指定します

図: 次に期間を選びます

図: 必要な期間を指定して，「CSVファイルをダウンロード」をクリックします

図: 気象データをcsv形式で入手できました!

ここでは事前に用意したcsv形式のファイルを使用しますので，データをダウンロードしていただく必要はありませんが，ご自身で，異なる地点や期間，異なる気象要素を各自でダウンロードして独自の分析に挑戦してみて下さい．

1.2.2 東京電力の電力消費量データ

また，ここでは、ビジネスデータとして東京電力の電力消費量データを使用します．気象庁AMeDASの気温データと同様に，2017年1月～12月の1時間毎の1年間分の東京電力の電力消費量を入手します．ファイル名は，tepco2017.csvとします．このデータを使う理由は，気象データと同様に連続的（1時間毎）で長期間（1年以上）にわたり入手可能なオープンデータであるだけでなく，気象との相関の高いパラメータであると考えられるからです．つまり，ここでは東京電力の電力消費量は東京の気温と関係しているのではないか？という仮説を立てたことになるわけです． （注意：以下の例では2016年のデータをダウンロードしています）

東京電力でんき予報のリンク： http://www.tepco.co.jp/forecast/index-j.html

図: 東京電力でんき予報のホームページの「過去の電力使用実績データ」をクリックします

図: 期間を選択します

図: 電力消費量データをcsv形式で入手できました!

ここでは事前に用意したcsv形式のファイルを使用しますので，データをダウンロードしていただく必要はありませんが，ご自身で，異なる地域の電力消費量データやその他のオープンデータを各自でダウンロードして分析してみて下さい．

1.3 Pythonプログラムの構成と実行

このプログラムは，本日配布した資料一式の中のpythonプログラムamedas_diagnostic_analytics.pyを，Jupyter Notebook用に再収録したものです．どちらも同じ内容のpythonプログラムとなっています．学習をする際にはこのJupyter Notebook環境を利用した方がプログラムと説明文の対応が分かりやすいと思います．ですが，pythonプログラミングに慣れてきたら上記のamedas_diagnostic_analytics.pyを直接利用して実行する方が便利です．上記のpythonプログラムを直接編集して実行する場合は，pythonの統合環境であるSpyderを利用することをお勧めします．Spyderは，Anaconda Promptで，

$ spyder

と打つことで起動します．また，左下のWindowsのスタートボタンから，「すべてのプログラム」→「Anaconda3 (64bit)」→「Spyder (Anaconda3)」でも起動できます．

図: Spyderの起動画面

上の図のように，Spyderの画面は3つに別れていて，左側の画面がプログラムエディタとなっています．文法のミスやエラーの場所などを色つきで教えてくれる優れものです．また，タブを押すことで候補となる変数やメソッドが表示されます．右上側の画面は変数エクスプローラーで，変数の型，サイズ，値を一覧で確認できます．右下側の画面はコンソールで，実行結果やエラーなどが表示されます．pythonプログラムを開いたら，上の緑色の三角形の「ファイルを実行（F5）」を押すことでプログラムが実行されます．Spyderでも，一度，pythonプログラムamedas_diagnostic_analytics.pyを動かしてみて下さい．

さて，Jupyter Notebookに戻りましょう．このプログラムは，以下のような構成となっています．

　1. ライブラリのインポート
　2. 気象庁AMeDASの気温データ（csv形式）を読む
　3. 東京電力の電力消費量データ（csv形式）を読む
　4. 2つのデータを結合する
　5. データを整理する
　6. 時系列図を作成する
　7. 散布図と回帰直線を作成する
　　(5～7の仮説と検証のプロセスを繰り返す)

早速，１つずつ手を動かしてみましょう．

1.3.1 ライブラリのインポート

以下のセルから，プログラムが始まります．上のプルダウンメニューが「Code」となっているセルがpythonプログラムに相当します．セルにカーソルを合わせて上にある「Run」ボタンを押して実行してみて下さい．以降のプログラムは上のプログラムがあることを前提にして組み立てられてますので，上から順番に実行する必要がありますので注意して下さい．

上の「View」タブの中の「Toggle Line Numbers」をクリックすると，セルに行番号が表示されるので表示させて読み進めて下さい．

まずは，必要なライブラリをインポートします．2行目で数値計算ライブラリnumpyをインポートします．import numpy as npとすることで，プログラム中ではnpと省略してインポートできます．作図のための，matplotlibとseabornも4行目と6行目でインポートします．また，8行目のdatetimeは日付（日時や時間）の処理ができる便利なライブラリで気象データ分析には必須となります．また，10行目でデータの前処理や整理に不可欠なpandasもインポートします．最後に，16行目に線形回帰直線を推定するのに必要な機械学習のライブラリsklearn (scikit-learn)の中のlinear_modelもインポートします．

Jupyter Notebookの中だけで有効なマジックコマンド%precison 3（小数点3桁まで）と%matplotlib inline（インラインで画像を表示）も忘れずにつけましょう．

# numpyをnpという別名でインポートします．
import numpy as np
# matplotlibをpltという別名でインポートします．
import matplotlib.pyplot as plt
# Seabornをインポートします．
import seaborn as sns
# datetimeは日時データを処理する際に便利なメソッドです．インポートします．
from datetime import date, datetime, timedelta
# pandasをpdという別名でインポートします．
import pandas as pd
# matplotlibで時系列図を作成するときには以下をインポートします
from pandas.plotting import register_matplotlib_converters
# これを登録しておきます．
register_matplotlib_converters()
# sklearn(scikit-learn)は機械学習関連のライブラリーです．インポートします．
from sklearn import linear_model
%precision 3
%matplotlib inline

1.3.2 気象庁AMeDASの気温データ（csv形式）を読む

さて次に，気象庁AMeDASの気温データ（amedas2017.csv）を読み込んで，pandasのデータフレームとしてメモリしましょう．表計算ソフトのExcelで例えるならば，Excelを起動してcsvファイルを開いた状態です．

このデータを見ると，上から5行目まではデータとは関係のないヘッダー情報となっているようです．6行目から最終行目（8765行目）までのデータ部が必要な情報です．また，A列は日時と時刻，B列が気温，C列とD列は品質や均質に関する参考情報となっています．データ分析に必要なデータは前の2列分（A列とB列）のみです．こういったデータの特性を考慮しながら必要な情報だけをメモリする必要があります．

インポートしたpandasには，csvファイルを読み込むメソッドpd.read_csv があります．これを使うとcsvファイルをそのままデータフレームに割り当てることができます．後のデータ整理ではpandasのデータフレームを使うと便利なので好都合です．ここでは，以下のセルのようにreadamedasという関数（3行目）の中でデータの読み込みを行います．この関数は分かりやすく複数行で書かれていますが，メソッドpd.read_csvを使えば，実質1行分のプログラミングでデータを読み込むことができます．このpd.read_csv（14行目）の引数として，15～20行目の6種類を入力する必要があります．

pd.read_csv(

filename：ファイル名（ここでは，amedas2017.csv）
encoding：文字エーコーディング（ここでは，Shift_JIS）
skiprows：上から何行分を読み飛ばすか（ここでは，5行）
names：列のラベル（ここでは，A列～D列を['date_time','temp','dummy1','dummy2']とつけます）
parse_dates：datetime型で読み込む列を構文解析する（ここでは，'date_time'の列をdatetime型で読み込む）
index_col：インデックスとする列（ここでは，'datetime'）

)

図: ここがポイント！

この様にして読み込んだデータは，pandasのデータフレーム型で割り当てられて，関数の中のdfという変数から，メインプログラムではamedasという変数に戻されて格納されます（30行目）．その直後の32行目では，

amedas=amedas.drop(['dummy1','dummy2'],axis=1)

のメソッドdropにより，ラベルが'dummy1'と'dummy2'となっている列が削除されています．axis=1とすると列の方向で削除されるという意味です．axis=0とすれば行の方向で削除されます．

それでは実行してみましょう．

# 気象庁アメダスの気温の時系列データを読み込んで，
# DataFrameに割り当てる関数
def readamedas(filename,skipline): 
    # ここがポイント！
    # pandasのread_csvというメソッドでcsvファイルを読み込みます．
    # 引数として，
    # [0]入力ファイル名
    # [1]エンコーディング
    # [2]読み飛ばす行数，
    # [3]column名
    # [4]datetime型で読み込むcolumn名
    # [5]indexとするcolumn名
    # を与える．
    df=pd.read_csv(
       filename,
       encoding='Shift_JIS',
       skiprows=skipline,
       names=['date_time','temp','dummy1','dummy2'],
       parse_dates={'datetime':['date_time']}, 
       index_col='datetime'
       )
    return df

# 気象庁AMeDAS（東京）の気温（℃）のcsvファイルの名前を指定します．
filename1='amedas2017.csv'
# csvファイルの上から5行目まではデータではないため呼び飛ばします．
skipline1=5
# 気象庁AMeDAS（東京）の気温（℃）のcsvファイルを読み込んで，
# pandasのDataFrame（tepco）に割り当てる関数を呼び出します．
amedas=readamedas(filename1,skipline1)
# DataFrame(amedas)の中のdummy1とdummy2の列を削除する．
amedas=amedas.drop(['dummy1','dummy2'],axis=1)
amedas

	temp
datetime
2017-01-01 01:00:00	5.1
2017-01-01 02:00:00	4.1
2017-01-01 03:00:00	4.0
2017-01-01 04:00:00	3.0
2017-01-01 05:00:00	3.6
...	...
2017-12-31 20:00:00	2.8
2017-12-31 21:00:00	3.3
2017-12-31 22:00:00	2.7
2017-12-31 23:00:00	2.7
2018-01-01 00:00:00	2.3

8760 rows × 1 columns

日付（'datetime'）と気温（'temp'）からなる表がプログラムの下に現れれば成功です．dummy1とdummy2も消えていることを確認して下さい．これは，amedasというデータフレームの内容を省略して表示しています．2017年1月1日01時00分00秒～2018年1月1日0時0分0秒までの1時間毎の気温データ（合計8760行）を読み込めました．

1.3.3 東京電力の電力消費量データ（csv形式）を読む

次に，同様に，東京電力でんき予報の電力消費量データ（tepco2017.csv）を読み込んで，pandasのデータフレームとしてメモリしてみましょう．これも，表計算ソフトのExcelで例えるならば，Excelを起動してcsvファイルを開いた状態となります．

このデータを見ると，上から3行目まではデータとは関係のないヘッダー情報となっているようです．4行目から最終行（8763行目）までのデータ部が必要な情報です．また，A列は日時，B列が時刻，C列が電力消費量となっています．気温データのときのように不必要な列（'dummy1'や'dummy2'）はありませんが，気温データと違って日時と時刻が2つの列に別れています．こういったデータの特性を考慮しながら情報をメモリする必要があります．

同じく，csvファイルを読み込むメソッドpd.read_csvを使って，csvファイルをそのままデータフレームに割り当てます．ここでは，以下のセルのようにreadtepcoという関数（3行目）の中でデータの読み込みを行います．この関数も分かりやすく複数行で書かれていますが，メソッドpd.read_csvを使えば，実質1行分のプログラミングでデータを読み込むことができます．このpd.read_csv（14行目）の引数として，次の6種類を設定する必要があります．

pd.read_csv(

filename：ファイル名（ここでは，tepco2017.csv）
encoding：文字エーコーディング（ここでは，Shift_JIS）
skiprows：上から何行分を読み飛ばすか（ここでは，3行）
names：列のラベル（ここでは，A列～C列を['date','time','power']とする）
parse_dates：datetime型で読み込む列を構文解析する（ここでは，'date'と'time'の2つの列をまとめてdatetime型で読み込む）
index_col：インデックスとする列（ここでは，'datetime'）

)

この様にして読み込んだデータは，pandasのデータフレーム型で割り当てられて，関数の中のdfという変数から，メインプログラムではtepcoという変数に戻されて格納されます（30行目）．

それでは実行してみましょう．

# 東京電力の電力消費量の時系列データを読み込んで，
# DataFrameに割り当てる関数
def readtepco(filename,skipline): 
    # ここがポイント！
    # pandasのread_csvというメソッドでcsvファイルを読み込みます．
    # 引数として，
    # [0]入力ファイル名
    # [1]エンコーディング
    # [2]読み飛ばす行数，
    # [3]column名
    # [4]datetime型で読み込むcolumn名
    # [5]indexとするcolumn名
    # を与える．
    df=pd.read_csv(
       filename,
       encoding='Shift_JIS',
       skiprows=skipline,
       names=['date','time','power'],
       parse_dates={'datetime':['date','time']},
       index_col='datetime'
       )
    return df

# 東京電力の電力消費量（10^6 kW）のcsvファイルの名前を指定します．
filename2='tepco2017.csv'
# csvファイルの上から3行目まではデータではないため呼び飛ばします．
skipline2=3
# 東京電力の電力消費量（10^6 kW）のcsvファイルを読み込んで，
# pandasのDataFrame（tepco）に割り当てる関数を呼び出します．
tepco=readtepco(filename2,skipline2)
tepco

	power
datetime
2017-01-01 00:00:00	2783
2017-01-01 01:00:00	2634
2017-01-01 02:00:00	2520
2017-01-01 03:00:00	2438
2017-01-01 04:00:00	2389
...	...
2017-12-31 19:00:00	3527
2017-12-31 20:00:00	3443
2017-12-31 21:00:00	3335
2017-12-31 22:00:00	3200
2017-12-31 23:00:00	3141

8760 rows × 1 columns

日付（'datetime'）と電力消費量（'power'）からなる表が現れれば成功です．これは，tepcoというデータフレームの内容を省略して表示しています．2017年1月1日0時0分0秒～2017年12月31日23時0分0秒までの1時間毎の電力消費量データ（合計8760行分）を読み込めたことが確認できました．

1.3.4 2つのデータを結合する

これまで，気温データの時系列amedasと電力消費量データの時系列tepcoを作ることができましたが，今のままでは別々の表となっていて扱いづらいです．データ分析をする際には，2つの時系列データを並べて表示するようなデータの集約が必要となってきます．例えば，表計算ソフトExcelで例えると，1つの表からもう1つの表へ時間軸が合うようにコピーした状態であると言えます．

pandasには2つのデータフレームを結合するメソッドpd.merge()があります．次のセルにあるように，

data=pd.merge(amedas, tepco, how="outer", left_index=True, right_index=True)

と実行することで，amedasのインデックスdatetimeとtepcoのインデックスdatetimeと紐付けながら，2つのデータを横方向に結合することができます．howは，outerとすれば全てのインデックスが残るように結合することになり（完全外部結合），innerとすれば両方のデータに含まれるインデックスだけが残るように結合することになります（内部結合）．left_index=Trueと設定すると，左側のデータフレームのインデックスを結合のキーとして用います．right_index=Trueと設定すると，右側のデータフレームのインデックスを結合のキーとして用います．

図: ここがポイント！

それでは，以下のセルを実行してみましょう．

# ここがポイント！
# 2つのDataFrame（amedasとtepco）を結合します．
data=pd.merge(amedas, tepco, how="outer", left_index=True, right_index=True)
data

	temp	power
datetime
2017-01-01 00:00:00	NaN	2783.0
2017-01-01 01:00:00	5.1	2634.0
2017-01-01 02:00:00	4.1	2520.0
2017-01-01 03:00:00	4.0	2438.0
2017-01-01 04:00:00	3.0	2389.0
...	...	...
2017-12-31 20:00:00	2.8	3443.0
2017-12-31 21:00:00	3.3	3335.0
2017-12-31 22:00:00	2.7	3200.0
2017-12-31 23:00:00	2.7	3141.0
2018-01-01 00:00:00	2.3	NaN

8761 rows × 2 columns

このような形にデータを集約することで，時系列図や散布図を作成したり，互いの相関関係を調べることが可能になります．ここでは，how="outer"という設定にしましたが，how="inner"としたらどうなるか確認してみてください．

1.3.5 データの整理（NaNを取り除く）

上の2つのデータの結合の際には，実は，わざと問題が生じるように完全外部結合（how="outer"）としています．どのような問題が生じたかと言えば，よく表を見てみると，データの最初の行（2017年1月1日0時0分0秒）と終わりの行（2018年1月1日0時0分0秒）にそれぞれNaNという文字が入ってしまっています．これはデータが存在しないことを意味するものです．データ分析の際にはこのような欠測値を扱わなくてはいけないことが，しばしばありますが，このようにNaNが入ったままだとこの先の分析に支障がありますので，どちらか一方にでもNaNがある場合にはその行を削除するようなデータの整形をしたいと思います．表計算ソフトExcelで例えるならば，フィルターでNaNを探してその行を削除するプロセスに該当します．

次のセルは，NaNを取り除くための関数dataprocess0（3行目）を定義し，それを呼び出すことでNaNを除外したデータフレームdata0を作成するプログラムとなっています（16行目）．関数dataprocess0の中では，6行目でNaNを除去しています．

df1=df.dropna(how='any')

pandasのdropnaメソッド で NaNを除去することができます．how='any'の場合には，'temp'と'power'のどちらか一方に欠測 がある時刻（行）は削除されることになります．一方，how=allの場合には，'temp'と'power'の両方ともに欠測がある時刻（行）が削除されることになります．

図: ここがポイント！

また，関数dataprocess0の戻り値は，データフレームdata0ともう１つcorrがあり（12行目），作成されたデータフレームdata0の中の'temp'と'power'の時系列データの間の相関係数（corr）が表示されます（18行目）．関数dataprocess0の中の8行目で相関係数を計算しています．

corr=df1[name1].corr(df1[name2])

pandasのメソッドcorr()を呼び出して，データフレームdf1（ここでは先に作成したdata）中のname1（ここでは'temp'）の列データdf1['temp']とname2（ここでは'power'）の列データdf1['power']の間の相関係数を計算しています．相関係数は，pythonのpandasの機能で簡単に計算できます（ExcelでもCORREL関数で簡単に計算できます）．

図: ここがポイント！

相関係数とは2つのデータ$x$と$y$の間の関係性の強さを表す統計学的指標です．相関係数は2つのデータの共分散を2つのデータの標準偏差の積で割ってものさしを揃えたようなものです．$x$の平均が$\overline{x}$，$y$の平均が$\overline{y}$とすると，相関係数$r$を式で表すと，

$$ r=\displaystyle{\frac{\frac{1}{n} \sum_{i=1}^{n} (x_i - \overline{x}) (y_i - \overline{y}) }{ \sqrt{\frac{1}{n} \sum_{i=1}^{n} (x_i - \overline{x})^2} \sqrt{\frac{1}{n} \sum_{i=1}^{n} (y_i - \overline{y})^2} }} $$

となります．相関係数は，データの単位に関係なく，-1から1までの範囲の値を取ります．相関係数が1に近いほど正の相関が強く（一方が上がれば，もう一方も上がる），相関係数が-1に近いほど負の相関が強くなります（一方が上がれば，もう一方は下がる）．分野や論者によって，相関係数の数値の見方は異なりますが，目安として，絶対値が0.2から0.4なら弱い相関，0.4から0.7なら中程度，0.7以上ならば強い相関という人が多いようです．ただし，相関係数が高いからといって，直接の因果関係があるということでは必ずしもありません．「気温が高いとアイスクリームが売れる」ことには因果関係があります．また，「気温が高いと水難事故が増える」ことにも因果関係があります．しかしながら，「アイスクリームが売れれば水難事故が増える」ことには因果関係があるとは言えません．因果関係がなくても高い相関が生まれることがあり，このことを 疑似相関 と言います．相関係数が高ければ必ず因果関係があるわけではないので注意してください．

図: 相関係数の値の目安

それでは以下のセルを実行してみましょう．

# ２つの時系列データからNaNを取り除いて，
# DataFrameに割り当てる関数
def dataprocess0(df,name1,name2):
    # ここがポイント！
    # 2つのデータのどちらかにNaNがあれば，その行を削除します．
    df1=df.dropna(how='any')
    # ここがポイント！
    # 処理した2つのname1とname2のデータを取り出し，
    # 2つのデータ間の相関係数を計算します．
    corr=df1[name1].corr(df1[name2])
    # 処理したデータと相関係数を返します．
    return df1,corr

# データ前処理0：　dataprocess0関数を呼び出して，
#             どちらかにNaNがある行を取り除きます．
data0,corr=dataprocess0(data,'temp','power')
# 処理された2つのデータの間の相関係数を表示します．
print("Dataprocess0... Corr=",corr)
data0

Dataprocess0... Corr= 0.005663728756090282

	temp	power
datetime
2017-01-01 01:00:00	5.1	2634.0
2017-01-01 02:00:00	4.1	2520.0
2017-01-01 03:00:00	4.0	2438.0
2017-01-01 04:00:00	3.0	2389.0
2017-01-01 05:00:00	3.6	2394.0
...	...	...
2017-12-31 19:00:00	4.0	3527.0
2017-12-31 20:00:00	2.8	3443.0
2017-12-31 21:00:00	3.3	3335.0
2017-12-31 22:00:00	2.7	3200.0
2017-12-31 23:00:00	2.7	3141.0

8759 rows × 2 columns

結果を確認してみましょう．まず，NaNを取り除く処理がなされたデータフレームdata0の中味を見てみると，確かにNaNの最初と最後のNaNを含む行が抜けていることが分かります．データ整理がうまくいきました．

また，計算された相関係数はCorr= 0.00566となりました．ちゃんと計算できたようですが，気温と電力消費量との間は無相関のようです．この結果から，気温と電力消費量との間には全く関係性がないと言ってしまってよいのでしょうか？何が間違っているのでしょうか？もう少し分析を進めましょう．

1.3.6 時系列図を作成する

次に，時系列図を作成して調べてみましょう．データを適切な形にグラフ化することで，その特徴がみえてきます．数値だけでは分かりにくい例外的な特徴も，グラフにすることによって抽出することできます．ここでは，pythonのライブラリmatplotlibを使って2つの時系列データを一枚の時系列図に収める方法を見ていきましょう．表計算ソフトExcelで例えるならば，3つの列データ（'datetime'，'temp'，'power'）を選択して，ツールバーの中の折れ線グラフボタンを押すことで時系列グラフを作成できます．

次のセルは，時系列図を作成する関数timeseries（2行目）を定義し，それを呼び出すことでデータフレームdata0の中の2つ列データから作図しています（34行目）．関数timeseriesの中では，まず，データフレームdfから横軸のデータX（4行目）と2つの縦軸（気温と電力消費量）のデータY1（6行目）とY2（8行目）を取り出しています．2つの時系列図を1枚の図の中に重ねて表示するため，メソッドplt.subplot()では図の配列は1行1列が指定されています（12行目）．2つ目の図もx軸は共有することになるため，ここでメソッドtwinx()を呼び込みます（14行目）．ax1から始まるオブジェクトは，ここでは，気温の時系列図に関わる設定をしています．一方で，ax2から始まるオブジェクトは，電力消費量の時系列図に関わる設定をしています．15行目～18行目で時系列図を描画し，19行目～26行目で図のタイトルや軸ラベルの設定をし，27行目～30行目で図の凡例を設定しています．

それでは以下のセルを実行してみましょう．

# ２つの時系列データから時系列図を作成する関数
def timeseries(df,name1,name2): 
    # dfのインデックス（時間）をXとする
    X=df.index
    # dfのname1列を指定してデータを取り出し，numpy配列で値をY1に与える．
    Y1=df.loc[:,[name1]].values
    # dfのname1列を指定してデータを取り出し，numpy配列で値をY2に与える．
    Y2=df.loc[:,[name2]].values
    # 時系列図の大きさを指定
    plt.figure(figsize=(20, 10))
    # 1つ目(name1)のグラフを1行1列の1つ目に
    ax1=plt.subplot(1,1,1)
    # 2つ目(name2)のグラフのx軸を共有する
    ax2=ax1.twinx()
    # 1つ目(name1)の時系列 
    ax1.plot(X,Y1,color='blue',label=name1)
    # 2つ目(name2)の時系列 
    ax2.plot(X,Y2,color='red',label=name2)
    # グラフのタイトル
    ax1.set_title("Timeseries:"+name1+" and "+name2)
    # x軸のラベル
    ax1.set_xlabel('Time')
    # y軸（左側の第1軸）のラベル
    ax1.set_ylabel('Temperature [degree C]')
    # y軸（右側の第2軸）のラベル
    ax2.set_ylabel('Electric Power [$10^6$kW]')
    # 1つ目(name1)の凡例（左上に置く） 
    ax1.legend(loc='upper left')
    # 2つ目(name1)の凡例（右上に置く）
    ax2.legend(loc='upper right')
    return

# 処理されたデータを用いて時系列図を作成します．
timeseries(data0,'temp','power')

青い線がax1の気温データ，赤い線がax2の電力消費量のデータの時系列図が出てきました．気温の縦軸Y1は左側に，電力消費量の縦軸Y2は第2軸として右側に現れています．さて，完成した時系列図をみて考察してみましょう．次のようなことに気づきます．

　・気温の高い夏季に電力消費量も上がる傾向にある（季節変化）
　・気温の低い冬季になると電力消費量が再び上昇する傾向にある（季節変化）
　・気温の高い昼間と気温の低い夜間の電力消費量の差が大きい（特に夏季に）（日変化）
　・気温では明瞭でないが，電力少量は1週間スケールでも周期的に変化している（週変化）

他にも何か読み取れることはあるでしょうか？ここでの「直感」がデータ分析の仮説を立てるうえで重要です．

1.3.7 散布図と回帰直線を作成する

時系列図で様々な特徴が読み取れました．次に，横軸に気温，縦軸に電力消費量を取ることで，各時刻のデータをプロットした散布図をmatplotlibで作成してみましょう．

2つのデータの間に関係性があるときには，散布図上では，一方が上がれば一方が上がる関係（正の相関）や，一方が下がれば一方が上がると関係（負の相関）が見られるはずです．2つのデータの間に関係性がないときには（相関係数が0に近いときには），このように一本の線に乗らず全体的に散らばったような分布（無相関）となるはずです．相関係数を計算するときには一緒に散布図も確認することが重要です．

図: 正の相関と負の相関の散布図

また，正の相関や負の相関のあるデータに対しては，散布図上で回帰直線（線形回帰モデル）を引きます．これは，散布図の縦軸（目的変数）を$y$，横軸（説明変数）を$x$として，最も当てはまりのよい直線$y=ax+b$の係数（傾き$a$と切片$b$）を最小二乗法により推定することを指します．当てはまりの良い回帰直線が作れたならば，気温$x$から電力消費量$y$を予測することも可能となります．表計算ソフトExcelの場合には，散布図上で右クリックして「近似曲線の追加」を選択することで線形近似線を推定することができます．また，ExcelのSLOPE関数（傾きaの推定）やINTERCEPT関数（切片bの推定）でも回帰直線を推定できます．

次のセルは，散布図を作成する関数scatter（2行目）を定義し，それを呼び出すことでデータフレームdata0の中の2つの列データから散布図を作成しています（42行目）．関数scatterの中では，まず，データフレームdfから横軸の気温データX（7行目）と縦軸の電力消費量のデータY（9行目）を取り出しています．

4行目では，scikit-learnライブラリから，線形回帰モデルを作るクラスlinear_model.LinearRegression()からclfというインスタンスを作成しています．12行目のメソッドclf.fit()で， 説明変数Xと目的変数Yから当てはまりの良い直線を最小二乗法で推定します．回帰直線の傾きaはclf.coef_に（15行目），切片bはclf.intercept_で取り出すことができます（18行目）．また，回帰直線の当てはまりの良さを表す決定係数$R^2$はclf.score()で取り出すことができます（21行目）．

図: ここがポイント！

24行目以降で，散布図と回帰直線を作図しています．27行目で$x$と$y$の散布図を作成し，30行目で回帰直線を引いています．clf.predict(X)は回帰直線Y=aX+bにより推定されたYを計算しています．32行目では作成された回帰直線を文字型equation（23行目）で置いています．33行目から38行目は，グラフのタイトル，横軸Xのラベル，縦軸Yのラベルを設定しています．

図: ここがポイント！

それでは以下のセルを実行してみましょう．

# ２つの時系列データから散布図を作成する関数
def scatter(df,name1,name2):
    # ここがポイント！
    # scikit-learnの線形回帰モデルのクラスを呼び出す
    clf = linear_model.LinearRegression()
    # 説明変数Xにはname1を割り当てる（numpy配列）
    X=df.loc[:,[name1]].values
    # 説明変数Yにはname2を割り当てる（numpy配列）
    Y=df.loc[:,[name2]].values
    # ここがポイント！
    # Y=aX+bの予測モデルを作成する
    clf.fit(X,Y)
    # ここがポイント！
    # 回帰係数a
    slope=clf.coef_
    # ここがポイント！
    # 切片b
    intercept=clf.intercept_
    # ここがポイント！
    # 決定係数R2（回帰直線の当てはまりの良さ）
    r2=clf.score(X,Y)
    # 文字列"Y=aX+b (R2=r2)"
    equation="  y = "+str('{:.1f}'.format(slope[0][0]))+" x +"+str('{:.0f}'.format(intercept[0]))+" (R2="+str('{:.2f}'.format(r2))+")"
    # 散布図の大きさを指定
    plt.figure(figsize=(8, 8))
    # 散布図のプロット
    plt.plot(X, Y, 'o')
    # ここがポイント！
    # 散布図上に回帰直線を引く
    plt.plot(X, clf.predict(X))
    # 文字列"Y=aX+b (R2=r2)"を図の左上に置く
    plt.text(np.nanmin(X), np.nanmax(Y), equation)
    # グラフのタイトル
    plt.title("Scatter diagram:"+name1+" and "+name2)
    # x軸のラベル
    plt.xlabel('Temperature [degree C]')
    # y軸のラベル
    plt.ylabel('Electric Power [$10^6$kW]')    
    return

# 処理されたデータを用いて散布図を作成します．
scatter(data0,'temp','power')

横軸Xを気温，縦軸Yを電力消費量とする散布図を作成することができました．この散布図を見ると，全体にばらついた分布をしていて，確かに相関係数は0に近くなりそうです．ただ，XとYの間に全く関係性がないわけではなさそうです．特徴的なのは18℃を中心としたV字の分布をしていることが分かります．18℃よりも高いところでは右上がりの変化をしていますが（正の相関），18℃よりも低いところでは右下がりの変化をしています（負の相関）．

線形回帰モデルはY=0.4X+3259となりました．例えば気温が20℃の時は，電力消費量は（0.4 x 20 + 3259）＝ 3,267 [万kW] と予測できます．また，気温が10℃と30℃の時は，それぞれ，3,263 [万kW] と3,271 [万kW] となります．気温が何度であっても予測結果に大きな違いは生じていません．決定係数$R^2$は0.00と低く，予測式としては使い物にならなさそうです．

これらの結果から，気温と電力消費量との間で 無相関 となった理由は夏季と冬季を一緒に分析したからだと考えられます．よって，更なる特徴を把握するためには，データを更に加工して分析する必要がありそうです．

1.3.8 データの整理と分析（仮説１：夏の期間のデータの抽出）

そこで，新たな仮説を立ててみましょう．先の時系列図や散布図の分析で，夏季と冬季では異なる傾向を示していることが分かりました．新たな仮説１として，夏季と冬季に区別して分析すればより明瞭な関係性が得られるものと期待されます．ここでは，データフレームdata0の中から，夏季（7月1日～9月30日まで）のみのデータを抽出して，再度，同様の分析を行ってみましょう．表計算ソフトExcelで同様の処理をする場合には，不要な期間のデータ（行）を見つけて手動で表から削除するか，フィルターをかけることで，データを抽出することができます．

次のセルは，夏季の期間のデータを抽出する関数dataprocess1（3行目）を定義し，それを呼び出すことでデータフレームdata0から夏季のみのデータからなるデータフレームdata1を作成するプログラムとなっています（15行目）．関数dataprocess1の中では，6行目で夏季のみのデータを取り出しています．

df1=df[(df.index.date >= date(2017,7,1) ) & (df.index.date <= date(2017,9,30))]

datetimeのdateというライブラリを用いることで日付処理を行うことができます．この式を見ると，データフレームdfの中のインデックスの日付（ここではdatetimeの中のdate）が2017年7月1日よりも後であり，かつ（&），日付が2017年9月30日より前であれば，夏季のデータであるとしてデータを残すと判断しています．この条件に当てはまらなかったデータは除外されます．

図: ここがポイント！

また，関数dataprocess1の戻り値は，データフレームdata1ともう１つcorrがあり（15行目），作成されたデータフレームdata1の中の'temp'と'power'の時系列データの間の相関係数（corr）が表示されます（17行目）．関数dataprocess1の中の9行目で相関係数を計算しています．

corr=df1[name1].corr(df1[name2])

pandasのcorr()メソッドを呼び出してで，データフレームdf1（ここでは先に作成したdata）中のname1（ここでは'temp'）の列データdf1['temp']とname2（ここでは'power'）の列データdf1['power']の間の相関係数を計算しています．

それでは以下のセルを実行してみましょう．

# ２つの時系列データから夏の期間（7/1から9/30まで）を取り出して，
# DataFrameに割り当てる関数
def dataprocess1(df,name1,name2):
    # ここがポイント！
    # 夏の期間（7/1から9/30まで）のデータを取り出します．
    df1=df[(df.index.date >= date(2017,7,1) ) & (df.index.date <= date(2017,9,30))]
    # 処理した2つのname1とname2のデータを取り出し，
    # 2つのデータ間の相関係数を計算します．
    corr=df1[name1].corr(df1[name2])
    # 処理したデータと相関係数を返します．．    
    return df1,corr

# データ前処理１：　dataprocess1関数を呼び出して，
#             夏の期間（7/1から9/30まで）だけを取り出します．
data1,corr=dataprocess1(data0,'temp','power')
# 処理された2つのデータの間の相関係数を表示します．
print("Dataprocess1... Corr=",corr)
data1

Dataprocess1... Corr= 0.7510366253230192

	temp	power
datetime
2017-07-01 00:00:00	23.9	2730.0
2017-07-01 01:00:00	23.0	2570.0
2017-07-01 02:00:00	23.3	2491.0
2017-07-01 03:00:00	22.7	2488.0
2017-07-01 04:00:00	23.1	2478.0
...	...	...
2017-09-30 19:00:00	19.2	3057.0
2017-09-30 20:00:00	18.7	2929.0
2017-09-30 21:00:00	18.7	2765.0
2017-09-30 22:00:00	18.2	2626.0
2017-09-30 23:00:00	18.6	2464.0

2208 rows × 2 columns

新たに作成されたデータフレームdata1の出力結果を見ると，確かに夏季の7月1日から9月30日までのデータを抽出することができています．データの行数も2208行と減っています．また，相関係数は0.7510とかなり高くなりました．夏季のデータのみを取り出すことで「やや強い相関」となりました．ここで立てた仮説１「夏季と冬季に区別して分析すればより明瞭な関係性が得られる」は正しかったということになります．

更に，時系列図と散布図を作成してみましょう．先ほど，時系列図を作成する関数timeseriesと散布図を作成する関数scatterを作っていますので，それをそのまま流用すればOKです．ここでは2つの関数を呼び出すだけとなります．先との違いは，引数として入力するデータフレームがdata0からdata1に変わったことだけです．

それでは実行してみましょう．

# 処理されたデータを用いて時系列図を作成します．
timeseries(data1,'temp','power')
# 処理されたデータを用いて散布図を作成します．
scatter(data1,'temp','power')

夏季のデータのみを抽出したデータフレームdata1から時系列図と散布図が作成できたと思います．

データフレームdata1の時系列を見ると，気温の高い8月に電力消費量が大きい傾向がよりはっきりと見えてきました．また，気温も電力消費量も明瞭な日変化があることがはっきりと見て取れます．また，データフレームdata1の散布図を見ると，先のV字のような分布ではなく，右上がりの分布となっていることが分かります．

推定された線形回帰モデルはY=156.6X-610となりました．決定係数$R^2$は0.56となり，先の推定式と比べるとかなり精度が高くなったと言えます．しかしながら，まだ，全体的に回帰直線の周りで広くばらついた分布をしているため，より精度の高いモデルを作るためには更なる仮説を立ててデータを整理する必要があると考えられます．

1.3.9 データの整理と分析（仮説２：12～15時のデータの抽出）

さらに仮説を立ててみましょう．先の時系列図や散布図の分析で，気温も電力消費量も明瞭な日変化を示すことが分かりました．気温の高い日中には冷房需要が高くなると想定されます．そこで，新たな仮説２として，気温が最も高くなる12～15時のみのデータのみを抽出して分析すればより明瞭な関係性が得られると考えられます．ここでは，データフレームdata1の中から，12時，13時，14時，15時のみのデータを抽出して，再度，同様の分析を行ってみましょう．表計算ソフトExcelで同様の処理をする場合には，フィルターで該当する時刻のデータのみを抽出して，それ以外を削除することになります．

次のセルは，夏季の期間のデータを抽出する関数dataprocess2（3行目）を定義し，それを呼び出すことでデータフレームdata1から12～15時のみのデータからなるデータフレームdata2を作成するプログラムとなっています（15行目）．関数dataprocess1の中では，6行目で12～15時のみのデータを取り出しています．

df1=df[(df.index.hour >= 12)&((df.index.hour <= 15))]

データフレームdfの中のインデックス（ここではdatetime）から時間のみを取り出し(.hour)，それが12～15の数字であれば，12～15時のデータであると判断してデータを残します．それ以外のデータは除外されます．

図: ここがポイント！

また，関数dataprocess2の戻り値は，データフレームdata2ともう１つcorrがあり（15行目），作成されたデータフレームdata2の中の'temp'と'power'の時系列データの間の相関係数（corr）が表示されます（17行目）．関数dataprocess2の中の9行目で相関係数を計算しています．

corr=df1[name1].corr(df1[name2])

pandasのcorr()メソッドを呼び出して，データフレームdf1（ここでは先に作成したdata）中のname1（ここでは'temp'）の列データdf1['temp']とname2（ここでは'power'）の列データdf1['power']の間の相関係数を計算しています．

それでは以下のセルを実行してみましょう．

# ２つの時系列データから12～15時だけを取り出して，
# DataFrameに割り当てる関数
def dataprocess2(df,name1,name2):
    # ここがポイント！
    # 14時のデータだけを取り出します．
    df1=df[(df.index.hour >= 12)&((df.index.hour <= 15))]
    # 処理した2つのname1とname2のデータを取り出し，
    # 2つのデータ間の相関係数を計算します．
    corr=df1[name1].corr(df1[name2])
    # 処理したデータと相関係数を返します．
    return df1,corr

# データ前処理2：　dataprocess2関数を呼び出して，
#             12～15時だけを取り出します．
data2,corr=dataprocess2(data1,'temp','power')
# 処理された2つのデータの間の相関係数を表示します．
print("Dataprocess2... Corr=",corr)
data2

Dataprocess2... Corr= 0.7821019235290686

	temp	power
datetime
2017-07-01 12:00:00	24.1	3382.0
2017-07-01 13:00:00	23.0	3340.0
2017-07-01 14:00:00	22.9	3297.0
2017-07-01 15:00:00	22.7	3220.0
2017-07-02 12:00:00	28.7	3363.0
...	...	...
2017-09-29 15:00:00	23.5	3458.0
2017-09-30 12:00:00	22.2	2981.0
2017-09-30 13:00:00	23.6	2995.0
2017-09-30 14:00:00	23.8	2979.0
2017-09-30 15:00:00	23.1	2948.0

368 rows × 2 columns

新たに作成されたデータフレームdata2の出力結果を見ると，データの行数も368行とさらに減り，確かに12～15時のみのデータが抽出されていることが分かります．相関係数は0.7821と相関係数が上昇しました．ここで立てた仮説２「気温が最も高くなる12～15時のみのデータを抽出して分析すればより明瞭な関係性が得られる」も正しかったといえるでしょう．

更に，時系列図と散布図を作成してみましょう．先ほど，時系列図を作成する関数timeseriesと散布図を作成する関数scatterを作っていますので，それをそのまま流用すればOKです．ここではこれらの関数を呼び出すだけとなります．先との違いは，引数として入力するデータフレームがdata2に変わったことだけです．

それでは実行してみましょう．

# 処理されたデータを用いて時系列図を作成します．
timeseries(data2,'temp','power')
# 処理されたデータを用いて散布図を作成します．
scatter(data2,'temp','power')

12時～15時のデータのみ抽出されたデータフレームdata2から時系列図と散布図が作成できたと思います．

データフレームdata2の時系列を見ると，気温と電力消費量は大まかな変化については一致していますが，一部で細かな変化を追えていないように見えます．また，データフレームdata2の散布図を見ると，先のdata1と同様に右上がりの分布で正の相関があることが分かりますが，data1と同様にばらついた分布のまま解消されていないようです．

また，推定された線形回帰モデルはY=132.4X+296となりました．この推定式の決定係数$R^2$は0.61となり，先の推定式の決定係数（0.56）に比べて改善が見られます．全体的に回帰直線の周りで広くばらついた分布をしているため，より精度の高いモデルを作るためには，別の仮説を立ててデータを整理する必要があると言えるでしょう．

1.3.10 データの整理と分析（仮説３：休日のデータの削除）

さらに別の仮説を立ててみましょう．先の時系列図や散布図の分析で，電力消費量には1週間毎の変化があることが分かりました．時系列データを注意深く見てみると，電力消費量は5日間は高く，2日間は低いといった変化を繰り返しているように見えます．これは，平日と休日（土曜日・日曜日）の違いを反映した変化であると推測されます．そこで，新たな仮説３として，土曜日と日曜日のデータを取り除いて平日のみを分析すればより明瞭な関係性が得られると考えました．ここでは，データフレームdata2の中から，平日のみのデータを抽出して，再度，同様の分析を行ってみましょう．表計算ソフトExcelで同様の処理をする場合には，WEEKDAY関数で曜日を判定して，フィルターをかけて土曜日と日曜日のデータを除外するプロセスとなります．

次のセルは，夏季の期間のデータを抽出する関数dataprocess3（3行目）を定義し，それを呼び出すことでデータフレームdata2から平日のみのデータからなるデータフレームdata3を作成するプログラムとなっています（15行目）．関数dataprocess3の中では，6行目で平日のみのデータを取り出しています．

df1=df[(df.index.weekday != 5) & (df.index.weekday != 6) ]

データフレームdfの中のインデックス（ここではdatetime）から曜日情報を取り出し(.weekday)，「5」であれば土曜日と，「6」であれば日曜日と判定します．ちなみに，0は月曜日です．最終的に，土曜日でもなく，かつ，日曜日でもないデータが抽出されます．

図: ここがポイント！

また，関数dataprocess3の戻り値は，データフレームdata3ともう１つcorrがあり（15行目），作成されたデータフレームdata3の中の'temp'と'power'の時系列データの間の相関係数（corr）が表示されます（17行目）．関数dataprocess3の中の9行目で相関係数を計算しています．

corr=df1[name1].corr(df1[name2])

pandasのcorr()メソッドを呼び出して，データフレームdf1（ここでは先に作成したdata）中のname1（ここでは'temp'）の列データdf1['temp']とname2（ここでは'power'）の列データdf1['power']の間の相関係数を計算しています．

それでは以下のセルを実行してみましょう．

# ２つの時系列データから土曜日・日曜日のデータを取り除いて，
# DataFrameに割り当てる関数
def dataprocess3(df,name1,name2):
    # ここがポイント！
    # 土曜日（weekday==5）と日曜日（weekday==6）を取り除く．
    df1=df[(df.index.weekday != 5) & (df.index.weekday != 6) ]
    # 処理した2つのname1とname2のデータを取り出し，
    # 2つのデータ間の相関係数を計算します．
    corr=df1[name1].corr(df1[name2])
    # 処理したデータと相関係数を返します．
    return df1,corr

# データ前処理3：　dataprocess3関数を呼び出して，
#             土曜日・日曜日を取り除きます．
data3,corr=dataprocess3(data2,'temp','power')
# 処理された2つのデータの間の相関係数を表示します．
print("Dataprocess3... Corr=",corr)
data3

Dataprocess3... Corr= 0.8804545257470545

	temp	power
datetime
2017-07-03 12:00:00	31.9	4375.0
2017-07-03 13:00:00	31.3	4578.0
2017-07-03 14:00:00	31.0	4664.0
2017-07-03 15:00:00	30.4	4611.0
2017-07-04 12:00:00	28.3	4150.0
...	...	...
2017-09-28 15:00:00	23.7	3579.0
2017-09-29 12:00:00	22.7	3342.0
2017-09-29 13:00:00	23.2	3463.0
2017-09-29 14:00:00	23.5	3485.0
2017-09-29 15:00:00	23.5	3458.0

260 rows × 2 columns

新たに作成されたデータフレームdata3の出力結果を見ると，土曜日と日曜日のデータが確かに削除されていることが分かります．データの行数も260行とさらに減っています．相関係数は0.880となりより強い相関となりました．ここで立てた仮説３「土曜日と日曜日のデータを取り除いて分析すればより明瞭な関係性が得られる」は正しかった言えます．仮説の通り，土曜日と日曜日は休日であるため都市部での電力需要量が減り，結果として平日とは異なる傾向を示すと言えます．data1やdata2の散布図のばらつきの主因となっていたと考えられます．より当てはまりのよい線形回帰モデルY=aX+bを作るためには，土日（あるいは逆に平日）を取り除いたサンプルを用いる必要があるわけです．

更に，時系列図と散布図を作成してみましょう．既に時系列図を作成する関数timeseriesと散布図を作成する関数scatterを作っていますので，それをそのまま流用すればOKです．ここでは2つの関数を呼び出すだけとなります．先との違いは，引数として入力するデータフレームがdata3に変わったことだけです．

それでは実行してみましょう．

# 処理されたデータを用いて時系列図を作成します．
timeseries(data3,'temp','power')
# 処理されたデータを用いて散布図を作成します．
scatter(data3,'temp','power')

平日のデータのみを抽出したデータフレームdata3から時系列図と散布図が作成できたと思います．

データフレームdata3の時系列を見ると，気温の変化と電力消費量の変化がほとんど重なっていることが分かります．ただ，8月中旬については気温の変化との対応が悪い期間が存在しているようです．また，データフレームdata3の散布図を見ると，先のdata1やdata2に比べて，ばらつきの少ない右上がりの分布となっていることが分かります．「夏季」「12～15時」「平日」の気温と電力消費量の間には明瞭な関係性があると言えるでしょう．

また，推定された線形回帰モデルはY=136.1X+372となりました．推定式自体にはあまり変化が見られませんが，決定係数$R^2$は0.78となり，より精度の高い線形回帰モデルを作ることができました．

1.3.11 データの整理と分析（仮説４：お盆・祝日のデータの削除）

またまた仮説をたててみましょう．これまで見てきたように，気温と電力消費量との関係性には平日と休日とで明瞭な違いがあることが分かりました．それならば，夏季休暇（お盆休み）や祝日も取り除いて評価した方が良さそうです．先の時系列図からも，8月中頃のお盆頃において，気温と電力消費量の変化に大きな隔たりがあることが確認されていますので，そのようなデータを更に除外することによって，より当てはまりの良い線形回帰モデルが作れそうです．そこで，新たな仮説４として，お盆（8/7～8/20）と祝日（7/17，8/11, 9/18，9/23）のデータを取り除いて分析すればより明瞭な関係性が得られるのではないかと考えました．ここでは，データフレームdata3の中から，お盆と祝日のデータを削除して，再度，同様の分析を行ってみましょう．表計算ソフトExcelで同様の処理をする場合には，お盆や祝日の真偽を判定にする関数はありませんので，手動で行を削除することになります．

次のセルは，お盆と祝日以外の期間のデータを抽出する関数dataprocess4（3行目）を定義し，それを呼び出すことでデータフレームdata3からお盆と祝日以外のデータからなるデータフレームdata4を作成するプログラムとなっています（27行目）．関数dataprocess4の中では，4～18行目で平日のみのデータを取り出しています．

df1=df[(df.index.date < date(2017,8,7) ) | (df.index.date > date(2017,8,20))]
df1=df1[(df1.index.date != date(2017,7,17))]
df1=df1[(df1.index.date != date(2017,8,11))]
df1=df1[(df1.index.date != date(2017,9,18))]
df1=df1[(df1.index.date != date(2017,9,23))]

図: ここがポイント！

データフレームdfの中のインデックス（ここではdatetime）から，6行目では2017年8月7日よりも前と，2017年8月20日よりも後のデータのみを抽出しています（つまり，お盆の2週間は削除されます）．また，9行目では7月17日（海の日）以外のデータを抽出しています．また，12行目では8月11日（山の日）以外のデータを抽出しています（6行目のお盆で既に削除されています）．また，15行目では9月18日（敬老の日）以外のデータを抽出しています．最後に，18行目では9月23日（秋分の日）以外のデータを抽出しています．

また，関数dataprocess4の戻り値は，データフレームdata4ともう１つcorrがあり（23行目），作成されたデータフレームdata4の中の'temp'と'power'の時系列データの間の相関係数（corr）が表示されます（29行目）．関数dataprocess4の中の21行目で相関係数を計算しています．

corr=df1[name1].corr(df1[name2])

pandasのcorr()メソッドを呼び出して，データフレームdf1（ここでは先に作成したdata）中のname1（ここでは'temp'）の列データdf1['temp']とname2（ここでは'power'）の列データdf1['power']の間の相関係数を計算しています．

それでは以下のセルを実行してみましょう．

# ２つの時系列データからお盆や祝日のデータを取り除いて，
# DataFrameに割り当てる関数
def dataprocess4(df,name1,name2):
    # ここがポイント！
    # お盆休暇（8/7から8/20まで）以外残す
    df1=df[(df.index.date < date(2017,8,7) ) | (df.index.date > date(2017,8,20))] 
    # ここがポイント！
    # 海の日以外を残す
    df1=df1[(df1.index.date != date(2017,7,17))] 
    # ここがポイント！
    # 山の日以外を残す
    df1=df1[(df1.index.date != date(2017,8,11))] 
    # ここがポイント！
    # 敬老の日以外を残す
    df1=df1[(df1.index.date != date(2017,9,18))] 
    # ここがポイント！
    # 秋分の日以外を残す
    df1=df1[(df1.index.date != date(2017,9,23))] 
    # 処理した2つのname1とname2のデータを取り出し，
    # 2つのデータ間の相関係数を計算します．
    corr=df1[name1].corr(df1[name2])
    # 処理したデータと相関係数を返します．
    return df1,corr

# データ前処理4：　dataprocess4関数を呼び出して，
#             お盆・祝日を取り除きます．
data4,corr=dataprocess4(data3,'temp','power')
# 処理された2つのデータの間の相関係数を表示します．
print("Dataprocess4... Corr=",corr)
data4

Dataprocess4... Corr= 0.9133959244892919

	temp	power
datetime
2017-07-03 12:00:00	31.9	4375.0
2017-07-03 13:00:00	31.3	4578.0
2017-07-03 14:00:00	31.0	4664.0
2017-07-03 15:00:00	30.4	4611.0
2017-07-04 12:00:00	28.3	4150.0
...	...	...
2017-09-28 15:00:00	23.7	3579.0
2017-09-29 12:00:00	22.7	3342.0
2017-09-29 13:00:00	23.2	3463.0
2017-09-29 14:00:00	23.5	3485.0
2017-09-29 15:00:00	23.5	3458.0

212 rows × 2 columns

新たに作成されたデータフレームdata4の出力結果を見ると，表示に制限があるため祝日やお盆のデータが削除されているかは直接確認できませんが．行数212行に減っていますのでdataprocess4は上手く機能していると考えられます．また，相関係数は更に上昇し，0.913となり最高値を記録しています．ここで立てた仮説４「お盆（8/7～8/20）と祝日（7/17，8/11, 9/18，9/23）のデータを取り除いて分析すればより明瞭な関係性が得られる」は正しかったと言えるでしょう．data3の時系列に見られた8月中旬の電力消費量の落ち込みはお盆における需要低下であると言え，それを適切に除去することで気温との関係性はより明瞭に現れると考えられます．

更に，時系列図と散布図を作成してみましょう．既に時系列図を作成する関数timeseriesと散布図を作成する関数scatterを作っていますので，それをそのまま流用すればOKです．ここでは2つの関数を呼び出すだけとなります．先との違いは，入力するデータフレームがdata4に変わったことだけです．

それでは実行してみましょう．

# 処理されたデータを用いて時系列図を作成します．
timeseries(data4,'temp','power')
# 処理されたデータを用いて散布図を作成します．
scatter(data4,'temp','power')

お盆・祝日のデータを除去したデータフレームdata4から時系列図と散布図が作成できたと思います．

まず，データフレームdata4の時系列を見ると，気温の変化と電力消費量の変化がほぼ重なっていることが分かります．data3で見られた，夏季休暇（お盆休み）の急激に電力消費量が低下している期間を除去することで，より一層，気温の変化との対応が良くなったように見えます．また，データフレームdata4の散布図を見ると，data3に比べて，より線形回帰モデルの近くにプロットが集まっていることが分かります．

また，推定された線形回帰モデルはY=129.4X+624となりました．推定式自体にはあまり変化が見られませんが，決定係数$R^2$は更に上昇して0.83となり，実用に耐える精度を有する線形回帰モデルを作ることができたと言えるでしょう．ただし，この回帰式が使えるのは，「夏季」「12～15時」「平日（祝日を除く）」ときだけである点に注意してください．

このような分析の結果から，第一次近似的に気温と電力消費量との間には関係があり，正確に電力消費量を予測するためには，さらに，季節，時間，曜日等を考慮する必要があることがわかりました．次に，機械学習により電力消費量の予測モデルを作る予測的分析を行ってみましょう．

1.4 まとめ

●　データ分析とは，仮説と検証の繰り返しにより，データの中に潜む本質にせまること．
●　まずは，仮説に従って，データを適切に要約・整理します．
●　データを時系列図で示すことにより，特徴や傾向を読み取ります．
●　散布図により，2つのデータの間の関係性や異常値を読み取ります．
●　相関係数により，2つのデータの間の相関関係の強さを数値で表すことができます．
●　以上のプロセスを経て，仮説を検証することができます．

データ分析により，誰も気づかなかった真実を発見できると楽しいですね．そこから「新たなビジネス」が生まれます．

1.5 著作権について

Copyright (c) 気象データ×IT勉強会，吉野　純，岐阜大学工学部附属応用気象研究センター 2022 All rights reserved.
＜利用条件＞
本書は、本書に記載した要件・技術・方式に関する内容が変更されないこと、および出典を明示いただくことを前提に、無償でその全部または一部を複製、翻案、翻訳、転記、引用、公衆送信等して利用できます。なお、全体または一部を複製、翻案、翻訳された場合は、本書にある著作権表示および利用条件を明示してください。
＜免責事項＞
本書の著作権者は、本書の記載内容に関して、その正確性、商品性、利用目的への適合性等に関して保証するものではなく、特許権、著作権、その他の権利を侵害していないことを保証するものでもありません。本書の利用により生じた損害について、本書の著作権者は、法律上のいかなる責任も負いません。