Implementing the Maths Garden Update Algorithm
In meow: Unified Framework for Computer Adaptive Testing Simulations

knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
library(meow)

The Maths Garden algorithm (Klinkenberg, Straatemeier, and van der Maas, 2011) is a gradient-based, Elo-style estimation method for computer adaptive practice systems. It updates person abilities and item difficulties on the fly, which makes it well suited to real-time educational applications.

Mathematical foundation

The algorithm updates both abilities ($\theta$) and difficulties ($b$) from prediction errors:

$$\theta_j^{new} = \theta_j + K_\theta \sum_{i \in I_j} (S_{ij} - E(S_{ij})), \qquad b_i^{new} = b_i + K_b \sum_{j \in J_i} (E(S_{ij}) - S_{ij}),$$

where $S_{ij} \in {0, 1}$ is the observed response, $I_j$ and $J_i$ are the items answered by person $j$ and the people who answered item $i$, and $K_\theta$, $K_b$ are learning rates. The expected response follows the Rasch (1PL) model,

$$E(S_{ij}) = \frac{1}{1 + e^{-(\theta_j - b_i)}}.$$

Intuitively, performing better than expected raises an ability estimate; an item answered correctly more often than expected becomes easier.

Implementation in `meow`

update_maths_garden() follows the parameter update contract (vignette("parameter-update")): it reads the administered responses from the matrix state and returns updated pers and item data frames. The prediction errors are aggregated per respondent and per item with tapply() rather than explicit loops:

update_maths_garden <- function(pers, item, R, admin, K_theta = 0.1, K_b = 0.1) {
  idx    <- which(admin != 0, arr.ind = TRUE)
  person <- idx[, 1]
  itm    <- idx[, 2]
  resp   <- R[idx]

  E_Sij <- stats::plogis(pers$theta[person] - item$b[itm])

  dtheta <- tapply(resp - E_Sij, person, sum)
  pers$theta[as.integer(names(dtheta))] <-
    pers$theta[as.integer(names(dtheta))] + K_theta * dtheta

  db <- tapply(E_Sij - resp, itm, sum)
  item$b[as.integer(names(db))] <-
    item$b[as.integer(names(db))] + K_b * db

  list(pers = pers, item = item)
}

Using it

Learning rates are passed through update_args:

sim <- meow(
  select_fun  = select_max_info,
  update_fun  = update_maths_garden,
  data_loader = data_simple_1pl,
  data_args   = list(N_persons = 100, N_items = 50),
  update_args = list(K_theta = 0.05, K_b = 0.05)
)
head(sim$results[, 1:3])

Extending the algorithm

Because update functions are ordinary R functions, variations are easy. The following adds adaptive learning rates (shrinking as a respondent answers more items) and bounds on the estimates, while staying within the matrix contract:

update_maths_garden_adaptive <- function(pers, item, R, admin,
                                         base_K = 0.1, decay = 0.05,
                                         bounds = c(-4, 4)) {
  idx    <- which(admin != 0, arr.ind = TRUE)
  person <- idx[, 1]
  itm    <- idx[, 2]
  resp   <- R[idx]
  E_Sij  <- stats::plogis(pers$theta[person] - item$b[itm])

  n_person <- tapply(resp, person, length)
  err_p    <- tapply(resp - E_Sij, person, sum)
  who_p    <- as.integer(names(err_p))
  K_p      <- base_K / (1 + n_person * decay)
  pers$theta[who_p] <- pers$theta[who_p] + K_p * err_p
  pers$theta <- pmin(pmax(pers$theta, bounds[1]), bounds[2])

  n_item <- tapply(resp, itm, length)
  err_i  <- tapply(E_Sij - resp, itm, sum)
  who_i  <- as.integer(names(err_i))
  K_i    <- base_K / (1 + n_item * decay)
  item$b[who_i] <- item$b[who_i] + K_i * err_i
  item$b <- pmin(pmax(item$b, bounds[1]), bounds[2])

  list(pers = pers, item = item)
}

Practical notes

Start with modest learning rates ($K_\theta = K_b = 0.1$) and check that estimates stabilize across iterations.
Large learning rates can make estimates oscillate; bounding the estimates helps.
The algorithm assumes a Rasch model; if your data need discrimination parameters, consider an MLE updater or a 2PL extension.