rmult.clm: Simulating Correlated Ordinal Responses Conditional on a Marginal Cumulative Link Model Specification

Description

Simulates correlated ordinal responses assuming a cumulative link model for the marginal probabilities.

Usage

rmult.clm(clsize = clsize, intercepts = intercepts, betas = betas,
  xformula = formula(xdata), xdata = parent.frame(), link = "logit",
  cor.matrix = cor.matrix, rlatent = NULL)

Value

Returns a list that has components:

Ysim: the simulated ordinal responses. Element ($i$,$t$) represents the realization of $Y_{it}$.
simdata: a data frame that includes the simulated response variables (y), the covariates specified by xformula, subjects' identities (id) and the corresponding measurement occasions (time).
rlatent: the latent random variables denoted by $e^{O1}_{it}$ in Touloumis (2016).

Arguments

clsize: integer indicating the common cluster size.
intercepts: numerical vector or matrix containing the intercepts of the marginal cumulative link model.
betas: numerical vector or matrix containing the value of the marginal regression parameter vector associated with the covariates (i.e., excluding intercepts).
xformula: formula expression as in other marginal regression models but without including a response variable.
xdata: optional data frame containing the variables provided in xformula.
link: character string indicating the link function in the marginal cumulative link model. Options include 'probit', 'logit', 'cloglog' or 'cauchit'. Required when rlatent = NULL.
cor.matrix: matrix indicating the correlation matrix of the multivariate normal distribution when the NORTA method is employed (rlatent = NULL).
rlatent: matrix with clsize columns containing realizations of the latent random vectors when the NORTA method is not preferred. See details for more info.

Author

Anestis Touloumis

Details

The formulae are easier to read from either the Vignette or the Reference Manual (both available here).

The assumed marginal cumulative link model is $$Pr(Y_{it}\le j |x_{it})=F(\beta_{tj0} +\beta^{'}_{t} x_{it})$$ where $F$ is the cumulative distribution function determined by link. For subject $i$, $Y_{it}$ is the $t$-th ordinal response and $x_{it}$ is the associated covariates vector. Finally, $\beta_{tj0}$ is the $j$-th category-specific intercept at the $t$-th measurement occasion and $\beta_{tj}$ is the $j$-th category-specific regression parameter vector at the $t$-th measurement occasion.

The ordinal response $Y_{it}$ is obtained by extending the approach of McCullagh (1980) as suggested in Touloumis (2016).

When $\beta_{tj0}=\beta_{j0}$ for all $t$, then intercepts should be provided as a numerical vector. Otherwise, intercepts must be a numerical matrix such that row $t$ contains the category-specific intercepts at the $t$-th measurement occasion.

betas should be provided as a numeric vector only when $\beta_{t}=\beta$ for all $t$. Otherwise, betas must be provided as a numeric matrix with clsize rows such that the $t$-th row contains the value of $\beta_{t}$. In either case, betas should reflect the order of the terms implied by xformula.

The appropriate use of xformula is xformula = ~ covariates, where covariates indicate the linear predictor as in other marginal regression models.

The optional argument xdata should be provided in ``long'' format.

The NORTA method is the default option for simulating the latent random vectors denoted by $e^{O1}_{it}$ in Touloumis (2016). To import simulated values for the latent random vectors without utilizing the NORTA method, the user can employ the rlatent argument. In this case, element ($i,t$) of rlatent represents the realization of $e^{O1}_{it}$.

References

Cario, M. C. and Nelson, B. L. (1997) Modeling and generating random vectors with arbitrary marginal distributions and correlation matrix. Technical Report, Department of Industrial Engineering and Management Sciences, Northwestern University, Evanston, Illinois.

Li, S. T. and Hammond, J. L. (1975) Generation of pseudorandom numbers with specified univariate distributions and correlation coefficients. IEEE Transactions on Systems, Man and Cybernetics 5, 557--561.

McCullagh, P. (1980) Regression models for ordinal data. Journal of the Royal Statistical Society B 42, 109--142.

Touloumis, A. (2016) Simulating Correlated Binary and Multinomial Responses under Marginal Model Specification: The SimCorMultRes Package. The R Journal 8, 79--91.

Touloumis, A., Agresti, A. and Kateri, M. (2013) GEE for multinomial responses using a local odds ratios parameterization. Biometrics 69, 633--640.

Examples

Run this code

## See Example 3.2 in the Vignette.
set.seed(12345)
sample_size <- 500
cluster_size <- 4
beta_intercepts <- c(-1.5, -0.5, 0.5, 1.5)
beta_coefficients <- matrix(c(1, 2, 3, 4), 4, 1)
x <- rep(rnorm(sample_size), each = cluster_size)
latent_correlation_matrix <- toeplitz(c(1, 0.85, 0.5, 0.15))
simulated_ordinal_dataset <- rmult.clm(clsize = cluster_size,
  intercepts = beta_intercepts, betas = beta_coefficients, xformula = ~x,
  cor.matrix = latent_correlation_matrix, link = "probit")
head(simulated_ordinal_dataset$simdata, n = 8)

## Same sampling scheme except that the parameter vector is time-stationary.
set.seed(12345)
simulated_ordinal_dataset <- rmult.clm(clsize = cluster_size, betas = 1,
  xformula = ~x, cor.matrix = latent_correlation_matrix,
  intercepts = beta_intercepts, link = "probit")
## Fit a GEE model (Touloumis et al., 2013) to estimate the regression
## coefficients.
library(multgee)
ordinal_gee_model <- ordLORgee(y ~ x, id = id, repeated = time,
  link = "probit", data = simulated_ordinal_dataset$simdata)
coef(ordinal_gee_model)

Run the code above in your browser using DataLab