目次

1.はじめに
 2.データセット
 3.データの中身を確認
 4.データの特徴を読み取る
 5.機械学習の実装
 6.ベストスコアと各モデルの結果を比較
 7.予測
 8.まとめ
 9.参考文献

1.はじめに

現在オンラインスクール受講中のデータ分析初学者が、
LightGBMや線形重回帰など様々なモデルを使い、カリフォルニア州の住宅価格予測を実装してみました。

環境

Macbook pro
python3.8
Google Colabo

2.データセット

まず、sklearn.datasets.fetch_california_housing()関数を呼び出します。戻り値として辞書オブジェクトを返します。その配下のdataにデータが、targetにラベルが、feature_namesに特徴名が、関数の引数にas_frame=Trueを含めた場合はframeに全9列のpandasデータフレーム（※target列は最後の列）が格納されています。

import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt


from sklearn.datasets import fetch_california_housing

data_housing = fetch_california_housing(as_frame=True)

3.データの中身を確認

#説明変数と目的変数をそれぞれ変数に格納
X = data_housing.data
y = data_housing.target
print("data_shape:{}".format(X.shape))
print("taget_shape:{}".format(y.shape))

説明変数と目的変数のデータを分けてみました。
Xは説明変数、yは目的変数です。

#表示結果
data_shape:(20640, 8)
taget_shape:(20640,)

データの中身
・MedInc（median income）：ブロックグループの世帯ごとの「所得」の中央値。
・HouseAge（median house age）：ブロックグループの「築年数」の中央値。
・AveRooms（avarage number of rooms）：各ブロックグループ内にある世帯ごとの「部屋数」の平均値。
・AveBedrms（avarage number of bedrooms）：各ブロックグループ内にある世帯ごとの「寝室数」の平均値。
・Population：ブロックグループの「人口」（＝居住する人の総数）
・AveOccup（average occupancy rate）：各ブロックグループ内にある世帯ごとの「世帯人数」の平均値。
。Latitude：ブロックグループの中心点の「緯度」。値が＋方向に大きいほど、そのブロックグループは北にあります。
・Longitude：ブロックグループの中心点の「経度」。値が－方向に大きいほど、そのブロックグループは西にあります。
・MedHouseVal（median house value）：「住宅価格」（100,000ドル＝10万ドル単位）の中央値。

MedIncからLongitudeまでを説明変数、
MedHouseVal（median house value）を目的変数とします。

pandasのdataframeに変換し、説明変数の中身を確認してみます。

df_housing = pd.DataFrame(X, columns=data_housing.feature_names)
df_housing.head()

さらに、目的変数を"Price"としてこのdataframeに入れます。

#目的変数の"price"を同じデータフレームに加える
df_housing["Price"] = y
df_housing.head()

次に、欠損値がないか確かめます。

df_housing.isnull().sum()

MedInc        0
HouseAge      0
AveRooms      0
AveBedrms     0
Population    0
AveOccup      0
Latitude      0
Longitude     0
dtype: int64

訓練データセット、検証データセット共に欠損値はありませんでした。

次にデータの総数、各属性に対する総数、データのタイプを確認します。

df_housing.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 20640 entries, 0 to 20639
Data columns (total 9 columns):
 #   Column      Non-Null Count  Dtype  
---  ------      --------------  -----  
 0   MedInc      20640 non-null  float64
 1   HouseAge    20640 non-null  float64
 2   AveRooms    20640 non-null  float64
 3   AveBedrms   20640 non-null  float64
 4   Population  20640 non-null  float64
 5   AveOccup    20640 non-null  float64
 6   Latitude    20640 non-null  float64
 7   Longitude   20640 non-null  float64
 8   Price       20640 non-null  float64
dtypes: float64(9)
memory usage: 1.4 MB

値は全てfloat型、文字型は一切ないのでこのままデータを使っていきます。

4.データの特徴を読み取る

各カラムをヒストグラムに表します。

#ヒストグラムで可視化
df_housing.hist(bins=50, figsize=(20, 15), label='TrainData')
plt.legend()
plt.savefig('train_hist.png')
plt.show()

"HouseAge"一番右の最大値に注目してみると、異常に数値が伸びていることが確認できます。
実際は築50年以上の住宅が存在しますが、「最大値を50年」としてまとめられているので、データ数が異常に伸びていると推測できます。"Price"に関しても同様です。

次に、相関関係をみてみます。

#相関関係
corr = df_housing.corr()
corr

相関関係は、下の図のように表せます。
一般的には、0.7を超えたら相関は強く、0.2以下だと相関は弱いとされてます。

相関係数のだけだと、全体の関係性が見にくいのでヒートマップで表示します。

#相関関係を可視化
sns.heatmap(corr)

このヒートマップは、色が白いほど正の相関が強く、黒いほど負の相関が強いことを表しています。
ここで、"MedInc"と"price"に注目します。相関係数は0.688と0.7には満たないものの他の特徴量と比べても、目的変数"Price"との相関は強いと捉えることができそうです。

そこで、"MedInc"と"price"の関係を可視化してみます。

#"MedInc","Price"の散布図
fig, ax = plt.subplots()
ax.scatter(x = df_housing['MedInc'], y = df_housing['Price'],s=5)
plt.ylabel('Price', fontsize=15)
plt.xlabel('MedInc', fontsize=15)
plt.show()

"MedInc"(ブロックグループの世帯ごとの「所得」の中央値)と"Price"(住宅価格)が関係性が高そうですね。

最後に、統計量を確認します。
.describe()で、データの統計量を見ることができます。
count（データ数）
mean（平均値）
std（標準偏差（平均値からの分散））
min（最小値）
25%（第１四分位数）
50%（第２四分位数（中央値））
75%（第３四分位数）
max（最大値）

#統計量を見る
df_housing.describe()

まず、平均値と第二四分位数(50%)の値が乖離度合いを見てみます。

mean_median = df_housing.describe().loc[["mean","50%"],:]
mean_median

"population"ブロックグループの人口数の乖離が気になりますが、上記の相関関係から目的変数に与える影響はほぼないと考えられます。
全体的に乖離は低いと確認できます。

次に、最小値、最大値、四分位数の各値との乖離度合いを確認します。

min_max = df_housing.describe().loc[["min","25%","50%","75%","max"],:]
min_max

"MedInc"に注目すると、最大値と第二四分位数に大きな開きがあります。
そこで、正規分布と照らし合わせて考察したいと思います。

正規分布と呼ばれる分布では，平均±標準偏差の範囲に全体のデータの約68%が，平均±2標準偏差の範囲に全体の約95%が，平均±3標準偏差の範囲には全体の約99.7%がそれぞれ含まれています。(下記の図参照）

"MedInc"の正規分布では、最大値の15万ドルは正規分布の平均±3標準偏差の範囲にあると推測でき、外れ値として処理しても良いかもしれません。

5.機械学習の実装

実際に機械モデルを実装し、予測していきます。

#説明変数と目的変数を切り分ける
X = df_housing.drop(['Price'], axis=1)
y = df_housing['Price'].copy()

訓練データと検証データを切り分けます。
今回全データの80%を訓練データセットとして使い、残りの20%をテストデータとします。

#学習用データセットとテスト用データセットに分ける
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y,test_size=0.2, random_state=0)

各データセットの型をデータサイズを確認しておきます。

# 訓練データと検証データのサイズを確認。
print ("X_train :", X_train.shape)
print ("y_train :", y_train.shape)
print ("X_test :", X_test.shape)
print ("y_test :", y_test.shape)

#表示
X_train : (16512, 8)
y_train : (16512,)
X_test : (4128, 8)
y_test : (4128,)

今回使用するモデルは、LightGBM、線形重回帰、ランダムフォレスト、XGBoostでそれぞれの結果を比較します。

from sklearn.linear_model import LinearRegression
from sklearn.ensemble import RandomForestRegressor
import xgboost as xgb
import lightgbm as lgb

max_score = 0 #最大スコアと最も良いモデル名を代入する入れ物を作成しておく（変数の初期化）
best_model = ""

#線形重回帰
model = LinearRegression()
model.fit(X_train,y_train)
score = model.score(X_test,y_test)
print(f"train score:{model.score(X_train,y_train)}")
print(f"test score:{model.score(X_test,y_test)}" )  
a_1 = model.score(X_train,y_train)
b_1 = model.score(X_test,y_test)
if max_score < score:
    max_score = score
    best_model = "線形重回帰"

#ランダムフォレスト
model = RandomForestRegressor()
model.fit(X_train,y_train)
score = model.score(X_test,y_test)
a_2 = model.score(X_train,y_train)
b_2 = model.score(X_test,y_test)
if max_score < score:
    max_score = score
    best_model = "ランダムフォレスト"

# XGBoost 
model = xgb.XGBRegressor()
model.fit(X_train,y_train)
score = model.score(X_test,y_test)
a_3 = model.score(X_train,y_train)
b_3 = model.score(X_test,y_test)
if max_score < score:
    max_score = score
    best_model = "XGBoost"

# LightGBM 
model = lgb.LGBMRegressor()
model.fit(X_train,y_train)
score = model.score(X_test,y_test)
a_4 = model.score(X_train,y_train)
b_4 = model.score(X_test,y_test)
if max_score < score:
    max_score = score
    best_model = "LightGBM"

6.ベストスコアと各モデルの結果を比較

上記でモデルのベストスコアの算出と、訓練データのスコア、テストデータのスコアを算出しました。
下記が結果になります。

ベストスコアを出したモデルと決定係数の表示


print("モデル:{}".format(best_model))

モデル:LightGBM

print("決定係数:{}".format(max_score))

決定係数:0.8386495586073639

各モデルの訓練データスコア、テストデータスコア

model_dic ={"線形重回帰": [a_1,b_1],"ランダムフォレスト":[a_2,b_2], "XGBoost":[a_3,b_3],"LightGBM":[a_4,b_4]}
df_model = pd.DataFrame(model_dic)
df_model

ベストスコアはLightGBMでした。
決定係数は0.838650（表の値は四捨五入されています）でした。

LightGBMは、XGboost(勾配ブースティング)の強化版のイメージです。
勾配ブースティングには、「予測精度は高いが、計算時間が長い」という特徴がありますが、LightGBMは「予測精度を保ったまま計算時間を大きく削減できる」という特徴があります。

7.予測

最後にベストモデルのLightGBMを用いて、価格の予測をします。

lgb_pred = model.predict(X_test)
price_pred = pd.DataFrame({
    'Price': y_test,
    'Prediction': lgb_pred
})

price_pred.head()

8.まとめ

カリフォルニア州の住宅価格の予測を様々なモデルで検証してみました。
初めてデータの読み込みから予測までやったので相当な時間を要しましたが、
予測までのワークフローなど非常に勉強になりました。

データの外れ値処理を施したり、モデルにハイパーパラメータを最適化することにより精度の高い予測ができると思いますので、今後も勉強に励みたいと思います。

ありがとうございました。

9.参考文献

この記事は以下の情報を参考に書きました。
-①カリフォルニアの住宅価格の表形式データセット
-②カリフォルニア住宅価格の予測
-③Kaggler がよく使う「LightGBM」とは？【機械学習】