--- title: Getting started with PP author: Manuel Reif output: rmarkdown::html_vignette vignette: > %\VignetteEngine{knitr::rmarkdown} %\VignetteIndexEntry{Getting started with PP} %\VignetteEncoding{UTF-8} --- ```{r setup, include = FALSE} knitr::opts_chunk$set(collapse = TRUE, comment = "#>") ``` PP estimates person parameters (abilities) when the item parameters are already known, for example from an earlier calibration or an item bank. It is the right package if 1. you have item responses of a group of persons, 2. you have item parameters for these items, which are treated as fixed, 3. the items follow the 1PL to 4PL model, the (generalized) partial credit model, or a mix of both, and you want ability estimates (MLE, WLE, MAP, EAP or a robust estimate), jackknife estimates, plausible values or person fit statistics. This vignette shows the estimation functions. Person fit is the topic of `vignette("intro_pfit")`. Estimation and person fit in one step and models fitted with the eRm package are shown in `vignette("ppass_erm")`. ## A 3PL test The data set `fourpl_df` contains the responses of 60 persons to 12 items, together with an id and the sex of each person. The item parameters are attributes of the data set. ```{r get} library(PP) data(fourpl_df) head(fourpl_df[, 1:7]) diff_par <- attr(fourpl_df, "diffpar") slope_par <- attr(fourpl_df, "slopes") low_par <- attr(fourpl_df, "lowerA") items <- as.matrix(fourpl_df[, -(1:2)]) ``` `PP_4pl()` covers all dichotomous models. Which model is used depends on the item parameters you supply: difficulties (`thres`) alone give the 1PL model, adding `slopes` gives the 2PL model, adding `lowerA` gives the 3PL model, and adding `upperA` as well gives the 4PL model. Here the difficulties, slopes and lower asymptotes are known, so the 3PL model is used. ```{r mle} res_mle <- PP_4pl(respm = items, thres = diff_par, slopes = slope_par, lowerA = low_par, type = "mle") head(res_mle$estimates, 8) ``` The result is a list. Its main part is the data.frame `$estimates` with one row per person, in the order of the response matrix: - `estimate` and `SE`: the ability estimate and its standard error, - `nsteps`: the number of Newton steps needed for this person, - `status`: `"ok"` if the estimate converged inside `range` (default `c(-6, 6)`); `"extreme"` for persons with all answers right or wrong, who have no finite maximum likelihood estimate; `"boundary"` if there is no solution inside `range` and the estimate is set to its limit; `"no_response"` for persons without any answer; `"not_converged"` if the estimation did not converge within `maxsteps` steps. ```{r status} summary(res_mle, nrowmax = 5) ``` Person 35 answered only item 8 correctly, the item with the highest guessing parameter, and all easier items wrong (item 7 is missing). Under the 3PL model the likelihood of this pattern keeps rising as the ability decreases, so there is no maximum inside `range`. The estimate is set to the lower limit, the status is `"boundary"`, and `PP_4pl()` gives a warning. `coef()` returns the estimates as a vector and `as.data.frame()` the table, so the estimates can easily be added to the person data: ```{r bind} dat <- data.frame(fourpl_df[, 1:2], as.data.frame(res_mle)[, c("estimate", "SE")]) head(dat) ``` Missing responses (`NA`) are allowed; the ability of a person is estimated from the items they answered. ## Choosing the estimation type `type` selects the estimator: - `"mle"`: maximum likelihood. Persons with all answers in the lowest or highest category get `-Inf`/`Inf`. - `"wle"` (default): weighted likelihood (Warm, 1989). It is less biased than the MLE and finite for extreme scores. - `"map"` and `"eap"`: maximum and mean of the posterior with a normal prior `N(mu, sigma2)` (default `N(0, 1)`). The prior may differ between persons. The EAP is computed on a grid of ability values; its SE is the posterior standard deviation. - `"robust"`: a Huber-type estimator (Schuster & Yuan, 2011), which downweights unexpected answers to items far from the person's ability. ```{r types} res_wle <- PP_4pl(respm = items, thres = diff_par, slopes = slope_par, lowerA = low_par) res_eap <- PP_4pl(respm = items, thres = diff_par, slopes = slope_par, lowerA = low_par, type = "eap", mu = ifelse(fourpl_df$sex == "w", 0.2, -0.2), sigma2 = 0.8) head(cbind(mle = coef(res_mle), wle = coef(res_wle), eap = coef(res_eap))) ``` ## Measurement precision `plot()` shows where the test measures accurately and whether the persons are located there. The upper panel is the distribution of the estimates, the lower panel the standard error of each person. Because the item parameters are known, the standard error of the whole test follows from the test information: $1/\sqrt{I(\theta)}$. The ticks at the bottom mark the items, at the ability where each item gives the most information. One response in `fourpl_df` is missing, so the curve of all 12 items is drawn on request (`curve = TRUE`); it holds for every person who answered all items. ```{r precision, fig.width = 6, fig.height = 4.5, fig.alt = "Histogram of the WLE estimates above, and below the standard error of each person against the estimate with the standard error curve of the test"} plot(res_wle, which = "precision", curve = TRUE) ``` ## Polytomous and mixed tests `PP_gpcm()` estimates abilities for the generalized partial credit model. With all slopes equal to 1 (the default) this is the partial credit model. `thres` has one column per item with its thresholds; items with fewer categories have `NA` below their last threshold. A first row of zeros is optional. ```{r gpcm} thres <- rbind(c(-1.5, -1, -0.5, 0, 0.5), c(0, 0.5, 1, 1.2, NA), c(1, 1.5, NA, 2, NA)) slopes <- c(1, 1.2, 0.8, 1.1, 0.9) set.seed(1234) x_gpcm <- sim_gpcm(thres = thres, alpha = slopes, theta = rnorm(200)) res_gpcm <- PP_gpcm(respm = x_gpcm, thres = thres, slopes = slopes, type = "wle") head(res_gpcm$estimates) ``` For a test with dichotomous and polytomous items, use `PPall()`. Each column of `thres` holds the thresholds of one item: a 4PL item has a single threshold (its difficulty), a GPCM item has several, and this decides the model of each item (`findmodel()` shows it). Lower and upper asymptotes apply to the 4PL items only. ```{r mixed} thres_mix <- rbind(c(-1, 0, 1, -0.5, 0.3), c(NA, NA, NA, 0.6, 1.1)) findmodel(thres_mix) x_mix <- sim_mixed(thres = thres_mix, alpha = c(1, 1.3, 0.9, 1, 1), lowerA = c(0.2, 0.1, 0.15, 0, 0), theta = rnorm(200)) res_mix <- PPall(respm = x_mix, thres = thres_mix, slopes = c(1, 1.3, 0.9, 1, 1), lowerA = c(0.2, 0.1, 0.15, 0, 0), type = "map") summary(res_mix$estimates$estimate) ``` `sim_4pl()`, `sim_gpcm()` and `sim_mixed()` simulate responses for given item and person parameters, as used above. ## Jackknife `JKpp()` leaves out each answered item once, estimates the ability again and combines the pseudovalues into a jackknife estimate and standard error (Wainer & Wright, 1980). `cmeth = "median"` and `cmeth = "AMT"` (Andrews' sine estimator) give robust estimates; the standard error is the same for all three. ```{r jk} jk <- JKpp(res_wle) head(jk$estimates) head(JKpp(res_wle, cmeth = "AMT")$estimates) ``` ## Plausible values Plausible values are draws from the posterior distribution of each person. They need a prior, so `PV()` takes an estimate of type `"eap"` or `"map"` and uses its prior (here different for women and men). ```{r pv} set.seed(1) pv <- PV(res_eap, npv = 5) head(pv$pvdraws) ``` By default (`approx = TRUE`) the values are drawn from a normal distribution with the posterior mean and standard deviation. With `approx = FALSE` they are drawn from the posterior itself. ## References Schuster, C., & Yuan, K. H. (2011). Robust estimation of latent ability in item response models. *Journal of Educational and Behavioral Statistics, 36*(6), 720–735. Wainer, H., & Wright, B. D. (1980). Robust estimation of ability in the Rasch model. *Psychometrika, 45*(3), 373–391. Warm, T. A. (1989). Weighted likelihood estimation of ability in item response theory. *Psychometrika, 54*, 427–450.