忍者ブログ
×

[PR]上記の広告は3ヶ月以上新規記事投稿のないブログに表示されています。新しい記事を書く事で広告が消えます。

分散分析の効果量
# サンプル
## Anova関数の分散分析表
                SS num Df Error SS den Df         F       Pr(>F)
(Intercept) 3444.5      1       77      7 313.13636 4.537639e-07
afact        217.5      3       71     21  21.44366 1.347309e-06

d family
# Cohenのf, f二乗

f二乗 = 要因の平方和/誤差の平方和 = F値 * (df1/df2)
21.444366*(3/21)
## つまり、F値に自由度をかけただけ。逆に言えば、F値は要因の"平均"平方和/誤差の"平均"平方和 だから

f はsqrt(f二乗)
##被験者間 一元配置分散分析のときはこれをつかう

r family
# 相関比と偏相関比 (相関比 = イータ二乗) がある
# 相関比 (イータ二乗)

要因の平方和/全体平方和
217.5/(77+217.5+71)
## 全体平方和は要因の平方和と誤差の平方和の総計。spssでは"修正総和"と表示される (spssで"総和"とされているものは切片の平方和も足したもので、相関比の計算には使わない (と思う)
### Anova関数では"SS"と"Error SS"を総計し、(Intercept) の"SS"は引くこと
### aov関数では要因、およびResidualsの"Sum Sq"で表示されている部分の総計

# 偏相関比 (偏イータ二乗)
要因の平方和/(要因の平方和+誤差平方和)
一要因の被験者間分散分析は偏相関比 = 相関比。要因の平方和+誤差平方和 = 全体平方和になるので。被験者内一要因だと切片の誤差平方和 (被験者要因の平方和) を加える必要がある
## Anova関数ではIntercept以外のSS, aov関数ではSum Sq
例だと
217.5/(217.5+71)
####基本は上なのだが、F値は自由度と平方和の関数なので、以下も可能
 [効果Aの自由度(df1)×効果AのF値]/[効果Aの自由度(df1)×効果AのF値+効果Aの検定に用いる誤差の自由度(df2)]
(3*21.44366)/((3*21.44366)+21)
別の計算手続きとしては,F値にdf1をかけてdf2で割り (上のf二乗と同じ) ,その値を「その値に1を足した値」で割る。
((21.44366*3)/21)/(((21.44366*3)/21)+1)
参考: http://home.hiroshima-u.ac.jp/nittono/QA.html

# 偏オメガ二乗
## 反復測定ありだとまた違う式らしい
偏相関比の母集団推定値
要因の平方和 - (要因の自由度*誤差の平均平方和) / 全体平方和
[効果Aの自由度(df1)×(効果AのF値-1)]/ [効果Aの自由度(df1)×(効果AのF値-1)+全データ数]

f二乗の母集団推定値 = オメガ二乗 /(1-オメガ二乗)
## 標本推定によるf二乗とは異なる


そのほかの効果量
# 多変量分散分析 (MANOVA)

1 - Wilks's lambda
これがmultivariate 相関比らしい。多変量偏相関比はまた違うらしい

# カイ二乗検定
2x2 ファイ係数
2x2以外 クラメールの連関係数

重回帰や相関ではそのままR二乗やrを見ればいい

参考文献
水本 篤・竹内 理 (2008). 研究論文における効果量報告のために――基礎概念と注意点―― 英語教育研究, 31, 57-66.
http://www.mizumot.com/files/EffectSize_KELES31.pdf


PR
効果量とは
d family (Cohenのd, f2) とr family (r, eta) がある。
d family は平均値差を標準化したもの、r family は変数間の関連度合いを表したもの
参考:
http://www.ec.kagawa-u.ac.jp/~hori/spss/tokidoki2.html#22
http://home.hiroshima-u.ac.jp/nittono/QA.html#Stat

t検定の効果量
# サンプル
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) # 女子
## 心理統計学の基礎、p. 9 変化量データ。
tres <- t.test(y1, y2, paired=F, var.equal=T)
t.v <- tres$statistic
df <- tres$parameter
tres
        Two Sample t-test

data:  y1 and y2
t = -2.4332, df = 38, p-value = 0.01978
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -4.7632071 -0.4367929
sample estimates:
mean of x mean of y
     1.95      4.55

d family
Cohenのd

平均値差 / sqrt(分散の和/2)
(mean(y1)-mean(y2))/(sqrt((var(y1)+var(y2))/2))
  -0.7694321
## 対応あり、なし、どちらでも同じ (?)
## t値を2倍し、Nの平方根で割るという方法もあるようだが、数値が一致しなくてよくわからない:
2*t.v /sqrt(length(y1)+length(y2))
## これはHedgeのgというらしい。分母を自由度 ( = 38) ではなくNの総計 ( = 40) にするとCohen's dに一致する。分母を自由度にした場合は、効果量の分母で標本分散を使ったときのCohen's dに一致する
2*t.v /sqrt(df)
fv1 <- var(na.omit(y1))*(length(na.omit(y1))-1)/length(na.omit(y1))
fv2 <- var(na.omit(y2))*(length(na.omit(y2))-1)/length(na.omit(y2))
(mean(y1)-mean(y2))/(sqrt((fv1+fv2)/2))
p.rep.t関数で算出されるのは標本の方みたい
## 対応なしで各群の人数が異なる場合は以下を分母にする
s1 <- (var(y1))
s2 <- var(y2)
n1 <- length(y1)
n2 <- length(y2)
sqrt(((n1-1)*s1+(n2-1)*s2)/(n1+n2))

r family
r

sqrt(t.v^2/(t.v^2 + df))
## 群を01としたときの点双列相関係数に等しい
gr <- c(rep(1, 20), rep(0, 20))
y12 <- c(y1, y2)
cor(gr, y12)

#aovの結果
Error: a:s
          Df  Sum Sq Mean Sq F value    Pr(>F)   
a          3 217.500  72.500 21.4437 1.347e-06 ***
s          7  77.000  11.000  3.2535   0.01682 * 
Residuals 21  71.000   3.381 

# Anovaの結果
                SS num Df Error SS den Df         F       Pr(>F)
(Intercept) 3444.5      1       77      7 313.13636 4.537639e-07
afact        217.5      3       71     21  21.44366 1.347309e-06

Anova関数の被験者要因は
平方和: (Intercept) のError SS = 77
平均平方和: (Intercept) のError SS/ (Intercept) のden DF 77/7 = 11
誤差平方和: 要因のError SS/要因のden Df = 71/21 = 3.380952
F値: 平均平方和/誤差平方和 = 11/3.380952 = 3.253521
要因自由度: Interceptのden Df = 7
誤差自由度: 要因のden Df = 21


## aov結果のとりだし
aovrestb <- summary(aov(result ~ a + s + Error(a:s), RB))
## aovrestb[[1]][[1]]["Sum Sq"]
library(ctv)
install.views("Psychometrics")
install.views("SocialSciences")
install.views("Multivariate")
dat <- read.table("https://blog.cnobi.jp/v1/blog/user/89d80905c7038b4121822249e9062fba/1256404631",header=T)

dat1 <- dat

mod1sd <- sd(dat1$mod1)
mod1ab <- dat1$mod1 - mod1sd
mod1bl <- dat1$mod1 - (-1*mod1sd)
mod2sd <- sd(dat1$mod2)
mod2ab <- dat1$mod2 - mod2sd
mod2bl <- dat1$mod2 - (-1*mod2sd)
dat2 <- cbind(dat1, mod1ab, mod1bl, mod2ab, mod2bl)

lmres.whl <- lm(dv ~ idv*mod1*mod2, data = dat2)
summary(lmres.whl)

lmres.1bl.2bl <- lm(dv ~ idv*mod1bl*mod2bl, data = dat2)
lmres.1bl.2ab <- lm(dv ~ idv*mod1bl*mod2ab, data = dat2)
lmres.1ab.2bl <- lm(dv ~ idv*mod1ab*mod2bl, data = dat2)
lmres.1ab.2ab <- lm(dv ~ idv*mod1ab*mod2ab, data = dat2)
summary(lmres.1bl.2bl) #LL
summary(lmres.1bl.2ab) #LH
summary(lmres.1ab.2bl) #HL
summary(lmres.1ab.2ab) #HH

# 標準化値を調べる
dat3 <- dat1
dat3$idv <- scale(dat1$idv)
dat3$mod1 <- scale(dat1$mod1)
dat3$mod2 <- scale(dat1$mod2)
dat3$dv <- scale(dat1$dv)
smod1sd <- sd(dat3$mod1)
smod1ab <- dat3$mod1 - smod1sd
smod1bl <- dat3$mod1 - (-1*smod1sd)
smod2sd <- sd(dat3$mod2)
smod2ab <- dat3$mod2 - smod2sd
smod2bl <- dat3$mod2 - (-1*smod2sd)
dat4 <- cbind(dat3, smod1ab, smod1bl, smod2ab, smod2bl)

slmres.whl <- lm(dv ~ idv*mod1*mod2, data = dat4)
summary(slmres.whl)
slmres.1bl.2bl <- lm(dv ~ idv*smod1bl*smod2bl, data = dat4)
slmres.1bl.2ab <- lm(dv ~ idv*smod1bl*smod2ab, data = dat4)
slmres.1ab.2bl <- lm(dv ~ idv*smod1ab*smod2bl, data = dat4)
slmres.1ab.2ab <- lm(dv ~ idv*smod1ab*smod2ab, data = dat4)
summary(slmres.1bl.2bl) # LL
summary(slmres.1bl.2ab) # LH
summary(slmres.1ab.2bl) # HL
summary(slmres.1ab.2ab) # HH


プロフィール
HN:
tao
HP:
性別:
非公開
職業:
会社員
趣味:
アウトドア、自転車、ジョギング、英語学習
自己紹介:
・千葉在住のサラリーマンです。データ分析っぽいことが仕事。
・今年英検1級取得。今はTOEIC高得点を目指して勉強中。
・興味のあることは野球、アウトドア、英語学習、統計、プログラミング、PC関係などなど。
ブログ内検索
freead
順位表
プロ野球データFreak
セリーグ順位表
パリーグ順位表