depower (version 2025.1.20)

glmm_poisson: GLMM for Poisson ratio of means

Description

Generalized linear mixed model for two dependent Poisson outcomes.

Usage

glmm_poisson(data, test = "wald", ci_level = NULL, ...)

Value

A list with the following elements:

Slot	Subslot	Name	Description
1		chisq	\(\chi^2\) test statistic for the ratio of means.
2		df	Degrees of freedom.
3		p	p-value.
4		ratio	Estimated ratio of means (sample 2 / sample 1).
4	1	estimate	Point estimate.
4	2	lower	Confidence interval lower bound.
4	3	upper	Confidence interval upper bound.
5		mean1	Estimated mean of sample 1.
5	1	estimate	Point estimate.
5	2	lower	Confidence interval lower bound.
5	3	upper	Confidence interval upper bound.
6		mean2	Estimated mean of sample 2.
6	1	estimate	Point estimate.
6	2	lower	Confidence interval lower bound.
6	3	upper	Confidence interval upper bound.
7		item_sd	Estimated standard deviation of the item (subject) random intercept.
7	1	estimate	Point estimate.
7	2	lower	Confidence interval lower bound.
7	3	upper	Confidence interval upper bound.
8		n1	Sample size of sample 1.
9		n2	Sample size of sample 2.
10		method	Method used for the results.
11		test	Type of hypothesis test.
12		alternative	The alternative hypothesis.
13		ci_level	Confidence level of the interval.
14		hessian	Information about the Hessian matrix.
15		convergence	Information about convergence.

Arguments

data

(list)
A list whose first element is the vector of Poisson values from sample 1 and the second element is the vector of Poisson values from sample 2. Each vector must be sorted by the subject/item index and must be the same sample size. NAs are silently excluded. The default output from sim_bnb().

test

(String: "wald"; "c("wald", "lrt"))
The statistical method used for the test results. test = "wald" performs a Wald test and optionally the Wald confidence intervals. test = "lrt" performs a likelihood ratio test and optionally the profile likelihood confidence intervals (means and ratio). The Wald interval is always used for the standard deviation of the item (subject) random intercept.

ci_level

(Scalar numeric: NULL; (0, 1))
If NULL, confidence intervals are set as NA. If in (0, 1), confidence intervals are calculated at the specified level. Profile likelihood intervals are computationally intensive, so intervals from test = "lrt" may be slow.

...

Optional arguments passed to glmmTMB::glmmTMB().

Details

Uses glmmTMB::glmmTMB() in the form

glmmTMB(
  formula = value ~ condition + (1 | item),
  data = data,
  family = stats::poisson
)

to model dependent Poisson outcomes \(X_1 \sim \text{Poisson}(\mu)\) and \(X_2 \sim \text{Poisson}(r \mu)\) where \(\mu\) is the mean of sample 1 and \(r\) is the ratio of the means of sample 2 with respect to sample 1.

The hypotheses for the LRT and Wald test of \(r\) are

$$ \begin{aligned} H_{null} &: log(r) = 0 \\ H_{alt} &: log(r) \neq 0 \end{aligned} $$

where \(r = \frac{\bar{X}_2}{\bar{X}_1}\) is the population ratio of arithmetic means for sample 2 with respect to sample 1 and \(log(r_{null}) = 0\) assumes the population means are identical.

When simulating data from sim_bnb(), the mean is a function of the item (subject) random effect which in turn is a function of the dispersion parameter. Thus, glmm_poisson() has biased mean estimates. The bias increases as the dispersion parameter gets smaller and decreases as the dispersion parameter gets larger. However, estimates of the ratio and standard deviation of the random intercept tend to be accurate. In summary, the Poisson mixed-effects model fit by glmm_poisson() is not recommended for the BNB data simulated by sim_bnb(). Instead, wald_test_bnb() or lrt_bnb() should typically be used instead.

References

hilbe_2011depower

hilbe_2014depower

See Also

glmmTMB::glmmTMB()

Examples

#----------------------------------------------------------------------------
# glmm_poisson() examples
#----------------------------------------------------------------------------
library(depower)

set.seed(1234)
d <- sim_bnb(
  n = 40,
  mean1 = 10,
  ratio = 1.2,
  dispersion = 2
)

lrt <- glmm_poisson(d, test = "lrt")
lrt

wald <- glmm_poisson(d, test = "wald", ci_level = 0.95)
wald

#----------------------------------------------------------------------------
# Compare results to manual calculation of chi-square statistic
#----------------------------------------------------------------------------
# Use the same data, but as a data frame instead of list
set.seed(1234)
d <- sim_bnb(
  n = 40,
  mean1 = 10,
  ratio = 1.2,
  dispersion = 2,
  return_type = "data.frame"
)

mod_alt <- glmmTMB::glmmTMB(
  formula = value ~ condition + (1 | item),
  data = d,
  family = stats::poisson,
)
mod_null <- glmmTMB::glmmTMB(
  formula = value ~ 1 + (1 | item),
  data = d,
  family = stats::poisson,
)

lrt_chisq <- as.numeric(-2 * (logLik(mod_null) - logLik(mod_alt)))
lrt_chisq
wald_chisq <- summary(mod_alt)$coefficients$cond["condition2", "z value"]^2
wald_chisq

anova(mod_null, mod_alt)