目的

シェッフェの方法による線形比較を行う

使用法

scheffe(n, m, u, g1, g2, conf.level=0.95)

引数

n           各群のデータ個数のベクトル
m           各群の平均値のベクトル
u           各群の不偏分散のベクトル
g1          第 1 グループの指定
g2          第 2 グループの指定
conf.level  線形比較の信頼区間を計算する信頼率

ソース

インストールは，以下の 1 行をコピーし，R コンソールにペーストする
source("http://aoki2.si.gunma-u.ac.jp/R/src/scheffe.R", encoding="euc-jp")

# シェッフェの方法による線形比較
scheffe <- function( n,                      # 各群のデータ個数のベクトル
                        m,                      # 各群の平均値のベクトル
                        u,                      # 各群の不偏分散のベクトル
                        g1,                     # 第一グループの指定
                        g2,                     # 第二グループの指定
                        conf.level=0.95)        # 線形比較の信頼区間を計算する信頼率
{
        stopifnot(length(n) == length(m), length(m) == length(u), n > 1, u > 0, floor(n) == n, floor(g1) == g1, floor(g2) == g2)
        method <- "シェッフェの方法による線形比較"
        data.name <- paste(deparse(substitute(g1)), "and", deparse(substitute(g2)))
        ng <- length(n)                              # 群の数
        k1 <- ng-1
        nc <- sum(n)                         # 全体のデータ数
        dfw <- nc-ng                         # 群内平方和の自由度
        Vw <- sum(u*(n-1))/dfw                       # 群内平均平方（群内不偏分散）
        n1 <- length(g1)                     # 第一グループにまとめる群数
        n2 <- length(g2)                     # 第二グループにまとめる群数
        g0 <- (1:ng)[-c(g1, g2)]             # どちらのグループにも含まれない群の番号
        n0 <- ng-n1-n2                               # どちらのグループにも含まれない群の数
        weight <- rep(c(1/n1, -1/n2, 0),     # 重み（合計すると 0 になる）
                        c(n1, n2, n0))[order(c(g1, g2, g0))]
        theta <- sum(weight*m)                       # θ推定量
        Vtheta <- Vw*sum(weight^2/n)         # θの分散
        conf.int <- theta-c(1, -1)*          # θの信頼区間
                sqrt(k1*qf(1-conf.level, k1, dfw, lower.tail=FALSE)*Vtheta)
        attr(conf.int, "conf.level") <- conf.level
        F0 <- theta^2/k1/Vtheta                      # F 値
        p <- pf(F0, k1, dfw, lower.tail=FALSE)       # P 値
        return(structure(list(statistic=c(theta=theta, "V(theta)"=Vtheta, F=F0),
                parameter=c(df1=k1, df2=dfw), p.value=p, conf.int=conf.int,
                method=method, data.name=data.name, contrast=list(g1, g2)), class="htest"))
}


使用例

> n <- c(8, 11, 22, 6)	# 例数ベクトル
> m <- c(135.83, 160.49, 178.35, 188.06)	# 平均値ベクトル
> sd <- c(19.59, 12.28, 15.01, 9.81)	# 標準偏差ベクトル

> scheffe(n, m, sd^2, 1:2, 3:4)	# 1，2群と3，4群の線形比較

	シェッフェの方法による線形比較

data:  1:2 and 3:4 
theta = -35.0450, V(theta) = 23.4094, F = 17.4880, df1 = 3, df2 = 43, p-value = 1.435e-07
95 percent confidence interval:
 -49.12185 -20.96815 

> scheffe(n, m, sd^2, 1:3, 4)# 1，2，3 群と4群の線形比較

	シェッフェの方法による線形比較

data:  1:3 and 4 
theta = -29.8367, V(theta) = 42.8136, F = 6.9310, df1 = 3, df2 = 43, p-value = 0.0006583
95 percent confidence interval:
 -48.87380 -10.79954 

> scheffe(n, m, sd^2, 2:3, 4) # 2，3群と4群の線形比較

	シェッフェの方法による線形比較

data:  2:3 and 4 
theta = -18.6400, V(theta) = 43.9185, F = 2.6371, df1 = 3, df2 = 43, p-value = 0.06169
95 percent confidence interval:
 -37.9212072   0.6412072 

 解説ページ