An End-to-End GLMM Case Study of Rural Clinic Visits

From Clinic Counts to Comparable Rates

A clinic that sees more visits may serve more patients or may have a higher visit rate. Counts alone cannot distinguish these explanations. In a rural-clinic utilization analysis, the denominator matters: a clinic observed over 1,000 patient-years has had more opportunity to record visits than one observed over 500 patient-years. Clinics are also clusters: repeated observations from one clinic may be more alike than observations from different clinics.

This case study uses hypothetical visit counts from eight rural clinics, each observed in two calendar years. The goal is to estimate annual change while accounting for patient-year exposure and persistent differences among clinics, then report clinic-level rates and comparisons with appropriate uncertainty. The data below are illustrative, not observations from an actual health system.

Assemble and Inspect the Data

Each row in the analysis data represents one clinic-year. Counts are visits recorded during that year; exposure is the corresponding time at risk, in patient-years. The compact table shows the eight clinic records in each year, with exposure of 1,000 patient-years for every clinic-year.

Before fitting a model, verify that counts are nonnegative integers, exposures are positive and expressed in the same units, and clinic identifiers are consistent across years. The total is 1,600 visits in 8,000 patient-years for 2024 and 1,930 visits in 8,000 patient-years for 2025. A crude rate is a useful descriptive check, but it does not account for the repeated observations within clinics.

Fit and Interpret the Exposure-Adjusted Model

The log link makes the mean visit rate multiplicative. Since \(\mu_{jt}/E_{jt}\) is the expected rate per patient-year, the model implies \(\log(\mu_{jt}/E_{jt})=\beta_0+\beta_1 I(t=2025)+b_j\). Thus, for a clinic with \(b_j=0\), \(\exp(\beta_0)\) is its modelled 2024 rate per patient-year, and \(\exp(\beta_1)\) is the 2025-to-2024 rate ratio. The random intercept permits clinics to have different underlying rates, while partially pooling estimates toward the across-clinic distribution.

For a transparent worked analysis, suppose a correctly converged fit to these hypothetical data gives rounded estimates \(\hat\beta_0=-1.609\), \(\hat\beta_1=0.190\), \(SE(\hat\beta_1)=0.080\), and \(\hat\tau^2=0.025\). The intercept corresponds to \(\exp(-1.609)=0.200\) visits per patient-year, or 200 visits per 1,000 patient-years, for a clinic with random effect zero in 2024. These figures are illustrative rounded model output; a fresh software fit can differ slightly because of estimation and rounding.

The Wald calculation uses a large-sample normal approximation for the fixed-effect estimate. It should not be treated as an exact test, particularly with very few clinics. Report the estimate, interval, and model context alongside the p-value; the p-value alone is not a measure of the size or clinical importance of the change.

Compare Clinic-Level Rates

A fitted clinic random effect is a deviation on the log-rate scale. For a specified year, a clinic’s conditional rate per 1,000 patient-years is \(1{,}000\exp(\hat\beta_0+\hat\beta_1 I(t=2025)+\hat b_j)\). A comparison of two clinics in the same year cancels the common fixed year effect: the conditional rate ratio for clinics A and B is \(\exp(\hat b_A-\hat b_B)\). These are conditional, model-based comparisons, not raw rates and not automatically population-average contrasts.

A clinic-effect estimate is typically shrunk toward zero compared with an estimate based only on that clinic’s observed rate. This partial pooling stabilizes estimates, especially for clinics with less exposure, but it means that a model-based clinic rate is not simply that clinic’s count divided by exposure. If many clinic pairs are presented as inferential comparisons, the report should also address multiplicity or clearly label the comparisons as descriptive.

Diagnostics, Assumptions, and Alternatives

The Poisson GLMM assumes that, conditional on the random intercept and fixed effects, the count variance equals the conditional mean. The normal random-intercept distribution is an assumption about between-clinic variation. The offset assumes exposure is measured adequately and that expected counts scale proportionally with exposure. Conditional independence also matters: after accounting for the clinic effect and year, there should not be unmodelled dependence that materially changes the uncertainty.

Check optimizer convergence and inspect warnings, gradients, and estimates near parameter boundaries before interpreting output. Examine observed and simulated counts or simulation-based residuals by clinic and year; investigate systematic residual patterns, unusually many zeros, and clinics with disproportionate influence. Compare observed rates with fitted rates, while remembering that fitted random effects are informed by the same observations. A conspicuous clinic can reflect data error, genuine heterogeneity, or model misspecification; it is not a reason to delete a clinic without justification.

If conditional variation is greater than the Poisson model can explain, a negative-binomial mixed model may be a useful alternative. If there are excess zeros for a substantive reason, consider whether a zero-inflated or hurdle model is appropriate rather than choosing one solely because it improves a fit statistic. A Poisson model without random effects ignores within-clinic dependence and may understate uncertainty. A model with clinic fixed effects instead estimates a separate parameter for each clinic, does not use the same partial pooling, and can be cumbersome when many clinics are present. Each alternative changes the assumptions and often the target of inference; compare fit and scientific purpose, not only a single score.

With only eight clinics, estimates of \(\tau^2\), clinic effects, and their uncertainty may be imprecise. A random-effects distribution can also be a poor description if clinic differences are strongly non-normal or if important clinic-level predictors are omitted. The random intercept does not establish why clinics differ, and the estimated annual association does not show that calendar year caused a change in utilization. The analysis also cannot correct errors in the visit definition, incomplete exposure measurement, or changes in patient composition that have not been represented in the model.

Make the Report Reproducible

A report should let a reader reconstruct the analysis decisions, not merely view a polished table. State the hypothetical or real data source, the unit of analysis, the outcome definition, the exposure unit, the time period, and how missing or excluded records were handled. Specify the model equation, distribution, link, offset, random-effects structure, estimation method, software and package versions, and any convergence warnings. Preserve the analysis data and executable code or a version-controlled source file, and identify the saved model output used for reported estimates.

Present crude totals and rates as descriptive context, followed by adjusted estimates with units, confidence intervals, and the conditioning interpretation. Label the clinic comparisons as model-based conditional rates; include the exposure unit in every rate table. Report diagnostic findings and limitations rather than implying that a successful optimizer proves the model is correct. A useful report also makes clear which values are raw summaries and which are fitted estimates, and avoids excessive decimal precision relative to the data.

Common Mistakes and Final Interpretation

• Comparing raw counts. A larger count can result from greater exposure. Compare rates with the appropriate exposure denominator and model the exposure offset.
• Calling \(\exp(\hat b_j)\) a causal clinic effect. It is a fitted clinic deviation under the model, not an explanation of why the clinic rate differs.
• Ignoring clustering. Treating repeated years from one clinic as independent without a clinic effect can misrepresent uncertainty.
• Reading a conditional estimate as population-average. State whether a rate or contrast conditions on clinic effects; the two targets need not coincide.
• Reporting a p-value as the finding. Give the rate ratio, interval, units, and assumptions, and describe the p-value as evidence relative to the specified null under the model.