試験の数理 その2（項目反応理論の数理モデル）

[試験の数理 その1（問題設定とデータの生成）](https://qiita.com/takuyakubo/items/43d56725952e67032b49) の続きです。

前回は、「問題設定とデータの生成について」でした。
今回は、「項目反応理論で使われる数理モデルについて」です。

用いた環境は、

* python 3.8
* numpy 1.19.2
* matplotlib 3.3.1

です。

## 項目特性曲線

前回書いたようにここでの問題意識は、試験の結果が与えられた時に**受験者の能力と問題の難易度を推定**することでした。これを推定するために、ある特定の問題に注目して、その問題をどれぐらいの能力の受験者が正答することができるのかをグラフにしてみます。極端な場合には次のようなグラフが得られます。

```python
import numpy as np
import matplotlib.pyplot as plt
x = np.linspace(-4, 4, 41)
y = x > 1.3
plt.step(x, y)
plt.xlabel("leaner's skill")
plt.ylabel("probability which learner can answer correctly")
plt.show()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/30977/72a2ac23-c2e1-c481-39c0-009a6d819158.png)

横軸は受験者の能力を表す数値、縦軸はその能力の学習者がその問題を正答する確率です。このような、グラフのことを**項目特性曲線**といいます。
この例は、学習者の能力が1.3を超えた時に常に正答し、そうでない時に常に誤答するというものになります。このような問題があれば、この問題を試験に入れることで、学習者の能力が1.3を超えているかどうかを測定することができます。

ここで、学習者の能力の数値とは何？という疑問が当然浮かぶことと思います。結論を言ってしまうと、この数値の絶対値に意味はありません。ただし、相対的に２つの問題があった時、例えば上を問1として、もう一つ問2が
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/30977/bc965f90-28cb-b137-aecb-16cd15b5adea.png)
という項目特性曲線を持つとしましょう。このとき、問1は問2より難しいと判断することができます。実際、

| |受験者 1|受験者 2|受験者 3|
|:---:|:---:|:---:|:---:|
|問 1|誤|誤|正|
|問 2|誤|正|正|

という結果が得られた時、能力の値は、受験者1 < 受験者2 < 受験者3となるでしょう。また、この状況に置いて、「問1(正) かつ　問2(誤)」は起こりません。

さて、この状況は正答できるか誤答するかの２者が決定的になっており、少し極端です。実際の統計処理を考えると、正答誤答はある程度確率的であることが予想されます。特に、正誤の境あたりの能力の受験者は問題によっては、ある程度の確率で正答することができるでしょう。そういう意味で、実際に項目特性曲線として使用されるのは、例えば、次のようなものになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/30977/c6fb770e-c2f7-e8ef-6b9d-b354d6190cd2.png)
性質の良い問題（難易度と正答率に逆転がおこらない問題）であれば、確率分布の累積密度関数が、モデルとしては良さそうです。
項目反応理論では、その中でも数理的にも扱いやすい関数として、ロジスティック分布をよく用います。ロジスティック分布を用いたモデルを**ロジスティックモデル**(logistic model)といい。問題毎のparameterの数に応じて、1 ~ 3 parameter logistic modelがよく知られています[^1]。
以下では、この1 ~ 3 parameter logistic modelについて述べます。

## 1 parameter logistic model (1PL model, Rasch model)

1 parameter logistic model は項目特性曲線としては扱われている中で最も簡単なもののうちの1つで、以下の式で表せます。

```math
\Pr\{u_{ij} = 1|\theta, a, b\} = \frac{1}{1 + \exp(-a(\theta_j - b_i))}
```
ここで、前回の記事と同様に問題に関するparameterは$i$で、受験者に関するparameterは$j$で添字づけます。$u_{ij}$は問題$i$を受験者$j$が正答できるかどうかを示す確率変数です。問題の総数を$I$、受験者の総数を$J$とすると、このモデルのparamterの数は問題に関わるものが$I$個($=b_i$)、受験者に関わるものが$J$個($=\theta_j$)、全体に関わるものが1個($=a$)の$I+J+1$個あります。
項目特性曲線を描くと次のようになります。

```python
a = 3
def L1P(b, x, a=a):
    return 1 / (1 + np.exp(-  a * (x - b)))

x = np.linspace(-4, 4, 41)
for b in np.linspace(-2, 2, 5):
    y = partial(L1P, b)(x)
    plt.plot(x, y, label=f"{a=}, {b=}")
plt.xlabel("leaner's skill")
plt.ylabel("probability which learner can answer correctly")
plt.legend()
plt.show()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/30977/06089fa2-b468-9967-d337-0923bbe45ed6.png)![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/30977/2ff95113-0e35-e7e5-b1fa-f1dec95e8f88.png)
見ての通り、傾きが一斉に変わるのが特徴です。つまり、問題毎の難易度の識別のしやすさが同じとなります。なお、このように$a$は識別のしやすさに関わる量なので、**識別力**と呼ばれます。識別力の範囲は正の実数で、識別力が大きい方が識別しやすいことがわかります。また、$b$は問題の難かしさを表すことから、**困難度**と呼ばれます。困難度のとりうる範囲は実数全体[^2]で、困難度が高いほど問題は難しいと判定されます。
このモデルは1960年代初頭にデンマークの数学者Raschが研究していたことにより、Raschモデルとも呼ばれます。

## 2 parameter logistic model (2PL model)

2 parameter logistic model は標準的なモデルであり、pythonのpackageである[pyirt](https://github.com/junchenfeng/pyirt)に搭載されている唯一の[^3]モデルです。Raschモデルとほとんど同じ数式であり、次のように表せます。

```math
\Pr\{u_{ij} = 1|\theta, a, b\} = \frac{1}{1 + \exp(-a_i(\theta_j - b_i))}
```

つまり、問題によらなかった識別力$a$が問題に依存するようになった($a \rightarrow a_j$)ということです。このモデルのparamterの数は問題に関わるものが$2I$個($=a_i, b_i$)、受験者に関わるものが$J$個($=\theta_j$)の$2I+J$個あります。
項目特性曲線を描くと次のようになります。

```python
def L2P(a, b, x):
    return 1 / (1 + np.exp(-  a * (x - b)))
x = np.linspace(-4, 4, 41)
for idx in range(5):
    a = 2 * (idx + 1) / 5
    b = -2.0 + idx
    y = partial(L2P, a, b)(x)
    plt.plot(x, y, label=f"{a=}, {b=}")
plt.xlabel("leaner's skill")
plt.ylabel("probability which learner can answer correctly")
plt.legend()
plt.show()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/30977/2adb7e04-cb15-6cd1-a007-28f12d07c6b3.png)
## 3 parameter logistic model (3PL model)
3 parameter logistic model は2PL modelに**当て推量**とよばれる量を加えたモデルです。TOEICなどの試験では問題が択一式（3択問題、4択問題など）となっています。この択一問題で受験者に正答する能力がなかった場合、受験者がとる行動とはなんでしょうか？それは、ランダムに選択することです。このようなときには、例えば4択問題であれば、最低でも25%の正答率が確保されることとなります。この25%の部分が当て推量です。数式で表すと、次のようになります。

```math
\Pr\{u_{ij} = 1|\theta, a, b, c\} = c_i + \frac{1 - c_i}{1 + \exp(-a_i(\theta_j - b_i))}
```
ここで、$c_i$が問題$i$における当て推量であり、とりうる値の範囲は$0 \leq c_i \leq 1$となります。このモデルのparamterの数は問題に関わるものが$3I$個($=a_i, b_i, c_i$)、受験者に関わるものが$J$個($=\theta_j$)の$3I+J$個あります。
項目特性曲線を描くと次のようになります。

```python
def L3P(a, b, c, x):
    return c + (1 - c) / (1 + np.exp(-  a * (x - b)))
x = np.linspace(-4, 4, 41)
for idx in range(5):
    a = 2 * (idx + 1) / 5
    b = -2.0 + idx
    c = (4 -  idx) / 10
    y = partial(L3P, a, b, c)(x)
    plt.plot(x, y, label=f"{a=}, {b=}, {c=}")
plt.xlabel("leaner's skill")
plt.ylabel("probability which learner can answer correctly")
plt.legend()
plt.show()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/30977/eed07786-c5a4-4f35-a2ba-a504dbddc3e5.png)

## 次回
3PL modelについて、parameter推定の方法について紹介します。
[試験の数理 その3（3PL modelの最適化）](https://qiita.com/takuyakubo/items/9ec1d3e3a1ae7f1f7594)

## 参考文献
* 書籍 [項目反応理論・理論編―テストの数理 (統計ライブラリー)(amazon.co.jp)](https://www.amazon.co.jp/dp/4254126697/ref=cm_sw_r_tw_dp_x_JuiyFb2VA01PX)
* 書籍 [項目反応理論 入門編 (第2版) (統計ライブラリー)(amazon.co.jp)](https://www.amazon.co.jp/dp/4254127952/ref=cm_sw_r_tw_dp_x_mrmyFb8XMJ9Z3)


[^1]: どれほど能力の高い受験者でも正答率が$d_i( < 1)$にしか届かないことを表現する4 parameter logistic modelもあるが、あまり見ないので割愛する。
[^2]: 数値計算など実際には適当な範囲に限定する。
[^3]: 2020年9月16日現在