Specification of a random slopes model with two grouping factors

Question

Suppose people take a 10-item exam. For each item, $k = 1,2,…,10$, exactly one rater assigns a score to exactly one person, with the constraint that no person sees the same rater twice. There are $i = 1,2,…,5000$ people and $j = 1,2,…,500$ raters. So people and raters are partially crossed, with most person-rater pairs having 0 data points, and the rest having exactly 1 data point.

Update: items were randomly assigned to positions.

Suppose I want to estimate a random slopes model where scores vary linearly with item position and the variability in scores is due to person intercept variability, person slope variability, rater intercept variability, rater slope variability, and residual error variability – in R,

lmer(rating ~ itemPosition + (itemPosition | personID) + (itemPosition | raterID))

Firstly, is there a better way to write the model down in non-matrix form than what I have below? Because there are two grouping factors, I am confused about the indices. If I let $q$ index data points, I want something like

$Y_q = \beta_0 + \beta_1 ItemPosition_q + p_{0i} + p_{1i} ItemPosition_q + r_{0j} + r_{1j} ItemPosition_q + e_q$

where $p_{0i} \sim N(0, \sigma_{p0}^2)$ denotes a person-level random intercept, $r_{0j} \sim N(0, \sigma_{r0}^2)$ denotes a rater-level random intercept, $p_{1i} \sim N(0, \sigma_{p1}^2)$ denotes a person-level random slope, $r_{1i} \sim N(0, \sigma_{r1}^2)$ denotes a rater-level randoms slope, $e_q \sim N(0, \sigma^2)$ is the residual error, the person random effects are uncorrelated with the rater random effects, the random effects are uncorrelated with the error term, $cov(p_{0i}, p_{1i}) = \sigma_{p01}$, and $cov(r_{0i}, r_{1i}) = \sigma_{r01}$. More specifically, can I replace the index $q$ with $ijk$? Secondly, do I have the correct expression for $$Var(p_{0i} + p_{1i} ItemPosition_q + {r_0j} + r_{1j} ItemPosition_q + e_q) ?$$ I get $(\sigma_{0p}^2 + 2\sigma_{p01} ItemPosition_q + \sigma_{1p}^2 ItemPosition_q^2)$ + $(\sigma_{0r}^2 + 2\sigma_{r01} ItemPosition_q + \sigma_{1r}^2 ItemPosition_q^2)$ + $\sigma^2.$ I want to use this as the denominator to calculate the proportion of variance due to persons, for example. Lastly, any general comments about the appropriateness of the model are also welcome but not at all necessary (as I haven’t provided any information about the data).

Answer 1

I believe both your model formulation and variance are fine for your purposes (and with the stated independences). Indeed, indexing with ijk is a good idea (note that you can then replace ItemPosition simply with k).

Word of warning though: your question does not explicitly hold this: this kind of variance is the variance between subjects for given question.

Second word of warning: is it reasonable to expect linearity of the rates with item position?

Answer 2

Let's think of a rating of test item $k$ as a function of a three latent variables: the intrinsic difficulty of the question, $u_k$, an item specific ability of test subject $i$, $u_{ik}$, a question specific rating 'leniency', $u_{jk}$.

Assuming a rating is a linear function of these factors we can say:

$$ \text{Rating}_{ijk} = u_k + u_{ik} + u_{jk} + \varepsilon_{ijk} $$

If $u_k$ is drawn at random conditional on $k$, we can get a valid estimate of item position specific ability learning and leniency from this model:

$$ \text{Rating}_{ijk} = \beta_0 + u_{ik} + u_{jk} + \varepsilon_{ijk} $$

Where $\beta_0$ is the average difficulty of the questions and variance of difficulties is aborbed into the error term. This model is estimable but it has a lot of parameters. If you are interested in detecting linear trends then you can simplify the model as follows.

$$ \begin{align} u_{ik} = p_{0i} + p_{1i}\text{Item Position}_k \\ u_{jk} =r_{0j} + r_{1j}\text{Item Position}_k \end{align} $$

And in the combined model:

$$ \text{Rating}_{ijk} = \beta_0 + p_{0i} + p_{1i}\text{Item Position}_k + r_{0j} + r_{1j}\text{Item Position}_k + \varepsilon_{ijk} $$

This is similar to your proposed model except it does not include a predictor for specific item difficulty, which makes sense because if item difficulties are random, your $\beta_1$ term should equal 0.

You might consider adding higher order polynomials of item position to allow for more flexibility in learning and leniency rates.