myscratchpad

# ちょっと復習
## 心理統計学の基礎、p. 9 変化量データ。
y1 <- c(4, 3, -3, 4, 1, 4, 0, 6, 2, -6, 2, 6, 6, -2, 1, 1, 6, -2, 7, -1) # 男子
y2 <- c(9, 0, 6, 5, 5, 5, 2, 11, 5, 3, 4, 7, 4, 10, -2, 6, 1, 2, 4, 4) # 女子
n1 <- 20
n2 <- 20
m1 <- mean(y1)
m2 <- mean(y2)
(tres <- t.test(y1, y2, var.equal=T))

# 標本分散の場合
s1 <- var(y1)*((n1-1)/(n1))
s2 <- var(y2)*((n2-1)/(n2))
sps <- sqrt(((s1*n1+s2*n2)/(n1+n2-2)))
(m1-m2)/sps

# 不偏分散の場合
v1 <- var(y1)
v2 <- var(y2)
spe <- sqrt(((v1*(n1-1)+v2*(n2-1))/(n1+n2-2)))
(m1-m2)/spe

# つまるところ、プールされた標準偏差の分子は2群それぞれの平方和を足したもので、それを(n1+n2-2) = 自由度　で割る
# Peasonのrは対応のないt検定だけ

## t.testの結果オブジェクト、もしくはt値とサンプルサイズを与えて効果量を算出する
tef <- function(tres=NULL, tv, n1, n2) {
  if (is.null(tres)) {
      tv <- tv
      df <- n1+n2-2
      }
      else {
       tv <- tres$statistic
       df <- tres$parameter
      }
   #dv1 <- (2*tv)/sqrt(df)
   dv2 <- (2*tv)/sqrt(df+2)
   r1 <- sqrt(tv^2/(tv^2+df))
   #print("# G*Powerと一致するのはdv2")
   #c(dv1=dv1, dv2=dv2, r=r1)
   res <- c(dv2, r1)
   names(res) <- c("Cohen's d", "Pearson's r")
   res
 }

tef(tres)
tef(tv=2.20, n1=50, n2=60)

# 芝・南風原 (1990). 行動科学における統計解析法　p251　より
n1v <- c(50,30,29,25,20,20,15,15)
n2v <- c(60,35,30,30,26,24,18,15)
tvv <- c(2.20,0.91,0.76,1.15,-0.06,0.93,1.64,0.38)

dvs <- list()
for (i in 1:8) {
  n1i <- n1v[i]
  n2i <- n2v[i]
  tvi <- tvv[i]
  x <- tef(tv=tvi, n1=n1i, n2=n2i)
  dvs[[i]] <- x
  }
dvs


# 効果量dを他の色んな効果量に変換するExcelシート
http://www.psychsystems.net/Manuals/StatsCalculators/Effect_Size_Calculator%2017.xls