In the mean time, here are some key routines written in PowerBasic.
This is a very powerful programming language with abosultely amazing support. Their on-line forums contain megabytes of reuseable code and answers to an amazing range of questions. If you are in the measurement field and program in PowerBasic, please drop me a line - LMR @ Xedres.org (no X and no spaces).
The PowerBasic code is very similar to the code in Active Server Pages, VBScript, and Visual Basic. I hope you find it useful.
' IRT CALIBRATION
' Maximum likelihood estimation via Newton Raphson.
' L.Rudner Aug 2003
'
' See Lord, F.M. (1980). Applications of Item Response Theory
' to Practical Testing Problems. Erlbaum - Page 61 and Chapter 12
FUNCTION probt( theta AS DOUBLE, a AS DOUBLE, b AS DOUBLE, c AS DOUBLE ) AS DOUBLE
' probability of a correct response on item j for examine of ability th
LOCAL dem AS DOUBLE, pr AS DOUBLE
pr = 0
dem = 1.0# + EXP( -1.700# * a * ( theta - b ))
pr = c + (( 1.0# - c ) / dem )
IF pr < .001 THEN pr = .001
IF pr > .999 THEN pr = .999
FUNCTION = pr
END FUNCTION
FUNCTION mse( obs() AS WORD, th() AS DOUBLE, a AS DOUBLE, b AS DOUBLE, c AS DOUBLE, N AS DWORD, nobs AS DWORD ) AS DOUBLE
'mean square error
LOCAL i AS DWORD, pr AS DOUBLE
LOCAL ms AS DOUBLE
ms = 0
nobs = 0
FOR i = 1 TO N
IF obs( i ) < > 9 THEN
INCR nobs
pr = probt( th( i ), a, b, c )
ms = ( obs( i ) - pr ) ^ 2
END IF
NEXT
FUNCTION = ms / Nobs
END FUNCTION
FUNCTION lik( obs() AS WORD, th() AS DOUBLE, a AS DOUBLE, b AS DOUBLE, c AS DOUBLE, N AS DWORD, nobs AS DWORD ) AS EXT
'log likelihood
LOCAL i AS DWORD, pr AS DOUBLE
LOCAL lk AS EXT
lk=1.0
nobs = 0
FOR i = 1 TO N
IF obs( i ) < > 9 THEN
INCR nobs
pr = probt( th( i ), a, b, c )
lk = lk * ( pr^obs(i)) * ((1 - pr ) ^ (1-obs(i)))
END IF
NEXT
FUNCTION = LOG(lk)
END FUNCTION
SUB calib_theta( seeth AS DOUBLE, theta AS DOUBLE, ii AS WORD, a () AS DOUBLE, b() AS DOUBLE, c() AS DOUBLE, obs() AS WORD )
' 3PL
LOCAL nl AS DWORD
LOCAL tdem AS DOUBLE, tnum AS DOUBLE, tho AS DOUBLE
LOCAL i AS DWORD, j AS DWORD
LOCAL pr AS DOUBLE
' see - standard error - returned
' th - theta in and updated
' ii - number of administered items
' ui - response (correct=1, wrong=0, missing/skip=9)
'
' formula is from Lord, p 61
nl = 0
DO
tdem = 0.0#
tnum = 0.0#
tho = theta
INCR nl
FOR i = 1 TO ii
j = i
pr = probt( theta, a ( i ), b( i ), c( i ))
tnum = tnum + ( obs( i ) - pr ) * 1.7 * ( a ( i ) / ( 1.0# - c( i ))) * ( pr - c( i )) / pr
tdem = tdem + 1.7# * 1.7# * a ( i ) * a ( i ) * (( 1.0# - pr ) / pr ) * (( pr - c( i )) / ( 1.0# - c( i ))) ^ 2
NEXT i
theta = theta + ( tnum / tdem )
IF theta > 3! THEN theta = 3: tho = 3
IF theta < - 3! THEN theta = - 3: tho = - 3
LOOP UNTIL ABS( theta - tho ) < .00001 OR nl > 200
seeth = 1 / SQR( tdem )
IF seeth > 99! THEN seeth = 99!
END SUB
SUB Calib_a( seea AS DOUBLE, obs() AS WORD, th() AS DOUBLE, a AS DOUBLE, b AS DOUBLE, c AS DOUBLE, N AS DWORD )
LOCAL nl AS DWORD
LOCAL tdem AS DOUBLE, tnum AS DOUBLE, da AS DOUBLE
LOCAL i AS DWORD, j AS DWORD
LOCAL pr AS DOUBLE, aold AS DOUBLE, qr AS DOUBLE
nl = 0
DO
aold = a
tdem = 0.0#
tnum = 0.0#
INCR nl
FOR i = 1 TO N
IF obs( i ) < > 9 THEN
pr = probt( th( i ), a, b, c )
qr = 1 - pr
da = 1.7# * ( th( i ) - b ) * qr * ( pr - c ) / ( 1 - c )
tnum = tnum + (( obs( i ) - pr ) / ( pr * qr )) * da
tdem = tdem + ( da * da ) / ( pr * qr )
END IF
NEXT
a = a + ( tnum / tdem )
IF a > 10! THEN a = 10: aold = a
IF a < 0.0 THEN a = .0001: aold = a
LOOP UNTIL ABS( a - aold ) < .001 OR nl > 200
seea = 1 / SQR( tdem )
IF seea > 99! THEN seea = 99!
END SUB
SUB Calib_b( seeb AS DOUBLE, obs() AS WORD, th() AS DOUBLE, a AS DOUBLE, b AS DOUBLE, c AS DOUBLE, N AS DWORD )
LOCAL nl AS DWORD
LOCAL tdem AS EXT, tnum AS EXT, db AS DOUBLE
LOCAL i AS DWORD, j AS DWORD
LOCAL pr AS DOUBLE, bold AS DOUBLE, qr AS DOUBLE
nl = 0
DO
bold = b
tdem = 0.0#
tnum = 0.0#
INCR nl
FOR i = 1 TO N
IF obs( i ) < > 9 THEN
pr = probt( th( i ), a, b, c )
qr = 1 - pr
db = - 1.7# * a * Qr * ( pr - c ) / ( 1 - c )
tnum = tnum + (( obs( i ) - pr ) / ( pr * qr )) * db
tdem = tdem + ( db * db ) / ( pr * qr )
END IF
NEXT
b = b + ( tnum / tdem )
IF b > 5! THEN b = 5.0000: bold = b
IF b < - 5! THEN b = - 5.0000: bold = b
LOOP UNTIL ABS( b - bold ) < .0001 OR nl > 200
seeb = 1 / SQR( tdem )
IF seeb > 99! THEN seeb = 99!
END SUB
SUB Calib_c( seec AS DOUBLE, obs() AS WORD, th() AS DOUBLE, a AS DOUBLE, b AS DOUBLE, c AS DOUBLE, N AS DWORD )
LOCAL nl AS DWORD
LOCAL tdem AS EXT, tnum AS EXT, dc AS DOUBLE
LOCAL i AS DWORD, j AS DWORD
LOCAL pr AS DOUBLE, cold AS DOUBLE, qr AS DOUBLE
nl = 0
DO
cold = c
tdem = 0.0#
tnum = 0.0#
INCR nl
FOR i = 1 TO N
IF th( i ) < 0.0 AND obs( i ) < > 9 THEN ' restrict c estimation
pr = probt( th( i ), a, b, c )
qr = 1 - pr
dc = Qr / ( 1 - c )
tnum = tnum + (( obs( i ) - pr ) / ( pr * qr )) * dc
tdem = tdem + ( dc * dc ) / ( pr * qr )
END IF
NEXT
IF tdem < > 0.0 THEN c = c + ( tnum / tdem ) ELSE c = c + .01
IF c > 1! THEN c = .999: cold = c
IF c < 0.0 THEN c = .0001: cold = c
LOOP UNTIL ABS( c - cold ) < .0001 OR nl > 200
seec = 1 / SQR( tdem )
IF seec > 99! THEN seec = 99!
END SUB
SUB calib_abc( seea AS DOUBLE, seeb AS DOUBLE, seec AS DOUBLE, obs() AS WORD, th() AS DOUBLE, a AS DOUBLE, b AS DOUBLE, c AS DOUBLE, N AS DWORD )
' calibrate a, c, c given theta
LOCAL nl AS DWORD
LOCAL ao AS DOUBLE, bo AS DOUBLE, co AS DOUBLE
DO
INCR nl
ao = a
bo = b
co = c
CALL Calib_b( seeb, obs(), th(), a, b, c, N )
CALL Calib_a( seea, obs(), th(), a, b, c, N )
CALL Calib_c( seec, obs(), th(), a, b, c, N ) ' comment out to use original
LOOP UNTIL nl > 199 OR ( ABS( bo - b ) < .0001 AND ABS( co - c ) < .0001 AND ABS( ao - a ) < .0001 )
' PRINT nl; " iterations"
END SUB
FUNCTION PBMAIN
DIM th( 20000 ) AS DOUBLE, obs( 20000 ) AS WORD
DIM obsstr( 20000 ) AS STRING
LOCAL N AS DWORD
LOCAL aa$
LOCAL a AS DOUBLE, b AS DOUBLE, c AS DOUBLE
LOCAL seea AS DOUBLE, seeb AS DOUBLE, seec AS DOUBLE
LOCAL i AS DWORD, j AS DWORD
LOCAL ab$
LOCAL mserr AS DOUBLE
LOCAL pool$
LOCAL ni AS WORD, nobs AS DWORD
END FUNCTION