忍者ブログ
1

×

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

# 2件法の01データの場合はKuder-Richardsonの公式というのを使うらしい。全く知らなかった。知らないことが多いものである。日本語でぐぐると99件しかヒットしない。皆知らないのだろう。
 # 最近4件法や5件法は順序尺度として扱うべしとの意見があるが、信頼性を出すときはどうするんだろう。クロンバックでいいのだろうか
 # このサイトに少し解説。以下適当な訳と覚書
  # 1. 多くの研究でそうしているように、クロンバックのアルファを使う。その仮定が満たされている限り
  # 2. ポリコリック相関行列からクロンバックのアルファ的なものを計算する (B D Zumbo, A M Gadermann, and C Zeisser. Ordinal versions of coefficients alpha and theta for likert rating scales. Journal of Modern Applied Statistical Methods, 6: 21-29, 2007.)
 # 3. 01データならKuder-Richardsonの公式を使う
 # 4. 他にもいろんな信頼性指標がある (R E Zinbarg, W Revelle, I Yovel, and W Li. Cronbach’s α, Revelle’s β, and McDonald’s ωh: Their relations with each other and two alternative conceptualizations of reliability. Psychometrika, 70(1): 123-133, 2005.)
 # A good review is Bruce Thompson. Score Reliability. Contemporary Thinking on Reliability issues. Sage Publications, 2003. だそうだ。まあとりあえずはクロンバックを使おう。気が向いたらRevelle先生のやつを読もう
# 関連記事:Rで級内相関

# Kuder-Richardsonの公式
 # 以下参考URL。感謝
 # http://www.eflclub.com/elvin/publications/2003/itemanalysis.html
 # http://www.gifted.uconn.edu/siegle/research/Instrument%20Reliability%20and%20Validity/Reliability.htm
 # http://www.jstage.jst.go.jp/article/rika/17/3/17_181/_article/-char/ja

 dat <- data.frame(
 q001 = c(1,1,0,0,1,0,0,0,0),
 q002 = c(1,1,0,1,1,1,0,0,0),
 q003 = c(1,1,1,0,1,0,0,0,0),
 q004 = c(1,1,1,1,1,0,0,0,0),
 q005 = c(1,1,1,1,1,1,0,1,0),
 q006 = c(1,0,1,1,1,0,0,1,0),
 q007 = c(1,1,1,0,0,0,1,0,0),
 q008 = c(1,1,1,1,0,1,1,0,0),
 q009 = c(1,1,1,1,0,0,0,0,0),
 q010 = c(1,1,0,0,0,0,0,0,0)
 )
 dat

 (ni <- ncol(dat)) # 項目数
 (ns <- nrow(dat)) # 参加者数
 (mt <- mean(rowSums(dat))) # 合計得点の平均点
 (sdt <- sd(rowSums(dat))) # 合計得点の標準偏差
 (sds <- sqrt(var(rowSums(dat))*((ns-1)/ns))) # 合計得点の標本標準偏差

 pc <- colSums(dat)/ns # 各項目ごとの正答率 (1の個数/全参加者数)
 pw <- 1-pc # 各項目ごとの誤答率 (0の個数/全参加者数)
 pcpw <- pc*pw
 spcpw <- sum(pcpw)


# Kuder-Richardson 20
 (ni/(ni-1)) * (((sds*sds)-spcpw)/(sdt*sdt)) # 不偏標準偏差を使ったとき
 (ni/(ni-1)) * (((sds*sds)-spcpw)/(sds*sds)) # 標本標準偏差を使ったとき

# Kuder-Richardson 21
# おそらく、ローデータがないときに使う簡易版
ni/(ni-1) * (1-(mt-mt*mt/ni)/(sdt*sdt)) # 不偏標準偏差を使ったとき
ni/(ni-1) * (1-(mt-mt*mt/ni)/(sds*sds)) # 標本標準偏差を使ったとき

# よくわからんけど標本標準偏差を使うのかな



# おまけ。カッパ係数
 # 参考:http://www.mizumot.com/stats/kappa.htm
a <- c(2,1,3,2,2,3,2,1,3,3,3,3,2,1,2,3,1,3,1,2)
b <- c(2,1,3,2,3,3,3,1,3,3,3,3,2,1,2,3,2,3,1,2)
x <- table(a,b) # データフレームではなくtableで集計表にする
x

library(psych)
print.digits(cohen.kappa(x), digits=5)
 # 普通はunweighted kappaのestimateを見る

# psychパッケージ以外。まあpsychパッケージが一番いいけど
library(vcd)
Kappa(x)

PR


cprp <- function(x, p=1) {
xn <- addmargins(x)
cp <- t(t(xn)/xn[nrow(xn),]*p)[-nrow(xn), -ncol(xn)]
rp <- (xn/xn[, ncol(xn)]*p)[-nrow(xn), -ncol(xn)]
pt <- prop.table(x)
return(list(adm=xn, cp=cp, rp=rp, pt=pt))
}

x <- matrix(c(10,10,20, 15, 5,20, 25,10,5, 5,5,10), nr=3,
dimnames=list(paste("r", 1:3, sep=""), paste("c", 1:4, sep="")))
x
cprp(x)

x <- rnorm(1000)

## 正規性の検定
hist(x)
qqnorm(x)
qqline(x)
shapiro.test(x) # p > .05なら正規分布

ks.test(x=x, y="pnorm", mean=mean(x), sd=sd(x))

## ks.testは「2群の分布が等しい」という帰無仮説を検定するもの。正規性の検定の場合、yに正規分布を指定している
x2 <- rexp(1000) # 正規分布じゃないやつ
ks.test(x=x, y=x2)


## 多変量正規性
library(mvnormtest)
data(EuStockMarkets)
C <- t(EuStockMarkets[15:29,1:4])
mshapiro.test(C)


下のspssのページをもとにした。
http://www.nyu.edu/its/statistics/Docs/intracls.html
下のページも参照。感謝
http://www.hs.hirosaki-u.ac.jp/~pteiki/research/stat/S/icc/

dat <- data.frame(
rater1=c(9,6,8,7,10,6),
rater2=c(2,1,4,1,5,2),
rater3=c(5,3,6,2,6,4),
rater4=c(8,2,8,6,9,7)
)

# psycパッケージのICC関数を使う。psychパッケージはマジ便利
library(psych)
ICC(dat)
## 勉強のため、ICCのコードの一部を貼り付けておく
x <- dat
n.obs <- dim(x)[1]
nj <- dim(x)[2]
x.s <- stack(x)
x.df <- data.frame(x.s, subs = rep(paste("S", 1:n.obs, sep = ""),
nj))
aov.x <- aov(values ~ subs + ind, data = x.df)
s.aov <- summary(aov.x)
stats <- matrix(unlist(s.aov), ncol = 3, byrow = TRUE)
MSB <- stats[3, 1]
MSW <- (stats[2, 2] + stats[2, 3])/(stats[1, 2] + stats[1,
3])
MSJ <- stats[3, 2]
MSE <- stats[3, 3]
(ICC1 <- (MSB - MSW)/(MSB + (nj - 1) * MSW)) # ICC1. Single raters absolute
(ICC2 <- (MSB - MSE)/(MSB + (nj - 1) * MSE + nj * (MSJ - MSE)/n.obs)) # ICC2. Single random raters
(ICC3 <- (MSB - MSE)/(MSB + (nj - 1) * MSE)) # ICC3. Single fixed raters
(ICC12 <- (MSB - MSW)/(MSB)) # ICC1k. Average random raters
(ICC22 <- (MSB - MSE)/(MSB + (MSJ - MSE)/n.obs)) # ICC2k. Average random raters
(ICC32 <- (MSB - MSE)/MSB) # ICC3k. Average fixed raters


head(x.df)
aov.x
s.aov
round(stats, 3)
alpha(dat) # 質問紙の信頼性係数を求めるときと同じだが、列には項目ではなく評定者がならぶ
gplotパッケージplotmeans関数を使ってみよう
こちらにも解説がある

もしくはpsychパッケージのerror.bars.by関数

凡例表示を覚えないとどうしようもないが

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