gpcm parametrization

[bea...@gmx.de]

hi phil and my fellow mirt-users,

i have a (maybe more general) question regarding the parametrization of the gpcm-model in mirt.

so far, when i used ltm to estimate my models, the gpcm-parameters were the points on the theta-scale were the probability functions of two adjacent categories met.

how can i transform the gpcm-parameters given by mirt into these?

see below an example with the science-data

library(ltm)
library(mirt)

dat = sapply(Science, as.numeric)[, c(-2,-5:-6)]

ltm_model = gpcm(dat)
mirt_model = mirt(dat, model = 1, itemtype="gpcm")

coef(ltm_model)
coef(mirt_model)

plot(ltm_model, item = 1)
abline(v = coef(ltm_model)[1,1:3]) # the threshold-params as vertical lines

thanks a lot!

[Phil Chalmers]

Hi Felix,

I've received lots of requests to convert the parametrization I use to
the more traditional metrics, but I just haven't put it high enough on
the priority list to do it (github pull request anyone?). Here's how
to convert the gpcm into what you see in ltm (the slopes should be the
same):

ds <- coef(mirt_model, digits = 10)[[2]][-1]
a <- coef(mirt_model, digits = 10)[[2]][1]

dprime <- ds/a
ltmds <- numeric(3)
for(i in 1:3) ltmds[i] <- -(dprime[i+1]-dprime[i])
ltmds

A similar argument can be used for the nominal model, since they share
similar parametrizations. Hope that helps!

Phil

> --
> You received this message because you are subscribed to the Google Groups "mirt-package" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to mirt-package...@googlegroups.com.
> For more options, visit https://groups.google.com/groups/opt_out.

[bea...@gmx.de]

Thanks a lot! Best, Felix

[Phil Chalmers]

All right, I caved. On the Github dev you can pass the option IRTpars = TRUE to the coef() function, and if appropriate the traditional IRT metric will be returned. Still needs a few more itemtypes, but it's mostly done. Cheers.

Phil

[bea...@gmx.de]

Hi Phil, could you provide similar code for the graded response model? Unfortunately, when i install the dev-version from github, my code breaks, so i like to stick with the 0.9.0 version. best, felix

[Matthew Sigal]

What breaks in your code when using the latest dev-version?

[bea...@gmx.de]

i think it was about the order/variables in the parameter-dataframe that i got with mirt(), changed and then gave back into mirt().

[Phil Chalmers]

Hey Felix,

From the github code for the graded model I have

for(i in 2:ncat)
par[i] <- -par[i]/par[1]
names(par) <- c('a', paste0('b', 1:(length(par)-1)))

which converts the metric to the traditional a/b format. Cheers.

Phil

On Fri, Oct 11, 2013 at 1:07 PM, <bea...@gmx.de> wrote:
> Hi Phil, could you provide similar code for the graded response model? Unfortunately, when i install the dev-version from github, my code breaks, so i like to stick with the 0.9.0 version. best, felix
>

[Phil Chalmers]

That makes sense, since *someone* in this thread suggested the
data.frame input should be a little safer :-). Redefining your old
data.frame with names from the new one should do the trick though.

Phil

On Sat, Oct 12, 2013 at 11:40 AM, <bea...@gmx.de> wrote:
> i think it was about the order/variables in the parameter-dataframe that i got with mirt(), changed and then gave back into mirt().
>

[bea...@gmx.de]

ah, great, thanks :)

[bea...@gmx.de]

Dear Phil,

i used the code mentioned above to write two functions to transform the GPCM-parameters and i would like to ask, if you could check the numbers, since the params are somewhat different to what coef(..., IRTpars = TRUE) reports.

Initially, i have the following parameters from a unidimensional GPCM:

slope step1 step2 step3 step4
2.35 0.41 1.27 2.13 3.02

i convert them using my function (see below) and get that:
a1 d0 d1 d2 d3 d4
2.35 0 -0.964 -3.948 -8.954 -16.05

I take those parameters to fix a model and estimate fscores. This seems to work flawlessly, the results make sense. however, using coef(..., IRTpars = TRUE) gives me some parameters, which are different from the ones i putted in the model.

a b1 b2 b3 b4 b5
par 2.35 0 0.836 2.531 5.087 8.532

Is this kind of a different parametrization than my first one (because there is b1 = 0)? Do the values make sense?

Thanks a lot!
Felix

If necessary, here are the functions i used for transformation:

transform_mirt_to_thresholds = function(slope, ds)
{
ltmds <- c(slope, -diff(ds/slope))
names(ltmds) = c("slope", paste0("threshold", 1:(length(ds)-1)))
return(ltmds)
}

transform_thresholds_to_mirt = function(slope, thresholds)
{
back = c(0, -thresholds)
for (i in 3:length(back)) back[i] = back[i-1] + back[i]
back = c(slope, back*slope)
names(back) = c("a1", paste0("d", 0:(length(back)-2)))
return(back)
}

[Phil Chalmers]

I'm not sure the transform_mirt_to_thresholds is working the same as the package. The mirt2traditional() internal function has this block for the gpcm (where par <- coef(mod)[[itemnumber]]): 

        ds <- par[-1]/par[1]

        newd <- numeric(length(ds)/2)

        for(i in 1:length(newd))

            newd[i] <- -(ds[i+length(newd)] - ds[i])

        par <- c(par[1], newd)

        names(par) <- c('a', paste0('b', 1:length(newd)))

So working this out backwards should be pretty straightforward. Cheers.

Phil

}

--
You received this message because you are subscribed to the Google Groups "mirt-package" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mirt-package...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

[bea...@gmx.de]

Mhm, it doesn't seem to work the same. If you consider my earlier example, the results from coef(mirt_model, IRTpars = TRUE) are different numbers than the ltm results. the models itself seem to be the same, when you look at the EAP estimates.

library(ltm)
library(mirt)
dat = sapply(Science, as.numeric)[, c(-2,-5:-6)]
ltm_model = gpcm(dat)
mirt_model = mirt(dat, model = 1, itemtype="gpcm")
coef(ltm_model)
coef(mirt_model)

coef(mirt_model, IRTpars = TRUE)
plot(fscores(mirt_model)[, "F1"],factor.scores(ltm_model,method="EAP")$score.dat$z1)

[Phil Chalmers]

You were completely right, the code was out of sync in the last update since now I include the fixed 'scoring' coefficients as well. It should looks like this:

        ds <- par[-1]/par[1]

        ds <- ds[-c(1L:ncat)]

        newd <- numeric(length(ds)-1L)

        for(i in 2:length(ds))

            newd[i-1L] <- -(ds[i] - ds[i-1L])

        par <- c(par[1], newd)

        names(par) <- c('a', paste0('b', 1:length(newd)))

Thanks for noticing this Felix, the update has been pushed.

Phil

[bea...@gmx.de]

Great, thanks a lot! works nicely! felix