BIKE SHARING DEMAND データ分析

Aidemy様の専門実践教育訓練講座を受講、これの最終成果物として記事作成。
Kaggleの上位2.9%という結果。
講座修了の合否に使用される記事は順序立て・考察等で冗長な為、詳細な考察を少々省いたものを作成。
この記事を見ると同じ分析ができて、且つなぜその分析を行うかという理由も理解できるように情報量は残した（つまり、依然膨大な記事）。
目次から気になる点へ飛んでいただけると良いと思う。

講座修了用の「冗長なバージョン」はこちら：

https://qiita.com/SumMae/items/f0f221735f94c06c609e

# <font color="Red">目次</font>
- [**はじめに～記事概略まとめ**](#はじめに記事概略まとめ)
- [**分析環境**](#分析環境)
- [**テーマについて**](#テーマについて)
    - [Bike Sharing Demand 概略](#bike-sharing-demand-概略) 
    - [分析用データ](#分析用データ)
    - [評価指標](#評価指標)
    - [目標値設定](#目標値設定)
        - [目標RMSLE値について掘り下げ](#目標rmsle値について掘り下げ)
- [**データ分析1：データ確認**](#データ分析1データ確認)
    - [各変数](#各変数)
    - [情報要約](#情報要約)
    - [年、月、日、時刻、曜日への分割](#年月日時刻曜日への分割)
    - [年月日等データ確認](#年月日等データ確認)
    - [天候データ確認](#天候データ確認)
    - [casual、registeredについて確認](#casualregisteredについて確認)
    - [量的変数の可視化](#量的変数の可視化)
- [**最初の機械学習：LinearRegression**](#最初の機械学習linearregression)
    - [結果](#結果)
    - [スコア確認1](#スコア確認1)
    - [課題](#課題)
- [**説明変数の編集**](#説明変数の編集) 
    - [atempの算出](#atempの算出)
- [**目的変数の扱い**](#目的変数の扱い) 
- [**LightGBMによる予測**](#lightgbmによる予測)
    - [LightGBM最適化での学習](#lightgbm最適化での学習)
    - [スコア確認2](#スコア確認2)
    - [考察](#考察)
- [**スコアを上げるための施策**](#スコアを上げるための施策)
    - [説明変数の作成](#説明変数の作成)
        - [通勤利用者、レジャー利用者](#通勤利用者レジャー利用者)
        - [1日雨が続く](#1日雨が続く)
        - [利用者数が多い条件](#利用者数が多い条件)
        - [クリスマス](#クリスマス)
        - [monthを削除](#monthを削除)
    - [作成した説明変数まとめ](#作成した説明変数まとめ)
        - [casual、registered での使用](#casualregistered-での使用)
    - [スコア確認3](#スコア確認3)
    - [目的変数の分布について](#目的変数の分布について)
    - [スコア確認4](#スコア確認4)
- [**まとめ**](#まとめ)

# <font color="Red">はじめに～記事概略まとめ</font>
- Kaggleコンペ、[Bike Sharing Demand](https://www.kaggle.com/competitions/bike-sharing-demand/overview)のデータ分析を実施
- スコアは上位3%に収まる
- Aidemy様のデータ分析講座を受けた最終成果物として記事作成

具体的な実施項目。
- 内容：ワシントンD.C.に於けるレンタル自転車の利用量（1時間当たり）
- 目的：利用者数を予測。天候、日時、休祝日情報より
- 分析：
    - <details><summary>各データについて確認(クリックで展開)</summary>

        - 事実上の外れ値あり　→　実際の天気予報を入手し補完
        - 理論式から「体感温度」を再計算し、外れ値を補完
            - →データに施されているノーマライズを解読し再計算


    </details>

    - <details><summary>機械学習モデル・・・LinearRegressionとLightGBM(クリックで展開)</summary>

        - 線形モデルで学習、不適と判断
        - LightGBM（決定木ベース）で学習、良好
            - oputunaによるハイパーパラメータの最適化

    </details>

    - <details><summary>特徴量の再考(クリックで展開)</summary>

        - 「ユーザ属性」情報がない　→　模した特徴量を作成
        - 時系列に並べると急峻に降下する点あり　→　この特徴量を作成
        - 逆に、ピークとなる箇所の特徴を作成
        - クリスマス以降を休日扱いにする
        - 目的変数の分布確認、log(x+1)で変換する

    </details>

    - LightGBMで再学習、評価結果は上位2.9%
    - まとめ

以下はこれらの詳細となる。

# <font color = "Red">分析環境</font>
|使用環境||
|:-:|:-:|
|Kaggle Notebook |2024年5月現在の仕様|
|言語|Python3|
|Google Chrome   |ウェブブラウザとして |
|Windows 11|使用OS|

# <font color="Red">テーマについて</font>
## <font color="HotPink">Bike Sharing Demand 概略</font>
Kaggleのコンペティション、[Bike Sharing Demand](https://www.kaggle.com/competitions/bike-sharing-demand/overview)。

ワシントンD.C.におけるレンタル自転車の利用者数の予測問題。
- データ：1時間ごとの利用者数と気象データ
- データ期間：2011年1月1日～2012年12月31日
- 目的変数：1時間当たりの利用者数
- 説明変数：気候と曜日等（詳細は後述）
- 訓練用データ：毎月1日～19日までのデータ
- テスト用データ：毎月20日以降のデータ

適当な図だが、予測対象を視覚的に説明。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/4315f346-0f60-ecde-9da6-ea45afc07045.png)


## <font color="HotPink">分析用データ</font>
Kaggle提供のデータ。
提供データは[こちら](https://www.kaggle.com/competitions/bike-sharing-demand/data)（Kaggleのデータタブのリンク）。データの出典は[こちら](https://archive.ics.uci.edu/dataset/275/bike+sharing+dataset)。
見かけ上欠損値がなく、扱いやすいデータ。
Bike Sharingサービスの仕様上自動的にデータ取得される為か。
但し、ところどころ記録のない時刻（行ごと無し）や、事実上の欠損値となるデータがある。

## <font color="HotPink">評価指標</font>
評価は**RMSLE値** (**R**oot **M**ean **S**quare **L**og **E**rror の略)で、レンタル数の予測値と実際のレンタル実績数との差を評価。
``` math
RMSLE \ = \ \sqrt{\frac{1}{n}\sum^n_{i=1}\Bigl(log (p_i+1) - log(a_i +1)\Bigr)^2}
```
- n　:　データ数。1時間を1データとして扱う
- $p_i$：　予測値
- $a_i$：　実際の値
- $log$：　自然対数

## <font color="Hotpink">目標値設定</font>
<font color="Red">**上位5%以上**</font>を目標に定める（Kaggleの「Leaderboard」上の順位）。
RMSLE値約0.39が上位5%。
「RMSLE値がxxxの場合、モデルの精度がこの程度だ」といった数値からの精度目標設定は難解なので、
目標はKaggleの順位にした（目標値は偏差値換算で約66）。

### <font color = "Orange">目標RMSLE値について掘り下げ</font>
少々横道に逸れる。
「RMS」はよく使用するが「RMSLE」は今回初めて触れるので少し詳しく、別記事にまとめる

https://qiita.com/SumMae/items/bb22533f9cf4b951c75b

結局、図的に理解すると、
 - RMSE  ：　切片で誤差の幅を見る
 - RMSLE ：　傾きで誤差の角度を見る

|RMSE|RMSLE|
|:-|:-|
|![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/a275f2ba-2031-4d86-3ba8-95614ffd865a.png)|![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/e030fe8c-f534-9d26-92e9-1a05c4dba4a4.png)|

# <font color="Red">データ分析1：データ確認</font>
## <font color="HotPink">各変数</font>
|変数名|内容|
|:--|:--|
|datetime|日時|
|season|1 =春、2 =夏、3 =秋、4 =冬|
|holiday|1=祝日、0=他（土日は0扱い）|
|workingday|1=平日（土日祝を除く日）、0=土日祝|
|weather|1=晴れ、2＝曇り・霧、3＝軽い雨雪、4＝強い雨雷雪|
|temp|気温(℃)。アメリカだがセルシウス温度|
|atemp|体感気温(℃)|
|humidity|湿度(%)|
|windspeed|風速（単位不明）|
|casual|利用者数、うち非会員|
|registered|利用者数、うち会員|
|count|利用者数|

countが目的変数。
casual、registeredも説明変数の内訳という意味で目的変数に属する。
目的変数の配置は前述の通り。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/4315f346-0f60-ecde-9da6-ea45afc07045.png)

## <font color ="HotPink">情報要約</font>
一部重複するがまとめておく。
- 欠損値なし（見かけ上。事実は一部欠損値は0または0相当に置換されている模様。後述）
- 説明変数は9種類
- 目的変数はcount、int型
- count、registered、casualはtrain.csvのみに含まれる
- データの日付は2011-1-1 00:00:00　～　2012-12-31 23:00:00の間
- train.csvは毎月1日～19日までのデータ
- test.csvは毎月20日以降のデータ
- datetime以外の説明変数は全て数値型（カテゴリカルなものはエンコード済み）
- datetime（日付に相当する説明変数）の型はobject型
- countは最小で1人、最大で977人と、かなり数値の開きがある

## <font color="HotPink">年、月、日、時刻、曜日への分割</font>
ほぼ前処理不要なデータが提供されているが、「datetime」だけ修正し、年・月・日・時刻にする。
また、曜日情報も参考として欲しい。

*曜日　：　経験として金曜昼はムスリムの礼拝で路上駐車が多くて通行に悩まされた経験があり、「金曜昼だけレンタル自転車で移動する集団がいるのでは？」と期待した。
（結果的には、その特徴は顕在化しなかった。残念。）*

併せて、これらデータは最終的にcategory変数属性にする。

## <font color="HotPink">年月日等データ確認</font>
コードが正しいか確認を兼ねて、count総数で表示している。勉強中故。

各月のcount数を表示。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/a65ee631-b905-ebb3-8770-762b514f9d4e.png)

seasonでも同様であった。

曜日、時刻、workingday、holiday関連。
棒グラフは総数カウントを、線グラフは平均値を表示した。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/7836b4f1-1f61-6f42-8d56-c6515c08f785.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/6a5dac4f-5a5b-4f1c-86c0-e8110096424f.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/45406e25-f245-f70f-b35f-07c54b008630.png)

特に平日と休日とでは分布が大きく異なる。

## <font color="HotPink">天候データ確認</font>
weatherについて。
countの総数表示。赤点はratio(%)。

各天候がいくつ出現したか表示。weather=4は3回しか出現していない。
よって、最終的にweather=4をweather=3に統合した。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/baeda453-89c6-de60-197d-806700c30068.png)

各時間に対するweatherの分布（countの平均値を表示）
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/8c1c9614-eeca-125a-c268-9c69a692371e.png)


各天候に於けるcount値の平均値を時刻に対してプロットすると、分布の形はほぼ同じで量が小さくなる。
この分布の形状は先に示した、全日付に対して合計した「count vs hour」の形状にほぼ一致している。天候による変化は単純な全数の変化のみで、どこかの時刻に対して特別な変化はしない。

## <font color="HotPink">casual、registeredについて確認</font>
count = casual + registered　の関係となる。
casual：非登録ユーザ、registered：登録ユーザ、という違い。
これらの確認。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/5de7a5c2-ebf8-21a5-e104-254f583289b8.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/249799e0-b09a-311f-0dad-0d8c49bdb510.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/b7b61ee6-5606-a274-43c4-8b0950851174.png)

特筆すべきは、ヒストグラム。
特にcasualで顕著。0過剰で右に下がる分布。λ=1のポアソン分布をイメージさせるように分布。
registeredも似ているが、100あたりにも緩いピークを持つような分布。

他、casualでは特に季節的に気候が良い（monthに於ける4月～9月）と利用者数が増えてcountに占める割合が上がる。

## <font color="HotPink">量的変数の可視化</font>
一度に描いた。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/5f9be2f1-15b5-098d-1898-d986b0ee9b10.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/2796d99e-0db2-f682-b29a-6bd40caf46cf.png)

得られる情報のうち一部を抜粋。
- casualはregisteredと比較してtemp、atempが低いと数が顕著に少なくなる
- humidity、windspeedは全体の数値から離れた位置に0となる値を示す。欠損（データ取得できなかったケース）を0で埋めている可能性が高い
- casualとregister（またはcount）間の分布をみると「ハート形」となっており、2種類の異なる系統の確率分布が存在しているように見て取れる

余談だが、atemp=30あたりでhumidity、casualの分布が割れており、またatemp=20あたりでwinspeedの分布が割れている。
atempとtempとは相関の強い散布図を示している中、いくつかの点でatemp=12あたりにtempに依らず一定値を取る点群がある。
外れ値、というより欠損値が12あたりにアサインされているような雰囲気。

# <font color="Red">最初の機械学習：LinearRegression</font>
最初にシンプルに機械学習モデルLinear Regressionでフィッティングして得られる結果を確認。
結果はKaggleのスコア1.33972。上位92.9%。

- 使用モデル：
    - scikit-learn提供のLinear Regerssionモデル
- データ分割：
    - train用データをKFold交差検証で5分割して「train」と「val」データに分ける
- 結果の評価：
    - RMSLE値と、予測値 vs 真値 のプロットを描く
- 使用データ
    - オリジナルのtrain.csvをベースとし、以下変更を加えた変数
        - datetimeを年、月、日、時刻、曜日に分けた
        - datetime、registered、casual を削除した
        
## <font color="HotPink">結果</font>
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/c79f2c61-4492-a735-221a-0fc8dd17218b.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/66f47395-11c3-4e1a-b940-5d8d1f8dea89.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/857eb018-d120-70b6-832a-39203aa3f64c.png)

本来直線に乗るべき、予測値 vs 真値　のプロットを描くと
直線でなく　$\sqrt{x}$　か　log（x）のような分布になっている。

yearの係数が課題。2011年：2012年比は約1.67であったがyearの係数は82。

## <font color="HotPink">スコア確認1</font>
Kaggleにsubmitしたときのスコア。
 - スコア：1.33972
 - 順位」3012位 / 3242人	
 - 上位%：92.9%
 
ただの線形モデルフィットでも下位7%には勝っているらしい。
これの逆だとすると、目標の「上位5%以上」はなかなか厳しそう。

## <font color="HotPink">課題</font>
因子を書き出す

- <details><summary>説明変数の編集（クリックで展開）</summary>

    - temp、atempが強い相関を示しており、共線性が問題となる
    - atempの分布に歪みがあるので正しいatempを算出する、またはtempのみの使用とする
    - 線形モデルなので、seasonや月、時刻等のカテゴリカル変数のエンコード順を昇順または降順にすることで変数の線形近似が良くなると考える
    - humidity、windspeedにおいて、他の値から飛んで0となっているものがあり、欠損値を0で埋めたものの可能性がある。0を適正に「欠損値処理」する
    - 他、場合により量的データに対して標準化や対数化などの常套手段を検討

</details>

- <details><summary>目的変数の扱い（クリックで展開）</summary>

    - casual、registeredで分布の特徴が異なる
    - countを直接求めるのではなく、casual、registeredをそれぞれ求めて和を取る
</details>

- <details><summary>L2正則化（Ridge回帰）（クリックで展開）</summary>
  
    - 説明変数yearの係数が82と大きい（次に大きい数値は9）。L2正則化を行うRidge回帰が望ましい
</details>

# <font color="Red">説明変数の編集</font>
## <font color="HotPink">atempの算出</font>
これは長くなるので別記事にした。

https://qiita.com/SumMae/items/d0548ee2156965bdde12

実施内容と得られた結果。

|実施内容|得られた結果|
|:-|:-|
|事実上の欠損値を確認|humidity=0、atemp=12.12。特定の日付でこの値を示しており、欠損値と判断|
|実際の天気予報を取得|humidity=0を正しく埋める<br>windspeed=0は正しいことを確認<br>tempのノーマライズ式を仮定、検証|
|tempのデノーマライズ|「atemp=12.12」はデノーマライズしたtempの値が0の時のatempの値（切片）であった|
|atemp算出|デノーマライズしたtemp、windspeed、humidityよりatempを算出|
|at_norm取得|atempに施されているノーマライズ式を仮定、検証<br>算出した正しいatempにもノーマライズを適用<br>「at_norm」と名付け、欠損無し一部不連続点を排除した正確なatempとした|

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/dd5c4e26-79f0-63b3-986a-888f840575d7.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/d88f74dc-bd92-220e-a97b-6a1926e4074a.png)

- humidity=0は事実上の欠損値
 → 実際の天気予報より取得し補完

- windspeed=0は正常な値と判断（天気予報情報との照合、一般の風速計の特徴、ヒストグラムのビン幅より）
　→ 微風が目的変数に与える影響はほぼないと考え、不要に値の編集をするのを避けるべく補完はしなかった

よって、humidity=0は正しく補完を完了、windspeed=0は補完の必要なしとした。

# <font color="Red">目的変数の扱い</font>
予測対象を個別に「casual」と「registered」にした。
count = registered + casual　という関係なので和を取って最終予測値にする。
resisteredとcasualとは異なる確率分布を持っているように見られる為。
registered、casualの分布は以下。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/7aac219c-3472-1b20-76f1-b60cfc56036c.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/c1de5e29-7057-b6e4-92e3-e6059b143e61.png)

散布図より、registered vs casualについてworkingdayで色分けするときれいに傾向が分離された。
registeredは勤務日にて高い数値を出し、casualではその逆である。
勤務日か否かが、「casual」と「registered」の特徴を分かつ項目のようである。


# <font color="Red">LightGBMによる予測</font>
目的変数の分布（ゼロ過剰右裾長）より、線形モデルでは予測が厳しい。
- 線形モデル：最小二乗法を使用
- 最小二乗法：誤差が正規分布に従うという前提で考えている
- 目的変数の分布：0過剰の右裾長分布で正規分布とは遠く離れている

→　LightGBM採用：決定木手法のモデルで、且つ動作の軽い手法

KFold交差分割検証でデータを5分割、分解したvalidation用データを使用してrmsle評価用の関数を作成して評価。
同時にKFoldのtrain用データにも同評価を行い、train、validationともに学習曲線を描けるようにした。
KFold交差検証後の本番のモデルでは全データを用い、best_iteration_数だけn_estimatorsを回した。
num_leavesは学習用データの母数に比例する値なので、最適化で得られた値を1.25倍して使用した。
（本番モデル作成ではvalidation用データは確保・使用せず、early_stoppingは使用しない）

## <font color="HotPink">LightGBM最適化での学習</font>
パラメータをチューニングしてLightGBMを使用する。
パラメータはLightGBMの公式ドキュメントと各パラメータの意味より検討した。
（以下のパラメータ名は「sklearn API」をベースにした名称で記載している）

- 最適化対象
    - num_leaves
    - min_data_in_leaf
- 最適化手法
    - optunaによるベイズ最適化
    - early_stoppingを使用。n_estimatorsの値はここから得る
- early_stopping
    - KFold交差検証でvalidationデータを使用
    - validationデータセットに対してrmsleを計算し5ラウンド以上変化しない場合early_stopping
- n_estimatorsの最適値について
    - n_estimators　はearly_stoppingから得られるbest_iteration_の値を最適値とする
- 評価関数
    - モデル学習時のeval_metricには自作したRMSLE値算出関数を設定した

最適化の一例を示す。
optunaの標準関数により、最適値を得るまでの結果をグラフプロットできる。
casualについてのこれを示す。

＜最適化＞
目的変数 = casual、一例として
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/3ea1b922-ae35-0f1e-c050-71ceb4782236.png)
|最適化History|両パラメータ間のコンター図|
|:-|:-|
|![Optuna_History_LGBM_Casual_Optimal01.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/defcbc58-79dc-9413-5bb1-d38ab69aae9b.png)|![Optuna_Contour_LGBM_Casual_Optimal01.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/f831192c-db17-4bbf-a152-7d0dc550fad6.png)|

min_data_in_leafが4とかなり小さい。


＜KFold交差検証＞

目的変数 = casual
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/bef3c6ea-9ee7-3e68-cfb4-28bd94e7c816.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/7a785a22-99e4-1983-bebd-53e6f9f19290.png)


目的変数 = registered
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/e2649bd7-f740-83fc-897e-cccc8e252052.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/d0b08b5b-6ded-1255-47ac-47483f2440a8.png)

＜最終的な学習結果＞
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/bdfec5c8-7f55-df97-00ce-f14e1f3f892a.png)


目的変数 = casualとregistered、合計してcountを得る
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/be6dbe0d-5db8-fa26-23e3-b1384844c541.png)

各RMSLE値
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/768699d8-b7cb-da07-a3a4-4075e3615767.png)

 - 目的変数=casual＋registeredで算出 　：　RMLE値=0.246

## <font color ="HotPink">スコア確認2</font>
Kaggle上でのスコア。
 - スコア：0.47609
 - 順位：1116位 / 3242人	
 - 上位%：34.4%

順位は上記30%に近づいた。目標は上位5%。

## <font color="HotPink">考察</font>

過学習が懸念されるが、他の施策を行うこととする。

まず、以下より一見して過学習と考えられる。
- casual、registeredの「min_data_in_leaf」値が4前後であり、かなり小さい
- KFoldのRMLSE値に於いて、trainとvalidの差が大きい
- KFoldの散布図を見ると、trainではデータのバラつきが小さいがvalidでは大きい
- 学習曲線ではvalidが早期に収束しているのに対しtrainでは下降を続けている

最も理解しやすい過学習要因は「min_data_in_leaf=4前後」と極端に小さいこと。
しかし、これはvalidation用データを用いて最適化したものなので、これを無暗に大きな値にするのは短絡的であり正しい行為ではないと考える。
実はmin_data_in_leafが小さくなるのは、目的変数の分布がゼロで最大を示し右裾長になっているので仕方がないと考えられる。
滅多に起こらない事象が多くを占める分布に対して決定木を作ろうとすると、かなり分岐しないとレアな事象を特定（分離）できないと思われる。
よって、min_data_in_leafは小さい値になって然るべきとイメージできる。

よって、「汎化のための施策」ではなく、他に考えられる施策を行う。

# <font color="Red">スコアを上げるための施策</font>
スコアを上げるための方針を2点実施する。
- 説明変数の作成
- 目的変数の分布考察

## <font color="HotPink">説明変数の作成</font>

以下2点を考慮したい。
- 「ユーザ属性」に相当する情報
- 利用者数の少ない特徴を抽出


<ins>■ユーザ属性について</ins>
　一般に、予測するうえで「ユーザ属性」情報が重要になる。
　例えば「ユーザアンケート」を想像する。
　「年齢」、「性別」、「利用目的」等のユーザの属性の情報を得て、これを基に傾向を予測することが多い。
　今回の説明変数では気象条件と日時のみ。

<ins>■利用者数の少ない特徴について</ins>
　Kaggleコンペ色が強い考え方だが、予測結果はRMLSEで評価されるので小さな値を示す箇所で正確に予想するとスコアが伸びやすい。

### <font color="Orange">通勤利用者、レジャー利用者</font>
workingday=0とworkingday=1とで分布が異なる。
registeredの分布をみると通勤時間、退社時間でピークを持っており「通勤利用者」の像がイメージできる。
これを基に特徴量を作成する。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/a135722e-545a-7094-3a21-dd5c587f88a3.png)

|特徴名|値|条件|
|:-|:-|:-|
|commute|0,1|workingday = 1 , hour = 7 or 8 or 17 or 18|
|leisure|0,1|workingday = 0 , hour = 11 to 17 |

「commute」が通勤利用者。付随して、「leisure」としてレジャー利用者というものも作成した。
本当は「祝日でも働いており、通勤利用している人」がいるように見える。
ユーザ属性で、「祝日関係ない業務に従事している」というユーザ層が分かれば精度が上がりそうなのだが、残念。

### <font color="Orange">1日雨が続く</font>

時系列でデータを見る。
データ量の都合上、日単位でグループ平均してプロットした。
workingday、holiday、weatherも併せてプロットする。
これらは上下動が視認できるように、適当にスケーリングして表示している。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/134d6813-11cd-8233-731e-73599c4fc266.png)

ところどころ急峻な落ち込みを示す箇所があり、これは「weather」と関係が強い（日単位で平均している値であることに注意）。
「weather」は瞬時値でありこれとは正確に対応していないようなので、長期間にわたる悪天候を特徴化したい。
確かに、「今日は一日中雨」という長期の系と、「昼1時間だけ雨」という瞬時値の系とではレンタルバイク利用に与える影響度は異なりそう。

「precipitation」（降雨量）としてweather値の24時間移動平均を作成する（重み無し平均）。

また、「precipitation >  2」 となるデータに「bad_weather」 という0/1の値を付与する。

日単位グループ平均でbad_weather=1に該当する個所を赤で塗りつぶして表示する。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/c22ee5b8-0ab9-53ed-d2db-a5875b8923c2.png)
急峻な落ち込みの多くをカバーしているように思う。

時刻単位の例は、代表2日を表示。
bad_weather=1に該当するエリアを赤色で塗りつぶして表示した。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/2133b552-d7bc-2a0b-488f-2c04f48b9f32.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/598f04e1-af68-7cd5-c5a1-a9410994e949.png)

休祝日ごとに落ち込みがあるので分かりにくいが、うまく特徴をとらえていると思う。
また、bad_weather（およびprecipitation）にかからない落ち込みはweather=3の瞬時値が捉えている。

作成した特徴量まとめ。

特徴名|値|条件
:-|:-|:-
bad_weather|0,1|precipitation>2

*precipitation : weatherの24時間移動平均

### <font color="Orange">利用者数が多い条件</font>
「precipitatoin」、「bad_weather」で利用者数が極小になる条件を特徴化した。
折角なので、これの逆も行う。方法は同じ。
日単位グループの時系列データにて、極大ピークを示す点（年月日）を15個ピックアップ。
その共通条件（条件範囲）を調べた。
特筆すべきはworkingday=0が多く特にcasualの値が大きいという条件が目立った。
これをcomfortと命名し、以下に特徴化する。
|特徴名|値|条件|
|:-|:-|:-|
|comfort|0,1|18 < at_norm < 35 &<br> 40 < humidity < 80 &<br> windspeed < 18 &<br> precipitation < 1.7 & <br> workingday = 0|

### <font color ="Orange">クリスマス</font>
クリスマスが「holiday」になっていなかった。
欧米では年始の正月休みがない代わりにクリスマス以降休暇に入る文化があると経験している。
ドイツ系企業で勤務していたが、12/23を最終日として、以降連絡がつかなくなる。
アメリカでも同様であろう。
思い切って12月24日以降は全てholiday=1にしてしまう。
workingdayも併せて0にする。

### <font color="Orange">monthを削除</font>
日単位グループ平均のデータを見ると、特にbad_weatherの関係で落ち込む月がちらほらある。
本記事の序盤で、「monthは2011年、2012年とで傾向は同じ」と述べたが、よく見ると多少異なる。
月情報は19日間分の訓練用データ情報で、うち5日間程度悪天候が続く月などがある。
2011年と2012年とでこの様子は同様ではないことに気が付く。


以下、monthの分布。各年のcount量に対する各月の占める割合を示す。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/1b28c106-0644-d3c9-675c-af837218291e.png)

seasonでは以下の通り。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/a351b8b8-9244-d14f-55af-cd3eaeba6735.png)

monthを削除してseasonを採用し、情報の平準化をしたほうが良いと判断した。


## <font color ="HotPink">作成した説明変数まとめ</font>

|変数名|表現する特徴|条件設定|
|:-|:-|:-|
|commute|通勤利用者|workingday=1の通勤ピーク時刻|
|leisure|レジャー利用者|workingday=0の日中|
|precipitation|利用者が極小となる条件|長時間悪天候を示す|
|comfort|極大となる箇所、casualをベースに作成|ピーク個所の気象条件を並べた|
|その他|クリスマスを休暇に|12月24日以降をholidayに設定|

### <font color ="Orange">casual、registered での使用</font>
|casualで使用|registeredで使用|
|:-|:-|
|precipitation|precipitation|
|---|commute
|leisure|leisure|
|bad_weather|---|
|comfort|---|

それぞれ、casual、registeredの特徴を記したものはもう一方の特徴量にあてはめが効かないようで、
実験的に交差検証したところ上記組み合わせで使用するとスコアが伸びた。
もともとの変数作成モチベーションと照らし合わせると合点がいく結果である。

## <font color="HotPink">スコア確認3</font>
Kaggle上でのスコア。

スコア：0.42761
順位」536位 / 3242人
上位%：16.5%

順位は前回の半分になった。
特徴量の追加は非常によく働いた。

## <font color="HotPink">目的変数の分布について</font>
目的変数をlog(x+1)で変換することにした。

目的変数の分布がゼロ過剰の右裾長である。
この分布では、誤差が正規分布を仮定する重回帰分析の前提条件から外れる。
（最小二乗法の観点から、正確には撹乱項(残差の確率分布項)が独立に同一の正規分布に従うことを仮定）

最尤度推定により特別な確率分布を求める一般化線形モデルという手法があるが、これは棄却した。
一般化成型モデルでは推定に使う確率分布の妥当性を検証することが重要（尤度関数作成がモデルの肝となるため）。

目的変数はポアソン分布に類似した形を示している。
ポアソン分布の特徴は**平均値と分散とが同じ値（$\lambda$）**となる。
describe()関数で手軽に平均meanと標準偏差stdが得られるのでこれで比較。
casualでは、平均（mean）と分散（std値を2乗した値）は大きく異なっており、ポアソン分布には従わないことが分かる。
他、いくつかの分布の和で表現できないか確認・検討したが難しく、一般化線形モデルの適用は棄却した。
よって、一般化線形モデルの適用は難しいと考え棄却。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/9b49a047-1c45-ad5b-f064-5a1179f7b1d5.png)

（※メモ：DataFrameのdescribe()関数は標本標準偏差（n-1で割る）なので、ここで議論したい値と比べると若干値は大きめに記される。
一様分布であるdayOfWeek等に於いてmeanとstdとの関係が合わないのはこのためと思われる）

替わりに、目的変数にlogを取る手法があることに思い当たる。
目的変数に0が多数あるため、真数条件を考慮しlog(x+1)で変換する。
これのヒストグラムを描くと分布は極端な右裾長ではなくなる。
この分布を以て機械学習を行うと比較的多くのエリアで自然な分布（正規分布）に近い状況に対する予測に落とし込むことができる。
予測能力向上が見込めそうである。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/2e9476a0-9d5f-6950-62c9-258006f4e4f1.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3800268/e45a61bb-a6a0-3231-9328-7eaead55babb.png)

## <font color="HotPink">スコア確認4</font>
Kaggle上でのスコア。

スコア：0.37525
順位：94位 / 3242人
<font color="Red">**上位%：2.9%**</font>

飛躍的に改善。上位3%以上になった。これで目標達成。
ちなみにこの時のLightGBMのパラメータを確認すると、以下の表の通り。
min_data_in_leafは最初の最適化時には2や4という極小さい値であり過学習が懸念された。
最終的には微量ながらこれも改善されており、casualで6、registeredで23となった。

LightGBM パラメータ

＜casual＞
|ハイパーパラメータ|値|
|:-|:-|
|n_estimators|74|
|num_leaves|130|
|min_data_in_leaf|6|

＜registered＞
|ハイパーパラメータ|値|
|:-|:-|
|n_estimators|107|
|num_leaves|128|
|min_data_in_leaf|23|

# <font color = "Red">まとめ</font>
KaggleのBike Sharing Demand に予測を行った。

本コンペの内容の特徴は以下の通り。
- 予測する期間は内挿：毎月20日以降の予測
- 評価値はRMSLE値
- 目的変数の分布は0が多い右裾長の分布
- 目的変数はcasual、registeredという2種類の属性の和となる
- 説明変数は気象条件と年月日、平日、休日情報のみで利用者に関する情報はない
- 見かけ上欠損値は無いが、事実上の欠損値がある

分析時の特徴は以下の通り。
- casual、registeredで特色の異なる分布を示す
- 目的変数はゼロ過剰右裾長の分布であり、正規分布からかけ離れている
- 目的変数をlog(x+1)とすると機械学習の精度を上げやすい
- 決定木をベースとしたLightGBMで比較的良い予測、線形モデルでは分布の形状に追従できない（バリアンスがかなり高くなる）
- 利用者に関する情報がない為、分布からそれに相当するものを作ると良い

<ins>**考察**</ins>
分布がゼロ過剰なうえ、評価指標がRMSLEなのでスコアを上げることが難しい問題であった。
Kaggleリーダーボード上の1位のスコアも0.33756であり、凡その目安でプラス側40%、マイナス側30%の誤差量に相当。
一見すると悪いスコアだが、RMSLEではゼロ付近の誤差が大きく響く特徴を持ち、今回の分布はゼロが多いことに依ると考えられる。

今回の自己のスコアは0.37525、上位約2.9%。
「機械学習」を学んで3か月程度で取得した順位としては悪くないと思う。
一方、モデルは簡易にLightGBMを使用し特徴量と目的変数について考慮したのみ。
データに関する解釈性については現状で良いが、さらに一歩踏み込んでスコアアップをすることもできそう。
アンサンブルにより弱学習機の結果を合わせることでバリアンスを抑えることを試してみたい。
中心極限定理（平均値の寄り集め集団と母平均との誤差は正規分布に収束する）のように良い値に収束していくのかと思われる。
今回のモデルでは予測値と真値の散布図をプロットするとバイアスはあまり感じられないようなので、そこそこ良い結果を得られそう。

___<ins>**追記**</ins>___
→後日、アンサンブルを行った。
　これまでのEDAより線形モデルでもR2スコア0.8程度のスコアを出すようになっており係数も$\pm2$以内に収まっていた。
 Linear Regression、Linear SVR、SVR、Random Forest、LightGBMの5個でアンサンブルを行ったがスコアは悪くなった。
 唯一、Random ForestとLightGBMの2つのみを2:8で平均するとRMSLE=3.7422、上位2.84%になった。
 感触的に、5モデルでの結果を受け、次はモデルを25個ほど用意したときの効果を見たい印象。
 そこまでの多様なモデルは容易に作れなかったため、ここで断念。
<br>


最後に、Kaggleコンペと実務について。
コンペではスコアを競うという特性であり、今回の問題設定は現実のタスクとは離れた問題設定となっている。
各月の20日以降の予測という問題設定だが、実際のタスクでは未来の予測をしたいと思われる。
この場合、2011年→2012年で利用者数が約1.67倍増加していることより、2013年の予測も利用者数は増加しつつ気象や休祝日の影響を受けて変動するものと考えられる。
LightGBMによる決定木ベースのモデルでは外挿には弱いと考えており、また、利用者数増加量については今回のモデルでは「year」特徴量がその解釈を担っていると考えられる。
未来を予測したい場合、year相当の特徴量は年間の統計がすべて出た後に定義可能なので使用できない。
時系列データでの分析から増加量の傾向を特徴化して、これに気象や休祝日による変動の効果を考慮して、といったモデルを作成できると良いように思われる。
こうした分析まで検討できると実用的。

最後に、離職中なので就職活動に勤しまなければいけない。
雇ってくれる方募集。