longitudinal_grmtree: Fit a Longitudinal Graded Response Model Tree for Response...

In grmtree: Recursive Partitioning for Graded Response Models

Fit a Longitudinal Graded Response Model Tree for Response Shift Detection

Description

This function implements a tree-based longitudinal graded response model using model-based recursive partitioning (MOB) to detect measurement heterogeneity in patient-reported outcome measures (PROMs) measured at two time points. The Longitudinal GRMTree extends the cross-sectional grmtree to the longitudinal setting by embedding a constrained two-factor GRM within the MOB framework. The resulting tree identifies patient subgroups where the longitudinal measurement model differs, providing the foundation for subgroup-specific response shift characterization via rs_characterize.

Usage

longitudinal_grmtree(
  formula,
  data,
  n_items = NULL,
  na.action = na.omit,
  control = grmtree.control(),
  mtry = NULL,
  ...
)

Arguments

formula	A formula specifying the model structure with a wide-format response matrix on the left-hand side and partitioning variables on the right-hand side (e.g., resp_wide ~ age + sex + comorbidity_count). The response matrix must have 2 * n_items columns: the first n_items columns are Time 1 (T1) responses and the next n_items columns are Time 2 (T2) responses for the same items.
data	A data frame containing the variables in the formula. Must include the response matrix as a matrix-valued column (created via data$resp_wide <- as.matrix(...)) and all partitioning variables as separate columns.
n_items	Integer specifying the number of items per time point. If NULL (default), automatically detected as half the number of columns in the response matrix. Must satisfy ncol(response) == 2 * n_items.
na.action	How to handle missing values. Default is na.omit, which removes rows with any missing values in the response matrix or partitioning variables.
control	A list of control parameters created by grmtree.control. Key parameters include: minbucket Minimum number of observations in a terminal node. Should be at least 10 times the number of free parameters per node to ensure stable estimation. For 8 items with 5 categories each, the constrained longitudinal GRM has approximately 40 item parameters plus 3 latent parameters, suggesting minbucket >= 200. alpha Significance level for the parameter instability tests. Default is 0.05. p_adjust Method for adjusting p-values across covariates at each split. Options include "bonferroni" (default, applied locally during tree construction), "BH" (Benjamini-Hochberg, applied post-hoc with pruning), and others. See grmtree.control for details.
mtry	Number of variables randomly sampled as candidates at each split. If NULL (default), all variables are considered. Can be used for random forest extensions.
...	Additional arguments passed to the internal fitting function longitudinal_grmfit and ultimately to mirt.

Details

Overview

The Longitudinal GRMTree is a unified framework for response shift detection that operates in two phases:

Phase 1 (this function): Identifies patient subgroups where the constrained longitudinal measurement model differs. The constrained model represents the null hypothesis of no response shift (item parameters equal across time). MOB tests whether this null model's item parameters are stable across patient covariates and recursively partitions the sample where instability is detected.

Phase 2 (rs_characterize): Within each terminal node, tests whether item parameters actually change from T1 to T2 by comparing the constrained model to an unconstrained model using likelihood ratio tests. This phase characterizes response shift at the item level, classifying changes as recalibration, reprioritization, or both.

Constrained Two-Factor Longitudinal GRM

Let Y_{im} denote the response of individual i to item m at time t. The longitudinal GRM specifies two correlated latent factors:

\theta_{i,T1} \sim N(0, 1)

\theta_{i,T2} \sim N(\mu_{T2}, \sigma^2_{T2})

Cov(\theta_{i,T1}, \theta_{i,T2}) = \sigma_{12}

For item m at time t, the graded response model is:

P(Y_{im,t} \geq j | \theta_{i,t}) = \frac{\exp(a_{m,t}(\theta_{i,t} - b_{mj,t}))}{1 + \exp(a_{m,t}(\theta_{i,t} - b_{mj,t}))}

where a_{m,t} is the discrimination parameter and b_{mj,t} are threshold parameters for item m at time t.

Under the no-response-shift constraint:

a_{m,T1} = a_{m,T2} \quad \text{and} \quad b_{mj,T1} = b_{mj,T2} \quad \forall m, j

This means item parameters are identical across time, so any observed changes in responses are attributed to true latent change (\mu_{T2}) rather than changes in measurement properties.

The Longitudinal GRMTree Algorithm

Step 1: Global Model Estimation. Fit the constrained longitudinal GRM to all individuals at the root node, estimating item parameters \hat{\beta} and latent parameters (\mu_{T2}, \sigma^2_{T2}, \sigma_{12}) via maximum likelihood.

Step 2: Parameter Stability Testing. For each covariate X_p, compute individual score function contributions s(\hat{\beta}; \mathbf{y}_i) and test whether the item-parameter scores fluctuate systematically with X_p using structural change tests; the structural-parameter scores (the follow-up mean \mu_{T2} and the between-occasion covariance \sigma_{12}) are held out, so the null hypothesis is that the item parameters are stable across all values of X_p.

Step 3: Recursive Partitioning. If significant instability is detected (after p-value adjustment):

• Select the covariate X_p^* with strongest instability

• Find the optimal split point c^* maximizing the partitioned log-likelihood

• Split the sample: X_p^* \leq c^* vs. X_p^* > c^*

Step 4: Recursion. Repeat Steps 1–3 within each child node until no significant instability remains or the minimum node size is reached.

Step 5 (Post-hoc): Apply rs_characterize to test for response shift within each terminal node.

Formal Model Structure

The fitted Longitudinal GRMTree provides a piecewise constrained longitudinal GRM:

P(Y_{im,t} = k | \theta_{i,t}, \mathbf{x}_i) = \sum_{b=1}^B I(\mathbf{x}_i \in \mathcal{X}_b) \cdot P_b(Y_{im,t} = k | \theta_{i,t})

where B is the number of terminal nodes, \mathcal{X}_b is the covariate subspace defining terminal node b, and P_b is the node-specific constrained longitudinal GRM. Each terminal node contains:

• Node-specific constrained item parameters (equal across T1 and T2)

• Node-specific latent trait parameters: \mu_{T2,b}, \sigma^2_{T2,b}, r_{T1,T2,b}

Interpretation

The tree structure identifies patient subgroups whose item parameters differ, i.e. subgroups with differential item functioning in the constrained longitudinal measurement model. A split reflects:

• Different item discrimination patterns across subgroups

• Different threshold locations across subgroups

The structural parameters – the true latent change \mu_{T2}, the follow-up latent variance \sigma^2_{T2}, and the test-retest correlation – are estimated and reported at every node but are held out of the split test, so a subgroup difference in how much the construct truly changed, or in its stability over time, cannot by itself produce a split.

Crucially, the tree does not directly detect response shift (temporal change in item parameters within a subgroup). Response shift is tested in Phase 2 using rs_characterize, which relaxes the equality constraints within each terminal node. A tree with only one terminal node (no split) indicates that no covariate moderates the measurement model; it does not imply the absence of response shift, since rs_characterize can still detect uniform RS in the unsplit root node.

Relationship to Existing Methods

The Longitudinal GRMTree extends several existing approaches:

• GRMTree (Arimoro et al., 2026): Cross-sectional DIF detection using tree-based GRM. The longitudinal extension embeds a two-factor model instead of a single-factor model.

• LIRTree (Ames & Leventhal, 2021): Longitudinal tree-based IRT using Rasch/2PL models. The GRM extension allows item-specific discrimination parameters, which are important for PROMs where items vary in discriminating ability.

• Oort SEM (Oort, 2005): Response shift detection via structural equation modeling. The tree-based approach removes the requirement for pre-specified subgroups.

Practical Recommendations

• Minimum node size: Use at least 10–25 times the number of free parameters per node. For M items with K response categories, the constrained model has approximately M(K-1) + M + 3 parameters (M discriminations, M(K-1) thresholds, and 3 latent parameters). For 8 items with 5 categories: 8 \times 4 + 8 + 3 = 43 parameters, suggesting minbucket = 200--400.

• Number of covariates: Bonferroni correction becomes increasingly conservative with more covariates. Consider using p_adjust = "BH" for exploratory analyses with many covariates.

• Sample size: With typical PROM instruments (5–15 items), a total sample of at least 500–1000 is recommended for adequate power to detect meaningful splits.

• Response matrix preparation: Use prepare_longitudinal_data to construct the wide-format response matrix from separate T1 and T2 item columns.

Value

An object of class c("longitudinal_grmtree", "grmtree", "modelparty", "party") containing the fitted tree structure. The object inherits from modelparty and includes:

Tree structure

Accessible via standard partykit methods such as nodeids, data_party, and indexing with [[.

Node models

Each terminal node contains a fitted constrained longitudinal GRM with item parameters (discrimination and thresholds), latent trait means (\mu_{T2}), and latent trait covariance matrix (\sigma^2_{T2}, r_{T1,T2}).

info$n_items

Number of items per time point.

info$model_type

"longitudinal_grm".

info$p_adjust

The p-value adjustment method used.

info$call

The original function call.

Author(s)

Olayinka Imisioluwa Arimoro [e-mail], Lisa M. Lix, Tolulope T. Sajobi

References

Arimoro, O. I., Lix, L. M., Patten, S. B., Sawatzky, R., Sebille, V., Liu, J., Wiebe, S., Josephson, C. B., & Sajobi, T. T. (2025). Tree-based latent variable model for assessing differential item functioning in patient-reported outcome measures: a simulation study. Quality of Life Research. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1007/s11136-025-04018-6")}

Ames, A. J., & Leventhal, B. C. (2021). Application of a longitudinal IRTree model: response style changes over time. Educational and Psychological Measurement, 81(3), 561–582.

Oort, F. J. (2005). Using structural equation modeling to detect response shifts and true change. Quality of Life Research, 14(3), 587–598.

Samejima, F. (1969). Estimation of latent ability using a response pattern of graded scores. Psychometrika Monograph Supplement, 34, 100–114.

Sprangers, M. A., & Schwartz, C. E. (1999). Integrating response shift into health-related quality of life research. Social Science & Medicine, 48(11), 1507–1515.

Zeileis, A., Hothorn, T., & Hornik, K. (2008). Model-based recursive partitioning. Journal of Computational and Graphical Statistics, 17(2), 492–514.

See Also

rs_characterize for Phase 2 response shift characterization, prepare_longitudinal_data for data preparation, grmtree for cross-sectional DIF detection, grmtree.control for control parameters, mob for the underlying MOB framework

Examples

library(grmtree)

  # Load the synthetic longitudinal data
  data("grmtree_long_data", package = "grmtree")

  # Prepare the wide-format response matrix
  items_t1 <- c("MOS_Listen", "MOS_Info", "MOS_Advice_Crisis", "MOS_Confide",
                "MOS_Advice_Want", "MOS_Fears", "MOS_Personal", "MOS_Understand")
  ld <- prepare_longitudinal_data(
    data = grmtree_long_data,
    items_t1 = items_t1,
    items_t2 = paste0(items_t1, "_year1"),
    covariates = c("sex", "age", "residency", "job",
                   "education", "comorbidity_count", "ever_smoker")
  )

  # Phase 1: fit the longitudinal GRM tree
  ltree <- longitudinal_grmtree(
    resp_wide ~ sex + age + residency + job +
      education + comorbidity_count + ever_smoker,
    data = ld, n_items = 8,
    control = grmtree.control(minbucket = 200)
  )

  # Print tree structure
  print(ltree)

  # Plot threshold regions
  plot(ltree, type = "regions", tnex = 2L)

  # Phase 2: characterize response shift within each subgroup
  rs <- rs_characterize(ltree, p_adjust = "fdr",
    global_p_adjust = "bonferroni")
  print(rs)

grmtree documentation built on Sept. 2, 2026, 1:07 a.m.