はじめに
前回の記事で 1 seabornにて公開されているpenguinsのデータセットを使って重回帰分析を行った。seabornにはpenguinsのデータセット以外にもいくつかのデータセットが公開されている。そこで今回はdiamondsのデータセットを使って一元配置分散分析を実施してみる。
その前に
まずは定義を確認しておきます。
一元配置分散分析：1つのカテゴリ変数の3つ以上のグループに基づいて、1つの連続従属変数の平均を比較。
必要なライブラリ、データのインポート
必要なパッケージを読み込み、データの確認を行います。
# パッケージのインポート
import pandas as pd
import seaborn as sns
# ダイヤモンドデータセットをseabonrから読み込む
diamonds = sns.load_dataset("diamonds", cache=False)
# 先頭10行を確認しておく。
print(diamonds.head(10))
carat cut color clarity depth table price x y z
0 0.23 Ideal E SI2 61.5 55.0 326 3.95 3.98 2.43
1 0.21 Premium E SI1 59.8 61.0 326 3.89 3.84 2.31
2 0.23 Good E VS1 56.9 65.0 327 4.05 4.07 2.31
3 0.29 Premium I VS2 62.4 58.0 334 4.20 4.23 2.63
4 0.31 Good J SI2 63.3 58.0 335 4.34 4.35 2.75
5 0.24 Very Good J VVS2 62.8 57.0 336 3.94 3.96 2.48
6 0.24 Very Good I VVS1 62.3 57.0 336 3.95 3.98 2.47
7 0.26 Very Good H SI1 61.9 55.0 337 4.07 4.11 2.53
8 0.22 Fair E VS2 65.1 61.0 337 3.87 3.78 2.49
9 0.23 Very Good H VS1 59.4 61.0 338 4.00 4.05 2.39
今回は各color（色）のグレードに対して、price（価格）の平均が統計的に有意な差があるかどうかを一元配置分散分析を使って確認していきたいと思います。
データの準備
color（色）とprice（価格）のデータを確認しておきます。
まずはcolor（色）から。
# 各色のグレードのダイヤモンドの数を確認
print(diamonds["color"].value_counts())
G 11292
E 9797
F 9542
H 8304
D 6775
I 5422
J 2808
Name: color, dtype: int64
次に、価格のヒストグラムを作成して分布を確認しておく。
import matplotlib.pyplot as plt
fig = sns.histplot(diamonds["price"])
fig.set_xlabel("Residual Value")
fig.set_title("Histogram of Residuals")
plt.show()
分布が大きく偏っていて高い値の長い裾が確認出来ます。そこで対数変換を行い歪んだ分布が対称的になるかを確認します。これは対数変換が正規性の仮定が必要な場合の統計分析とモデリングで特に役立つからですね。
# 価格の対数を取り、3列目に挿入します。
diamonds.insert(10, "log_price", [math.log(price) for price in diamonds["price"]])
# ヒストグラムを作成（省略）
正規分布に近づいてきました。よってprice（価格）ではなく、log_price（価格を対数化した値）を連続変数として扱うことにします。
探索的データ分析（EDA）
次に箱ひげ図を作成して基本的な探索的データ分析（EDA）を行います。
# 色のグレード別に価格の分布を示すボックス プロットを作成する
sns.boxplot(x = "color", y = "log_price", data = diamonds, hue="color", palette="Set1")
plt.show()
箱ひげ図より、ダイヤモンドの各色の log_price（価格を対数化した値） の分布に多くの重複していることが分かりました。しかし、統計的に有意に異なるかどうかはまだ分かりません。
分散分析（ANOVA：Analysis of Variance）を実行するには、回帰モデルを作成する必要があります。これを行う為に、statsmodels.api パッケージと ols() 関数を使用します。そして、変数Xがcolor（色）である、単純な線形回帰モデルを作成し、C() を使用してカテゴリとしてコード化します。次に、モデルをデータに適合させ、モデルの回帰分析の結果を取得します。
# statsmodels.api パッケージと ols() 関数をインポート
import statsmodels.api as sm
from statsmodels.formula.api import ols
# 単純な線形回帰モデルを構築し、モデルを適合する
model = ols(formula = "log_price ~ C(color)", data = diamonds).fit()
# 分析結果の概要を取得
print(model.summary())
OLS Regression Results
==============================================================================
Dep. Variable: log_price R-squared: 0.026
Model: OLS Adj. R-squared: 0.026
Method: Least Squares F-statistic: 237.8
Date: Thu, 09 Jan 2025 Prob (F-statistic): 3.77e-301
Time: 20:22:14 Log-Likelihood: -76617.
No. Observations: 53940 AIC: 1.532e+05
Df Residuals: 53933 BIC: 1.533e+05
Df Model: 6
Covariance Type: nonrobust
=================================================================================
coef std err t P>|t| [0.025 0.975]
---------------------------------------------------------------------------------
Intercept 7.6169 0.012 625.984 0.000 7.593 7.641
C(color)[T.E] -0.0375 0.016 -2.370 0.018 -0.069 -0.006
C(color)[T.F] 0.1455 0.016 9.146 0.000 0.114 0.177
C(color)[T.G] 0.1727 0.015 11.219 0.000 0.143 0.203
C(color)[T.H] 0.3015 0.016 18.390 0.000 0.269 0.334
C(color)[T.I] 0.4061 0.018 22.250 0.000 0.370 0.442
C(color)[T.J] 0.5291 0.022 23.537 0.000 0.485 0.573
==============================================================================
Omnibus: 11794.122 Durbin-Watson: 0.059
Prob(Omnibus): 0.000 Jarque-Bera (JB): 2240.596
Skew: 0.064 Prob(JB): 0.00
Kurtosis: 2.010 Cond. No. 8.56
==============================================================================
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
モデルの出力結果の P>|t| 列を確認すると、color（色）に関連付けられた係数βの p 値はすべて 0.05 未満であることが確認できます。つまり、color（色）は log_price（価格を対数化した値） に影響を及ぼす、と言えます。ただし、各color（色）間で価格の平均に大きな差（優位な差）がある、とは断定できません。そこで、一元配置分散分析を用います。
一元配置分散分析（One-way ANOVA）の実施
帰無仮説と対立仮説
まず、帰無仮説($H_0$)と対立仮説($H_1$)を確認しておきます。
帰無仮説：
H_0 :
𝑝𝑟𝑖𝑐𝑒_D = 𝑝𝑟𝑖𝑐𝑒_E = 𝑝𝑟𝑖𝑐𝑒_F = 𝑝𝑟𝑖𝑐𝑒_G = 𝑝𝑟𝑖𝑐𝑒_H = 𝑝𝑟𝑖𝑐𝑒_I = 𝑝𝑟𝑖𝑐𝑒_J
ダイヤモンドの価格の平均は色のグレードによって差はない。
対立仮説：
H_1 : Not (𝑝𝑟𝑖𝑐𝑒_D = 𝑝𝑟𝑖𝑐𝑒_E = 𝑝𝑟𝑖𝑐𝑒_F = 𝑝𝑟𝑖𝑐𝑒_G = 𝑝𝑟𝑖𝑐𝑒_H = 𝑝𝑟𝑖𝑐𝑒_I = 𝑝𝑟𝑖𝑐𝑒_J)
ダイヤモンドの価格の平均は色のグレードによって異なる。
Pythonで一元配置分散分析（One-way ANOVA）の実施
# 一元配置分散分析（One-way ANOVA）
print(sm.stats.anova_lm(model, typ = 1))
df sum_sq mean_sq F PR(>F)
C(color) 6.0 1431.255783 238.54263 237.807767 3.767555e-301
Residual 53933.0 54099.661516 1.00309 NaN NaN
PR(>F) 列の p 値 が非常に小さいため、帰無仮説を棄却し、対立仮説を採用します。つまり、「ダイヤモンドの価格の平均は色のグレードによって異なる」が結論となります。
技術的な補足。
anova_lm の typ についてはここ 2 とかここ 3 などを参照されたい。
-
Pythonを使って重回帰分析を行う https://qiita.com/tatsu_sekine/items/4c50a0b0e0257761cc2f ↩
-
Anova – Type I/II/III SS explained | R-bloggers https://www.r-bloggers.com/2011/03/anova-%e2%80%93-type-iiiiii-ss-explained/ ↩
-
r - How to interpret type I, type II, and type III ANOVA and MANOVA? - Cross Validated https://stats.stackexchange.com/questions/20452/how-to-interpret-type-i-type-ii-and-type-iii-anova-and-manova ↩