こんにちは！ぼりたそです！
今回はCSVを入力するだけで重回帰分析を実行できるPythonコードを作成しましたので、ご紹介します。
この記事は以下のポイントについてまとめています。
- 重回帰分析について
- 今回のコードの仕組み
- 実装コードと解説
それでは順に解説していきます。
重回帰分析について
まずは重回帰分析について簡単に説明します。
重回帰分析は、統計学や機械学習の中でもよく用いられる分析手法の一つです。主に、複数の説明変数（特徴量）と1つの目的変数（ターゲット）との関係を分析するのに使われます。
例を挙げながら説明します。ある店舗の売上を広告費と入店者数から予測すると仮定します。
■説明変数
- 広告費：どれくらい店舗の広告に費用をかけたか
- 入店者数：どれくらいの人が店舗に入店したか
■目的変数
- 店舗の売上：どれくらい店舗で売上があったか
重回帰分析では、これらのデータを使って、目的変数である「売上」を説明変数である「広告費」「入店社数」から予測することが目的です。
数学的に表現すると、重回帰分析では以下のようなモデルを考えます
売上 = $a_1$ × 広告費 + $a_2$ × 入店者数
ここで、$a_1$, $a_2$はそれぞれ回帰係数と呼ばれるもので、変数の相関関係や重要度を把握することができます。
詳細については以下の記事を参考にしていただければと思います。
実装コードの仕組み
次に今回の実装コードの仕組みについてご説明します。
今回の実装コードについて入力データと出力データは以下の通りとなっています。
■入力データ
- 解析データが格納されたCSVファイル
目的変数が1列目、説明変数が2列目以降になっているCSVファイルを用意してください。インデックスやカラム名意外でデータに文字列などは含まないようにしてください。
■出力データ
- 標準編回帰係数&t値が格納されたCSVファイル
変数の重要度を見積もることができる標準編回帰係数とその回帰係数が統計的に意味があるのかを示すt値をそれぞれの変数ごとにCSVとして出力されるようになっています。
- VIFが格納されたCSVファイル
説明変数同士の多重共線性を評価する指標として使われるVIFをそれぞれの変数ごとにCSVとして出力されるようになっています。
- 変数間の相関係数が格納されたCSVファイル
説明変数同士の相関係数について格納されたCSVファイルが出力されるようになっています。主にVIFが高い変数について、他のどの変数と相関が高いのかを確認するためのファイルになります。
- 自由度調整済み決定係数が格納されたファイル
入力データの回帰モデルについて自由度調整済み決定係数が格納されたファイル。基本的にFold数10のクロスバリデーションで計算され、10回分の決定係数とその平均が出力されます。
実装コードと解説
それでは実装コードの紹介とコードの解説をしていきます。
今回の実装コードは以下の通りです。
import pandas as pd
import numpy as np
from sklearn.model_selection import cross_val_score, LeaveOneOut, KFold
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import StandardScaler
from statsmodels.regression.linear_model import OLS
from statsmodels.tools import add_constant
from statsmodels.stats.outliers_influence import variance_inflation_factor
# CSVファイルの読み込み
data = pd.read_csv('your_data.csv')#解析したいCSVファイルのディレクトリを入力
# 説明変数と目的変数に分ける
X = data.iloc[:,1:]
y = data.iloc[:,:1]
# データの標準化
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
y_scaled = scaler.fit_transform(y.values.reshape(-1, 1)).flatten()
# 定数項を追加
X_scaled_with_const = add_constant(X_scaled)
# 重回帰分析の実行
model = OLS(y_scaled, X_scaled_with_const).fit()
# 標準偏回帰係数とt値を取得
coefficients = model.params
t_values = model.tvalues
# 標準偏回帰係数とt値をCSVとして出力
results_df = pd.DataFrame({
'Variable': ['Intercept'] + list(X.columns),
'Standardized Coefficients': coefficients,
't-value': t_values
})
results_df.to_csv('standardized_coefficients_and_tvalues.csv', index=False)
# VIFの計算のための関数
def calculate_vif(X):
vif = pd.DataFrame()
vif['Feature'] = X.columns
vif['VIF'] = [variance_inflation_factor(X.values, i) for i in range(X.shape[1])]
return vif
# 定数項を除いたVIFを計算
X_scaled_no_const = pd.DataFrame(X_scaled, columns=X.columns)
vif_data = calculate_vif(X_scaled_no_const)
# VIFをCSVとして出力
vif_data.to_csv('vif.csv', index=False)
# 説明変数間の相関係数を計算
correlation_matrix = pd.DataFrame(X_scaled_no_const).corr()
# 相関係数マトリクスをCSVとして出力
correlation_matrix.to_csv('correlation_matrix.csv')
# クロスバリデーションによる決定係数の評価
if len(y) <= 10:
cv = LeaveOneOut()
else:
cv = KFold(n_splits=10, shuffle=True, random_state=1)
# クロスバリデーションで決定係数を評価
model_cv = LinearRegression()
scores = cross_val_score(model_cv, X_scaled, y_scaled, cv=cv, scoring='r2')
# Adjusted R²の計算
n = len(y)
p = X.shape[1]
adjusted_r2 = 1 - (1 - scores) * (n - 1) / (n - p - 1)
mean_adjusted_r2 = np.mean(adjusted_r2)
# Adjusted R²をCSVとして出力
adjusted_r2_df = pd.DataFrame({
'Fold': range(1, len(adjusted_r2) + 1),
'Adjusted R²': adjusted_r2
})
adjusted_r2_df.loc[len(adjusted_r2_df)] = ['Mean', mean_adjusted_r2]
adjusted_r2_df.to_csv('adjusted_r2.csv', index=False)
上記コードでは大きく5つのフローに分かれています。
- データ前処理および重回帰分析の実行
- 標準編回帰係数とt値の出力
- 説明変数の相関係数の出力
- VIFの出力
- 決定係数の出力
順に説明していきます。
データの前処理及び重回帰分析の実行
データの前処理および重回帰分析の実行は以下の部分で実行しています。
# CSVファイルの読み込み
data = pd.read_csv('your_data.csv')#解析したいCSVファイルのディレクトリを入力
# 説明変数と目的変数に分ける
X = data.iloc[:,1:]
y = data.iloc[:,:1]
# データの標準化
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
y_scaled = scaler.fit_transform(y.values.reshape(-1, 1)).flatten()
# 定数項を追加
X_scaled_with_const = add_constant(X_scaled)
# 重回帰分析の実行
model = OLS(y_scaled, X_scaled_with_const).fit()まずはCSVファイルのディレクトリを入力してCSVファイルを読み込み、列指定で目的変数は1列目、説明変数は2列目以降としてわけています。
次にデータの標準化として目的変数、説明変数ともにデータを標準化（min: 0, max: 1）しています。
最後に一応定数項を追加（ほぼ0ですが）して重回帰分析を実行しています。
標準編回帰係数とt値の出力
標準編回帰係数とt値の出力ですが、以下の部分で処理しています。
# 標準偏回帰係数とt値を取得
coefficients = model.params
t_values = model.tvalues
# 標準偏回帰係数とt値をCSVとして出力
results_df = pd.DataFrame({
'Variable': ['Intercept'] + list(X.columns),
'Standardized Coefficients': coefficients,
't-value': t_values
})
results_df.to_csv('standardized_coefficients_and_tvalues.csv', index=False)まずはOLSをインスタンス化したmodelについてparamsメソッドを用いて標準編回帰係数を取得しています。またtvaluesメソッドを使用してt値を取得しています。
次にresults_dfとして先ほど取得した標準編回帰係数とt値について変数をインデックスとしてデータフレームを作成しています。
最後にCSVとして出力しています。
VIFの出力
VIFの出力については以下の部分で処理しています
# VIFの計算のための関数
def calculate_vif(X):
vif = pd.DataFrame()
vif['Feature'] = X.columns
vif['VIF'] = [variance_inflation_factor(X.values, i) for i in range(X.shape[1])]
return vif
# 定数項を除いたVIFを計算
X_scaled_no_const = pd.DataFrame(X_scaled, columns=X.columns)
vif_data = calculate_vif(X_scaled_no_const)
# VIFをCSVとして出力
vif_data.to_csv('vif.csv', index=False)VIFの計算についてはまず、calculate_vif()という計算用の関数を作成しています。引数を説明変数のデータが格納されたデータフレームとしており、それぞれの変数が他の変数でどれだけ説明されるかを計算し、VIFを計算するようになっています。
この関数を使用してVIFを計算し、CSVとして出力しています。
変数間の相関係数の出力
変数間の相関係数の出力については以下の部分で処理しています。
# 説明変数間の相関係数を計算
correlation_matrix = pd.DataFrame(X_scaled_no_const).corr()
# 相関係数マトリクスをCSVとして出力
correlation_matrix.to_csv('correlation_matrix.csv')相関係数は定数項を除いた変数で作成したデータフレームに対して.corr()メソッドを用いて計算しています。
最後に相関係数をマトリクスとしてCSVで出力しています。
決定係数の出力
決定係数の出力は以下の部分で処理しています。
# クロスバリデーションによる決定係数の評価
if len(y) <= 10:
cv = LeaveOneOut()
else:
cv = KFold(n_splits=10, shuffle=True, random_state=1)
# クロスバリデーションで決定係数を評価
model_cv = LinearRegression()
scores = cross_val_score(model_cv, X_scaled, y_scaled, cv=cv, scoring='r2')
# Adjusted R²の計算
n = len(y)
p = X.shape[1]
adjusted_r2 = 1 - (1 - scores) * (n - 1) / (n - p - 1)
mean_adjusted_r2 = np.mean(adjusted_r2)
# Adjusted R²をCSVとして出力
adjusted_r2_df = pd.DataFrame({
'Fold': range(1, len(adjusted_r2) + 1),
'Adjusted R²': adjusted_r2
})
adjusted_r2_df.loc[len(adjusted_r2_df)] = ['Mean', mean_adjusted_r2]
adjusted_r2_df.to_csv('adjusted_r2.csv', index=False)今回、決定係数の算出はクロスバリデーションにより行なっており、データ数が10より少なければLeaveOneOutクロスバリデーションで10以上であれば10Foldクロスバリデーションを行ない、平均化した決定係数を出力します。
クロスバリデーションについては以下の記事にまとめてありますので、興味ある方はご参考ください。
また、決定係数について基本的に説明変数の数が多くなると1に近づいてしまうため、自由度調整済み決定係数としてadjust_r2を計算しています。
最終的には、クロスバリデーションで求めたそれぞれのFoldにおける決定係数とその平均値がCSVとして出力されるようになっています。
終わりに
以上がCSVを入力するだけで重回帰分析を実行し、結果を出力するコードの説明になります。使用しているライブラリもメジャーなものですし、非常に使いやすいコードなのではないでしょうか。ぜひ仕事などで使用する機会があればご参照いただければと思います。
