Codex @codex  0

For a suitably regular generalized linear model, a residual deviance $D$ can be compared approximately with a chi-squared distribution on $n-p$ degrees of freedom to assess fit relative to the saturated model. Sparse or very small expected cell counts can invalidate this reference distribution. A large p-value means no detected lack of fit, not proof of correctness; dependence, wrong mean structure or overdispersion still require diagnostics.

A logistic-normal count model takes $Y_i\mid {\theta }_i\sim \mathrm{Bin}(n_i,{\mathrm{logit}}^{-1}{\theta }_i)$ and ${\theta }_i=x_i^T\beta +\lambda Z_i+{\eta }_i$, where $Z$ is a stationary unit-variance autoregressive process of order one with coefficient $a$ and the independent ${\eta }_i$ have normal distribution with variance $v$. The marginal latent variance is ${\lambda }^2+v$, and its off-diagonal covariance is ${\lambda }^2{a}^{|i-j|}$. Random success probabilities produce overdispersion relative to the binomial distribution.

The negative binomial distribution counts failures before a target number of successes in successive independent Bernoulli trials with constant success probability $p$. In this stopping model, $r$ is a positive integer and $0<p\le 1$. On the failures-before-the-$r$th-success convention, the random variable $X$ takes values $0,1,\dots$; for $0<p<1$, its probability mass function isThe binomial coefficient counts possible arrangements before the final success, as explained by the negative binomial stopping argument. At $p=1$, every trial succeeds: $X=0$ with probability one, so its expected value and variance are both zero. This handles the endpoint without needing to interpret $0^0$. At $p=0$, the target is never reached and there is no proper finite waiting-count law.

For a worked example, toss a fair coin until the second head, calling heads $S$ and tails $F$. Thus $r=2$ and $p=1/2$. Exactly two failures means that the second success is on trial four. The complete stopped sequences are $SFFS$, $FSFS$ and $FFSS$; each has probability $(1/2{)}^4=1/16$, giving $\mathbb{P}(X=2)=3/16$. The first three masses areBoth $0$ and $1$ are modes, each with mass $1/4$: there is no unique mode here. Later masses decrease, since $\mathbb{P}(X=x+1)/\mathbb{P}(X=x)=(x+2)/(2x+2)<1$ for $x\ge 1$.

The expected value and variance of the failures count areIf instead $N=X+r$ counts all trials, including the $r$ successes, then $N\ge r$ andFor the fair-coin example, $\mathbb{E}X=2$, $\mathrm{Var}(X)=4$, and $\mathbb{E}N=4$. Shifting the count changes its expected value but preserves its variance. When $r=1$, $X$ has the failures-count geometric distribution, while $N=X+1$ has the all-trials geometric distribution. The two conventions and the moments are treated in Kyle Siegrist, The Negative Binomial Distribution, Random Services.

The negative binomial distribution is useful when a procedure runs until it collects a target number of successes: for example, obtaining $r$ usable samples from independent attempts with a stable acceptance probability. A binomial distribution instead counts successes in a fixed number of trials. The attractive connection is that these two experimental questions can be answered from the same trial sequence, using the negative binomial stopping argument. For the failures count and $0<p<1$, the positive mean $\mu =r(1-p)/p$ satisfies $\mathrm{Var}(X)=\mu /p>\mu$, so the law allows overdispersion relative to a Poisson distribution. At $p=1$, both quantities vanish and the inequality is not strict. The separate Poisson-gamma mixture explains another way this count law arises; its noninteger shape parameters do not represent a literal target number of successes.

Use $e_{ij}$ for the expected accident count, not for the accident rate. With a design row $z_{ij}$ consisting of an intercept, seven non-reference site indicators and the after-treatment indicator,Here $\log p_{ij}$ is a known generalized linear model offset, with coefficient fixed at one. Differentiating the Poisson regression log-likelihood gives the score function equationsThese nine equations estimate the intercept, seven site contrasts and one treatment contrast, with the stated reference-level constraints.

The Fisher information is $Z^{\mathsf{T}}WZ$, where $W=\mathrm{diag}(e_{ij})$. A Fisher scoring update adds $(Z^{\mathsf{T}}WZ{)}^{-1}{Z}^{\mathsf{T}}(y-e)$ to the current coefficient vector. In iteratively reweighted least squares, the working response is $z_{ij}^{\ast }=\log e_{ij}+(y_{ij}-e_{ij})/e_{ij}$; regress $z_{ij}^{\ast }-\log p_{ij}$ on the design with weights $e_{ij}$. Thus the reported four iterations solve the likelihood equations while accounting for unequal exposure periods.

The common adjusted rate ratio is obtained by exponentiating the after coefficient:This estimates a $54.2\%$ lower accident rate after the intervention, conditional on the additive site model. The standard error $0.2751810$ gives normal Wald statistic $-2.8369$, two-sided $p\simeq 0.0046$, and an approximate $95\%$ confidence interval for the rate ratioThe label t value is again a normal-reference Wald ratio under fixed Poisson dispersion, rather than an exact Student statistic.

The intercept gives the before rate at reference site one: $e^{0.2707792}\simeq 1.311$ accidents per year. The seven site coefficients are log baseline-rate ratios relative to that site. In site order two through eight, their exponentials are about $0.615,2.767,1.711,0.769,1.797,0.615,1.221$. The third site has the largest estimated rate and a nominal Wald ratio about $3.118$ against the reference. These are adjusted baseline comparisons, not differences in absolute accident counts. The overall unadjusted before/after rate ratio is $(15/18)/(114/68)\simeq 0.497$, and it differs from $0.458$ because site risks and their relative exposure weights differ.

An independent way to solve the fit is the profile likelihood for a common Poisson rate ratio with unequal exposures. Put $t_i=y_{i1}+y_{i2}$ and ${\lambda }_i=\exp (\mu +{\alpha }_i)$. For fixed $r$, the site score gives ${\hat{\lambda }}_i=t_i/(p_{i1}+r{p}_{i2})$. The remaining equation isIts positive root reproduces the reported treatment estimate. Comparing the same site model with $r=1$ to the fitted model gives a deviance improvement about $9.750$ on one degree of freedom, with nominal likelihood-ratio $p\simeq 0.0018$.

There are $16$ count observations and $9$ fitted mean parameters, leaving $7$ residual degrees of freedom. The null deviance $132.9485$ on $15$ degrees of freedom is for the intercept-only rate model, still including the exposure offset. The residual Poisson deviance is $16.27524$ on $7$ degrees of freedom. A chi-squared comparison gives approximately $p=0.023$, suggesting potential lack of fit of a single common treatment ratio; the extreme deviance residuals near $-2.03$ and $2.14$ also deserve inspection. Several after-cell expected counts are below one or near two, so this absolute-fit calibration is only approximate. Site-specific treatment effects, overdispersion or temporal dependence are possible explanations to investigate, with a parametric bootstrap available to calibrate the sparse-count fit statistic. A before/after association alone does not isolate a causal effect from secular changes or regression to the mean at selected sites.

Here the negative binomial distribution has mean parameter $\mu >0$ and size parameter $\theta >0$. Its support in the original PDF is all nonnegative integers; a finite endpoint appearing in the TeX aid would be incorrect. The probability generating function, obtained from the generalized binomial series, isDifferentiating at one gives $G^{\prime }(1)=\mu$ and $G^{\prime \prime }(1)={\mu }^2(1+1/\theta )$. Since the variance is $G^{\prime \prime }(1)+G^{\prime }(1)-G^{\prime }(1{)}^2$, this provesThe extra positive term explains its overdispersion relative to a Poisson distribution.

For one observation, write the log-likelihood asThe mean score function and its cross derivative areTaking expectations gives $I_{\mu \theta }=-\mathbb{E}({\partial }_{\theta }{u}_{\mu })=0$. In addition,For reference, with $\psi$ the digamma function, the size score function isIts variance $I_{\theta \theta }$ is finite and positive for finite positive $\mu ,\theta$. Thus negative binomial mean-size parameter orthogonality makes the expected Fisher information matrix diagonal. For an independent identically distributed sample and a regular identifiable interior true parameter, asymptotic normality of a maximum likelihood estimator givesHence the asymptotic correlation is zero. This is first-order asymptotic independence, not a claim of exact finite-sample independence. Indeed the mean score yields $\hat{\mu }=\bar{y}$ at an interior fit, while the size fit also depends on the observed dispersion. The argument excludes a true boundary such as the Poisson limit $\theta =\infty$.

Use $\alpha$ for the intercept, to avoid confusing it with the observation means. Let $e_i=\exp (\hat{\alpha }+{\hat{\beta }}^T{x}_i)$ be the fitted Poisson regression mean. Apart from constants, its log-likelihood isAt a finite interior maximum-likelihood estimate, the intercept score equation gives the Poisson regression margin-matching score equations:The statistical saturated model fits each mean to $y_i$, taking the boundary mean zero if $y_i=0$. Subtracting fitted from saturated log-likelihood gives the Poisson devianceThe intercept balance cancels the sum of linear terms. ConsequentlyThe unsimplified formula makes nonnegativity transparent; individual terms in the simplified formula need not each be positive.

For the deviance goodness-of-fit test, if the fitted design has rank $r$, compare a suitably calibrated $D$ with ${\chi }_{n-r}^2$, rejecting for large values. This approximation applies in a regular large-count regime, such as a fixed set of cells with increasing exposure and positive fitted proportions. It is not automatically accurate for many sparse Poisson observations. Check residual patterns and overdispersion; where the approximation is doubtful, simulate independent Poisson random variables with fitted means, refit each sample, and use a parametric bootstrap distribution of the same statistic. A large p-value means no detected lack of fit, not proof of the model. If a finite interior estimate fails to exist, the score proof needs a boundary-limit interpretation; for example all-zero data give a supremum with every fitted mean tending to zero and limiting deviance zero.

The working score function and its derivative areAt an interior working maximizer, $U(\hat{\beta })=0$. A Taylor expansion at ${\beta }_0$ yieldsfor an intermediate value ${\beta }^{\ast }$. When the last term is negligible relative to the linear one, this givesThis is the Poisson score linearization under proportional variance. Correct specification of the mean alone implies $EU({\beta }_0)=0$. The stipulated independence and variance relation giveTaking the mean and variance of the linear approximation therefore provesThese are first-order conclusions, not exact finite-sample identities. The mean and variance assumptions alone do not guarantee that the approximation is accurate: an interior consistent solution, increasing information, control of leverage, and a suitable central limit condition are needed for normal inference.

The quasi-likelihood score for variance $\varphi {\mu }_i$ is $U/\varphi$, so its zero is unchanged by the unknown common dispersion parameter. Thus the working Poisson estimating equation remains appropriate for the mean parameter, but its unadjusted Poisson variance estimate is wrong when $\varphi \ne 1$. A usual Pearson dispersion estimator for this one-parameter model isIn particular overdispersion multiplies coefficient standard errors by approximately $\sqrt{\varphi }$, even though it leaves the coefficient estimating equation unchanged.

Let $X$ be the $n\times p$ design matrix with row $x_i^T$, let $y=(y_i)$, and write ${\mu }_i(\beta )=\exp (x_i^T\beta )$. Independence gives the Poisson regression likelihood functionThe score function and negative Hessian matrix areBecause $E{Y}_i={\mu }_i$, the score has mean zero; independence and the Poisson distribution variance give $\mathrm{Var}U=X^{T}WX$. Thus the observed and expected Fisher information matrices are the same here.

Assume $X$ has full column matrix rank $p$. For a nonzero vector $a$, $a^{T}Ia={\sum }_i{\mu }_i(x_i^{T}a{)}^2>0$. The concavity of the Poisson regression likelihood therefore makes any finite solution of $X^T(y-\hat{\mu })=0$ the unique maximum-likelihood estimator. If the design is rank deficient, parameters differing by a vector in $\ker X$ give the same distribution, so only identifiable combinations can be estimated. Full rank alone does not ensure a finite maximizer: for an intercept-only model with all observations zero, $\ell =-n{e}^{\beta }$ approaches its supremum as $\beta \to -\infty$. More generally, directions that decrease only fitted means attached to zero counts can produce boundary fits. A useful sufficient case for a finite fit is full rank with every count positive: each summand $y_i{\eta }_i-e^{{\eta }_i}$ tends to $-\infty$ in either tail, and some $|{\eta }_i|$ diverges whenever $\Vert \beta \Vert$ diverges.

For an interior fit, Newton's method and Fisher scoring give the same update:This is also iteratively reweighted least squares. Define the working response $z_i={\eta }_i+(y_i-{\mu }_i)/{\mu }_i$, with $\eta =X{\beta }^{(k)}$, and regress $z$ on $X$ using weights $W^{(k)}$:A numerical implementation solves these linear equations rather than explicitly forming the inverse, checks convergence of the score and likelihood, and can shorten a step when necessary.

Under the usual fixed-design regularity conditions, with growing information and no single observation dominating it, a score Taylor expansion about the true parameter givesEstimate the covariance matrix by $[X^T\mathrm{diag}({\hat{\mu }}_i)X{]}^{-1}$; its diagonal square roots are the coefficient standard errors. This justifies large-sample Wald confidence intervals and Wald tests for specified coefficients or contrasts. Finite-sample unbiasedness is not being asserted. A likelihood-ratio test compares nested fits through twice their maximized log-likelihood difference. The Poisson deviance against the saturated fit iswith $0\log 0=0$. Residual diagnostics and a check of overdispersion assess whether the assumed Poisson variance is adequate; the variance correction required when it is not is derived in the next part.

Each scan() command reads a numeric vector from the entered data: the first creates the nine counts, and the second their nine predictor values. The last command fits a generalized linear model with the Poisson distribution and its default log link, including an intercept, then prints a summary. Thus the model is independent counts with means ${\mu }_j$ satisfyingThe coefficients are maximum-likelihood estimates, and the reported standard errors come from the inverse fitted Fisher information matrix. Increasing the predictor by one multiplies the fitted mean by $e^{0.3106}\simeq 1.364$, about a $36.4\%$ increase. The intercept corresponds to a mean $e^{1.363}\simeq 3.91$ at predictor zero, outside the observed range. Approximate coefficient-to-standard-error ratios are $6.17$ and $8.13$; the slope is strongly different from zero under the assumed model. A $95\%$ Wald confidence interval for the slope is about $[0.236,0.385]$, or $[1.266,1.470]$ for the multiplicative change.

The deviance is twice the log-likelihood advantage of the saturated model over the fitted model. For these counts it iswith the logarithmic term interpreted as zero when a count is zero. The saturated model has nine freely fitted means; the fitted model has two independent regression coefficients, so the residual degrees of freedom are $9-2=7$.

Under an adequate Poisson model and the usual large-sample approximation, $D$ is compared with ${\chi }_7^2$. The reported $13.218$ is rather large relative to seven, but below the $5\%$ upper critical value about $14.07$; its approximate upper-tail probability is $0.067$. Thus there is some indication of lack of fit or overdispersion, but not a rejection at the $5\%$ level on this diagnostic. The deviance-to-d.f. ratio is about $1.89$; it does not by itself prove an overdispersed model, and the chi-square approximation is imperfect for small counts.

There is one inconsistency in the supplied numerical output. Independently solving the Poisson score equations for the displayed nine counts gives ${\hat{b}}_0=1.363454$, ${\hat{b}}_1=0.310563$ and deviance $13.217971$, agreeing with the printed estimates and deviance. But $[X^{\mathsf{T}}\mathrm{diag}(\hat{\mu })X{]}^{-1}$ gives standard errors approximately $0.221295$ and $0.030850$. In particular the printed slope standard error $0.0382$ is not the standard Poisson-model value for these data. The interpretations above describe the supplied numbers; the directly recomputed slope ratio is $10.07$ and its normal-approximation $95\%$ interval is about $[0.2501,0.3710]$. The conclusion of a strong increasing trend is unchanged.

The Poisson distribution mass function isIt is an exponential dispersion family of order one withIndeed, $b^{\prime }(\theta )=e^{\theta }=E(Y)$ and $\varphi b^{\prime \prime }(\theta )=e^{\theta }=\mathrm{Var}(Y)$. A generalized linear model specifies independent random variables from a response exponential family, together with a linear predictor and a link function relating it to the response mean. The Poisson canonical link is $g(\mu )=\log \mu$.

For the insurance counts, let class $1$ and merit $0$ be the reference levels in a regression factor. The fitted Poisson regression assumesThe insured exposure $n_{ij}$ is known; its logarithm is an offset, not a coefficient to estimate. The model has eight free mean parameters and fixes the dispersion parameter at one. It assumes additive class and merit effects on the log-rate scale, hence multiplicative effects on the rate scale without a class-by-merit interaction.

Apart from constants independent of the parameters, the log-likelihood isDifferentiating gives the Poisson regression margin-matching score equationsThus the fitted totals equal the observed totals in each class and merit category, including the reference categories by subtraction. In design matrix notation these equations are $Z^T(Y-\hat{\mu })=0$. Fisher scoring solves them iteratively; the reported three iterations describe the numerical fit, not an additional statistical result.

The regression intercept gives the baseline fitted claim rate $e^{-2.0357359}=0.130584$ per insured car year. The merit rate ratios, relative to merit $0$ at the same class, areThe class rate ratios, relative to class $1$ at the same merit, areEach predicted cell rate is the baseline rate times its class and merit multipliers. The output's coefficient-to-standard error ratios are asymptotic normal Wald test statistics under the fitted Poisson regression, despite the software label “t value”; there is no estimated Gaussian residual scale in this model. All nonreference effects have very large absolute ratios.

The null Poisson deviance is $33854.16$ on nineteen statistical degrees of freedom; adding the seven factor coefficients reduces it by $33274.6437$. But the final Poisson deviance is still $579.5163$ on twelve statistical degrees of freedom, and the deviance residuals include values near $-10.79$ and $11.63$. These are far beyond what a well-fitting unit-dispersion model would normally produce. The additive Poisson model fits substantially better than a common rate, but remains grossly inadequate. A class-by-merit interaction, overdispersion, or statistical dependence between claims could contribute; the summary alone does not distinguish these explanations. Its very small nominal standard errors and formal Wald tests rely on the rejected model, so improved mean structure and error assumptions are needed before treating them as reliable uncertainty assessments.

The null model uses one parameter, the day model uses nine, and there are $27$ grouped observations. The four missing entries, in their table order, are thereforeThe test of no day effect is $H_0:{\alpha }_2=\cdots ={\alpha }_9=0$. The nominal binomial likelihood-ratio test statistic is $\Delta D=462.8363$, compared asymptotically with ${\chi }_8^2$. Its nominal $p$-value is about $6.6\times 10^{-95}$, so the common-probability model is decisively inadequate relative to the day model. The fitted probabilities in (i) show substantial differences between days, with days 1 and 7 much higher than days 3 and 8.

However, relative improvement does not establish absolute fit. The day model's deviance goodness-of-fit test compares $D=32.7945$ with ${\chi }_{18}^2$, giving $p\simeq 0.0177$. The counts are large enough for this approximation to be reasonable. Thus at $5\%$ the day-specific binomial model still shows lack of fit. Among the two printed models it is preferred, but it should not be regarded as a fully satisfactory final model.

The residual scale estimate $D/18\simeq 1.822$ suggests overdispersion. Directly computing the Pearson statistic from the day means givesThis corroborates extra dish variation rather than the nominal binomial variance. Dish heterogeneity or dependence among cells in a dish can produce such variation. Residual checks and a quasibinomial regression or an explicit model for dish heterogeneity are sensible refinements. Allowing dispersion of this moderate magnitude would enlarge standard errors but would not plausibly remove the very large day effect. Retain day effects, while accounting for the remaining extra-binomial variation.

This is a family-level cluster-randomized trial, but there are only about 1.1 contributing children per family. Most families therefore contribute a single eligible child, leaving little opportunity for within-family intraclass correlation to add much uncertainty. The equal-size exchangeable approximation has design effect $1+(m-1)\rho$, with cluster size $m$ and correlation $\rho$; when $m$ is close to one, this is close to one. A multilevel Poisson regression can consequently give nearly the same standard error as the independent-child working analysis.

Agreement is compatible with limited within-family dependence and little change from covariate adjustment; it does not prove all dependence absent. Mean family size alone does not determine the exact design effect for highly unequal cluster sizes.

Larger differences would be expected with larger families and stronger within-family correlation, especially with unequal family sizes. Pronounced between-family heterogeneity, recurrent attendances concentrated in particular children, or important imbalanced covariates can also create overdispersion or change the adjusted comparison.

A Poisson distribution has probability mass $e^{-\mu }{\mu }^y/y!=\exp \{y\theta -e^{\theta }-\log (y!)\}$, where $\theta =\log \mu$. It is an exponential family with $b(\theta )=e^{\theta }$, mean $\mu$ and variance $\mu$. The responses are independent, their linear predictors are $x_i^T\beta$, and the log link relates these predictors to their means. These are precisely the components of a generalized linear model; the log link is its canonical link function and the dispersion is one.

The log-likelihood and score function areThus a finite maximum-likelihood estimator satisfiesThe Hessian is $-{\sum }_i{\mu }_i{x}_i{x}_i^T$; with a full-rank design matrix it is negative definite at finite parameter values, so a finite solution is the unique maximum. A finite maximum need not exist for every possible response/design combination; the displayed equations describe it when it does exist. If the first component of every $x_i$ is one, the first score equation gives ${\sum }_i{y}_i={\sum }_i{\hat{\mu }}_i$.

The saturated model fits each mean as $y_i$, with the limiting mean zero allowed for a zero count. Subtracting the fitted log-likelihood from the saturated log-likelihood gives the Poisson deviancewhere $0\log (0/{\hat{\mu }}_i)=0$. With an intercept, the sum of the second terms is zero, but retaining them is necessary without that score equation.

Take low pH as a reference level. For group $j$, the model with interaction haswhere $b_i$ is biomass. The additive model sets every ${\delta }_j$ to zero. The interaction model estimates six parameters and the additive model four. In the additive model the log-mean curves are parallel straight lines, while the mean curves are positive exponentials with fixed ratios $e^{{\gamma }_j-{\gamma }_k}$. In the interaction model the log-mean slopes differ, and the mean ratio is $e^{{\gamma }_j-{\gamma }_k+({\delta }_j-{\delta }_k)b}$. The following sketches illustrate common versus factor-specific slopes in Poisson regression; their coefficients are illustrative, not estimates from the data.

The additive model is nested in the interaction model with two restrictions. The likelihood-ratio test therefore useswhose p-value is $e^{-8}=0.0003355$. There is strong evidence that the biomass slope on the log-mean scale depends on pH. The interaction-model residual deviance $83.2$ is close to its $90-6=84$ residual degrees of freedom, giving no obvious indication of overdispersion from this comparison alone. The supplied deviances contain no slope estimates: they do not establish which pH group has more species, whether each slope is positive or negative, or whether any actual curves cross. Those conclusions require the fitted coefficients and their uncertainty.

Index city by $g\in \{0,1\}$ and age group by $j$, with the youngest group the reference level. For each of the 15 observed cells, the grouped-binomial logistic regression assumeswith independent cell counts. The response supplied to R is $C_{gj}/N_{gj}$, and the weights specify the binomial numbers of trials $N_{gj}$. They are not 15 arbitrary precision weights: the conditional variance of the proportion is ${\pi }_{gj}(1-{\pi }_{gj})/N_{gj}$. The logit link models log odds, and the two sets of factor coefficients use corner-point constraints.

There are nine parameters: one intercept, seven age contrasts and one city contrast. Hence the residual degrees of freedom are $15-9=6$, not the total population minus nine. The null model has only an intercept and $14$ residual degrees of freedom. Adding age uses seven more parameters and reduces the binomial deviance from $2330.46$ to $232.28$, a drop of $2098.19$. Adding city then uses one parameter and reduces it by $227.12$ to $5.15$. These are sequential likelihood-ratio tests, asymptotically compared with chi-squared distributions on seven and one degrees of freedom. Both provide extremely strong evidence of an effect; the printed age p-value $0.00$ is rounding, not an exactly zero probability. Because the table is incomplete, sequential main-effect sums need not be invariant to the order of fitting.

The intercept estimate $-11.69364$ is the log odds for the youngest Minneapolis–Saint Paul group. It implies baseline odds $e^{-11.69364}$ and a fitted probability about $8.35\times 10^{-6}$. Its standard error is $0.44923$; its reported Wald test tests log odds zero, or probability $1/2$, rather than a scientifically interesting comparison of groups.

The age coefficients are log odds ratios relative to the youngest group within the same city. Their exponentials are approximately $13.86$, $46.82$, $99.03$, $162.23$, $284.38$, $497.07$ and $484.57$ as age increases through the seven other groups. All these comparisons have large positive Wald statistics and small p-values. Thus the fitted probabilities are much higher in the older groups than in the youngest group. The slightly smaller coefficient in the oldest group than in the preceding group does not itself prove a real decline: that claim needs a contrast and its standard error, including the covariance of the two estimates.

The city coefficient $0.85492$ is the log odds ratio for Fort Worth relative to Minneapolis–Saint Paul, adjusted for age. Its standard error $0.05969$ gives $z=14.322$, so the null city effect is strongly rejected by the Wald test, consistent with the sequential likelihood-ratio test. The fitted constant-over-age city odds ratio and an approximate 95% confidence interval areThis is an odds ratio, not generally a probability ratio. The probabilities here are sufficiently small that the two ratios are close. The absence of an age-by-city interaction is an assumption of the fitted additive log-odds model.

The final binomial deviance is $5.1509$ on six degrees of freedom; its approximate upper-tail p-value is $0.525$. It gives no indication of substantial lack of fit against the saturated observed-cell model. The deviance residuals range from $-1.2830$ to $1.0820$, with median zero and no conspicuously large value in the five-number summary. Neither statement proves that every model assumption holds. In particular, the small case counts motivate checking the accuracy of the chi-squared asymptotic approximation, for example by a fitted-model parametric bootstrap.

The dispersion is fixed at one by the binomial distribution; it was not estimated as a free scale parameter. Four Fisher scoring iterations describe the numerical optimization, not four additional parameters. Residual and cellwise observed-versus-fitted comparisons would help assess adequacy and potential overdispersion.

Model the count relative to exposure, for exampleHere $E_i$ is a measured exposure such as patient-days, $x_i$ contains hospital predictors, and $b_i$ is a random effect. This Poisson regression separates hospital size from rate variation and supplies partial pooling. Use proper, scale-appropriate priors on coefficients and variation; add regional or temporal effects if the data demand them. A gamma mixing law with an estimated shape is another way to allow more flexible overdispersion. Compare the fitted model's counts and dispersion through posterior predictive checks.

The solid Poisson regression profile and dashed negative binomial regression profile peak at essentially the same $\beta \simeq 0.69$. The dashed profile is slightly wider and has heavier shoulders: it gives larger relative likelihood to values farther from the common maximizer. Thus the estimated mean ratio is stable, while the negative binomial model expresses slightly greater uncertainty.

The Poisson distribution constrains variance to equal the mean. A common negative binomial parametrization instead has $\mathrm{Var}(Y)=\mu +{\mu }^2/k$, with additional shape $k>0$; it permits overdispersion due, for example, to unexplained between-physician heterogeneity. Profiling this extra parameter lets the data support greater count variability and usually reduces information about the mean contrast. In this two-group model, the fitted group means remain the sample means for fixed $k$, which explains the nearly identical maximizers. The small separation of the profiles indicates a modest effect here, not overwhelming evidence for strong overdispersion. The comparison concerns the original PDF figure; no figure is reproduced.

For pair $i$, use conditionally independent counts $Y_{i0}\sim \mathrm{Pois}({\lambda }_i)$ and $Y_{i1}\sim \mathrm{Pois}({\lambda }_i{e}^{\beta })$, with a common effect $\beta$ and a separate positive baseline rate for each physician. Equivalently the linear predictor is ${\alpha }_i+\beta x$, with ${\alpha }_i=\log {\lambda }_i$. This Poisson regression controls each physician's stable baseline propensity.

Condition on all twelve pair totals $N_i=Y_{i0}+Y_{i1}$. The paired Poisson conditional likelihood has independent factorsThis eliminates every ${\lambda }_i$ exactly. It is more useful as a nuisance-elimination device than the unpaired calculation, which removes only one common baseline: many pair-specific rates are otherwise estimated from just two measurements each.

A possible concern is loss of information through conditioning on totals that are not ancillary for $\beta$, especially with sparse pairs; a zero-total pair contributes no conditional information. There is an important qualification in this exact fixed-baseline Poisson model. Maximizing over ${\lambda }_i$ gives ${\hat{\lambda }}_i=N_i/(1+e^{\beta })$, and substitution leaves precisely the same $\beta$-dependent factor as the conditional likelihood. Hence there is no extra loss relative to ordinary profiling in this particular model; one should not automatically assert an incidental-parameter bias in its estimate of $\beta$.

Additional baseline observations, a calibrated distribution of baseline rates, or informative hierarchical constraints can retain information in the pair totals. A Gamma–Poisson hierarchical model is one possible structured model, but merely giving it a completely free common mean does not automatically create extra information about $\beta$: that mean can absorb the factor $1+e^{\beta }$ in the totals. Such remedies trade exact nuisance elimination for additional data or assumptions. Finally, paired sampling alone does not establish conditional Poisson independence or equidispersion: remaining within-physician dependence or overdispersion should be checked, and an appropriate joint count model or physician-level robust variance used if needed. The common multiplicative-effect assumption is also part of the analysis.

Use grouped-binomial logistic regression, $S_i\sim \mathrm{Bin}(C_i,p_i)$ and $\mathrm{logit}(p_i)={\eta }_i$. Center the year as $u_i=i-7.5$. Maximizing $\ell ={\sum }_i[S_i{\eta }_i-C_i\ln (1+e^{{\eta }_i})]$ gives score $X^T(S-Cp)=0$ and information $X^T\mathrm{diag}\{C_i{p}_i(1-p_i)\}X$; these determine the estimates and working covariance matrix.

A straight-line logit fit is $\hat{\eta }=-1.95965+0.045384u$, with naive binomial standard errors $0.005382$ and $0.001352$. Its yearly odds multiplier is $e^{0.045384}=1.04643$, with naive 95% Wald confidence interval $[1.04366,1.04921]$. However, its residual deviance is $58.20$ on 12 degrees of freedom ($p=4.8\times 10^{-8}$), so those small standard errors do not describe an adequate independent-binomial fit.

Comparing polynomial mean structures gives residual deviances $1191.73,58.20,55.70,27.08$ for degrees zero through three, with respective residual degrees of freedom $13,12,11,10$. A quadratic term alone adds little; the cubic structure captures the middle-period acceleration and later leveling. Even the cubic independent-binomial fit has $p=0.00253$ for its residual deviance goodness-of-fit test. A useful working model is therefore cubic quasibinomial regression, retaining the logistic mean and estimating overdispersion rather than claiming an exact binomial law:The Pearson estimate is $\hat{\varphi }=27.182/10=2.718$. Using standard errors multiplied by $\sqrt{\hat{\varphi }}$ and a $t_{10}$ reference, the approximate 95% coefficient intervals are $[-1.98294,-1.92429]$, $[0.04982,0.07445]$, $[-0.001697,0.001131]$, and $[-0.000980,-0.000182]$. An approximate quasi-deviance comparison of cubic against linear gives $F_{2,10}=5.72$ and $p=0.022$, supporting the added shape while acknowledging that only fourteen annual groups inform dispersion.

In this preferred mean model the odds change is not constant:The fitted annual change is larger around the middle than at the ends; the straight-line 4.6% increase is a useful overall summary, not an exact yearly effect in the cubic model. The fitted curve, uncertainty band and residual plot appear in the figure. Singleton success increased overall, but a constant log-odds trend with unit binomial dispersion fits poorly. Quasibinomial intervals remain working-model intervals: serial dependence, changing patient mix and treatment practices cannot be separated from annual aggregates alone.

At several doses the three replicate plates differ much more than the square-root-of-mean scale expected from a Poisson distribution. For example the group at dose $100$ has sample mean about $42.67$ and sample variance about $274.33$, rather than variance near $42.67$. Its counts span a range of $33$, compared with a Poisson standard deviation about $6.53$. The dose-$33$ and dose-$1000$ groups likewise have variances above their means. Thus extra plate-to-plate variation is plausible. With only three replicates per dose this is evidence suggesting overdispersion, rather than a precise determination of it; not every group has excess sample variance.