雪国の主婦が積雪を予想してみました

2021年12月6日 22:43

はじめに

初めまして。note初投稿・エンジニア転職を目指し、プログラミングを勉強しはじめた主婦です。出産をきっかけに自分の人生や働き方について深く考えるようになりました。一度きりの人生、自分のやりたいことをやってみたいという気持ちが強くなりプログラミングを学ぶことから始めました。まずは行動ということでこのブログが私の第一歩となります。初心者ですがどうぞ宜しくお願い致します。

自己紹介

30歳主婦子ども一人

今までパソコンを所持していなかったほどのIT初心者です

アイデミーデータ分析講座にてプログラミングを勉強中

目的

私は雪国新潟の出身で毎年必ず雪が降ります。しかし積雪量というのは毎年異なり、腰の高さまで降る年もあれば、除雪が必要のない位少ない年もあります。最近の積雪量は異常気象により予測が難しいものとなっておりますが、今年の積雪も気になるところではありますので積雪量が予測できるようにしていきたいと思います。

データの中身

【期間】2000年～2021年の11月～３月までの降雪量合計

【場所】新潟

実行環境

Google Colaboratory

最初のデータの状態

データは気象庁のホームページからダウンロードしました。

1.データの読み込み

ダウンロードしたデータが.csvの為、pandasを用いて読み込みます。

import pandas as pd
df = pd.raed_csv("./C:\Users\waiwa\OneDrive\data (2).csv")

データファイルのパスの取得方法につきましてはこちらを参考にしました。

https://ossan-tech.work/windows-filepath-copy/

パスの取得方法すら分からなかったIT初心者に優しく教えてくださったアイデミーの先生方には感謝です。

２．データの整理

予測したい値の時系列データの１列だけ(積雪量1列)のDataFrameを用意し、indexが集計時期となるようなデータフレームを作成します。またそのデータを時系列分析に適用して期間を予測していきます。

まず先頭の２行をカラムとします。

df = pd.read_csv("data.csv",encoding="SHIFT-JIS",header=2)
df = df[['年月日.1','降雪量合計(cm)']]
df = df.dropna()

カラム名を分かりやすいように年月日→'Date'、降雪量合計→'Snowfall'に変更します。

df = df.rename(columns={'年月日.1': 'Date', '降雪量合計(cm)': 'Snowfall'})

Dateをdatetimeに変換し、indexとします。

df['Date'] = pd.to_datetime(df['Date'])
df = df.set_index('Date')

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

import warnings
import itertools
import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
# %matplotlib inline
warnings.simplefilter('ignore')

orderの最適化関数

ここでパラメータを選択するために、時系列データ:DATA, パラメータs(周期):sを入力すると、最も良いパラメーターとそのBICを出力するselectparameteという関数を定義していきます。

def selectparameter(DATA, s):
   p = d = q = range(0, 2)
   pdq = list(itertools.product(p, d, q))
   seasonal_pdq = [(x[0], x[1], x[2], s)
                   for x in list(itertools.product(p, d, q))]
   parameters = []
   BICs = np.array([])
   for param in pdq:
       for param_seasonal in seasonal_pdq:
           try:
               mod = sm.tsa.statespace.SARIMAX(DATA,
                                               order=param,
                                               seasonal_order=param_seasonal)
               results = mod.fit()
               parameters.append([param, param_seasonal, results.bic])
               BICs = np.append(BICs, results.bic)
           except:
               continue
   return parameters[np.argmin(BICs)]

モデルの構築

いよいよモデルの構築です。

今回積雪量予測は時系列データなのでSARIMAモデルを用いて時系列解析をしていきます。

SARIMAモデルは、ARIMAモデルをさらに季節周期を持つ時系列データにも拡張できるようにしたモデルの為今回の予測に適したモデルとなっています。

best_params = selectparameter(df, 12)
SARIMA_model = sm.tsa.statespace.SARIMAX(df.Snowfall, order=best_params[0],
                                            seasonal_order=best_params[1],
                                            enforce_stationarity=False, enforce_invertibility=False).fit()

周期は月ごとのデータであることも考慮して1年は12カ月なので、s=12となります。

orderはselectparameter関数の0インデックス, seasonal_orderは1インデックスに格納していきます。

予測

pred = SARIMA_model.predict('2016-01-20', '2021-11-30	')

predに予測期間での予測値を代入していきます。

グラフの可視化

plt.plot(df)
plt.plot(pred, "r")
plt.show()

予測値は赤色でプロットします。

結果

こちらが実行結果となります。

このグラフを見ると、赤の部分の予測が正解(青の部分)に比べてある程度正しく予測することができていることが分かります。
ただ最大積雪量が正しく予測できていなかったので、もう少し正確に予測するためには例年の最初の月の最高気温や最低気温などもデータに入れると精度を上げることができるのかなと思いました。

まとめ

プログラミング初心者で右も左も分からず何度も挫折しそうになりましたが、先生方に教えていただきながらなんとか時系列分析を完成させることができました。プログラミングは本当に奥が深くたくさんの知識を必要とするものだと思いますが、その入り口に触れることができ勉強することの楽しさを思い出すことができました。今後も少しずつですがステップアップしていきたいと思います。