Yasuharu Okamoto, 2017.7
Analysis of Polytomous Items by Item Response Theory
A Stan script with Python
The well-known model for polytomous items is the graded response model (Samejima, 1969). For binary items, visit this website.
Let the number of categories of items be K, where categories are denoted by integers 1 to K. Denote a response (category) on item by , then
Set probabilities of as follows:
where
Parameters s and s appear in equations in the form , and for any value , we have
That is, the model is not identifiable. To make the model identifiable, set
For the above model, the following Stan script is prepared.
========================= irt_pol.sta ===============================
data {
int <lower = 3>
K;
//
Number of categories
int Npsn;
// Number of
persons
int Nitm;
//
Number of items
int Ntot;
//
Number of data, Ntot = Nprn
* Nitm
int<lower = 1, upper =
Npsn> IDpsn[Ntot]; // Person ID
int<lower = 1, upper =
Nitm> IDitm[Ntot]; // Item ID
int<lower = 1, upper =
K> Res[Ntot];
//
Response-> an integer from 1 to K
}
parameters {
real<lower = 0.0> a[Nitm];
real
b[Nitm - 1];
ordered[K - 1] c;
real
F[Npsn];
}
transformed parameters {
vector[K] p[Ntot];
vector[K] shft[Nitm];
for (k in 1:(K-1)){
shft[1][k] = c[k];
for
(j in 2:(Nitm)){
shft[j][k] = shft[1][k] + b[j - 1];
}
}
for (i in 1:Ntot){
p[i][K] = inv_logit(1.7 * a[IDitm[i]] * (F[IDpsn[i]] - shft[IDitm[i]][K
- 1]));
p[i][1] = 1.0 - inv_logit(1.7 * a[IDitm[i]] * (F[IDpsn[i]] -
shft[IDitm[i]][1]));
for
(k in 2:(K - 1)){
p[i][k] = inv_logit(1.7 * a[IDitm[i]] * (F[IDpsn[i]] -
shft[IDitm[i]][k-1]))
-inv_logit(1.7 * a[IDitm[i]] * (F[IDpsn[i]] - shft[IDitm[i]][k]));
}
}
}
model {
for (i in 1:Npsn)
F[i]
~ normal(0.0, 1.0);
for (i in 1:Ntot)
Res[i] ~ categorical(p[i]);
}
=================================================================
A Python script using the above Stan script is this. The script files and an example data file were archived into a file PrgIRTpol.zip, which can be down loaded by clicking the name PrgIRTpol.zip.
=========================================================
import pystan
from pystan import
StanModel
import pickle
def calc_mean_med( d
):
d.sort()
med = d[len(d) // 2]
sum = 0
for v in d:
sum
+= v
mean = sum / len(d)
return mean, med
#
# Prepare
the input data file
#
fn_in =
input("Input data file = ")
f_in = open(fn_in,
"r")
#
# Set the
contents of the input data file in the object data_f
#
data_f =
f_in.readlines()
f_in.close()
#
# Prepare
the output file
#
fn_out =
input("Output data file = ")
f_out = open(fn_out,
"w")
f_out.write("Input
data file = " + fn_in + "\n\n")
pos = 0
while True:
if len(data_f[pos]) > 0:
if
data_f[pos][0] == "/":
break
pos += 1
pos += 1
N = int(data_f[pos])
pos += 1
M = int(data_f[pos])
pos += 1
K = int(data_f[pos])
print('N = ', N,
' M = ', M, ' K = ', K)
X = []
for i in range(N):
X.append([])
pos += 1
temp_strs =
data_f[pos].split()
for j in range(M):
X[i].append(int(temp_strs[j + 1]))
IDpsn = []
IDitm = []
Res = []
for i in range(N):
for j in range(M):
IDpsn.append(i + 1)
IDitm.append(j + 1)
Res.append(X[i][j])
f_out.write('\nData...\n')
for i in range(N):
f_out.write('{0:>5}:
'.format(i+1))
for j in range(M):
f_out.write(' {0}'.format(X[i][j]))
f_out.write('\n')
Data = {'K': K,
'Npsn': N, 'Nitm': M, 'Ntot': N * M, 'IDpsn': IDpsn, 'IDitm': IDitm, 'Res':
Res}
fit =
pystan.stan(file = 'irt_pol.stan', data = Data, seed = 999,
pars = ['a', 'b', 'c'], n_jobs = 1)
print(fit)
a = fit['a']
b = fit['b']
c = fit['c']
a_smpls = []
b_smpls = []
c_smpls = []
for i in range(M):
a_smpls.append([])
for i in range(M-1):
b_smpls.append([])
for k in range(K-1):
c_smpls.append([])
for v in a:
for i in range(M):
a_smpls[i].append(v[i])
for v in b:
for i in range(M-1):
b_smpls[i].append(v[i])
for v in c:
for i in range(K-1):
c_smpls[i].append(v[i])
print('\nlen(a) = ',
len(a), ' len(a_smpls) = ',
len(a_smpls),
' len(a_smpls[0] = ',
len(a_smpls[0]))
print('\nlen(b) = ',
len(b), ' len(b_smpls) = ',
len(b_smpls),
' len(b_smpls[0] = ',
len(b_smpls[0]))
print('\nlen(c) = ',
len(c), ' len(c_smpls) = ',
len(c_smpls),
' len(c_smpls[0] = ',
len(c_smpls[0]))
f_out.write('\n')
for i in range(M):
mean_a, med_a =
calc_mean_med(a_smpls[i])
f_out.write('a[{0}] =
{1:<.5}(mean)
{2:<.5}(med)\n'.
format(i+1, mean_a, med_a))
f_out.write('\n')
f_out.write('b[1] =
0.0 (indentification
condition)\n')
for i in range(M-1):
mean_b, med_b =
calc_mean_med(b_smpls[i])
f_out.write('b[{0}] =
{1:<.5}(mean)
{2:<.5}(med)\n'.
format(i+2, mean_b, med_b))
f_out.write('\n')
for i in range(K-1):
mean_c, med_c =
calc_mean_med(c_smpls[i])
f_out.write('c[{0}] =
{1:<.5}(mean)
{2:<.5}(med)\n'.
format(i+1, mean_c, med_c))
f_out.close()
=======================================================
Run this Python script, then file names of input data and output data are asked (Figure 1).
Figure 1
A file name of the output data is any text file name. The input data file should be prepared in the following format shown in Figure 2.
・
・
・
Figure 2
On the line next
to the line with slash / at the head, the number of data (persons) is written (“200” in
Figure 2). On the next line to this line, the number of items is written (“10” in Figure 2), and on the following line,
the number of categories is written (“6” in Figure 2). After these 3 lines, the responses of persons are
written, one person’s data per line. Each line starts
with data ID, then responses are written with one or more than one spaces as
seperators.
After the output
file name is set, calculation starts. The numbers of persons, items, and
categories are displayed, then calculation by Stan starts (Figure 3).
Figure 3
Figure 4 shows that Stan took about 5 minutes to sample in a chain for this condition. Total time of sampling is the sum of these times for chains when the chains are run sequentially.
Figure 4
When the program ends, statistics by Stan is displayed (Figure 5).
Figure 5
Parts of the content of the output file are shown in Figures 6 and 7.
Figure 6
On the first line, the input data file name is printed out. Then, the read in data are printed out, so that you can check whether the intended data are correctly read in or not. After the read in data, the point estimates of parameters, mean and medians, are printed out (Figure 7).
Figure 7
As shown by Equation (1), any shifted values and for arbitrary constant can also be accepted. If we want to put the center category boundary at the origin 0, set and transform s and s as follows:
Reference
Samejima, F. (1969). Estimation of latent ability using a response pattern of graded scores. Psychometric Monograph, No. 17, Richmond: Psychometric Society.