dscore estimation under different transformations

[stefvanbuuren]

I would expect that ability estimation is insensitive to a linear transformation of the ability scale, but such turns out to be the case only approximately. The following script compares two scales (one D-score and one logit), that produces (slightly) different estimates.

transform <- c(41.10, 2.23)

# ability
data <- data.frame(
  age = rep(round(21/365.25, 4), 10),
  GSFIXEYE = c(NA, NA, 0, 0, 0, 1, 0, 1, 1, 1),
  GSRSPCH =  c(NA, NA, 0, 0, 1, 0, 1, 0, 1, 1),
  GSMLEG =   c(NA,  0, 0, 1, 0, 0, 1, 1, 0, 1))
items <- c("GSFIXEYE", "GSRSPCH", "GSMLEG")

keyd <- data.frame(item = items,
                   delta = gettau(items = items),
                   stringsAsFactors = FALSE)

zd <- ability(data, items = items, dec = 4, metric = "dscore", 
              key = keyd)$b

qpl <- ((-10:100) - transform[1]) / transform[2]
keyl <- data.frame(item = items,
                   delta = gettau(items = items),
                   stringsAsFactors = FALSE)
keyl$delta <- (keyl$delta - transform[1]) / transform[2]
zl <- ability(data, items = items, dec = 4, transform = transform, 
              qp = qpl, metric = "logit", key = keyl)$b

test_that("logit and dscore are identical", {
  expect_identical(zl, (zd - transform[1])/transform[2])
})

When tracking down differences between the two methods, I found that taking out the division (qp[2] - qp[1]) in normalize() will produce the same prior. After that, the next divergence appears in cpc <- t(exp(outer(0:m, qp) + c(0, -cumsum(delta)))) in posterior(). This suggest that the exponential transform here introduces instability. I have no time to further dive in and smooth out differences, and have put back (qp[2] - qp[1]) into normalise(), but evidently this is somewhat fishy.

Some options to pursue:

1. Perhaps we can bypass quadrature methods altogether, and use a normal approximation everywhere. The prior is normal anyway, and I remember seeing that a prior normal in combination with a logistic model produces a normal posterior (Albert and Chibb? Gelman's BDA book?). If so, this would considerably speed up and simplify calculations.
2. Study what happens in packages ltm, sirt or similar packages that can estimate EAP.
3. Choose one scale, derive the other by a linear transform. I would then choose for the D-score scale, and derive the logic form from that
4. Ot just live with the difference? It's not big, and in practice it may not matter.

[stefvanbuuren]

I have now found out that option 1 cannot work. The posteriors are not normal. See https://stefvanbuuren.name/dbook1/sec-dscoreestimation.html#numerical-example for an example that shows that it is skewed.
Currently on option 4.

[stefvanbuuren]

The above code will not work anymore in dscore 0.38.0 and above. Will update to account for the new function arguments.

[stefvanbuuren]

Solved in dscore 1.4.1