IRTで尺度のバリデーションを行う - 井出草平の研究ノート

データを用意

bfiを使う。NEO性格特性のN（神経症傾向）項目だけを抽出する。

データの説明はこちら

ides.hatenablog.com

library(psych)
data(bfi)
BFI_N <- bfi[, 1:5]  # N項目（1～5列目）を抽出
colnames(BFI_N) <- c("N1", "N2", "N3", "N4", "N5") # 変数名の変更
head(BFI_N)

      N1 N2 N3 N4 N5
61617  2  4  3  4  4
61618  2  4  5  2  5
61620  5  4  5  4  4
61621  4  4  6  5  5
61622  2  3  3  4  5
61623  6  6  5  6  5

段階反応モデル（GRM）の適合

library(mirt)

# GRMモデルの適合
mod_grm <- mirt(BFI_N, 1, itemtype = "graded")

# モデル概要の表示
summary(mod_grm)

> summary(mod_grm)
       F1    h2
E1 -0.660 0.436
E2 -0.788 0.621
E3  0.632 0.399
E4  0.759 0.576
E5  0.569 0.323

SS loadings:  2.355
Proportion Var:  0.471

Factor correlations:

   F1
F1  1

結果の確認

各項目のパラメータ（識別力と閾値）

coef_grm <- coef(mod_grm, simplify = TRUE)
print(coef_grm)

$items
       N1    d1     d2     d3     d4     d5
N1 -0.862 0.780 -0.641 -1.426 -2.391 -3.843
A2  1.840 5.575  3.936  3.027  1.214 -1.195
A3  2.526 5.748  4.052  2.956  1.020 -1.845
A4  1.047 3.511  2.337  1.748  0.742 -0.433
A5  1.701 5.109  3.326  2.244  0.628 -1.612

この部分は、各項目の識別力パラメータ (N1) と困難度/閾値パラメータ (d1-d5) を表している。

識別力パラメータ (N1):

各項目が、潜在特性（この場合は神経症傾向）をどれだけ識別できるかを表す指標。
値が大きいほど、その項目は潜在特性の高い人と低い人をよく区別できることを意味する。
一般的に、1以上であれば良好な識別力を持つとされている。
N1 (N1)の識別力パラメータは-0.862と負の値になっている。これは、この項目が他の項目とは逆の方向に機能していることを意味する。

困難度/閾値パラメータ (d1 〜 d5):

各項目に反応するために必要な潜在特性のレベルを表す。
値が大きいほど、その水準で反応するのが難しいことを意味する。
N2 (A2)のd1は5.575です。これは、N2でカテゴリ1からカテゴリ2に上がるためには、神経症傾向のレベルが5.575以上である必要があることを意味する。
一般的に、閾値パラメータは、d1 < d2 < d3 < d4 < d5 の順に大きくなる。

モデル適合度の確認

# モデルの適合度指標
print(mod_grm)

Log-likelihood = -19604.66
Estimated parameters: 30
AIC = 39269.33
BIC = 39447.45; SABIC = 39352.13

視覚化カテゴリ反応曲線（CRC）の描画

# カテゴリ反応曲線の描画
plot(mod_grm, type = "trace")

この図は、神経症傾向という潜在特性のレベル（θ）と、各項目・各カテゴリに反応する確率（P(θ)）の関係を示している。

横軸 (θ): 潜在特性である神経症傾向のレベルを示す。右に行くほど神経症傾向が高い。
縦軸 (P(θ)): 各項目・各カテゴリに反応する確率を示す。
曲線 (P1～P6): 各項目について、神経症傾向のレベルごとに、各カテゴリ（1～6）に反応する確率を示している。

各項目の特徴

N1: 明らかに他の項目と異なる挙動を示している。潜在特性が低いときにカテゴリ1に反応する確率が高く、高くなるにつれてその確率が低下している。逆転項目である可能性が極めて高い。また、真ん中の選択肢（P3）が他の項目に比べて選ばれにくい傾向が見られる。
A2, A3, A5: 潜在特性が高くなるにつれて、高いカテゴリ（特にP5）に反応する確率が高くなる、典型的なGRMのパターンを示している。しかし、真ん中の選択肢（P3, P4）が突出して選ばれやすい傾向が見られる。
A3: 他の項目に比べて、潜在特性が高いレベルでP5に反応する確率が最も高くなる。識別力が高いことを示唆している。
A4: 潜在特性が非常に高いときに、P6に反応する確率が他の項目よりも高くなる。非常に高い神経症傾向を測定するのに適した項目である可能性が高い。A2, A3, A5 と同様、真ん中が選択されやすい傾向にある。

閾値パラメータとの関連

カテゴリ反応曲線上で、各カテゴリの曲線が交差する点が、閾値パラメータに相当する。A3で閾値パラメータd2からd5の値が減少していたことは、カテゴリ反応曲線でP3, P4, P5, P6の曲線が左側（低い潜在特性レベル）で密集していることと整合する。他の項目でも、d3以降の閾値パラメータの減少傾向が見られ、カテゴリ反応曲線でもP4, P5, P6の曲線が左側に寄っていることが確認できる。

まとめ

N1は逆転項目である可能性が極めて高い。得点を反転して再分析する必要があるだろう。
全体的に、高いカテゴリ（特にカテゴリ5）に反応する確率が高く、神経症傾向を過大評価している可能性がある。
多くの項目で、真ん中の選択肢（P3, P4）が他のカテゴリに比べて高くなっており、回答が集中する傾向がある。これは、回答の偏り、または、項目の特性に起因する可能性がある。－A3は識別力が高く、A4は非常に高い神経症傾向を測定するのに適している可能性が示唆された。

テスト情報関数（TIF）の描画

# テスト情報関数の描画
plot(mod_grm, type = "info")

この図は、神経症傾向を測定する5項目（N1～A5）から構成されるテストが、潜在特性である神経症傾向（θ）のどのレベルにおいて、どれだけの情報を提供できるかを示している。

横軸 (θ)：神経症傾向のレベルを表す。右に行くほど神経症傾向が高い。
縦軸 (I(θ))：テスト情報量を表す。この値が大きいほど、その潜在特性レベルにおける測定の精度が高いことを意味する。
曲線：神経症傾向のレベルごとに、テスト情報量がどのように変化するかを示している。

このテスト情報関数から読み取れる情報

最も情報量が高い領域：
- このテストは、神経症傾向が-2から0付近（特に-1.5付近と-0.5付近にピークがある）の人々を最も正確に測定できることが分かる。 - つまり、神経症傾向が平均よりやや低いレベルから平均レベルの人々の能力を識別するのに最も適している。

情報量が低い領域：
- 神経症傾向が非常に低い（-4以下）または非常に高い（2以上）人々については、測定の精度が低いことが分かる。
- 特に、神経症傾向が非常に高い人々に対する情報量が極めて少ない。

曲線の形状：
曲線は2つのピークを持つ非対称な形状を示している。これは、項目によって測定の精度が高い潜在特性レベルが異なることを示唆している。特に、以前の分析で示唆されたA4の高い潜在特性レベルでの情報提供能力が、この形状に影響を与えている可能性がある。

まとめ

このテストは、神経症傾向が平均付近の人々を最も正確に測定できる一方で、神経症傾向が非常に高い、または非常に低い人々を正確に測定するには適していない可能性が高い。
特に、神経症傾向が非常に高い人々に対する情報量が極めて少ないため、高いレベルの神経症傾向を識別するためには、項目の追加や修正が必要かもしれない。
2つのピークは、項目によって測定精度の高いレベルが異なることを示唆している。各項目の項目情報関数を確認することで（plot(mod_grm, type = "infotrace", which.items = 1)のように項目を指定）、どの項目がどのレベルで高い情報を提供しているかを把握できる。

項目情報関数

plot(mod_grm, type = "infotrace", which.items = 1)

横軸 (θ)：神経症傾向のレベルを表す。右に行くほど神経症傾向が高い。
縦軸 (I(θ))：項目情報量を表す。この値が大きいほど、その潜在特性レベルにおけるN1の測定精度が高いことを意味する。
曲線：神経症傾向のレベルごとに、N1が提供する情報量の大きさを示している。

この項目情報関数から読み取れるN1の情報

最も情報量が高い領域：

N1は、神経症傾向が-2から0付近（特に-1付近）の人々について最も多くの情報を提供することが分かる。
つまり、N1は、神経症傾向が平均よりやや低いレベルから平均レベルの人々を識別するのに最も適した項目であると言える。

情報量が低い領域：

神経症傾向が非常に低い（-4以下）または非常に高い（2以上）人々については、N1はあまり情報を提供していない。
特に、神経症傾向が高い人々に対する情報量が非常に少ないことが見て取れる。

曲線の形状：

曲線は、左側に偏った非対称な形状をしている。これは、N1が、神経症傾向が低いレベルでより高い識別力を持ち、高いレベルでは識別力が低いことを示している。

考察

N1は、神経症傾向が平均付近の人々を識別するのに適した項目である一方、神経症傾向が非常に高い、または非常に低い人々を識別するには適していない可能性が高い。
特に、神経症傾向が高い人々に対する情報量が非常に少ないため、このレベルの神経症傾向を正確に測定するためには、N1以外の項目や、N1の修正が必要かもしれない。
これまでの分析で示唆されたN1の逆転項目の可能性を踏まえると、この項目情報関数の形状は、N1が他の項目とは異なる特性を持つことを裏付けている。

N2以降の項目情報関数の確認

plot(mod_grm, type = "infotrace", which.items = 2) # N2の項目情報関数
plot(mod_grm, type = "infotrace", which.items = 3) # N3の項目情報関数
plot(mod_grm, type = "infotrace", which.items = 4) # N4の項目情報関数
plot(mod_grm, type = "infotrace", which.items = 5) # N5の項目情報関数

残差の確認

residuals()関数を用いて、モデルの残差を確認し、局所依存性などの問題がないかを確認する必要がある。

residuals(mod_grm) # 項目間の残差の相関行列
residuals(mod_grm, type = "LD") # 局所依存性の検定統計量

出力。

LD matrix (lower triangle) and standardized residual correlations (upper triangle)

Upper triangle summary:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 -0.114  -0.083   0.006   0.004   0.084   0.117

        N1      N2      N3     N4     N5
N1          -0.114   0.088  0.071  0.112
N2 178.444          -0.091  0.064 -0.093
N3 106.804 114.158         -0.057  0.117
N4  69.523  56.609  44.637        -0.053
N5 174.981 119.044 188.730 38.707

この出力は、項目間の局所依存性、つまり、特定の項目ペアの間に想定される1因子GRMモデルでは説明できない関連性があるかどうかを示している。結果は、行列形式で表示され、下三角部分が行列形式で出力され、上三角部分に標準化残差の相関 が示されている。

行列の構造

下三角部分 (LD matrix): mirtパッケージ独自の局所依存性指標。
上三角部分 (standardized residual correlations): 標準化残差の相関。通常の残差の相関とほぼ同等と考えて問題ない。
対角成分: 値なし（-または空白）。

解釈

局所依存性指標 (下三角部分):
- N1とA2間（178.444）、N1とA5間（174.981）、A2とA5間（119.044）、A3とA5間（188.730）で、比較的大きな値を示している。
- これらの値は、N1, A2, A3, A5の間に、想定した1因子GRMモデルでは説明できない関連性（局所依存性）が存在する可能性を示唆している。
- 特に、A5は他の3項目との間で局所依存性指標が大きく、A5が何らかの特異な項目である、または、別の潜在特性を測定している可能性 が考えられる。
標準化残差の相関 (上三角部分):
- 全体的に絶対値が0.1を超えるものは少なく、最も大きいものでN1とA2の間の-0.114、次にA3とA5の間の0.117、N1とA5の間の0.112である。
- これらの値は、項目間に強い関連性がないことを示唆している。

考察

下三角部分の局所依存性指標と、上三角部分の標準化残差の相関とで、結果に乖離が見られる。これは、mirtの局所依存性指標が、通常の残差の相関とは異なる計算方法に基づいているためと考えられる。
下三角部分の局所依存性指標の方が、N1、A2、A3、A5の間に潜在的な問題を示唆しているため、こちらを重視すべき である。特に、A5は他の3項目との間で局所依存性指標が大きいため、注意が必要だ。
先に述べたように、真ん中の値が選択されやすいことも、この数値に影響を与えている可能性はある。

対処

N1の扱い:
- N1が逆転項目である可能性を改めて検討し、得点を反転させた場合の影響を確認する。
- N1をモデルから除外することも検討する。
A5の検討:
- A5が他の項目と異なる特性を持つ可能性（例えば、神経症傾向とは別の側面を測定している、A4とA5はともに「不安」に関する項目であるため、他の項目と比較して似たような回答がされやすく、そのために局所独立性が満たされていないなど）を検討する。
- A5をモデルから除外することも検討する。

経験的信頼性(empirical reliability)

empirical_rxx()で計算できるはずだが、このデータでは無理だったので、手動で計算した。

# 潜在特性スコアと標準誤差の取得
fscores_result <- fscores(mod_grm, full.scores = TRUE, full.scores.SE = TRUE)

# スコアと標準誤差を確認
head(fscores_result)

# 潜在特性スコアと標準誤差を分離
Theta_scores <- fscores_result[, "F1"]        # 潜在特性スコア
Theta_SE <- fscores_result[, "SE_F1"]        # 標準誤差

# 経験的信頼性を計算
SE_mean_square <- mean(Theta_SE^2)           # 標準誤差の二乗の平均
Theta_variance <- var(Theta_scores)          # 潜在特性スコアの分散
empirical_rel <- 1 - (SE_mean_square / Theta_variance)

# 結果の表示
print(empirical_rel)

経験的信頼性: 0.7157662

信頼性指標

alpha(BFI_N)

Some items ( N1 ) were negatively correlated with the first principal component and 
probably should be reversed.
To do this, run the function again with the 'check.keys=TRUE' option
Reliability analysis
Call: alpha(x = BFI_N)

  raw_alpha std.alpha G6(smc) average_r  S/N   ase mean   sd median_r
      0.43      0.46    0.53      0.14 0.85 0.016  4.2 0.74     0.32

    95% confidence boundaries
         lower alpha upper
Feldt      0.4  0.43  0.46
Duhachek   0.4  0.43  0.46

 Reliability if an item is dropped:
   raw_alpha std.alpha G6(smc) average_r  S/N alpha se  var.r med.r
N1      0.72      0.72    0.67     0.397 2.63   0.0087 0.0065 0.375
N2      0.28      0.30    0.39     0.096 0.43   0.0219 0.1095 0.081
N3      0.18      0.21    0.31     0.061 0.26   0.0249 0.1014 0.081
N4      0.25      0.30    0.44     0.099 0.44   0.0229 0.1604 0.104
N5      0.21      0.24    0.36     0.071 0.31   0.0238 0.1309 0.094

 Item statistics
      n raw.r std.r r.cor r.drop mean  sd
N1 2784 0.066 0.024 -0.39  -0.31  2.4 1.4
N2 2773 0.630 0.665  0.58   0.37  4.8 1.2
N3 2774 0.724 0.742  0.71   0.48  4.6 1.3
N4 2781 0.686 0.661  0.50   0.37  4.7 1.5
N5 2784 0.700 0.719  0.64   0.45  4.6 1.3

Non missing response frequency for each item
      1    2    3    4    5    6 miss
N1 0.33 0.29 0.14 0.12 0.08 0.03 0.01
N2 0.02 0.05 0.05 0.20 0.37 0.31 0.01
N3 0.03 0.06 0.07 0.20 0.36 0.27 0.01
N4 0.05 0.08 0.07 0.16 0.24 0.41 0.01
N5 0.02 0.07 0.09 0.22 0.35 0.25 0.01
警告メッセージ:
alpha(BFI_N) で:
  Some items were negatively correlated with the first principal component and probably
should be reversed.
To do this, run the function again with the 'check.keys=TRUE' option

alpha()関数を使用した際の結果から、次のような問題点が指摘されています：

N1項目が負の相関を持つ:

項目N1が、最初の主成分と負の相関を持っている。
これは、N1がスコア全体と反対の方向に寄与している可能性があることを示している。

全体の信頼性（α係数）が低い:

Cronbachのα係数は0.43（標準化後は0.46）で、一般的に望ましい水準（0.7以上）を下回っている。。