![[RjpWiki]](image/../Rlogo.png "[RjpWiki]")
α係数，信頼性係数の算出
青木先生のページ 解説 もあり
> library(psych)
> example(alpha)
alpha> r4 <- sim.congeneric()
alpha> alpha(r4)
Reliability analysis
Call: alpha(x = r4)
raw_alpha std.alpha G6(smc) average_r
0.74 0.74 0.69 0.42
Reliability if an item is dropped:
raw_alpha std.alpha G6(smc) average_r
V1 0.62 0.62 0.53 0.36
V2 0.66 0.66 0.57 0.39
V3 0.70 0.70 0.62 0.44
V4 0.74 0.74 0.66 0.49
Item statistics
r r.cor
V1 0.81 0.74
V2 0.78 0.67
V3 0.73 0.59
V4 0.68 0.50
alpha> r9 <- sim.hierarchical()
alpha> alpha(r9)
Reliability analysis
Call: alpha(x = r9)
raw_alpha std.alpha G6(smc) average_r
0.76 0.76 0.76 0.26
Reliability if an item is dropped:
raw_alpha std.alpha G6(smc) average_r
V1 0.71 0.71 0.70 0.24
V2 0.72 0.72 0.71 0.25
V3 0.74 0.74 0.73 0.26
V4 0.73 0.73 0.72 0.25
V5 0.74 0.74 0.73 0.26
V6 0.75 0.75 0.74 0.27
V7 0.75 0.75 0.74 0.27
V8 0.76 0.76 0.75 0.28
V9 0.77 0.77 0.76 0.29
Item statistics
r r.cor
V1 0.72 0.71
V2 0.67 0.63
V3 0.61 0.55
V4 0.65 0.59
V5 0.59 0.52
V6 0.53 0.43
V7 0.56 0.46
V8 0.50 0.39
V9 0.45 0.32
alpha> #an example of two independent factors that produce reasonable alphas
alpha> #this is a case where alpha is a poor indicator of unidimensionality
alpha> two.f <- sim.item(8)
alpha> alpha(two.f,keys=c(rep(1,4),rep(-1,4)))
Reliability analysis
Call: alpha(x = two.f, keys = c(rep(1, 4), rep(-1, 4)))
raw_alpha std.alpha G6(smc) average_r mean sd
0.59 0.59 0.63 0.15 0.011 0.3
Reliability if an item is dropped:
raw_alpha std.alpha G6(smc) average_r
V1 0.54 0.54 0.57 0.14
V2 0.55 0.55 0.58 0.15
V3 0.57 0.57 0.60 0.16
V4 0.55 0.55 0.59 0.15
V5 0.57 0.57 0.61 0.16
V6 0.55 0.55 0.59 0.15
V7 0.56 0.56 0.60 0.15
V8 0.54 0.54 0.58 0.14
Item statistics
n r r.cor r.drop mean sd
V1 500 0.55 0.47 0.34 0.0717 1.07
V2 500 0.52 0.42 0.30 0.0423 1.01
V3 500 0.47 0.35 0.25 -0.0055 0.98
V4 500 0.52 0.42 0.30 0.0556 0.97
V5 500 0.46 0.33 0.24 -0.0433 1.00
V6 500 0.50 0.40 0.29 0.0100 0.99
V7 500 0.49 0.37 0.26 -0.0203 0.98
V8 500 0.55 0.46 0.33 -0.0186 1.05
alpha> #an example with discrete item responses -- show the frequencies
alpha> items <- sim.congeneric(N=500,short=FALSE,low=-2,high=2,categorical=TRUE) #500 responses to 4 discrete items
alpha> a4 <- alpha(items$observed) #item response analysis of congeneric measures
alpha> a4
Reliability analysis
Call: alpha(x = items$observed)
raw_alpha std.alpha G6(smc) average_r mean sd
0.7 0.7 0.66 0.37 0.0045 0.74
Reliability if an item is dropped:
raw_alpha std.alpha G6(smc) average_r
V1 0.55 0.55 0.46 0.29
V2 0.62 0.62 0.55 0.36
V3 0.65 0.65 0.57 0.38
V4 0.71 0.72 0.63 0.46
Item statistics
n r r.cor r.drop mean sd
V1 500 0.81 0.74 0.62 0.004 0.99
V2 500 0.74 0.61 0.51 0.000 1.01
V3 500 0.71 0.57 0.47 0.062 1.07
V4 500 0.64 0.43 0.36 -0.048 0.98
Non missing response frequency for each item
-2 -1 0 1 2 miss
V1 0.05 0.27 0.35 0.27 0.06 0
V2 0.06 0.25 0.37 0.25 0.07 0
V3 0.06 0.25 0.34 0.25 0.10 0
V4 0.07 0.24 0.42 0.22 0.05 0
alpha> #summary just gives Alpha
alpha> summary(a4)
Reliability analysis
raw_alpha std.alpha G6(smc) average_r mean sd
0.7 0.7 0.66 0.37 0.0045 0.74
> library(ltm) > example(cronbach.alpha) crnbc.> # Cronbach's alpha for the LSAT data-set crnbc.> # with a Bootstrap 95% CI crnbc.> cronbach.alpha(LSAT, CI = TRUE, B = 500) Cronbach's alpha for the 'LSAT' data-set Items: 5 Sample units: 1000 alpha: 0.295 Bootstrap 95% CI based on 500 samples 2.5% 97.5% 0.210 0.356
> library(psy)
> example(cronbach)
crnbch> data(expsy)
crnbch> cronbach(expsy[,1:10]) ## not good because item 2 is reversed (1 is high and 4 is low)
$sample.size
[1] 27
$number.of.items
[1] 10
$alpha
[1] 0.1762655
crnbch> cronbach(cbind(expsy[,c(1,3:10)],-1*expsy[,2])) ## better
$sample.size
[1] 27
$number.of.items
[1] 10
$alpha
[1] 0.3752657
crnbch> datafile <- cbind(expsy[,c(1,3:10)],-1*expsy[,2])
crnbch> library(boot)
Attaching package: 'boot'
The following object(s) are masked _by_ '.GlobalEnv':
simplex
crnbch> cronbach.boot <- function(data,x) {cronbach(data[x,])[[3]]}
crnbch> res <- boot(datafile,cronbach.boot,1000)
crnbch> quantile(res$t,c(0.025,0.975)) ## two-sided bootstrapped confidence interval of Cronbach's alpha
2.5% 97.5%
-0.2710921 0.6082468
crnbch> boot.ci(res,type="bca") ## adjusted bootstrap percentile (BCa) confidence interval (better)
BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 1000 bootstrap replicates
CALL :
boot.ci(boot.out = res, type = "bca")
Intervals :
Level BCa
95% (-0.1680, 0.6404 )
Calculations and Intervals on Original Scale
> library(multilevel) > example(cronbach) crnbch> data(bhr2000) crnbch> cronbach(bhr2000[,2:11]) $Alpha [1] 0.8798196 $N [1] 5400
CRAN のパッケージ置場にあるので
install.packages( 'Rcmdr' )
で入手できる．
reliability 関数を用いて Cronbach's alpha を求められる．
> scale01 <- as.data.frame( cbind( Q1, Q2, Q3, Q4, Q5 ) ) > reliability( cov( scale01 ) )
reliability 関数の含まれる Rcmdr パッケージは R Commander とよばれるウィンドウを開いて操作するので，library() 関数や require() 関数で呼び出すとウィンドウが開いてしまう．単に reliability() を使いたいだけである場合非常に欝陶しい．
対策：R Commander を終了させ，save(reliability, file="reliability.R") とでもして保存し，次回以降は使いたいときに load("reliability.R") とすればよいでしょう。（←Rcmdr::reliability(x)でよいのでは？）
http://home.earthlink.net/~bmagill/MyMisc.html 個人の作成した関数集に含まれているが，出力はあまり丁寧ではない．
