最終更新日 : 2007.12.3
2本の回帰直線の比較(傾きの検定)
帰無仮説「2本の回帰直線は平行である(傾きは同じである)」を検定するための関数を定義します。共分散分析(連続型変数とカテゴリカル変数を混在させるを参照)は2つの直線が平行であると仮定して行われるものですから、この検定は共分散分析を行うための事前検定として重要となります。
# 関数two.slopeを定義する two.slope <- function(dat){ attach(dat) # 引数に指定されたデータフレームにパスを通す # データフレームより各群のデータを抽出 Y1 <- Y[group==1] X1 <- X[group==1] Y2 <- Y[group==2] X2 <- X[group==2] # 2本の回帰直線を求める model1 <- lm(Y1 ~ X1) model2 <- lm(Y2 ~ X2) # それぞれの残差を求め、残差平方和を求める m1.resid <- sum(resid(model1) ^ 2) m2.resid <- sum(resid(model2) ^ 2) SSr <- m1.resid + m2.resid # 各群の平均値 y1.M <- mean(Y1) y2.M <- mean(Y2) x1.M <- mean(X1) x2.M <- mean(X2) # 2つの回帰直線が平行であるとしたときのパラメータ推定値と残差平方和 sum1 <- sum((X1 - x1.M) * (Y1 - y1.M)) sum2 <- sum((X2 - x2.M) * (Y2 - y2.M)) sum3 <- sum((X1 - x1.M) ^ 2) sum4 <- sum((X2 - x2.M) ^ 2) b <- (sum1 + sum2) / (sum3 + sum4) a1 <- y1.M - b * x1.M a2 <- y2.M - b * x2.M sig1 <- sum((Y1 - (a1 + b * X1)) ^ 2) sig2 <- sum((Y2 - (a2 + b * X2)) ^ 2) SSrb <- sig1 + sig2 m <- length(Y1) # グループ1のデータ数 n <- length(Y2) # グループ2のデータ数 DF2 <- m + n - 4 # 自由度2 Fb <- (SSrb - SSr) / (SSr / (DF2)) # F統計量の計算 pb <- 1 - pf(Fb, df1 = 1, df2 = DF2) # p値の計算 cat("F = ", Fb, "\n") cat("p = ", pb, "\n") detach(dat) }
■ 使用方法
引数にはデータフレームを指定してください。またデータフレームの列名は必ず目的変数をY、説明変数をX、グループをあらわす変数をgroupとして作成してください。そうでないと計算できません。
> g.dat <- as.factor(rep(c(1, 2), c(5, 6))) # グループを表す変数 > x.dat <- c(1.3, 1.5, 1.5, 1.8, 2.0, 1.5, 1.5, 1.8, 1.8, 2.0, 2.2) # 説明変数 > y.dat <- c(17.3, 15.6, 14.1, 13.8, 10.2, 17.8, 15.8, 16.0, 13.6, 11.2, 10.2) # 目的変数 > mydata <- data.frame(Y = y.dat, X = x.dat, group = g.dat) # データフレームの作成 # 使用例 > two.slope(mydata) F = 0.1161944 # F値 p = 0.7432014 # p値
参考文献