Bayesian IRT Models: Unidimensional Binary Models

Overview

This SAS Web Example demonstrates how to fit one-, two-, three-, and four-parameter (1PL, 2PL, 3PL, and 4PL) unidimensional binary item response theory (IRT) models by using the MCMC procedure. The Analysis section briefly presents mathematical representations of the models. The Examples section presents one example for each of the four models. The discussion focuses on Bayesian model specification and PROC MCMC syntax. There is also sample SAS code for producing plots of item characteristic curves (ICC), item information curves (IIC), test characteristic curves (TCC), and test information curves (TIC).

Analysis

In unidimensional binary item response theory models, an instrument (test) consists of a number of items (questions) that require a binary response (for example, yes/no, true/false, or agree/disagree). The purpose of the instrument is to measure a single latent trait (or factor) of the subjects. The latent trait is assumed to be measurable and to have a range that encompasses the real line. An individual’s location within this range, θ, is assumed to be a continuous random variable. IRT models enable you to estimate the probability of a correct response for each item, to estimate the levels of the latent traits of the subjects, and to evaluate how well the items, individually and collectively, measure the latent trait.

In the models that follow, the probability that a subject will answer an item correctly (or affirmatively) is assumed to be a function of the subject’s location (θ) on the range of possible values of the latent trait:

p(x = 1 | θ ) = ƒ( θ )

It is also assumed that a subject’s response to an item is independent of the subject’s response to any other item, conditional on θ. The function ƒ( θ ) is known as an item response function (IRF) and is assumed to be a sigmoid, resembling a normal or logistic cumulative distribution. Choosing an IRT model amounts to fully specifying the functional form for ƒ( θ ). The models that are presented in this web example are all variations of a logistic function.

Note: The parameter estimates from a logistic model and a normal model that are fitted to the same data set differ by an approximately linear transformation; the parameter estimates from a normal model are approximately 1.7 times the estimates from a logistic model. You need to modify only one statement in a PROC MCMC model specification to switch between a normal model and a logistic IRT model.

One-Parameter Logistic (1PL) Model

In a 1PL model, the probability that a subject will respond correctly or affirmatively to an item is a function of θ_i and an item parameter, δ_j, which measures an item’s location on a latent difficulty scale. An item’s difficulty is defined as the level of θ_i that a random test subject must possess, such that the probability of a correct response is 0.5. The basic premise of the 1PL model is that the greater the distance between θ_i and δ_j, the greater the certainty in how a subject is expected to respond to an item. As the distance approaches zero, the more likely that there is a 50:50 chance that a subject will correctly respond to an item (De Ayala 2009). Mathematically, the 1PL model is written as

Screenshot 2023-12-15 at 11.15.24 AM.png

where i indexes subjects and j indexes items. The parameter α is known as the item discrimination parameter. It is related to the IRF’s slope and reflects how well an item discriminates among individuals located at different points along the continuum (De Ayala 2009). In the 1PL model,α is a singular parameter (it has no subscript). The fact that is common to all the items implies that all the IRFs are parallel. If you remove α from the model, implying that it has a value of 1, you get what is known as the Rasch model.

The item information function for the 1PL model is

Two-Parameter Logistic (2PL) Model

The 2PL model has the following form:

It is identical to the 1PL model, except that α now has a subscript, indicating that each item now has its own discrimination parameter. This means that the IRFs of the items are no longer constrained to be parallel. The cost of this flexibility is that if you have J items, you now have to estimate J – 1 additional parameters compared to the 1PL model.

The item information function for the 2PL model is

Three-Parameter Logistic (3PL) Model

In both the 1PL and 2PL models, the IRFs have a lower asymptote of 0, implying that as θ approaches –∞, the probability of a correct or affirmative response approaches 0. However, if a subject randomly selects a response (guesses), the probability of a correct or affirmative response would be greater than 0 regardless of the subject’s θ or the difficulty of the item. The 3PL model includes an additional parameter, 9, for each item, enabling you to estimate the lower asymptote. The 3PL model has the following form:

For some items, it might be that guessing behavior is unlikely. In such cases, you can easily constrain selected g_jparameters to be 0.

The item information function for the 3PL model is

Four-Parameter Logistic (4PL) Model

Although the 3PL model allows for the possibility that the lower asymptote is nonzero, the upper asymptote is still 1. This means that as θ approaches +∞, the probability of a correct or affirmative response approaches 1. However, subjects with a very large can still make clerical errors. The 4PL model includes a ceiling parameter, , for each item, enabling you to estimate the upper asymptote. The 4PL model has the following form:

The item information function for the 4PL model is

ICC, TIC, IIC, and TIC Plots

After you choose a particular model and sample the posterior distributions of the model parameters, you can produce plots of ƒ( θ ) versus θ for each item. Such plots are known as item characteristic curves. If you compute the sum of ƒ( θ ) over all the items, you get what is known as the test characteristic function. A test characteristic function provides the expected number of correct or affirmative answers as a function of θ. A plot of the test characteristic function versus θ is known as a test characteristic curve. You can also assess the informational value of the items, individually and collectively. You do this by first computing the individual item information functions, which are defined as the negative of the expected value of the second derivative of the log-likelihood function. A plot of an item information function versus θ is known as an item information curve. If you then compute the sum of all the item information functions, you get what is known as the test information function. A plot of the test information function versus θ is known as a test information curve.

Examples

There are several different ways to specify IRT models in PROC MCMC. The method that the following examples use requires the data to be in long form (also known as hierarchical form, panel data form, or longitudinal form). When the data are in long form, the rows of the data set are indexed by both person and item, and there is a single response variable. Having the data in this form enables you to specify most of the model parameters as random effects by using the RANDOM statement in PROC MCMC. Treating the parameters as random effects is possible because of the conditional independence assumption stated in the Analysis section. This approach simplifies model specification in PROC MCMC, and experimentation shows that using this approach can reduce execution times significantly. In the following examples, the reduction in execution times ranges from 18% for the 4PL model to 37% for the 1PL model compared to alternative model specifications.

Preparing the Data

The following examples use the IrtBinary data set, which is used in the Getting Started example in the IRT procedure chapter in the SAS/STAT User's Guide. There are 50 subjects, and each subject responds to 10 items. These 10 items have binary responses: 1 indicates correct and 0 indicates incorrect. The following DATA step reads the data and creates the variable Person, which indexes the rows. The data set is in wide form.

data IrtBinary;
   input item1-item10 @@;
   person=_N_;
   datalines;
1 0 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1
0 0 0 0 1 0 1 0 0 1 1 1 1 1 1 0 0 0 0 0
1 1 1 1 1 0 0 0 1 1 1 1 1 1 1 0 0 1 0 1
1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 1 1 1

   ... more lines ...

0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1
0 1 0 1 0 0 0 0 0 0 1 1 1 0 0 1 0 1 1 1
1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 0 1 1
;

The following SAS statements use PROC TRANSPOSE to reshape IrtBinary into long form. PROC TRANSPOSE saves the output data set IrtBinaryRE. IrtBinaryRE contains three variables: Person, which indexes the subjects; Item, which indexes the items; and Response, which records whether the subjects’ responses to each item are correct. IrtBinaryRE contains 10 rows of data for each person; IrtBinary has only 1 row per person.

proc transpose data=IrtBinary out=IrtBinaryRE(rename=(col1=response _NAME_=item));
  by person;
run;

proc print data=IrtBinaryRE(obs=10);
run;

Output 1 shows the first 10 observations of the output data set IrtBinaryRE.

Output 1: First 10 Observations of IrtBinaryRE

Obs	person	item	response
1	1	item1	1
2	1	item2	0
3	1	item3	1
4	1	item4	1
5	1	item5	1
6	1	item6	1
7	1	item7	1
8	1	item8	1
9	1	item9	1
10	1	item10	0

Example 1: 1PL Model

Recall that the probability of a correct or affirmative response in the 1PL model is

To perform Bayesian estimation of this model, you must specify prior distributions for θ_i, δ_{j ,}and α. The prior for θ_i is usually specified as a standard normal distribution; δ_jis also often assigned a normal prior. However, unless you actually have prior information that you want to incorporate into the model, it is common practice to assign a large variance to the prior for δ_j so that it is an uninformative or diffuse prior. The theoretical range of α encompasses the real line, so you might consider assigning a diffuse normal distribution as the prior for α. However, because of the indeterminacy of metric that is inherent in IRT models (De Ayala 2009, p. 41), you sometimes get a negative estimate for α, which would be interpreted to mean that subjects with lower θ_is have a higher probability of a correct response than subjects with higher θ_is. As a result, it is common practice for modelers to assign prior distributions such as a lognormal or a truncated normal distribution that restrict α to be nonnegative.

You can also treat the parameters of the prior distributions for θ_i, δ_{j ,}or α as hyperparameters and assign hyperprior distributions. This example assigns a standard normal distribution to θ_i and a normal distribution with a large variance to δ_j. Both θ_i and δ_jare treated as random effects, and their prior distributions are specified by using RANDOM statements in PROC MCMC; is assigned a truncated normal prior with a mean of 1 (the value implied by the Rasch model), a large variance, and a lower truncation boundary of 0.

The following statements specify the 1PL model by using PROC MCMC:

ods graphics on;
proc mcmc data=IrtBinaryRE nmc=20000 outpost=out1p  seed=1000 nthreads=-1 ;
  parms alpha;
  prior alpha ~normal(1, var=10, lower=0);
  random theta~normal(0, var=1) subject=person namesuffix=position nooutpost;
  random delta~normal(0, var=10) subject=item monitor=(delta) namesuffix=position;
  p=logistic(alpha*(theta-delta));
  model response ~ binary(p);
run;

The OUTPOST= option in the PROC MCMC statement saves the MCMC samples in a data set named Out1p. The NTHREADS=–1 option sets the number of available threads to the number of hyperthreaded cores available on your system.

The PARMS statement declares the parameter α, and the PRIOR statement specifies the prior distribution for α as a truncated normal distribution with a mean of 1, a variance of 10, and a lower truncation boundary of 0.

The first RANDOM statement specifies the prior distribution for θ_i as a standard normal distribution. The SUBJECT= option specifies that the variable Person identifies the subjects. The NOOUTPOST option suppresses the output of the posterior samples of the θ_i, random-effects parameters to the OUTPOST= data set; this reduces the execution time. However, if you want to perform analysis on the posterior samples of θ_i,, you can omit this option. The NAMESUFFIX= option specifies how the names of the random-effects parameters are internally created from the SUBJECT= variable; when you specify NAMESUFFIX=POSITION, the MCMC procedure constructs the parameter names by appending the numbers 1, 2, 3, and so on, where the number indicates the order in which the SUBJECT= variable appears in the data set.

The second RANDOM statement specifies the prior distribution for δ_j as a normal distribution with a mean of 0 and a variance of 10. The SUBJECT= option specifies that the variable Item identifies the subjects. The MONITOR= option requests that the random-effects parameters for δ_j be monitored and reported in the output.

The following programming statement generates the variable P, which contains the probability of a correct response. If you prefer a normal model instead of a logistic model, just replace the LOGISTIC function with the PROBNORM function. That is, specify the following:

p=probnorm(alpha*(theta-delta));

The MODEL statement specifies the conditional distribution of the data, given the parameters (the likelihood function). In this example, it specifies that Response has a binary (Bernoulli) distribution with a probability of success equal to P.

Output 2 shows the posterior summaries and intervals table for the 1PL model. Difficulty parameters measure the difficulties of the items. As the value of the difficulty parameter increases, the item becomes more difficult. You can observe that all the difficulty parameters are less than 0, which suggests that all the items in this example are relatively easy.

Output 2: Posterior Summaries and Intervals

The MCMC Procedure

Posterior Summaries and Intervals
Parameter	N	Mean	Standard Deviation	95% HPD Interval
alpha	20000	1.2610	0.1550	0.9659	1.5725
delta_1	20000	-1.1553	0.2555	-1.6669	-0.6742
delta_2	20000	-1.3856	0.2818	-1.9369	-0.8512
delta_3	20000	-1.2137	0.2628	-1.7489	-0.7182
delta_4	20000	-1.1694	0.2604	-1.6824	-0.6653
delta_5	20000	-1.1555	0.2618	-1.6446	-0.6316
delta_6	20000	-0.4676	0.2223	-0.8955	-0.0275
delta_7	20000	-0.5179	0.2179	-0.9402	-0.0925
delta_8	20000	-0.4244	0.2222	-0.8802	-0.0187
delta_9	20000	-0.3768	0.2223	-0.8371	0.0450
delta_10	20000	-0.5646	0.2221	-0.9727	-0.1226

Producing ICC, IIC, TCC, and TIC Plots

To produce ICC, IIC, TCC, and TIC plots, you use the means of the posterior samples that are saved in the OUTPOST= data set Out1p to compute the following information: < /p>

the probability of a correct or affirmative response for each item over a range of values of θ (ICC)
the sum of the probabilities across all items (TCC)
the item information functions for each item over a range of values of θ (IIC)
the sum of the item information functions across all items (TIC)

The following statements create the data set Plots1PL, from which you can produce ICC, TCC, IIC, and TIC plots by using PROC SGPLOT:

data plots;
  set out1p;
  array d[10] delta_1-delta_10;
  array p[10] p1-p10;
  array i[10] i1-i10;
  do theta=-6 to 6 by .25;
    do j=1 to 10;
      p[j]=logistic(alpha*(theta-d[j]));
	  i[j]=alpha**2*p[j]*(1-p[j]);
	end;
	tcf=sum(of p1-p10);
	tif=sum(of i1-i10);
	output;
  end;
run;

proc sort data=plots;
  by theta;
run;

proc means data=plots noprint;
  var p1-p10 i1-i10 tcf tif;
  by theta;
  output out=plots1pl
    mean=
    p5=lb_p1-lb_p10 lb_i1-lb_i10 lb_tcf lb_tif
    p95=ub_p1-ub_p10 ub_i1-ub_i10 ub_tcf ub_tif;
run;

The following statements use the data set Plots1PL and PROC SGPLOT to produce a plot of the item characteristic curve and a 95% credible interval for the first item:

proc sgplot data=plots1pl;
  title "Item Characteristic Curve";
  title2 "item=1";
  band x=theta upper=ub_p1 lower=lb_p1 / modelname="icc" transparency=.5
       legendlabel="95% Credible Interval";
  series x=theta y=p1 / name="icc";
  refline .5 / axis=y;
  xaxis label="Trait ((*ESC*){unicode theta})";
  yaxis label="Probability";
run;
title;

Figure 1 displays a plot of the ICC and a 95% credible interval from the 1PL model for the first item. The curve shows the how the probability of a correct answer varies as the value of the latent trait varies. The difficult parameter for the first item is –1.1553; the ICC curve intersects the horizontal line at this point. This indicates that a person with a below average trait value of –1.1553 has a 50% chance of answering the question correctly.

Figure 1: Item Characteristic Curve

The following statements produce a plot of the test characteristic curve that shows you the expected number of correct responses as a function of θ:

proc sgplot data=plots1pl;
  title "Test Characteristic Curve";
  band x=theta upper=ub_tcf lower=lb_tcf / modelname="tcc" transparency=.5
    legendlabel="95% Credible Interval";
  series x=theta y=tcf / name="tcc";
  xaxis label="Trait ((*ESC*){unicode theta})";
  yaxis label="Expected Number of Correct Answers";
run;
title;

Figure 2 displays a plot of the TCC and a 95% credible interval from the 1PL model. The TCC curve shows how the expected number of correct answers varies as the value of the latent trait varies. The graph indicates that a person with an average trait value of 0 is expected to answer approximately eight of the ten questions correctly.

Figure 2: Test Characteristic Curve

The following statements produce a plot of the item information curve for the first item from the 1PL model:

proc sgplot data=plots1pl;
  title "Item Information Curve";
  title2 "item=1";
  band x=theta upper=ub_i1 lower=lb_i1 / modelname="iic" transparency=.5
       legendlabel="95% Credible Interval";
  series x=theta y=i1 / name="iic";
  xaxis label="Trait ((*ESC*){unicode theta})";
  yaxis label="Information";
run;
title;

Figure 3 displays a plot of the IIC and a 95% credible interval of the first item from the 1PL model. The curve shows that an item provides its maximum information at the point where the latent trait value equals the item’s difficulty.

Figure 3: Item Information Curve

The following statements produce a plot of the test information curve from the 1PL model:

proc sgplot data=plots1pl;
  title "Test Information Curve";
  band x=theta upper=ub_tif lower=lb_tif / modelname="tic" transparency=.5
       legendlabel="95% Credible Interval";
  series x=theta y=tif / name="tic";
  xaxis label="Trait ((*ESC*){unicode theta})";
  yaxis label="Information";
run;
title;

Figure 4 displays a plot of the TIC and a 95% credible interval from the 1PL model. The curve shows that the instrument provides its maximum information for estimating the latent trait at a value of approximately –1.

Figure 4: Test Information Curve

Example 2: 2PL Model

The 2PL model extends the 1PL model by relaxing the assumption that the IRFs are parallel. This is accomplished by introducing the discrimination parameter α_j for each item. The probability of a correct or affirmative response is

The following statements specify the 2PL model by using PROC MCMC:

proc mcmc data=IrtBinaryRE nmc=20000 outpost=out2pl  seed=1000 nthreads=-1;
  random theta~normal(0, var=1) subject=person namesuffix=position nooutpost;
  random delta~normal(0, var=10) subject=item monitor=(delta) namesuffix=position;
  random alpha~normal(1, var=10) subject=item monitor=(alpha) namesuffix=position;
  p=logistic(alpha*(theta-delta));
  model response ~ binary(p);
run;

The only change from the 1PL model is that α_j is now treated as a random effect, so you no longer have a PARMS statement and a PRIOR statement for α_j ; instead, you now have a RANDOM statement. In the new RANDOM statement, the prior distribution for α_j is specified as normal with a mean of 1 and a variance of 10. The SUBJECT= option specifies that the variable Item identifies the subjects, and the MONITOR= option requests that the random-effects parameters for α_j be monitored and reported in the output. The statement that generates the variable P is changed to reflect the 2PL model’s probability of a correct response.

Output 3 shows the posterior summaries and intervals table for the 2PL model.

Output 3: Posterior Summaries and Intervals

The MCMC Procedure

Posterior Summaries and Intervals
Parameter	N	Mean	Standard Deviation	95% HPD Interval
delta_1	20000	-0.8835	0.2244	-1.3201	-0.4620
delta_2	20000	-1.0329	0.2354	-1.4805	-0.5908
delta_3	20000	-0.9305	0.2397	-1.4207	-0.5080
delta_4	20000	-0.9486	0.2498	-1.4387	-0.5075
delta_5	20000	-1.1769	0.3760	-1.9535	-0.5719
delta_6	20000	-0.5351	0.3314	-1.1239	-0.00139
delta_7	20000	-0.7262	0.4301	-1.6151	-0.0580
delta_8	20000	-0.5977	0.4177	-1.3986	0.00567
delta_9	20000	-0.4594	0.3216	-1.0210	0.1094
delta_10	20000	-0.7162	0.4297	-1.5048	-0.0539
alpha_1	20000	2.3510	0.7812	0.9933	3.9078
alpha_2	20000	2.5123	0.8854	1.0055	4.1709
alpha_3	20000	2.3206	0.8125	0.9553	3.9300
alpha_4	20000	1.9653	0.6567	0.8032	3.2152
alpha_5	20000	1.3349	0.4434	0.5179	2.2018
alpha_6	20000	1.1622	0.3918	0.4383	1.9667
alpha_7	20000	0.9063	0.3398	0.3241	1.5864
alpha_8	20000	0.9158	0.3437	0.2909	1.6251
alpha_9	20000	1.1062	0.3871	0.4154	1.9247
alpha_10	20000	1.0332	0.3801	0.2493	1.7199

The estimates for α_j are quite variable, indicating that the parallel IRF assumption of the 1PL model is likely to be too restrictive.

The SAS program that produces the ICC, IIC, TCC, and TIC plots is essentially the same as for the 1PL model except for slight modifications to the code that creates the Plots2PL data set. The modifications are needed to accommodate the new α_j parameters, and the probability and information function formulas are updated for the 2PL model. Only the DATA step that creates the data set Plots2PL is shown here; everything else is unchanged and is not shown to save space. The complete SAS program is included (for all the models) in the downloadable sample SAS program file that accompanies this document.

data plots;
  set out2pl;
  array alpha[10] alpha_1-alpha_10;
  array d[10] delta_1-delta_10;
  array p[10] p1-p10;
  array i[10] i1-i10;
  do theta=-6 to 6 by .25;
    do j=1 to 10;
      p[j]=logistic(alpha[j]*(theta-d[j]));
      i[j]=alpha[j]**2*p[j]*(1-p[j]);
    end;
    tcf=sum(of p1-p10);
    tif=sum(of i1-i10);
    output;
  end;
run;

Example 3: 3PL Model

The 3PL model extends the 2PL model by adding the parameters g_j, which are used to estimate the lower asymptote of the IRF. This addresses the effect that guessing has on estimating the probability of a correct or affirmative response. The probability of a correct response for the 3PL model is

The following statements specify the 3PL model by using PROC MCMC:

proc mcmc data=IrtBinaryRE nmc=20000 outpost=out3pl  seed=1000 nthreads=-1;
  random theta~normal(0, var=1) subject=person namesuffix=position nooutpost;
  random delta~normal(0, var=10) subject=item monitor=(delta) namesuffix=position;
  random alpha~normal(1, var=10) subject=item monitor=(alpha) namesuffix=position;
  random g~beta(2,2)  subject=item monitor=(g) namesuffix=position;
  p=g + (1 - g)*logistic(alpha*(theta-delta));
  model response ~ binary(p);
run;

There are only two changes from the 2PL model. The first is a new RANDOM statement for the g_j parameters. The prior distribution for g_j is specified as a beta distribution with both shape parameters equal to 2. The SUBJECT= option specifies that the variable Item identifies the subjects, and the MONITOR= option requests that the random-effects parameters for g_j be monitored and reported in the output. The second change is in the statement that generates the variable P. It is changed to reflect the 3PL model’s probability of a correct response.

Output 4 shows the posterior summaries and intervals table for the 3PL model.

Output 4: Posterior Summaries and Intervals

The MCMC Procedure

Posterior Summaries and Intervals
Parameter	N	Mean	Standard Deviation	95% HPD Interval
delta_1	20000	-0.4897	0.2780	-1.0604	0.0423
delta_2	20000	-0.6071	0.2911	-1.1538	-0.0436
delta_3	20000	-0.6204	0.3250	-1.2540	-0.00025
delta_4	20000	-0.4028	0.3426	-1.0915	0.2306
delta_5	20000	-0.3040	0.4829	-1.2224	0.6641
delta_6	20000	0.2840	0.4623	-0.6075	1.1311
delta_7	20000	0.4890	0.3453	-0.1676	1.1299
delta_8	20000	0.4414	0.6911	-0.7109	1.6356
delta_9	20000	0.4442	0.5806	-0.5662	1.5225
delta_10	20000	0.3211	0.4154	-0.4698	1.0784
alpha_1	20000	3.3154	1.5304	1.0233	6.3126
alpha_2	20000	3.7251	1.6654	1.1349	7.2304
alpha_3	20000	2.4956	1.1757	0.7381	4.9297
alpha_4	20000	3.0525	1.5377	0.8125	6.2402
alpha_5	20000	2.3538	1.3245	0.5325	5.1481
alpha_6	20000	2.3686	1.4587	0.3157	5.3921
alpha_7	20000	4.0087	1.9609	0.7411	7.7422
alpha_8	20000	1.9140	1.3366	0.0539	4.5239
alpha_9	20000	2.1372	1.4044	0.1464	4.9903
alpha_10	20000	2.9880	1.7142	0.4812	6.4742
g_1	20000	0.2888	0.1242	0.0543	0.5155
g_2	20000	0.3214	0.1344	0.0691	0.5781
g_3	20000	0.2674	0.1283	0.0330	0.5079
g_4	20000	0.3358	0.1362	0.0767	0.5900
g_5	20000	0.3905	0.1451	0.1247	0.6752
g_6	20000	0.3231	0.1239	0.0793	0.5504
g_7	20000	0.4161	0.0915	0.2308	0.5895
g_8	20000	0.3361	0.1251	0.0960	0.5670
g_9	20000	0.3270	0.1256	0.0778	0.5533
g_10	20000	0.3747	0.1130	0.1502	0.5869

The strictly positive 95% HPD intervals for the g_j parameters indicate that the lower asymptote is not likely to be 0, as implied by the 1PL and 2PL models. The inclusion of g_j in the model has a significant effect on the estimates for the δ_j and α_jparameters.

To create the ICC, IIC, TCC, and TIC plots, you must modify the DATA step that creates the Plots3PL data set to accommodate the new parameters and to update the probability and information function formulas for the 3PL model. The modified statements are as follows:

data plots;
  set out3pl;
  array alpha[10] alpha_1-alpha_10;
  array d[10] delta_1-delta_10;
  array p[10] p1-p10;
  array i[10] i1-i10;
  array g[10] g_1-g_10;
  do theta=-6 to 6 by .25;
    do j=1 to 10;
      p[j]=g[j] + (1 - g[j])*logistic(alpha[j]*(theta-d[j]));
      i[j]=alpha[j]**2*((p[j]-g[j])**2/(1-g[j])**2)*((1-p[j])/p[j]);
    end;
    tcf=sum(of p1-p10);
    tif=sum(of i1-i10);
    output;
  end;
run;

Example 4: Four-Parameter Logistic IRT Model

The 4PL model extends the 3PL model by adding the parameters c_j, which are used to estimate the upper asymptote of the IRF. This addresses the effect that clerical errors have on estimating the probability of a correct or affirmative response. The probability of a correct response for the 4PL model is

The following statements specify the 4PL model by using PROC MCMC:

proc mcmc data=IrtBinaryRE nmc=20000 outpost=out4pl  seed=1000 nthreads=-1 ;
  random theta~normal(0, var=1) subject=person namesuffix=position nooutpost;
  random delta~normal(0, var=10) subject=item monitor=(delta) namesuffix=position;
  random alpha~normal(1, var=10) subject=item monitor=(alpha) namesuffix=position;
  random g~beta(2,2)  subject=item monitor=(g) namesuffix=position;
  llc=lpdfbeta(c,2,2,g);
  random c~general(llc)  subject=item monitor=(c) initial=.999 namesuffix=position;
  p=g + (c - g)*logistic(alpha*(theta-delta));
  model response ~ binary(p);
run;

There are three changes from the 3PL model. Like the g_j parameters of the 3PL model, the c_j parameters are assigned a beta prior distribution. However, you need to restrict the c_j parameters, which are the upper asymptotes of the IRF, to be greater than the corresponding g_j parameters, which are the lower asymptotes of the IRF. You can accomplish this by using a truncated beta distribution and specifying g_j as the lower boundary for the distribution of c_j. However, the RANDOM statement does not directly support the truncated beta distribution. Therefore, you must use an indirect method that is documented in the section "Using Density Functions in the Programming Statements" of the chapter "The MCMC Procedure" in the SAS/STAT User's Guide. The LPDFBETA function computes the density of a beta distribution, and it permits both upper and lower truncation boundaries. So you use a programming statement to create a variable Llc and set it equal to the LPDFBETA function with appropriate parameters as follows:

llc=lpdfbeta(c,2,2,g);

You then specify a RANDOM statement for c_j, specify the prior distribution as a general distribution, and specify the variable Llc as the general distribution’s argument as follows:

random c~general(llc)  subject=item monitor=(c) initial=.999 namesuffix=position;

When you use the general distribution in a RANDOM statement, PROC MCMC requires you to specify the INITIAL= option to provide initial values for the random-effects parameters.

Finally, you need to update the programming statement that specifies the formula for P to reflect the appropriate probability for the 4PL model.

Output 5 shows the posterior summaries and intervals table for the 4PL model.

Output 5: Posterior Summaries and Intervals

The MCMC Procedure

Posterior Summaries and Intervals
Parameter	N	Mean	Standard Deviation	95% HPD Interval
delta_1	20000	-0.6864	0.2934	-1.2772	-0.1235
delta_2	20000	-0.8145	0.3476	-1.4489	-0.1897
delta_3	20000	-0.8283	0.2863	-1.3502	-0.2612
delta_4	20000	-0.7009	0.3720	-1.4358	-0.0342
delta_5	20000	-0.6294	0.7838	-1.7101	0.5643
delta_6	20000	-0.2232	0.5148	-1.0763	0.6580
delta_7	20000	0.0288	0.8732	-1.6972	1.4492
delta_8	20000	-0.2553	0.8495	-1.5033	0.9210
delta_9	20000	-0.3324	0.8700	-1.7061	1.0800
delta_10	20000	-0.1409	0.7486	-1.4970	1.1418
alpha_1	20000	4.5986	1.8234	1.4773	8.1507
alpha_2	20000	4.5088	1.8829	1.4037	8.1701
alpha_3	20000	4.9369	1.8818	1.5300	8.4674
alpha_4	20000	4.2284	1.8461	1.1504	7.8987
alpha_5	20000	3.8816	2.0848	0.6875	8.2821
alpha_6	20000	3.8972	1.8955	0.7710	7.6917
alpha_7	20000	3.8045	2.1270	0.5872	8.2902
alpha_8	20000	3.7436	2.0594	0.5369	7.9271
alpha_9	20000	3.4237	1.9482	0.3702	7.2261
alpha_10	20000	3.7657	2.0476	0.4096	7.8540
g_1	20000	0.3007	0.1277	0.0564	0.5327
g_2	20000	0.3389	0.1442	0.0825	0.6243
g_3	20000	0.2636	0.1322	0.0255	0.5040
g_4	20000	0.3255	0.1411	0.0626	0.5841
g_5	20000	0.4286	0.1440	0.1518	0.7060
g_6	20000	0.3069	0.1101	0.0879	0.5157
g_7	20000	0.4122	0.1054	0.2091	0.6278
g_8	20000	0.3342	0.1156	0.1169	0.5652
g_9	20000	0.3103	0.1272	0.0644	0.5527
g_10	20000	0.3824	0.1180	0.1349	0.6067
c_1	20000	0.9174	0.0416	0.8352	0.9903
c_2	20000	0.9270	0.0378	0.8536	0.9906
c_3	20000	0.9113	0.0404	0.8338	0.9827
c_4	20000	0.9047	0.0486	0.8149	0.9936
c_5	20000	0.8857	0.0563	0.7662	0.9777
c_6	20000	0.8316	0.0806	0.6826	0.9777
c_7	20000	0.8286	0.0899	0.6575	0.9837
c_8	20000	0.7883	0.0858	0.6326	0.9576
c_9	20000	0.7706	0.0893	0.5962	0.9398
c_10	20000	0.8191	0.0838	0.6671	0.9784

The estimated upper asymptotes are consistently less than 1. The inclusion of the c_j parameters also significantly affects the estimates of the δ_j and α_j parameters.

To create the ICC, IIC, TCC, and TIC plots, you must modify the DATA step that creates the Plots4PL data set to accommodate the new c_j parameters and to update the probability and information function formulas for the 4PL model. The modified DATA step is as follows:

data plots;
  set out4pl;
  array alpha[10] alpha_1-alpha_10;
  array d[10] delta_1-delta_10;
  array p[10] p1-p10;
  array i[10] i1-i10;
  array g[10] g_1-g_10;
  array c[10] c_1-c_10;
  do theta=-6 to 6 by .25;
    do j=1 to 10;
      p[j]=g[j] + (c[j] - g[j])*logistic(alpha[j]*(theta-d[j]));
      i[j]=(alpha[j]**2*(p[j]-c[j])**2*(g[j]-p[j])**2)/((g[j]-c[j])**2*p[j]*(1-p[j]));
    end;
    tcf=sum(of p1-p10);
    tif=sum(of i1-i10);
    output;
  end;
run;

References

De Ayala, R.J. (2009). The Theory and Practice of Item Response Theory. New York: Guilford Press.