Item Response Theory of Stage Theory
A simple sample script in Stan
Yasuharu Okamoto
Item response theory (IRT) usually assumes that the ability is continuous variable, but in psychology, it is not rare to assume that ability is represented by stages. A famous stage theory is one by Piaget. Okamoto (2013) presented a model of stage theory by extension of item response theory. In this website, we present a simple script in Stan of stage theory based on Okamoto (2013)fs model. For general introduction, see Levy and Mislevy (2016).
We assume three stages, stage-0, stage-1,
and stage-2, for simplicity. Stages cannot be observed, so may be called latent
classes. Probability of correct answer on an item is dependent on the stage in
which the person is, and increases from stage-0, stage-1 to stage-2. In Okamoto (2013), this increase in
probability of correct response is represented using inverse logit function
(cumulative logistic function). Also in the following script, increasing
probabilities are represented by inverse logistic function, inv_logit function, with stacking
simplex variables rescaled. Running this script, we get estimations of the
probabilities of a person in stages, 
, and the probabilities of correct response of a person in a stage, 
. Stan cannot use categorical parameters and stages are categorical
variables. Hence, the script in the following uses marginalization to treat
categorical parameters, that is,
A sample script for stage theory is given as follows:
"""
Yasuharu
Okamoto, 2018.11, 2019.02
"""
import pystan
import pickle
import matplotlib.pyplot as plt
import seaborn as sb
import numpy as np
X = ([[1]*5 + [1]*5] * 7) + ([[2]*5 +
[1]*5] * 7) + \
([[2]*5 + [2]*5] * 7)
print('X =')
for v in X:
print(v)
N = len(X)
print('N = ', N)
M = len(X[0])
print('M = ', M)
NStg = 3
print('NStg = ', NStg)
stan_data = {'N': N, 'M': M, 'NStg':
NStg, 'X': X}
stan_code = """
//
// Binary Responses
for Three Latent Classes
//
Y. Okamoto
//
data {
int N;
// Number of subjects
int M;
// Number of items
int NStg; // Number of stages
int X[N, M]; // Responses: 1 -> failure, 2 -> success
}
parameters {
simplex[NStg+1]
Cintvl[M]; // Intervals of criterions
simplex[NStg] theta[N];
}
transformed parameters {
vector[NStg] C[M];
real p[M, NStg];
real p_person[N, M];
vector[2] p_person_cat[N,
M];
for (i in 1:M){
for
(stg in 1:NStg) {
if ( stg == 1 ) {
C[i][1] = -5 + Cintvl[i][1] * 10.0;
}
else {
C[i][stg] = C[i][stg - 1] + Cintvl[i][stg] * 10.0;
};
}
}
for (i in 1:M) {
for
(stg in 1:NStg) {
p[i][stg] = inv_logit(C[i][stg]);
}
}
for (j in 1:N) {
for
(i in 1:M){
p_person[j][i] = 0.0;
for (stg in 1:NStg){
p_person[j][i] = p_person[j][i] + p[i][stg] * theta[j][stg];
}
}
}
for (j in 1:N) {
for
(i in 1:M) {
p_person_cat[j][i][2] = p_person[j][i];
p_person_cat[j][i][1] = 1.0 - p_person_cat[j][i][2];
}
}
}
model {
for (j in 1:N) {
for
(i in 1:M){
X[j][i] ~ categorical(p_person_cat[j][i]);
}
}
}
"""
sm = pystan.StanModel(model_code =
stan_code)
with open('model.pkl', 'wb') as g:
pickle.dump(sm, g)
fit = sm.sampling(data = stan_data,
n_jobs = 1)
p_samples = fit['p']
p_item_stage = np.empty((M, NStg))
for i in range(M):
for s in range(NStg):
p_item_stage[i][s] = p_samples[:,i,s].mean()
print('\n\nP(correct|item,stage)\n')
print('{0:>7s}{1:>9s}{2:>9s}{3:>9s}'.format('
', 'Stage-0', 'Stage-1', 'Stage-2'))
for i in range(M):
print('Item-{0} '.format(i),
end = '')
for s in range(NStg):
print('{0:>9.5f}'.format(p_item_stage[i][s]), end = '')
print(' ')
theta = fit['theta']
p_stage = np.empty((N, NStg))
for j in range(N):
for s in range(NStg):
p_stage[j][s] = theta[:,j,s].mean()
print('\n\nP(stage|person)\n')
print('{0:>10s}{1:>9s}{2:>9s}{3:>9s}'.format('
', 'stage-0', 'stage-1', 'stage-2'))
for j in range(N):
print('Person-{0:
<3d}'.format(j), end = '')
for s in range(NStg):
print('{0:>9.5f}'.format(p_stage[j][s]), end = '')
print(' ')
The data is set by this code of the above script:
X = ([[1]*5 + [1]*5] * 7) + ([[2]*5 +
[1]*5] * 7) + \
([[2]*5 + [2]*5] * 7)
Run the script, we get the following:
X =
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 1, 1, 1, 1, 1]
[2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
[2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
[2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
[2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
[2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
[2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
[2, 2, 2, 2, 2, 2, 2, 2, 2, 2]
The above output (X =) shows the data generated by the script. Numeral 1 indicates error response, and 2 correct one.
The data consists of 
persons, the first 7
person is supposed to be in stage-0, the next 7 person in stage-1, and the rest
in stage-2.
The number of items is 10, the first 5 items is easy, so stage-1 or stage-2 person can answer correctly, but stage-0 person cannot answer correctly on all items. The last 5 items can be answered correctly only by a person in stage-2.
In the Stan script, probabilities 
are represented by
parameter theta , probabilities 
by parameter p , and probabilities
by parameter p_person .
Estimation of parameters by Bayesian analysis with the Stan script gives these results:
Probabilties of correct responses on items
by a person in each stage are
P(correct|item,stage)
Stage-0 Stage-1 Stage-2
Item-0 0.06546 0.89806 0.97390
Item-1 0.06758 0.89594 0.97358
Item-2 0.06639 0.89750 0.97394
Item-3 0.06653 0.89628 0.97293
Item-4 0.06806 0.89767 0.97323
Item-5 0.02634 0.10168 0.93298
Item-6 0.02651 0.10224 0.93154
Item-7 0.02668 0.10318 0.93133
Item-8 0.02670 0.10106 0.93228
Item-9 0.02697 0.10188 0.93122
The probabilities of stage-0 person are nearly 0, probabilities of stage-1 person is high on item-0 to item-4, but low on item-5 to item-9. Persons in stage-2 can answer correctly for all items.
Mean probabilities that a person is in each of the stages are shown as follows:
P(stage|person)
stage-0 stage-1 stage-2
Person-0 0.77684 0.14212 0.08104
Person-1 0.77793 0.14119 0.08088
Person-2 0.77805 0.14182 0.08013
Person-3 0.77776 0.14293 0.07930
Person-4 0.77680 0.14301 0.08019
Person-5 0.77623 0.14279 0.08098
Person-6 0.77871 0.14103 0.08026
Person-7 0.16397 0.66866 0.16738
Person-8 0.16714 0.66856 0.16431
Person-9 0.16475 0.66813 0.16712
Person-10 0.16460 0.67035 0.16505
Person-11 0.16773 0.66397 0.16830
Person-12 0.16798 0.66553 0.16648
Person-13 0.16112 0.67423 0.16465
Person-14 0.08208 0.14140 0.77652
Person-15 0.08016 0.13770 0.78214
Person-16 0.08116 0.14124 0.77759
Person-17 0.08011 0.13994 0.77995
Person-18 0.07980 0.13745 0.78274
Person-19 0.08032 0.14088 0.77880
Person-20 0.08126 0.13986 0.77888
The probabilities for each person are adequately estimated. That is, person-0 to person-6 are estimated more probably to be in stage-0, person-7 to person-13 in stage-1, and person-14 to person-20 in stage-2.
The results in this website suggests that stage theory can be easily applied to real data, and statistically analyzed using Bayesian method
Visit this website for a sample script in PyMC3.
Reference
Levy, R. &
Mislevy, R. J. (2016). Bayesian
psychometric modeling. CRC Press
Okamoto, Y gItem Response Theory with Latent Classesh, Faculty of Integrated Arts and Social Sciences journal, 24, 103-107, 2013. Japan Womenfs University. (Click on http://id.nii.ac.jp/1133/00001684/)