回帰分析概説

　回帰分析とは，複数のデータ間の関係性を調べる統計的な手法のことである。 お互いに影響を与え合う値の関係性を調べる相関分析とは異なり， 回帰分析では，従属変数Y（目的変数）は，独立変数X（説明変数）から どの影響を与えるかの関係性を調べます。

　身長と握力の関係を考える場合，握力（目的変数）は身長（説明変数） の1要因のみではなく，性別や年齢など複数の要因にも関わります。 しかし，問題を簡単化するために，１つの主な要因のみ考える場合， 単回帰分析といい，複数要因の場合，重回帰分析という。 ここで単回帰分析のR言語の実装について説明します。

相関と回帰分析

相関：2変数のデータ間に相関があるかどうか相関検定を行い， 相関係数を求め，相関の強さを判断します。 \(r>0.7\)の場合，かなり強い相関と認められます。

回帰分析：2変数のデータ（独立変数Xと従属変数Y）において， 独立変数から対応する従属変数を回帰式 （線形方程式\(Y = a+bX\)）で予測することができます。 また，予測の精度を評価するために， 決定係数\(R^2(r^2, rはピアソン係数)\) (寄与率)が利用されます。\(R^2>0.5 (r > 0.7)\)の場合， 適用がよい（当てはまる率が50%以上）と判断されます。

　例えば，独立変数(身長)\(X\)と従属変数(握力)\(Y\)の 回帰式\(Y= -155.6+1.19X\)と求まったとすると，\(X\)に 身長158(cm)を代入したら，身長158の握力は，\(Y= -155.6+1.19×158=32.42\) で 予測できます。

回帰分析の手順

①帰無仮説：回帰式が成立しない（回帰係数(傾き)bが∞や0である）。

②尺度：独立変数Xと従属変数Yとも連続変数で，間隔と比率尺度のデータが 適用されます。

③正規分布：XとYは正規分布である必要性がないですが，正規分布であれば， 残差（Q-Q残差プロットで確認）が正規分布になる可能性が高まり， XとYの線形性を有することを補助的に検証されます。 正規分布ではない場合，有意な回帰分析ではない可能性があり，外れ値の処理 や変数変換を行うこともあります。

④ピアソン相関検定：XとYの間に１次線形関数（直線）で表すため， ピアソン相関係数\(r\)を求め，有意な相関であることを確認します。 一般に\(r>0.7\)なら，決定係数\(R^2(=r^2)>0.5\)， 有意な回帰分析になる可能性があります。XとYが正規分布ではない場合， 相関検定を行いません。

⑤回帰分析：回帰式 \(Y=a+bX\)を求め，\(p<0.05\)なら， 帰無仮説を棄却でき，有意な回帰分析となります。また， 回帰式の当てはまりの良さを評価するため，決定係数が0.5以上ならば， 当てはまりが良いと判断します。従って，決定係数を予測精度として使われます。

⑥信頼区間の推定：回帰式\(Y = a+bX\)の\(a\)と\(b\)が 良く当てはまる範囲を95%信頼区間に示されます。

　以上の③と④は，有意な回帰分析であるかどうかを判断するための補助的な手段で， 必ず行う必要はありません。

回帰分析の実例

　20名の大学生の身長と体重のデータを表２1に示します。身長（説明変数）から 体重（目的変数）を回帰分析を行い，身長が165.5(cm)の体重を回帰式で予測します。

①詳細な回帰分析

　下記のRプログラムについて説明します。データをベクトルとして，Xに独立変数（説明変数）， Yに従属変数（目的変数）を入力します。まず，XとYの正規性を シャピロ検定で行い，両方とも正規分布ならば，ピアソン相関検定を行い， ピアソン係数とp値が出力されます。 XまたはYが正規分布ではないときにピアソン係数が利用できないから， 有意な回帰分析にならない可能性があると表示します。

　回帰分析の結果について，回帰式の切片\(a\)と傾き\(b\)が求められ， それぞれのp値が示されます。両方とも\(p<0.05\)であれば，有意な回帰式となります。 回帰式の95%信頼区間や回帰式，決定係数などが出力されます。 また，散布図に回帰直線，回帰式と決定係数も示され，回帰分析の 結果を確認しやすくなります。

#回帰分析で散布図に回帰直線と回帰式，決定係数
# 例１データ(Xは説明変数，Yは目的変数)
X <- c(171.3,168.1,171.3,168.6,170.6,162.3,170.7,175.1,168.9,172.4,176.5,168.5,172.6,167.9,161.1,167.4,174.9,166.3,174.3,171.1)
Y <- c(71.5,66,75.5,75.5,79.5,65.5,78.5,79,70.5,80.5,94,71,85,70,61,73.5,78,68,80.5,76)

#データの正規性をシャピロ・ウィルク検定で確認
shapiro_x <- shapiro.test(X)
shapiro_y <- shapiro.test(Y)
cat("Xのシャピロ検定 p値: ", shapiro_x$p.value, "\n")
cat("Yのシャピロ検定 p値: ", shapiro_y$p.value, "\n")

# 正規性の判定 (p値 > 0.05なら正規分布とみなす)
if (shapiro_x$p.value > 0.05 & shapiro_y$p.value > 0.05) {
    # ピアソン相関係数の計算
  result <- cor.test(X, Y, method = "pearson")
  cat("\nピアソン相関係数: ", result$estimate, "\n")
  cat("相関検定p値: ", result$p.value, "\n")
  } else {
  cat("\nXまたはY正規分布でないため、ピアソン係数が利用できない。
よって，下記の回帰分析は有意ではない可能性がある。\n")
}

# 回帰モデルの作成
kaiki <- lm(Y ~ X)
cat("\n回帰分析の結果とｐ値 \n")
print(summary(kaiki)$coefficients)
# 回帰係数を取得
intercept <- coef(kaiki)[1]
slope <- coef(kaiki)[2]
# 決定係数 (R-squared) を取得
r_squared <- summary(kaiki)$r.squared
# 回帰直線の数式とR-squaredを作成
equation <- paste("\nY = ", round(intercept, 2), " + ", round(slope, 2), "X", 
                  "\nR-squared = ", round(r_squared, 2), sep = "")
cat("\n回帰式と決定係数: ", equation, "\n")

# 散布図の作成
plot(X, Y, main = "散布図における回帰直線，回帰式，決定係数",
     xlab = "X", ylab = "Y", pch = 19)
# 回帰直線を追加
abline(kaiki, col = "blue", lwd = 2)
# 数式とR-squaredをプロットに追加
text(x = min(X), y = max(Y)-1, labels = equation, pos = 4, col = "red")
# 残差を取得
residuals <- residuals(kaiki)
# Q-Qプロットの作成
qqnorm(residuals, main = "Q-Q 残差プロット")
qqline(residuals, col = "green", lwd = 2)  #正規分布の確認

実行した結果は次のようになる。

## Xのシャピロ検定 p値:  0.5006022
## Yのシャピロ検定 p値:  0.7806714
## 
## ピアソン相関係数:  0.8608575 
## 相関検定p値:  1.108128e-06
## 
## 回帰分析の結果とｐ値
##                Estimate Std. Error   t value     Pr(>|t|)
## (Intercept) -201.230938 38.4880397 -5.228402 5.678162e-05
## X              1.624642  0.2263479  7.177629 1.108128e-06
## 
## 回帰式と決定係数:  
## Y = -201.23 + 1.62X
## R-squared = 0.74

　残差Q-Qプロットでは，残差がプロット上で直線に沿って並ぶ場合, 残差が正規分布に従っていることを認められます。曲線の場合は， 正規性の仮定（XとYは線形性をもつ）に違反している可能性があり， 有意な回帰分析ではないと判断します。

　プログラムを実行した結果，回帰式 \(Y=-201.23 + 1.62X\)と 決定係数 \(R^2=0.74\)が表示されます。回帰式の\(X\)に165.5を代入して，Yを計算すると， 表２の計測結果になかった身長165.5の体重は\(-201.23+\) \(1.62×165.5\) \(\fallingdotseq66.9\)と予測できます。この予測値の当てはまる確率は74%以上になります。

②回帰式のみ出力

回帰分析の途中のデータとグラフを示されず，回帰式と決定係数のみを 表示するプログラムは次のようになる。

#回帰分析で散布図に回帰直線と回帰式，決定係数
# 例１データ(Xは説明変数，Yは目的変数)
X <- c(171.3,168.1,171.3,168.6,170.6,162.3,170.7,175.1,168.9,172.4,176.5,168.5,172.6,167.9,161.1,167.4,174.9,166.3,174.3,171.1)
Y <- c(71.5,66,75.5,75.5,79.5,65.5,78.5,79,70.5,80.5,94,71,85,70,61,73.5,78,68,80.5,76)

# 回帰モデルの作成
kaiki <- lm(Y ~ X)
# cat("\n回帰分析の結果とp値 \n")
# print(summary(kaiki)$coefficients)
# 回帰係数を取得
intercept <- coef(kaiki)[1]
slope <- coef(kaiki)[2]
# 決定係数 (R-squared) を取得
r_squared <- summary(kaiki)$r.squared
# 回帰直線の数式とR-squaredを作成
equation <- paste("\nY = ", round(intercept, 2), " + ", round(slope, 2), "X", "\nR-squared = ", round(r_squared, 2), sep = "")
cat("\n回帰式と決定係数: ", equation, "\n")
## 
## 回帰式と決定係数:  
## Y = -201.23 + 1.62X
## R-squared = 0.74

結論

　単回帰分析の場合，この２つのプログラムが使えます。データが多い場合， データフレームやCSVファイルを読み込むように修正すればよい。 また，詳細な統計的解析を行う場合，①のプログラムを使い， 回帰分析の結果とp値を確認し，有意でない回帰分析になった場合， その結果（回帰式と決定係数）を慎重に解析して取捨する。単純に回帰式と 決定係数\(R^2\)を求める場合，②のプログラムを使ったほうがよい。

  1. 石川朗，対馬栄輝，「リハビリテーション統計学」第2版， 中山書店，2024年3月。↩︎