一般化可能性理論

広島大学大学院 教育における統計分析の応用

1 一般化可能性理論

一般化可能性理論とは，測定誤差に焦点を当て，その誤差の成分と大きさを検討することで，尺度や評定の信頼性を検討する方法である。

古典的テスト理論の場合には，

\[ X_{pf} = \tau_{p} + E_{pf} \]

すなわち，個人\(p\)のテスト冊子\(f\)で得られた得点(観測値)\(X_{pf}\)は，個人の能力の真値\(\tau_{p}\)と個人およびテスト形式等にともなう誤差\(E_{pf}\)から成ると考える。\(E_{pf}\)が何によってもたらされるかは議論の対象にはならず，基本的にはテスト冊子を構成する項目の修正や加除によって内的整合性を高めていく。

一方，一般化可能性理論では，個人\(p\)のテスト冊子\(f\)で得られた得点(観測値)\(X_{pf}\)のばらつきをもたらす要因に着目する。個人\(p\)のテスト項目\(i\)に対する得点\(X_{pi}\)は，

\[ \begin{aligned} X_{pi} = {} & \mu && \text{[全体平均]} \\ & + \mu_p - \mu && \text{[受検者要因]} \\ & + \mu_i - \mu && \text{[項目要因]} \\ & + X_{pi} - \mu_p - \mu_i + \mu && \text{[受検者×項目要因]} \end{aligned} \]

からなり，\(X_{pi}\)の分散である\(\sigma^2 (X_{pi})\)は，以下の分散成分からなる。

\[ \begin{aligned} \sigma^2 (X_{pi}) = {} & \sigma_p^2 && \text{[受検者要因]} \\ & + \sigma_i^2 && \text{[項目要因]} \\ & + \sigma_{pi,e}^2 && \text{[受検者×項目要因]} \end{aligned} \]

一般化可能性理論は，古典的テスト理論に基づく一般的に用いられている内的整合性の指標である信頼性係数（Cronbach’s \(\alpha\)）とは異なる，以下のような特徴がある。

• 尺度や評定の結果がばらつく要因は，受検者だけでなく，尺度や評定の項目，評定者も含まれる。さらに，これらの要因の交互作用によってもばらつく。一般化可能性理論は，これらの要因のうち評定結果のばらつきをもたらすものを特定することができる。これを一般化可能性研究という。

• 一般化可能性研究で求めた各要因と交互作用の分散成分の推定値を用いると，信頼性係数に相当する一般化可能性係数を求めることができる。さらに，項目数や評定者数を変化させるシミュレーションを行うことで，一般化可能性係数を高めるためのシナリオを検討することができる。これを決定研究という。

2 分析例

• コードは以下のwebpageを参照した。

  • https://mizumot.com/handbook/?page_id=978

2.1 データ

2項目からなる評価を，評定者2名が，21人の受検者に対して行った場合。

受検者	評定者1_項目1	評定者1_項目2	評定者2_項目1	評定者2_項目2
1	1	2	1	2
2	5	5	4	4
3	5	5	4	4
4	4	5	4	4
5	2	2	2	2
6	3	3	4	3
7	1	1	1	1
8	2	1	2	2
9	2	3	3	3
10	2	2	2	2
11	3	2	3	2
12	5	5	5	5
13	4	5	5	5
14	4	5	4	5
15	5	4	4	5
16	5	5	4	5
17	4	4	5	5
18	2	2	2	2
19	4	5	3	4
20	4	3	4	3
21	4	4	3	4

2.2 一般化可能性研究

2.2.1 データをLong型にする

i1r1 <- data.frame(c(1, 5, 5, 4, 2, 3, 1, 2, 2, 2, 3, 5, 4, 4, 5, 5, 4, 2, 4, 4, 4))
i2r1 <- data.frame(c(2, 5, 5, 5, 2, 3, 1, 1, 3, 2, 2, 5, 5, 5, 4, 5, 4, 2, 5, 3, 4))
i1r2 <- data.frame(c(1, 4, 4, 4, 2, 4, 1, 2, 3, 2, 3, 5, 5, 4, 4, 4, 5, 2, 3, 4, 3))
i2r2 <- data.frame(c(2, 4, 4, 4, 2, 3, 1, 2, 3, 2, 2, 5, 5, 5, 5, 5, 5, 2, 4, 3, 4))
colnames(i1r1) <- c("score"); colnames(i1r2) <- c("score")
colnames(i2r1) <- c("score"); colnames(i2r2) <- c("score")
p <- c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21)

i1r1$person <- p; i1r1$rater <- 1; i1r1$item <- 1
i1r2$person <- p; i1r2$rater <- 2; i1r2$item <- 1
i2r1$person <- p; i2r1$rater <- 1; i2r1$item <- 2
i2r2$person <- p; i2r2$rater <- 2; i2r2$item <- 2

i1r1 <- i1r1[c(2,3,4,1)]; i1r2 <- i1r2[c(2,3,4,1)]
i2r1 <- i2r1[c(2,3,4,1)]; i2r2 <- i2r2[c(2,3,4,1)]
gdata <- rbind(i1r1, i2r1, i1r2, i2r2)

# person, rater, item は要因なので因子型にする
gdata$person <- factor(gdata$person)
gdata$rater  <- factor(gdata$rater)
gdata$item   <- factor(gdata$item)

gdata

   person rater item score
1       1     1    1     1
2       2     1    1     5
3       3     1    1     5
4       4     1    1     4
5       5     1    1     2
6       6     1    1     3
7       7     1    1     1
8       8     1    1     2
9       9     1    1     2
10     10     1    1     2
11     11     1    1     3
12     12     1    1     5
13     13     1    1     4
14     14     1    1     4
15     15     1    1     5
16     16     1    1     5
17     17     1    1     4
18     18     1    1     2
19     19     1    1     4
20     20     1    1     4
21     21     1    1     4
22      1     1    2     2
23      2     1    2     5
24      3     1    2     5
25      4     1    2     5
26      5     1    2     2
27      6     1    2     3
28      7     1    2     1
29      8     1    2     1
30      9     1    2     3
31     10     1    2     2
32     11     1    2     2
33     12     1    2     5
34     13     1    2     5
35     14     1    2     5
36     15     1    2     4
37     16     1    2     5
38     17     1    2     4
39     18     1    2     2
40     19     1    2     5
41     20     1    2     3
42     21     1    2     4
43      1     2    1     1
44      2     2    1     4
45      3     2    1     4
46      4     2    1     4
47      5     2    1     2
48      6     2    1     4
49      7     2    1     1
50      8     2    1     2
51      9     2    1     3
52     10     2    1     2
53     11     2    1     3
54     12     2    1     5
55     13     2    1     5
56     14     2    1     4
57     15     2    1     4
58     16     2    1     4
59     17     2    1     5
60     18     2    1     2
61     19     2    1     3
62     20     2    1     4
63     21     2    1     3
64      1     2    2     2
65      2     2    2     4
66      3     2    2     4
67      4     2    2     4
68      5     2    2     2
69      6     2    2     3
70      7     2    2     1
71      8     2    2     2
72      9     2    2     3
73     10     2    2     2
74     11     2    2     2
75     12     2    2     5
76     13     2    2     5
77     14     2    2     5
78     15     2    2     5
79     16     2    2     5
80     17     2    2     5
81     18     2    2     2
82     19     2    2     4
83     20     2    2     3
84     21     2    2     4

2.2.2 分散成分を推定し検討する

library(lme4)  # lme4パッケージを使う

Loading required package: Matrix

gstudy <- lmer(score ~ 1 + (1|person) + (1|rater) + (1|item) + 
               (1|person:rater) + (1|person:item) + (1|rater:item),
               data = gdata)

boundary (singular) fit: see help('isSingular')

gstudy.sum <- summary(gstudy) #結果を格納する

options(scipen = 100) # 結果の指数表記を避ける
print(summary(gstudy), digits = 4)

Linear mixed model fit by REML ['lmerMod']
Formula: 
score ~ 1 + (1 | person) + (1 | rater) + (1 | item) + (1 | person:rater) +  
    (1 | person:item) + (1 | rater:item)
   Data: gdata

REML criterion at convergence: 183.8

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.39362 -0.59574 -0.07708  0.61079  1.44297 

Random effects:
 Groups       Name        Variance         Std.Dev.   
 person:item  (Intercept) 0.09524129964838 0.308611892
 person:rater (Intercept) 0.07142912115816 0.267262270
 person       (Intercept) 1.51875187567108 1.232376515
 rater:item   (Intercept) 0.00000000003872 0.000006223
 item         (Intercept) 0.00000000000000 0.000000000
 rater        (Intercept) 0.00000000000000 0.000000000
 Residual                 0.13095058381007 0.361870949
Number of obs: 84, groups:  
person:item, 42; person:rater, 42; person, 21; rater:item, 4; item, 2; rater, 2

Fixed effects:
            Estimate Std. Error t value
(Intercept)    3.393      0.279   12.16
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')

なお，boundary (singular) fit というメッセージは，分散が0と推定された成分があることを示すものであり，この分析では問題ない。

varcomp <- VarCorr(gstudy)  

output <- rep(NA, length(varcomp)+1)
for(i in 1:length(varcomp)){
output[i] <- varcomp[[i]][,1]
}
names(output) <- c(names(varcomp), "Residual")
output <- sort(output, na.last=T, decreasing=T)
output["Residual"] <- attr(varcomp, "sc")[[1]]^2
round(output,3)  # 分散成分の大きい順に表示される

      person  person:item person:rater   rater:item         item        rater 
       1.519        0.095        0.071        0.000        0.000        0.000 
    Residual 
       0.131 

# 百分率で表示
round(output/sum(output),3)*100

      person  person:item person:rater   rater:item         item        rater 
        83.6          5.2          3.9          0.0          0.0          0.0 
    Residual 
         7.2 

2.3 一般化可能性研究で分かること

要因	分散成分	割合
受検者(\(p\))	1.519	83.6
評定者(\(j\))	0.000	0.0
項目(\(i\))	0.000	0.0
\(p \times j\)	0.071	3.9
\(p \times i\)	0.095	5.2
\(j \times i\)	0.000	0.0
\(p \times i \times j\)	0.131	7.2

各分散成分は，次のように読む。

• 得点が受検者によってばらつく→受検者の間に実力の違いがある。評価で捉えたいのはこのばらつきである。
• 得点が評定者によってばらつく→評定者によって評価が優しい/厳しいがある。
• 得点が項目によってばらつく→項目によって難しい/易しいがある。
• 受検者と評定者の組み合わせによってばらつく→評定者による受検者のえこひいきがある。
• 受検者と項目の組み合わせによってばらつく→受検者によって項目に対する得意・不得意がある。

この分析例では，受検者の分散成分が全体の83.6%を占めており，このばらつきが一番大きい。評定者と項目の分散成分は0であり，評定者による優しい/厳しいの違いや，項目による難易の違いは見られない。受検者と評定者の組み合わせ（3.9%），受検者と項目の組み合わせ（5.2%）によるばらつきも小さい。

2.4 決定研究

分散成分を\(p\)（受検者），\(j\)（評定者），\(i\)（項目），\(pj\)（\(p \times j\)），\(pi\)（\(p \times i\)），\(ji\)（\(j \times i\)） ，\(pji\)（\(p \times i \times j\)），\(N_j\)（評定者数），\(N_i\)（項目数）としたときに，以下の式に代入すると一般化可能性係数\(E\rho^2\)を求めることができる。

\[ E\rho^2 = \frac{\sigma_p^2}{\sigma_p^2 + \frac{\sigma_{pj}^2}{N_j} + \frac{\sigma_{pi}^2}{N_i} + \frac{\sigma_{pji}^2}{N_j N_i}} \]

Rでは以下のように求める。

# G研究でのoutputの出力を使う
var_p   <- output[["person"]]        # 受検者
var_pj  <- output[["person:rater"]]  # 受検者×評定者
var_pi  <- output[["person:item"]]   # 受検者×項目
var_pji <- output[["Residual"]]      # 受検者×評定者×項目（残差）

# 関数
G.coeff <- function(nj, ni){
  denom <- var_p + var_pj / nj + var_pi / ni + var_pji / (ni * nj)
  G  <- var_p / denom
  return(G)
}

# G研究を行ったのと同じ条件（評定者数2，項目数2）
G.coeff(nj = 2, ni = 2)

[1] 0.9289998

この場合は一般化可能性係数が高く信頼性が高いといえる。

さらに，\(N_j\)と\(N_i\)を変化させて，一般化可能性係数の変化をシミュレーションすることができる。

# 評定者1-4名，1-5項目での変化

# --- ここを先に実行 ---
item <- 1:5   # 項目数1から5まで
ni <- item
reliability <- function(nj, ni){
  denom <- var_p + var_pj / nj + var_pi / ni + var_pji / (nj * ni)
  G  <- var_p / denom
  return(G)
}
# --- ここまで ---

# nj = 数字 で人数を変えて分析
rater1 <- reliability(nj=1, ni)    # 評定者1名の場合
rater2 <- reliability(nj=2, ni)    # 評定者2名の場合
rater3 <- reliability(nj=3, ni)    # 評定者3名の場合
rater4 <- reliability(nj=4, ni)    # 評定者4名の場合

# 結果
rater1   # 評定者1名の場合

[1] 0.8361454 0.8916647 0.9118466 0.9222841 0.9286620

rater2   # 評定者2名の場合

[1] 0.8854751 0.9289998 0.9444747 0.9524072 0.9572309

rater3   # 評定者3名の場合

[1] 0.9032377 0.9421495 0.9558760 0.9628903 0.9671485

rater4   # 評定者4名の場合

[1] 0.9123890 0.9488649 0.9616804 0.9682189 0.9721848

2.4.1 グラフにする

グラフにすると分かりやすい。

# Macで図に日本語を表示するための指定（Windowsでは不要なので，Macのときだけ実行する）
if (Sys.info()["sysname"] == "Darwin") par(family = "HiraKakuPro-W3")
plot(c(0, 0),xlim=c(1, 5),ylim=c(0, 1),
     type = "n", xlab = "項目数", ylab = "一般化可能性係数")

# 軸の目盛り幅を指定
axis(side=2, at=c(0.2, 0.4, 0.6, 0.8, 1))
axis(side=1, at=1:5)

 # 凡例の追加
legend("bottomright", cex=0.7, 
       legend = c("評定者1名","評定者2名","評定者3名","評定者4名"), 
       lty = c(1,2,3,4), 
       pch = c(1,2,3,4), 
       col = c(1,2,3,4),
       bty = "n")

# y軸の0.6で横線を引いて，実線（lty=1）を0.5の太さ（lwd=0.5）で引く
abline(h=0.6, lwd=0.5, lty=1)

# 評定者数ごとにプロットの点や線の指定（点や線，色は凡例に合わせている）
points(ni, rater1, pch=1, col=1); lines(ni, rater1, col=1, lty=1)
points(ni, rater2, pch=2, col=2); lines(ni, rater2, col=2, lty=2)
points(ni, rater3, pch=3, col=3); lines(ni, rater3, col=3, lty=3)
points(ni, rater4, pch=4, col=4); lines(ni, rater4, col=4, lty=4)

3 演習

1. 以下のデータに対して一般化可能性研究を行い，評価方法の改善策を考えなさい。

2. 以下のデータに対して決定研究を行い，評価方法の改善策を考えなさい。

3.1 データ

受検者	評定者1_項目1	評定者1_項目2	評定者2_項目1	評定者2_項目2
1	0	1	1	2
2	3	1	4	2
3	2	1	2	0
4	1	0	2	1
5	1	2	2	1
6	4	3	4	4
7	1	2	1	1
8	2	0	2	0
9	1	1	1	2
10	1	1	1	0
11	1	1	2	1
12	1	1	2	1
13	2	1	1	1

3.2 一般化可能性研究

3.2.1 データをLong型にする

i1r1 <- data.frame(c(0, 3, 2, 1, 1, 4, 1, 2, 1, 1, 1, 1, 2))
i2r1 <- data.frame(c(1, 1, 1, 0, 2, 3, 2, 0, 1, 1, 1, 1, 1))
i1r2 <- data.frame(c(1, 4, 2, 2, 2, 4, 1, 2, 1, 1, 2, 2, 1))
i2r2 <- data.frame(c(2, 2, 0, 1, 1, 4, 1, 0, 2, 0, 1, 1, 1))
colnames(i1r1) <- c("score"); colnames(i1r2) <- c("score")
colnames(i2r1) <- c("score"); colnames(i2r2) <- c("score")
p <- c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13)

i1r1$person <- p; i1r1$rater <- 1; i1r1$item <- 1
i1r2$person <- p; i1r2$rater <- 2; i1r2$item <- 1
i2r1$person <- p; i2r1$rater <- 1; i2r1$item <- 2
i2r2$person <- p; i2r2$rater <- 2; i2r2$item <- 2

i1r1 <- i1r1[c(2,3,4,1)]; i1r2 <- i1r2[c(2,3,4,1)]
i2r1 <- i2r1[c(2,3,4,1)]; i2r2 <- i2r2[c(2,3,4,1)]
data.g2 <- rbind(i1r1, i2r1, i1r2, i2r2)

# person, rater, item は要因なので因子型にする
data.g2$person <- factor(data.g2$person)
data.g2$rater  <- factor(data.g2$rater)
data.g2$item   <- factor(data.g2$item)

data.g2

   person rater item score
1       1     1    1     0
2       2     1    1     3
3       3     1    1     2
4       4     1    1     1
5       5     1    1     1
6       6     1    1     4
7       7     1    1     1
8       8     1    1     2
9       9     1    1     1
10     10     1    1     1
11     11     1    1     1
12     12     1    1     1
13     13     1    1     2
14      1     1    2     1
15      2     1    2     1
16      3     1    2     1
17      4     1    2     0
18      5     1    2     2
19      6     1    2     3
20      7     1    2     2
21      8     1    2     0
22      9     1    2     1
23     10     1    2     1
24     11     1    2     1
25     12     1    2     1
26     13     1    2     1
27      1     2    1     1
28      2     2    1     4
29      3     2    1     2
30      4     2    1     2
31      5     2    1     2
32      6     2    1     4
33      7     2    1     1
34      8     2    1     2
35      9     2    1     1
36     10     2    1     1
37     11     2    1     2
38     12     2    1     2
39     13     2    1     1
40      1     2    2     2
41      2     2    2     2
42      3     2    2     0
43      4     2    2     1
44      5     2    2     1
45      6     2    2     4
46      7     2    2     1
47      8     2    2     0
48      9     2    2     2
49     10     2    2     0
50     11     2    2     1
51     12     2    2     1
52     13     2    2     1

3.2.2 分散成分を推定し検討する

library(lme4)  # lme4パッケージを使う
gstudy <- lmer(score ~ 1 + (1|person) + (1|rater) + (1|item) + 
               (1|person:rater) + (1|person:item) + (1|rater:item),
               data = data.g2)

gstudy.sum <- summary(gstudy) #結果を格納する

options(scipen = 100) # 結果の指数表記を避ける
print(summary(gstudy), digits = 4)

Linear mixed model fit by REML ['lmerMod']
Formula: 
score ~ 1 + (1 | person) + (1 | rater) + (1 | item) + (1 | person:rater) +  
    (1 | person:item) + (1 | rater:item)
   Data: data.g2

REML criterion at convergence: 127.7

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.2268 -0.5677  0.0429  0.5653  1.2259 

Random effects:
 Groups       Name        Variance Std.Dev.
 person:item  (Intercept) 0.314108 0.56045 
 person:rater (Intercept) 0.067306 0.25943 
 person       (Intercept) 0.397435 0.63042 
 rater:item   (Intercept) 0.006411 0.08007 
 item         (Intercept) 0.108970 0.33011 
 rater        (Intercept) 0.009616 0.09806 
 Residual                 0.224358 0.47366 
Number of obs: 52, groups:  
person:item, 26; person:rater, 26; person, 13; rater:item, 4; item, 2; rater, 2

Fixed effects:
            Estimate Std. Error t value
(Intercept)   1.4615     0.3323   4.398

varcomp <- VarCorr(gstudy)  

output <- rep(NA, length(varcomp)+1)
for(i in 1:length(varcomp)){
output[i] <- varcomp[[i]][,1]
}
names(output) <- c(names(varcomp), "Residual")
output <- sort(output, na.last=T, decreasing=T)
output["Residual"] <- attr(varcomp, "sc")[[1]]^2
round(output,3)  # 分散成分の大きい順に表示される

      person  person:item         item person:rater        rater   rater:item 
       0.397        0.314        0.109        0.067        0.010        0.006 
    Residual 
       0.224 

# 百分率で表示
round(output/sum(output),3)*100

      person  person:item         item person:rater        rater   rater:item 
        35.2         27.8          9.7          6.0          0.9          0.6 
    Residual 
        19.9 

3.3 一般化可能性研究で分かること

• 以下の表の空欄を埋めましょう（解答）

要因	分散成分	割合
受検者(\(p\))	0.397	35.2
評定者(\(j\))	0.010	0.9
項目(\(i\))	0.109	9.7
\(p \times j\)	0.067	6.0
\(p \times i\)	0.314	27.8
\(j \times i\)	0.006	0.6
\(p \times i \times j\)	0.224	19.9

• 以下の内容に沿って評価しましょう（解答例）
  • 得点は受検者によってばらつくか→ばらつくが，受検者の分散成分は全体の35.2%にとどまる。分析例（83.6%）と比べると，受検者の実力の違い以外の要因によるばらつきが大きい。
  • 得点が評定者によってばらつくか→ほとんどばらつかない（0.9%）。評定者による優しい/厳しいの違いは小さい。
  • 受検者と評定者の組み合わせによってばらつくか→小さい（6.0%）。評定者による受検者のえこひいきは目立たない。
  • 受検者と項目の組み合わせによってばらつくか→大きい（27.8%）。受検者によって項目に対する得意・不得意がある。項目そのものの難易の違い（9.7%）もある。

改善策（解答例）：評定者に関わる分散成分（\(j\)，\(p \times j\)，\(j \times i\)）は合わせて7.5%と小さく，評定者を増やしても得られるものは少ない。一方，\(p \times i\) が大きいので，2項目だけでは，どの項目が出題されたかによって受検者の得点が左右されてしまう。項目数を増やすことが改善策となる。

3.4 決定研究

# G研究でのoutputの出力を使う
var_p   <- output[["person"]]        # 受検者
var_pj  <- output[["person:rater"]]  # 受検者×評定者
var_pi  <- output[["person:item"]]   # 受検者×項目
var_pji <- output[["Residual"]]      # 受検者×評定者×項目（残差）

# 関数
G.coeff <- function(nj, ni){
  denom <- var_p + var_pj / nj + var_pi / ni + var_pji / (ni * nj)
  G  <- var_p / denom
  return(G)
}

# G研究を行ったのと同じ条件（評定者数2，項目数2）
G.coeff(nj = 2, ni = 2)

[1] 0.6169134

# 評定者1-4名，1-10項目での変化

# --- ここを先に実行 ---
item <- 1:10   # 項目数1から10まで
ni <- item
reliability <- function(nj, ni){
  denom <- var_p + var_pj / nj + var_pi / ni + var_pji / (nj * ni)
  G  <- var_p / denom
  return(G)
}
# --- ここまで ---

# nj = 数字 で人数を変えて分析
rater1 <- reliability(nj=1, ni)    # 評定者1名の場合
rater2 <- reliability(nj=2, ni)    # 評定者2名の場合
rater3 <- reliability(nj=3, ni)    # 評定者3名の場合
rater4 <- reliability(nj=4, ni)    # 評定者4名の場合

# 結果
rater1   # 評定者1名の場合

 [1] 0.3961646 0.5414839 0.6169153 0.6631020 0.6942897 0.7167640 0.7337291
 [8] 0.7469894 0.7576391 0.7663800

rater2   # 評定者2名の場合

 [1] 0.4635487 0.6169134 0.6933817 0.7391944 0.7697078 0.7914892 0.8078176
 [8] 0.8205131 0.8306666 0.8389721

rater3   # 評定者3名の場合

 [1] 0.4914102 0.6469540 0.7232645 0.7685937 0.7986250 0.8199846 0.8359546
 [8] 0.8483463 0.8582413 0.8663251

rater4   # 評定者4名の場合

 [1] 0.5066358 0.6630988 0.7391931 0.7841881 0.8139141 0.8350159 0.8507711
 [8] 0.8629833 0.8727268 0.8806815

3.4.1 グラフにする

# Macで図に日本語を表示するための指定（Windowsでは不要なので，Macのときだけ実行する）
if (Sys.info()["sysname"] == "Darwin") par(family = "HiraKakuPro-W3")
plot(c(0, 0),xlim=c(1, 10),ylim=c(0, 1),
     type = "n", xlab = "項目数", ylab = "一般化可能性係数")

# 軸の目盛り幅を指定
axis(side=2, at=c(0.2, 0.4, 0.6, 0.8, 1))
axis(side=1, at=1:10)

 # 凡例の追加
legend("bottomright", cex=0.7, 
       legend = c("評定者1名","評定者2名","評定者3名","評定者4名"), 
       lty = c(1,2,3,4), 
       pch = c(1,2,3,4), 
       col = c(1,2,3,4),
       bty = "n")

# y軸の0.6で横線を引いて，実線（lty=1）を0.5の太さ（lwd=0.5）で引く
abline(h=0.6, lwd=0.5, lty=1)

# 評定者数ごとにプロットの点や線の指定（点や線，色は凡例に合わせている）
points(ni, rater1, pch=1, col=1); lines(ni, rater1, col=1, lty=1)
points(ni, rater2, pch=2, col=2); lines(ni, rater2, col=2, lty=2)
points(ni, rater3, pch=3, col=3); lines(ni, rater3, col=3, lty=3)
points(ni, rater4, pch=4, col=4); lines(ni, rater4, col=4, lty=4)

改善策（解答例）：現在の条件（評定者2名，2項目）での一般化可能性係数は0.617で，十分とは言えない。一般化可能性係数が0.8に達するのは，評定者2名であれば7項目（0.808），3名であれば6項目（0.820），4名であれば5項目（0.814）のときである。評定者1名では，10項目にしても0.766にとどまる。評定者を2名から4名に倍増させても，2項目のままでは0.663にしかならないのに対し，評定者2名のまま項目を4項目に倍増させれば0.739になる。したがって，評定者を増やすよりも，項目を増やす方が効果的である。