Poisson aggregate likelihood needs zero exposure-rate covariance
No, the covariance cancels only when it is common to both arms. The bias followed the identity log{(1 + rho1 CV1(T) CV1(lambda)) / (1 + rho0 CV0(T) CV0(lambda))} (slope 0.994); at plausible dispersions with a differential correlation of at least 0.2, 27 of 64 cells had material bias or undercoverage, with bias up to 0.14 on the log rate ratio and coverage down to 0.17. Integrating aggregate arms over exposure-weighted covariate moments removed the bias everywhere (largest 0.008). Even with zero correlation the unweighted likelihood undercovered at large dispersions (0.84), because the finite-sample covariance is not zero. Read the study below.
Statement
For an aggregate Poisson arm, component ML-NMR multiplies total exposure by the mean modeled rate under the supplied covariate integration distribution. If unweighted participant moments are supplied, that is exact only when exposure has zero covariance with the modeled individual rate; full independence of exposure and covariates is sufficient but not necessary, and treatment-dependent mean follow-up alone is already absorbed by each arm’s own total exposure. The problem is within-arm association between exposure and rate-relevant covariates, which is what differential follow-up by prognosis produces, and publications almost never report exposure-weighted covariate summaries. This is a property of the Poisson aggregate likelihood rather than of any one package: multinma also fits Poisson models with a log time-at-risk offset.
Why it is open
The integration over the aggregate covariate distribution uses reported unweighted moments while the Poisson rate is a person-time-weighted functional of covariates, and sicker patients typically contribute less follow-up, which is exactly the association that breaks the zero-covariance requirement. The correct input, exposure-weighted covariate summaries, is essentially never published, and the current interface neither labels which weighting a supplied moment set represents nor models the weighting itself. There is no diagnostic that tells an analyst whether the omitted covariance is material for their contrast.
What has been tried
| Work | What it contributes |
|---|---|
| multinma 0.8.0 NEWS, Poisson time-at-risk offset | fits Poisson models with a log time-at-risk offset, so the same exposure-weighting requirement applies outside component models |
| cpaic 0.1.0 (development version, choxos/cpaic@9d150e9) | unreleased development R package on GitHub, not on CRAN; computes each aggregate Poisson arm as total exposure times the mean modeled rate over unweighted covariate integration points, with no exposure-weighted moment interface and no joint exposure-covariate model |
Probable solution or research direction
Label the weighting convention in the interface, so an analyst who can obtain exposure-weighted covariate moments can supply them and make the existing likelihood correct, and request those summaries from data holders as a reporting item. Where only unweighted moments exist, model exposure and covariates jointly so the person-time weighting is generated by the model rather than assumed away. Add a sensitivity calculation for the magnitude of the exposure-rate covariance that would move a target contrast, since the severity is entirely conditional on that quantity.
Our study
The aggregate Poisson likelihood is biased by the exposure-rate covariance, and exposure-weighted covariate moments remove the bias
The question put to it. The aggregate Poisson likelihood multiplies an arm’s total exposure by its mean modeled rate, which is exact only when exposure and rate are uncorrelated. The error is the factor 1 + rho CV(T) CV(lambda) per arm. How much bias does that put on a transported rate ratio at plausible levels of differential follow-up, and do exposure-weighted covariate summaries remove it?
What it found. No, the covariance cancels only when it is common to both arms. The bias followed the identity log{(1 + rho1 CV1(T) CV1(lambda)) / (1 + rho0 CV0(T) CV0(lambda))} (slope 0.994); at plausible dispersions with a differential correlation of at least 0.2, 27 of 64 cells had material bias or undercoverage, with bias up to 0.14 on the log rate ratio and coverage down to 0.17. Integrating aggregate arms over exposure-weighted covariate moments removed the bias everywhere (largest 0.008). Even with zero correlation the unweighted likelihood undercovered at large dispersions (0.84), because the finite-sample covariance is not zero.
Results.
- Identity check over 208 cells: slope 0.994 (SE 0.002), largest deviation 2.9 MCSE.
- Exposure-weighted moments: bias at most 0.008, coverage 0.908 to 0.965; borrowed exposure slope: bias at most 0.015 but coverage down to 0.69 because the slope’s uncertainty is ignored.
- Null control failed on coverage (0.838 to 0.957 at zero correlation), not bias; positive control’s cell selection corrected after the run; sensitivity range did not contain cells with positive treated-arm correlation (coverage down to 0.71).
Design.
- Estimand
- Target-population log marginal rate ratio, B versus C
- Methods compared
- aggregate integration over unweighted moments, exposure-weighted moments, moments tilted by a borrowed exposure slope, sensitivity interval over a declared correlation range
- Data-generating mechanism
- A minimal anchored network, individual data on A versus C and arm totals on B versus C, with lognormal exposure solved so each arm’s exposure-rate correlation and dispersions take registered values exactly. 234 attainable cells crossing the control-arm correlation, its difference in the treated arm, CV of exposure, CV of the rate and arm size; unweighted, exposure-weighted and borrowed-slope integration and a sensitivity bound. 1000 replicates per cell.
- Replicates
- 1000 per cell, 208 main cells and 26 control cells
What this does not answer. One covariate, one lognormal exposure law, no overdispersion, the same exposure mechanism in both studies; the likelihood is the study’s own code, not multinma. Peer review has not been done.
Peer review. Not yet done. The result stands on the study alone until it is.
Read it: PDF · Markdown · OpenDocument · Code and data: studies/CMP-20-poisson-exposure-rate-covariance · Protocol: registered before the run
Verification
Verdict. Unverifiable. Auditors disagreed and the disagreement was evidence-bearing on both sides.
The source titled this. “Poisson AgD likelihood assumes exposure independent of covariates” The title states independence of exposure and covariates as the required condition, but the audit showed the necessary condition is zero covariance between exposure and the modeled rate.
The source said. “The Poisson aggregate-data likelihood in component ML-NMR assumes exposure time is independent of covariates. When exposure depends on prognosis or on treatment, the exposure-weighted covariate moments that actually drive the rate differ from the unweighted moments reported in publications.”
What was wrong with it. The source stated independence of exposure and covariates as the required condition and framed the issue as package-specific. GPT-5.6 Sol read the Poisson Stan code at choxos/cpaic@9d150e9 and showed the necessary condition is zero covariance between exposure and the modeled rate, not independence, and that supplying exposure-weighted covariate summaries through the existing moment fields makes the same likelihood correct; the literature auditor noted the same requirement applies to multinma’s Poisson models, so this is field-level rather than a defect of one implementation.
Dissent. codex dissented (overstated): C2 and C3 are supported, but C1 states a sufficient condition as though it were necessary. The AgD Stan likelihood computes total exposure multiplied by the mean modeled rate under the supplied covariate integration distribution. If that distribution represents unweighted participants, exactness requires mean(e_i exp(eta_i)) = mean(e_i) mean(exp(eta_i)), equivalently zero covariance between exposure and the modeled rate. Full independence of exposure and covariates is sufficient but not necessary. Treatment-dependent mean follow-up alone is handled by each arm’s own total E; the problem arises from within-arm association between exposure and rate-relevant covariates. Conversely, supplying an exposure-weighted covariate distribution through the existing moment fields can make the same likelihood correct even though the interface does not label or model that weighting. No joint exposure-covariate model is implemented. Severity is high only when the omitted covariance is material.
- Literature and prior-art check
- open / supported-with-caveat. The statistical point is correct and non-trivial: a Poisson rate is a person-time-weighted functional of covariates, published covariate moments are person-weighted, and the two coincide only under independence of exposure and covariates, which is exactly what differential follow-up by prognosis violates. I found no paper addressing exposure-weighted covariate moments in an aggregate-data integrated likelihood. Caveat for the editor: this is a field-level property of the Poisson aggregate likelihood, not a cpaic defect. multinma also fits Poisson models with a log time-at-risk offset (its 0.8.0 NEWS records a bug fix to that offset), and the same independence assumption applies there, so framing it as ‘cpaic’s Poisson aggregate-data likelihood assumes…’ understates the reach of the problem. The claim that joint exposure and covariate models are not implemented in cpaic is a code fact I cannot verify from public sources.
- GPT-5.6 Sol, technical and source-code lens
- open / overstated. C2 and C3 are supported, but C1 states a sufficient condition as though it were necessary. The AgD Stan likelihood computes total exposure multiplied by the mean modeled rate under the supplied covariate integration distribution. If that distribution represents unweighted participants, exactness requires mean(e_i exp(eta_i)) = mean(e_i) mean(exp(eta_i)), equivalently zero covariance between exposure and the modeled rate. Full independence of exposure and covariates is sufficient but not necessary. Treatment-dependent mean follow-up alone is handled by each arm’s own total E; the problem arises from within-arm association between exposure and rate-relevant covariates. Conversely, supplying an exposure-weighted covariate distribution through the existing moment fields can make the same likelihood correct even though the interface does not label or model that weighting. No joint exposure-covariate model is implemented. Severity is high only when the omitted covariance is material.
- Grok 4.5, recency lens
- open / supported. cpaic man/cmlnmr.Rd and NEWS state that Poisson aggregate arms use E * mean(exp(eta)) over unweighted covariate integration points and assume exposure independent of covariates; there is no interface for exposure-weighted moments or a joint exposure–covariate model. Literature search for Poisson AgD ML-NMR with person-time depending on covariates did not surface an evaluated method that closes this for component ML-NMR. The statistical claim that exposure-weighted moments differ from unweighted moments under dependence is standard rate epidemiology, not a solved software gap. Queries: Poisson aggregate exposure covariate independent ML-NMR; cpaic poisson exposure; person-time weighted covariate moments NMA. All three claims remain supported; problem open.
- Decision path
-
R5:default-open → R7:counter-and-support-both-present
Source
Derived from C:L1378, scope:package-specific implementation resting on a field-level rate-outcome problem of the reviewed corpus.