Posterior Predictive Checks For Count Noise

Why Check a Fitted Count Model?

A Poisson or negative binomial model can produce a plausible fitted mean and still fail to reproduce important features of single-cell counts. It might underestimate the number of zeros, smooth away a small group of very large counts, or reproduce the mean while missing the spread. A posterior predictive check makes these potential failures concrete: it asks whether data generated from the fitted model resemble the observed data in ways relevant to the analysis.

The check is not a test that identifies the biological cause of a discrepancy. For example, a model producing too few zeros signals a mismatch in zero frequency, but that mismatch alone does not prove structural zero inflation or dropout. The interpretation must stay tied to the statistic being checked and to the model’s assumptions.

Definition: Let $y$ denote the observed count data, $\theta$ the model parameters, and $p(\theta\mid y)$ the posterior distribution of those parameters. A posterior predictive replicate $Y^{\mathrm{rep}}$ is generated by first drawing $\theta$ from $p(\theta\mid y)$, then drawing new data from the sampling model $p(Y^{\mathrm{rep}}\mid\theta)$. Its distribution is $$ p(Y^{\mathrm{rep}}\mid y)=\int p(Y^{\mathrm{rep}}\mid\theta)p(\theta\mid y)\,d\theta. $$

For count checks, $y_i$ can denote the count for a specified gene in cell $i$, with $i=1,\ldots,n$. The replicate $Y_i^{\mathrm{rep}}$ represents a new count under the same model and the same design information, such as the same cells’ exposures or covariates. The observed counts are compared with replicate counts, not replaced by them. If cells have different fitted means, each replicate should retain those cell-specific means rather than treating the pooled counts as identically distributed.

A common conditional model is $Y_i\mid\theta\sim\operatorname{Poisson}(\lambda_i)$, independently across cells conditional on the model parameters. Its conditional mean and variance are both $\lambda_i$. A negative binomial model may instead be parameterized by mean $\mu_i$ and size $r$, with conditional variance $\mu_i+\mu_i^2/r$. Under this parameterization, its zero probability is $(r/(r+\mu_i))^r$. These are model assumptions, not universal properties of count data: independence, the mean structure, and any shared dispersion assumptions need to be reasonable for the application.

Choose Discrepancies That Target Count Features

A discrepancy statistic $T(y)$ reduces a data feature to a number that can be computed both for the observed data and for each replicate. The direction of concern should be stated before comparing: for example, “the model may underpredict the number of zeros” suggests looking at the upper tail of a zero-count statistic. Useful checks for count noise include:

Mean: the average count, or a mean calculated within groups of cells having similar fitted means.
Variance: the sample variance or a variance-to-mean comparison, again preferably accounting for variation in fitted means.
Upper tail: the number or proportion of counts at least as large as a specified threshold, or the maximum count.
Zero frequency: the number or proportion of cells with count zero.

For each statistic, generate $B$ posterior predictive replicates and compute $T(Y^{\mathrm{rep},b})$ for $b=1,\ldots,B$. A simple graphical check places the observed value alongside the distribution of replicated values. A one-sided posterior predictive tail area can be estimated by the fraction of replicates at least as extreme as the observed value in the chosen direction. This fraction summarizes the comparison; it is not, in general, a uniformly calibrated frequentist $p$-value. The data were used to form the posterior and are then compared with data generated from that posterior.

Specify the conditional model.
State the count distribution, its mean structure, dispersion assumptions, and any conditional independence assumptions.

Choose a targeted statistic.
Define a count feature and the direction of mismatch that would be concerning, before interpreting its comparison.

Generate posterior predictive replicates.
For each draw from the posterior, simulate counts under the fitted model while retaining the observed design information.

Compare and diagnose.
Locate the observed statistic in the replicated distribution, examine relevant plots, and describe the discrepancy without attributing a cause the check cannot identify.

What Posterior Predictive Variation Contains

Replicates vary for two reasons: counts vary even at fixed parameters, and the parameters themselves are uncertain after observing the data. The posterior predictive distribution includes both sources. For a scalar statistic $T(Y^{\mathrm{rep}})$, the law of total expectation gives

E[T(Y^{\mathrm{rep}})\mid y] = E_{\theta\mid y}\!\left[E\{T(Y^{\mathrm{rep}})\mid\theta\}\right].

The corresponding variance decomposition is

\operatorname{Var}(T(Y^{\mathrm{rep}})\mid y) = E_{\theta\mid y}\!\left[\operatorname{Var}\{T(Y^{\mathrm{rep}})\mid\theta\}\right] + \operatorname{Var}_{\theta\mid y}\!\left[E\{T(Y^{\mathrm{rep}})\mid\theta\}\right].

The first term is the average sampling variation at fixed parameters. The second is variation in the statistic’s conditional expectation across plausible parameter values. Omitting posterior uncertainty by simulating only at a single point estimate can therefore make the reference distribution too narrow, especially when the data provide limited information about dispersion or a tail feature.

Worked Example: Checking the Poisson Zero Frequency

Suppose, hypothetically, that 40 cells are modeled as independent Poisson counts, each with fixed mean $1.5$. The observed data contain 12 zero counts. Under this model, the probability of a zero in one cell is $$ p_0=P(Y_i=0)=e^{-1.5}\approx 0.223130. $$

Let $Z$ be the number of zeros in a replicate. Under the stated independence and common-mean assumptions, $Z\sim\operatorname{Binomial}(40,p_0)$. Thus $$ E(Z)=40p_0=40(0.223130)\approx 8.925, $$ and $$ \operatorname{SD}(Z)=\sqrt{40p_0(1-p_0)} =\sqrt{40(0.223130)(0.776870)} \approx 2.633. $$

The observed value, 12, is about $(12-8.925)/2.633\approx 1.168$ standard deviations above the model’s expected zero count. This calculation is a useful scale comparison, not a normal-theory test or a proof of model failure. A direct posterior predictive check would simulate the parameters if they were uncertain, generate replicated counts, and compare the observed zero count with the resulting distribution. If cell means differ, the zero indicators have different probabilities $e^{-\lambda_i}$; the common-probability binomial calculation above would then not apply.

Worked Example: Checking a Rare Upper-Tail Count

Consider a hypothetical fitted count model and a prespecified threshold for a cell-level upper-tail event. Suppose the model-based posterior predictive probability of that event for each of 20 conditionally independent cells is $q=5/256=0.01953125$, treated here as fixed for the arithmetic. The observed data have two or more such events. Let $T$ count these events in a replicate. Then $T\sim\operatorname{Binomial}(20,q)$, and the probability of at least two is $$ P(T\geq 2)=1-P(T=0)-P(T=1) =1-(1-q)^{20}-20q(1-q)^{19}. $$

Substitution gives $$ P(T\geq 2)=1-(251/256)^{20}-20(5/256)(251/256)^{19}\approx 0.0574. $$

So, under these fixed-probability assumptions, about 5.74% of replicate datasets would contain at least two threshold exceedances. The value describes the reference model’s predictive tail for this statistic; it is not the probability that the model is true. In an actual posterior predictive analysis, $q$ can vary across posterior draws and cells, so replicating the counts under those draws is preferable to substituting one fixed $q$.

Worked Example: Comparing Mean, Variance, and Upper Tail Under a Negative Binomial Model

Suppose, hypothetically, five counts for one gene are $0,1,2,4,8$. For illustration, consider a negative binomial model with common mean $\mu=3$ and size $r=2$, treating these parameters as fixed. The observed sample mean is $$ \bar{x}=(0+1+2+4+8)/5=15/5=3. $$ The sample variance is $$ s^2=\frac{(0-3)^2+(1-3)^2+(2-3)^2+(4-3)^2+(8-3)^2}{5-1} =\frac{40}{4}=10. $$

The model variance for one count is $\mu+\mu^2/r=3+9/2=7.5$. Thus the observed sample variance exceeds the model’s single-count variance, while the observed mean equals the model mean. That is a descriptive warning to check the replicate distribution, not a conclusion based on comparing two numbers alone.

The model’s zero probability is $$ P(Y_i=0)=\left(\frac{r}{r+\mu}\right)^r =\left(\frac{2}{5}\right)^2=0.16. $$ The expected zero count among five cells is $5(0.16)=0.8$. For the upper tail, the negative binomial probability mass function under this parameterization is $P(Y=k)=(k+1)(0.4)^2(0.6)^k$. Therefore, $$ P(Y\geq 8) =(0.4)^2\sum_{k=8}^{\infty}(k+1)(0.6)^k =(0.4)^2(0.6)^8\left(\frac{9}{0.4}+\frac{0.6}{0.4^2}\right) \approx 0.07054. $$

The expected number of counts at least 8 in five cells is $5(0.07054)\approx 0.3527$. The observed count of 8 is therefore worth examining in replicated datasets. It does not by itself establish that the model underestimates the upper tail: a single maximum in five observations is highly variable, and the model parameters would ordinarily also be uncertain. Posterior predictive replicates provide the joint reference distributions for the mean, variance, zero count, and threshold exceedances.

Reading Discrepancies and Diagnosing Failure

A useful check compares like with like. For a heterogeneous mean structure, an overall variance can be large simply because fitted means differ across cells. A pooled variance-to-mean ratio can then suggest excess dispersion even when each cell’s conditional distribution is appropriate. In that setting, compare counts to cell-specific predictive distributions, or examine summaries within groups of similar fitted means. Likewise, a single pooled zero percentage may obscure the fact that low-mean cells should have more zeros than high-mean cells.

Patterns across several checks can help localize a mismatch. An observed mean near the replicated range but a variance above nearly all replicated variances suggests that the model captures the average level better than the spread. A zero count above the replicated range indicates underprediction of zeros for the checked data and design. A high upper-tail count can point to insufficient tail variation or a misspecified mean pattern. These observations are model diagnostics, not unique explanations. Different mechanisms can create similar discrepancies, and a check on one summary can miss failures in another part of the distribution.

Diagnostics should also consider the simulation itself. Too few posterior draws can make estimated tails unstable; repeated runs or more draws can reveal Monte Carlo noise. Inspect whether posterior draws adequately represent uncertain parameters, especially dispersion. Check that replicated counts are nonnegative integers, use the intended exposures and covariates, and follow the stated dependence assumptions. For posterior draws with different fitted means, calculate the discrepancy on every replicate using the same definition used for the observed data.

Alternatives, Limitations, and a Common Mistake

Posterior predictive checks can be presented as overlays of observed and replicated count distributions, or as replicated distributions of targeted statistics. They differ from residual checks, which summarize deviations from fitted values on a chosen scale and can be sensitive to that scale. They also differ from a parametric bootstrap check, which commonly estimates parameters from the observed data, simulates datasets at those estimates, and may refit the model for each simulated dataset. Posterior predictive checks average over a posterior distribution of parameters; bootstrap procedures have a different construction and inferential interpretation. Neither approach removes the need to choose diagnostics that address the scientific question.

A limitation of posterior predictive checking is that the same data inform the posterior and the comparison. As a result, a model can appear adequate on checks closely related to the statistics it was fitted to, particularly when parameters are flexible. A check can also have low sensitivity to a mismatch it was not designed to detect. Use several interpretable discrepancies, and treat their results as evidence about predictive adequacy, not a final model-selection rule.

Mistake: “The observed zero count is larger than the model’s expected zero count, so the data prove dropout or zero inflation.”
Correction: A difference from an expectation is not enough; compare the observed statistic with its full predictive distribution, including sampling and parameter uncertainty. Even a clear excess-zero discrepancy says that the fitted model underpredicts zeros for that check. It does not identify the mechanism responsible.
Mistake: Checking only the mean because the model’s purpose is to estimate average expression.
Correction: A matching mean does not imply a matching variance, tail, or zero frequency. Check the count features relevant to the intended use.

Key takeaway: A posterior predictive check compares observed count summaries with summaries generated from the fitted model’s posterior predictive distribution. Its conclusion should identify the feature that is or is not reproduced, while avoiding claims about causes that the check cannot distinguish.

Check Your Understanding

Answer each question, then compare with the solutions below.

In one sentence, distinguish a posterior predictive replicate from a new draw at a single fixed parameter estimate.
A Poisson model assigns each of 30 independent cells mean 2. What are the expected zero count and its standard deviation?
A negative binomial model has $\mu=4$ and $r=2$. Calculate its variance and zero probability.
If a fitted model reproduces the observed mean but nearly always generates fewer zeros than observed, what does the check support, and what does it not establish?
Why can a pooled variance be misleading when the cells have different fitted means?

Solutions. (1) A posterior predictive replicate first draws parameters from their posterior given the observed data, then draws new counts conditional on those parameters; a plug-in replicate fixes parameters at one estimate and omits their posterior uncertainty. (2) The zero probability is $e^{-2}\approx0.1353$. Thus the expected zero count is $30e^{-2}\approx4.060$, and its standard deviation is $\sqrt{30e^{-2}(1-e^{-2})}\approx1.874$. (3) The variance is $4+4^2/2=12$, and the zero probability is $(2/(2+4))^2=1/9\approx0.1111$. (4) The check supports the conclusion that the fitted model underpredicts zero frequency for the chosen data and design. It does not prove dropout, structural zeros, or any other particular mechanism. (5) Different fitted means create between-cell variation in addition to conditional count variation, so a pooled variance can be large even if the conditional model is adequate within cells of comparable means.