From 13a7dcf86bfb03a6871f4767a52fb33472247807 Mon Sep 17 00:00:00 2001
From: "seb.dalgarno@gmail.com" <seb@poissonconsulting.ca>
Date: Fri, 8 May 2026 11:26:15 -0700
Subject: [PATCH 1/3] Add zero-truncated gamma-Poisson functions

Adds log_lik_gamma_pois_zt(), dev_gamma_pois_zt(), ran_gamma_pois_zt(),
and res_gamma_pois_zt() for the zero-truncated gamma-Poisson (negative
binomial) distribution, mirroring the existing _zi family.

The likelihood is the standard gamma-Poisson density divided by
1 - P(0; lambda, theta). Random samples use rejection sampling (zero
draws are redrawn until accepted). Residuals centre on the truncated
mean lambda / (1 - P(0)) and the standardised form scales by the
truncated variance derived analytically. Deviances reuse the
untruncated NB deviance approximation since the saturated MLE for
the truncated likelihood has no closed form.

Tests follow the patterns established for _zi and cover known values,
NA handling, vectorisation, and the simulate workflow.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
---
 NAMESPACE                     |  4 +++
 NEWS.md                       |  1 +
 R/dev.R                       | 29 +++++++++++++++++++++
 R/log-lik.R                   | 20 ++++++++++++++
 R/ran.R                       | 36 +++++++++++++++++++++++++
 R/res.R                       | 33 +++++++++++++++++++++++
 man/dev_gamma_pois_zi.Rd      |  4 +++
 man/dev_gamma_pois_zt.Rd      | 39 ++++++++++++++++++++++++++++
 man/log_lik_bern.Rd           |  1 +
 man/log_lik_beta.Rd           |  1 +
 man/log_lik_beta_binom.Rd     |  1 +
 man/log_lik_binom.Rd          |  1 +
 man/log_lik_exp.Rd            |  1 +
 man/log_lik_gamma.Rd          |  1 +
 man/log_lik_gamma_pois.Rd     |  1 +
 man/log_lik_gamma_pois_zi.Rd  |  1 +
 man/log_lik_gamma_pois_zt.Rd  | 46 ++++++++++++++++++++++++++++++++
 man/log_lik_lnorm.Rd          |  1 +
 man/log_lik_neg_binom.Rd      |  1 +
 man/log_lik_norm.Rd           |  1 +
 man/log_lik_pois.Rd           |  1 +
 man/log_lik_pois_zi.Rd        |  1 +
 man/log_lik_skewnorm.Rd       |  1 +
 man/log_lik_student.Rd        |  1 +
 man/log_lik_unif.Rd           |  1 +
 man/ran_bern.Rd               |  1 +
 man/ran_beta_binom.Rd         |  1 +
 man/ran_binom.Rd              |  1 +
 man/ran_gamma.Rd              |  1 +
 man/ran_gamma_pois.Rd         |  1 +
 man/ran_gamma_pois_zi.Rd      |  1 +
 man/ran_gamma_pois_zt.Rd      | 45 ++++++++++++++++++++++++++++++++
 man/ran_lnorm.Rd              |  1 +
 man/ran_neg_binom.Rd          |  1 +
 man/ran_norm.Rd               |  1 +
 man/ran_pois.Rd               |  1 +
 man/ran_pois_zi.Rd            |  1 +
 man/ran_skewnorm.Rd           |  1 +
 man/ran_student.Rd            |  1 +
 man/res_bern.Rd               |  1 +
 man/res_beta_binom.Rd         |  1 +
 man/res_binom.Rd              |  1 +
 man/res_gamma.Rd              |  1 +
 man/res_gamma_pois.Rd         |  1 +
 man/res_gamma_pois_zi.Rd      |  1 +
 man/res_gamma_pois_zt.Rd      | 49 +++++++++++++++++++++++++++++++++++
 man/res_lnorm.Rd              |  1 +
 man/res_neg_binom.Rd          |  1 +
 man/res_norm.Rd               |  1 +
 man/res_pois.Rd               |  1 +
 man/res_pois_zi.Rd            |  1 +
 man/res_skewnorm.Rd           |  1 +
 man/res_student.Rd            |  1 +
 tests/testthat/test-dev.R     | 28 ++++++++++++++++++++
 tests/testthat/test-log-lik.R | 29 +++++++++++++++++++++
 tests/testthat/test-ran.R     | 33 +++++++++++++++++++++++
 tests/testthat/test-res.R     | 44 +++++++++++++++++++++++++++++++
 57 files changed, 482 insertions(+)
 create mode 100644 man/dev_gamma_pois_zt.Rd
 create mode 100644 man/log_lik_gamma_pois_zt.Rd
 create mode 100644 man/ran_gamma_pois_zt.Rd
 create mode 100644 man/res_gamma_pois_zt.Rd

diff --git a/NAMESPACE b/NAMESPACE
index 680601a1..bf0d7dc7 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -38,6 +38,7 @@ export(dev_binom)
 export(dev_gamma)
 export(dev_gamma_pois)
 export(dev_gamma_pois_zi)
+export(dev_gamma_pois_zt)
 export(dev_lnorm)
 export(dev_neg_binom)
 export(dev_norm)
@@ -67,6 +68,7 @@ export(log_lik_exp)
 export(log_lik_gamma)
 export(log_lik_gamma_pois)
 export(log_lik_gamma_pois_zi)
+export(log_lik_gamma_pois_zt)
 export(log_lik_lnorm)
 export(log_lik_neg_binom)
 export(log_lik_norm)
@@ -106,6 +108,7 @@ export(ran_binom)
 export(ran_gamma)
 export(ran_gamma_pois)
 export(ran_gamma_pois_zi)
+export(ran_gamma_pois_zt)
 export(ran_lnorm)
 export(ran_neg_binom)
 export(ran_norm)
@@ -120,6 +123,7 @@ export(res_binom)
 export(res_gamma)
 export(res_gamma_pois)
 export(res_gamma_pois_zi)
+export(res_gamma_pois_zt)
 export(res_lnorm)
 export(res_neg_binom)
 export(res_norm)
diff --git a/NEWS.md b/NEWS.md
index f34816fb..52a4a3aa 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -2,6 +2,7 @@
 
 # extras 0.8.0.9001
 
+- Added `log_lik_gamma_pois_zt()`, `dev_gamma_pois_zt()`, `ran_gamma_pois_zt()`, and `res_gamma_pois_zt()` for the zero-truncated gamma-Poisson (negative binomial) distribution.
 - Added `log_lik_exp()`.
 - Added `log_lik_beta()`.
 - Added 'log_lik_unif()'.
diff --git a/R/dev.R b/R/dev.R
index 4c2d17b1..44bc07a0 100644
--- a/R/dev.R
+++ b/R/dev.R
@@ -244,6 +244,35 @@ dev_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0, res = FALSE) {
   dev_res(x, lambda * (1 - prob), dev)
 }
 
+#' Zero-Truncated Gamma-Poisson Deviances
+#'
+#' Deviances for the zero-truncated gamma-Poisson (negative binomial)
+#' distribution. The deviance is approximated by the untruncated
+#' gamma-Poisson deviance; deviance residuals are centred on the
+#' truncated mean `lambda / (1 - P(0; lambda, theta))`. The saturated
+#' MLE for the zero-truncated likelihood differs from the observation
+#' and has no closed form; the approximation is accurate for typical
+#' `x` and small zero-mass `P(0; lambda, theta)`.
+#'
+#' @inheritParams params
+#' @param x A whole numeric vector of values greater than or equal to 1.
+#'
+#' @return An numeric vector of the corresponding deviances or deviance residuals.
+#' @family dev_dist # make live when complete
+#' @export
+#'
+#' @examples
+#' dev_gamma_pois_zt(c(1, 3, 4), 3, 2)
+dev_gamma_pois_zt <- function(x, lambda = 1, theta = 0, res = FALSE) {
+  dev <- dev_gamma_pois(x, lambda = lambda, theta = theta, res = FALSE)
+  if (vld_false(res)) {
+    return(dev)
+  }
+  log_p0 <- dnbinom(0, mu = lambda, size = 1 / theta, log = TRUE)
+  trunc_mean <- lambda / -expm1(log_p0)
+  dev_res(x, trunc_mean, dev)
+}
+
 #' Log-Normal Deviances
 #'
 #' @inheritParams params
diff --git a/R/log-lik.R b/R/log-lik.R
index 9b2cfe72..a3ffbdc6 100644
--- a/R/log-lik.R
+++ b/R/log-lik.R
@@ -205,6 +205,26 @@ log_lik_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0) {
   log(lpois)
 }
 
+#' Zero-Truncated Gamma-Poisson Log-Likelihood
+#'
+#' Density of the gamma-Poisson (negative binomial) distribution conditional
+#' on `x >= 1`.
+#'
+#' @inheritParams params
+#' @param x A whole numeric vector of values greater than or equal to 1.
+#'
+#' @return An numeric vector of the corresponding log-likelihoods.
+#' @family log_lik_dist
+#' @export
+#'
+#' @examples
+#' log_lik_gamma_pois_zt(c(1, 3, 4), 3, 1)
+log_lik_gamma_pois_zt <- function(x, lambda = 1, theta = 0) {
+  log_lik <- dnbinom(x, mu = lambda, size = 1 / theta, log = TRUE)
+  log_p0 <- dnbinom(0, mu = lambda, size = 1 / theta, log = TRUE)
+  log_lik - log1p(-exp(log_p0))
+}
+
 #' Log-Normal Log-Likelihood
 #'
 #' @inheritParams params
diff --git a/R/ran.R b/R/ran.R
index 8c659585..aac41de2 100644
--- a/R/ran.R
+++ b/R/ran.R
@@ -97,6 +97,42 @@ ran_gamma_pois_zi <- function(n = 1, lambda = 1, theta = 0, prob = 0) {
   ran_neg_binom(n = n, lambda = lambda, theta = theta) * ran_bern(n, prob = 1 - prob)
 }
 
+#' Zero-Truncated Gamma-Poisson Random Samples
+#'
+#' Draws from the gamma-Poisson (negative binomial) distribution conditional
+#' on the value being greater than or equal to 1. Implemented via rejection
+#' sampling: zero draws are redrawn until accepted.
+#'
+#' @inheritParams params
+#' @return A numeric vector of the random samples.
+#' @family ran_dist
+#' @export
+#'
+#' @examples
+#' ran_gamma_pois_zt(10, lambda = 3, theta = 1)
+ran_gamma_pois_zt <- function(n = 1, lambda = 1, theta = 0) {
+  chk_whole_number(n)
+  chk_gte(n)
+  if (n == 0) {
+    return(integer(0))
+  }
+  lambda <- rep_len(lambda, n)
+  theta <- rep_len(theta, n)
+  out <- integer(n)
+  remaining <- seq_len(n)
+  while (length(remaining) > 0) {
+    draws <- as.integer(stats::rnbinom(
+      length(remaining),
+      mu = lambda[remaining],
+      size = 1 / theta[remaining]
+    ))
+    accepted <- draws > 0
+    out[remaining[accepted]] <- draws[accepted]
+    remaining <- remaining[!accepted]
+  }
+  out
+}
+
 #' Log-Normal Random Samples
 #'
 #' @inheritParams params
diff --git a/R/res.R b/R/res.R
index a1f1496a..aa9361a7 100644
--- a/R/res.R
+++ b/R/res.R
@@ -173,6 +173,39 @@ res_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0, type = "dev",
   )
 }
 
+#' Zero-Truncated Gamma-Poisson Residuals
+#'
+#' Residuals for the zero-truncated gamma-Poisson (negative binomial)
+#' distribution. Centring and standardisation use the truncated mean
+#' `lambda / (1 - P(0; lambda, theta))` and the truncated variance.
+#'
+#' @inheritParams params
+#' @param x A whole numeric vector of values greater than or equal to 1.
+#'
+#' @return An numeric vector of the corresponding residuals.
+#' @family res_dist
+#' @export
+#'
+#' @examples
+#' res_gamma_pois_zt(c(1, 2, 5), 2, 1)
+res_gamma_pois_zt <- function(x, lambda = 1, theta = 0, type = "dev", simulate = FALSE) {
+  chk_string(type)
+  if (!vld_false(simulate)) {
+    x <- ran_gamma_pois_zt(length(x), lambda = lambda, theta = theta)
+  }
+  log_p0 <- dnbinom(0, mu = lambda, size = 1 / theta, log = TRUE)
+  one_minus_p0 <- -expm1(log_p0)
+  trunc_mean <- lambda / one_minus_p0
+  trunc_var <- (lambda + lambda^2 * (1 + theta)) / one_minus_p0 - trunc_mean^2
+  switch(type,
+    data = x,
+    raw = x - trunc_mean,
+    standardized = (x - trunc_mean) / sqrt(trunc_var),
+    dev = dev_gamma_pois_zt(x, lambda = lambda, theta = theta, res = TRUE),
+    chk_subset(x, c("data", "raw", "dev", "standardized"))
+  )
+}
+
 #' Log-Normal Residuals
 #'
 #' @inheritParams params
diff --git a/man/dev_gamma_pois_zi.Rd b/man/dev_gamma_pois_zi.Rd
index 790d700f..c67a9d00 100644
--- a/man/dev_gamma_pois_zi.Rd
+++ b/man/dev_gamma_pois_zi.Rd
@@ -29,4 +29,8 @@ Zero-Inflated Gamma-Poisson Deviances
 \examples{
 dev_gamma_pois_zi(c(1, 3, 4), 3, 2)
 }
+\seealso{
+Other dev_dist # make live when complete: 
+\code{\link{dev_gamma_pois_zt}()}
+}
 \concept{dev_dist # make live when complete}
diff --git a/man/dev_gamma_pois_zt.Rd b/man/dev_gamma_pois_zt.Rd
new file mode 100644
index 00000000..0e279209
--- /dev/null
+++ b/man/dev_gamma_pois_zt.Rd
@@ -0,0 +1,39 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/dev.R
+\name{dev_gamma_pois_zt}
+\alias{dev_gamma_pois_zt}
+\title{Zero-Truncated Gamma-Poisson Deviances}
+\usage{
+dev_gamma_pois_zt(x, lambda = 1, theta = 0, res = FALSE)
+}
+\arguments{
+\item{x}{A whole numeric vector of values greater than or equal to 1.}
+
+\item{lambda}{A non-negative numeric vector of means.}
+
+\item{theta}{A non-negative numeric vector of the dispersion for the mixture
+models (student, gamma-Poisson and beta-binomial).}
+
+\item{res}{A flag specifying whether to return the deviance residual as
+opposed to the deviance.}
+}
+\value{
+An numeric vector of the corresponding deviances or deviance residuals.
+}
+\description{
+Deviances for the zero-truncated gamma-Poisson (negative binomial)
+distribution. The deviance is approximated by the untruncated
+gamma-Poisson deviance; deviance residuals are centred on the
+truncated mean \verb{lambda / (1 - P(0; lambda, theta))}. The saturated
+MLE for the zero-truncated likelihood differs from the observation
+and has no closed form; the approximation is accurate for typical
+\code{x} and small zero-mass \verb{P(0; lambda, theta)}.
+}
+\examples{
+dev_gamma_pois_zt(c(1, 3, 4), 3, 2)
+}
+\seealso{
+Other dev_dist # make live when complete: 
+\code{\link{dev_gamma_pois_zi}()}
+}
+\concept{dev_dist # make live when complete}
diff --git a/man/log_lik_bern.Rd b/man/log_lik_bern.Rd
index cefc8642..440e1440 100644
--- a/man/log_lik_bern.Rd
+++ b/man/log_lik_bern.Rd
@@ -30,6 +30,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_beta.Rd b/man/log_lik_beta.Rd
index 565d05dc..7be6aa73 100644
--- a/man/log_lik_beta.Rd
+++ b/man/log_lik_beta.Rd
@@ -31,6 +31,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_beta_binom.Rd b/man/log_lik_beta_binom.Rd
index 4fc365e5..3e2e91ae 100644
--- a/man/log_lik_beta_binom.Rd
+++ b/man/log_lik_beta_binom.Rd
@@ -45,6 +45,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_binom.Rd b/man/log_lik_binom.Rd
index 53369af5..863dd5b6 100644
--- a/man/log_lik_binom.Rd
+++ b/man/log_lik_binom.Rd
@@ -32,6 +32,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_exp.Rd b/man/log_lik_exp.Rd
index 2dbbd29f..95187fa0 100644
--- a/man/log_lik_exp.Rd
+++ b/man/log_lik_exp.Rd
@@ -29,6 +29,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_gamma.Rd b/man/log_lik_gamma.Rd
index ae073723..7e1c3b47 100644
--- a/man/log_lik_gamma.Rd
+++ b/man/log_lik_gamma.Rd
@@ -31,6 +31,7 @@ Other log_lik_dist:
 \code{\link{log_lik_exp}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_gamma_pois.Rd b/man/log_lik_gamma_pois.Rd
index 9041c389..17059b56 100644
--- a/man/log_lik_gamma_pois.Rd
+++ b/man/log_lik_gamma_pois.Rd
@@ -32,6 +32,7 @@ Other log_lik_dist:
 \code{\link{log_lik_exp}()},
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_gamma_pois_zi.Rd b/man/log_lik_gamma_pois_zi.Rd
index fed48ce0..145f3772 100644
--- a/man/log_lik_gamma_pois_zi.Rd
+++ b/man/log_lik_gamma_pois_zi.Rd
@@ -35,6 +35,7 @@ Other log_lik_dist:
 \code{\link{log_lik_exp}()},
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_gamma_pois_zt.Rd b/man/log_lik_gamma_pois_zt.Rd
new file mode 100644
index 00000000..5e46fd9a
--- /dev/null
+++ b/man/log_lik_gamma_pois_zt.Rd
@@ -0,0 +1,46 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/log-lik.R
+\name{log_lik_gamma_pois_zt}
+\alias{log_lik_gamma_pois_zt}
+\title{Zero-Truncated Gamma-Poisson Log-Likelihood}
+\usage{
+log_lik_gamma_pois_zt(x, lambda = 1, theta = 0)
+}
+\arguments{
+\item{x}{A whole numeric vector of values greater than or equal to 1.}
+
+\item{lambda}{A non-negative numeric vector of means.}
+
+\item{theta}{A non-negative numeric vector of the dispersion for the mixture
+models (student, gamma-Poisson and beta-binomial).}
+}
+\value{
+An numeric vector of the corresponding log-likelihoods.
+}
+\description{
+Density of the gamma-Poisson (negative binomial) distribution conditional
+on \code{x >= 1}.
+}
+\examples{
+log_lik_gamma_pois_zt(c(1, 3, 4), 3, 1)
+}
+\seealso{
+Other log_lik_dist: 
+\code{\link{log_lik_bern}()},
+\code{\link{log_lik_beta}()},
+\code{\link{log_lik_beta_binom}()},
+\code{\link{log_lik_binom}()},
+\code{\link{log_lik_exp}()},
+\code{\link{log_lik_gamma}()},
+\code{\link{log_lik_gamma_pois}()},
+\code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_lnorm}()},
+\code{\link{log_lik_neg_binom}()},
+\code{\link{log_lik_norm}()},
+\code{\link{log_lik_pois}()},
+\code{\link{log_lik_pois_zi}()},
+\code{\link{log_lik_skewnorm}()},
+\code{\link{log_lik_student}()},
+\code{\link{log_lik_unif}()}
+}
+\concept{log_lik_dist}
diff --git a/man/log_lik_lnorm.Rd b/man/log_lik_lnorm.Rd
index 7d7d3e02..515f5724 100644
--- a/man/log_lik_lnorm.Rd
+++ b/man/log_lik_lnorm.Rd
@@ -33,6 +33,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
 \code{\link{log_lik_pois}()},
diff --git a/man/log_lik_neg_binom.Rd b/man/log_lik_neg_binom.Rd
index 8d47937e..4b7f7adc 100644
--- a/man/log_lik_neg_binom.Rd
+++ b/man/log_lik_neg_binom.Rd
@@ -33,6 +33,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_norm}()},
 \code{\link{log_lik_pois}()},
diff --git a/man/log_lik_norm.Rd b/man/log_lik_norm.Rd
index 8d1b3c34..c48f3ca2 100644
--- a/man/log_lik_norm.Rd
+++ b/man/log_lik_norm.Rd
@@ -32,6 +32,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_pois}()},
diff --git a/man/log_lik_pois.Rd b/man/log_lik_pois.Rd
index 26e0de54..d56f178e 100644
--- a/man/log_lik_pois.Rd
+++ b/man/log_lik_pois.Rd
@@ -30,6 +30,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_pois_zi.Rd b/man/log_lik_pois_zi.Rd
index 8f378a9f..10e45c06 100644
--- a/man/log_lik_pois_zi.Rd
+++ b/man/log_lik_pois_zi.Rd
@@ -33,6 +33,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_skewnorm.Rd b/man/log_lik_skewnorm.Rd
index 93e60f2a..3a1e25e1 100644
--- a/man/log_lik_skewnorm.Rd
+++ b/man/log_lik_skewnorm.Rd
@@ -38,6 +38,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_student.Rd b/man/log_lik_student.Rd
index d392447c..1bb261ca 100644
--- a/man/log_lik_student.Rd
+++ b/man/log_lik_student.Rd
@@ -35,6 +35,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/log_lik_unif.Rd b/man/log_lik_unif.Rd
index a2996430..27da9ef6 100644
--- a/man/log_lik_unif.Rd
+++ b/man/log_lik_unif.Rd
@@ -32,6 +32,7 @@ Other log_lik_dist:
 \code{\link{log_lik_gamma}()},
 \code{\link{log_lik_gamma_pois}()},
 \code{\link{log_lik_gamma_pois_zi}()},
+\code{\link{log_lik_gamma_pois_zt}()},
 \code{\link{log_lik_lnorm}()},
 \code{\link{log_lik_neg_binom}()},
 \code{\link{log_lik_norm}()},
diff --git a/man/ran_bern.Rd b/man/ran_bern.Rd
index 2ccc2555..3be44ce3 100644
--- a/man/ran_bern.Rd
+++ b/man/ran_bern.Rd
@@ -29,6 +29,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_beta_binom.Rd b/man/ran_beta_binom.Rd
index bfca892e..e80beb49 100644
--- a/man/ran_beta_binom.Rd
+++ b/man/ran_beta_binom.Rd
@@ -42,6 +42,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_binom.Rd b/man/ran_binom.Rd
index 7084843e..78823fc0 100644
--- a/man/ran_binom.Rd
+++ b/man/ran_binom.Rd
@@ -31,6 +31,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_gamma.Rd b/man/ran_gamma.Rd
index 923a88fe..edba7369 100644
--- a/man/ran_gamma.Rd
+++ b/man/ran_gamma.Rd
@@ -30,6 +30,7 @@ Other ran_dist:
 \code{\link{ran_binom}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_gamma_pois.Rd b/man/ran_gamma_pois.Rd
index 50e0c60b..57d39ef2 100644
--- a/man/ran_gamma_pois.Rd
+++ b/man/ran_gamma_pois.Rd
@@ -31,6 +31,7 @@ Other ran_dist:
 \code{\link{ran_binom}()},
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_gamma_pois_zi.Rd b/man/ran_gamma_pois_zi.Rd
index 5a70d40c..bb574037 100644
--- a/man/ran_gamma_pois_zi.Rd
+++ b/man/ran_gamma_pois_zi.Rd
@@ -34,6 +34,7 @@ Other ran_dist:
 \code{\link{ran_binom}()},
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_gamma_pois_zt.Rd b/man/ran_gamma_pois_zt.Rd
new file mode 100644
index 00000000..7d027d9f
--- /dev/null
+++ b/man/ran_gamma_pois_zt.Rd
@@ -0,0 +1,45 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/ran.R
+\name{ran_gamma_pois_zt}
+\alias{ran_gamma_pois_zt}
+\title{Zero-Truncated Gamma-Poisson Random Samples}
+\usage{
+ran_gamma_pois_zt(n = 1, lambda = 1, theta = 0)
+}
+\arguments{
+\item{n}{A non-negative whole number of the number of random samples to
+generate.}
+
+\item{lambda}{A non-negative numeric vector of means.}
+
+\item{theta}{A non-negative numeric vector of the dispersion for the mixture
+models (student, gamma-Poisson and beta-binomial).}
+}
+\value{
+A numeric vector of the random samples.
+}
+\description{
+Draws from the gamma-Poisson (negative binomial) distribution conditional
+on the value being greater than or equal to 1. Implemented via rejection
+sampling: zero draws are redrawn until accepted.
+}
+\examples{
+ran_gamma_pois_zt(10, lambda = 3, theta = 1)
+}
+\seealso{
+Other ran_dist: 
+\code{\link{ran_bern}()},
+\code{\link{ran_beta_binom}()},
+\code{\link{ran_binom}()},
+\code{\link{ran_gamma}()},
+\code{\link{ran_gamma_pois}()},
+\code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_lnorm}()},
+\code{\link{ran_neg_binom}()},
+\code{\link{ran_norm}()},
+\code{\link{ran_pois}()},
+\code{\link{ran_pois_zi}()},
+\code{\link{ran_skewnorm}()},
+\code{\link{ran_student}()}
+}
+\concept{ran_dist}
diff --git a/man/ran_lnorm.Rd b/man/ran_lnorm.Rd
index 03ac9612..d1cef892 100644
--- a/man/ran_lnorm.Rd
+++ b/man/ran_lnorm.Rd
@@ -32,6 +32,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
 \code{\link{ran_pois}()},
diff --git a/man/ran_neg_binom.Rd b/man/ran_neg_binom.Rd
index 19567ec2..062c1750 100644
--- a/man/ran_neg_binom.Rd
+++ b/man/ran_neg_binom.Rd
@@ -32,6 +32,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_norm}()},
 \code{\link{ran_pois}()},
diff --git a/man/ran_norm.Rd b/man/ran_norm.Rd
index 4bc7ce72..8ab610d8 100644
--- a/man/ran_norm.Rd
+++ b/man/ran_norm.Rd
@@ -31,6 +31,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_pois}()},
diff --git a/man/ran_pois.Rd b/man/ran_pois.Rd
index 92c70b62..716bbb10 100644
--- a/man/ran_pois.Rd
+++ b/man/ran_pois.Rd
@@ -29,6 +29,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_pois_zi.Rd b/man/ran_pois_zi.Rd
index cb707984..de49883a 100644
--- a/man/ran_pois_zi.Rd
+++ b/man/ran_pois_zi.Rd
@@ -32,6 +32,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_skewnorm.Rd b/man/ran_skewnorm.Rd
index 805c4d11..7d616456 100644
--- a/man/ran_skewnorm.Rd
+++ b/man/ran_skewnorm.Rd
@@ -37,6 +37,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/ran_student.Rd b/man/ran_student.Rd
index ffbdb0bc..4d5e8f66 100644
--- a/man/ran_student.Rd
+++ b/man/ran_student.Rd
@@ -34,6 +34,7 @@ Other ran_dist:
 \code{\link{ran_gamma}()},
 \code{\link{ran_gamma_pois}()},
 \code{\link{ran_gamma_pois_zi}()},
+\code{\link{ran_gamma_pois_zt}()},
 \code{\link{ran_lnorm}()},
 \code{\link{ran_neg_binom}()},
 \code{\link{ran_norm}()},
diff --git a/man/res_bern.Rd b/man/res_bern.Rd
index 13cf16c4..bff0965d 100644
--- a/man/res_bern.Rd
+++ b/man/res_bern.Rd
@@ -33,6 +33,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_beta_binom.Rd b/man/res_beta_binom.Rd
index b7fd0675..621cdf05 100644
--- a/man/res_beta_binom.Rd
+++ b/man/res_beta_binom.Rd
@@ -53,6 +53,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_binom.Rd b/man/res_binom.Rd
index ad2fc8e1..ec20d03c 100644
--- a/man/res_binom.Rd
+++ b/man/res_binom.Rd
@@ -35,6 +35,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_gamma.Rd b/man/res_gamma.Rd
index 6b98bcfd..9f1bd8cb 100644
--- a/man/res_gamma.Rd
+++ b/man/res_gamma.Rd
@@ -34,6 +34,7 @@ Other res_dist:
 \code{\link{res_binom}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_gamma_pois.Rd b/man/res_gamma_pois.Rd
index c8d276cd..f36e9b03 100644
--- a/man/res_gamma_pois.Rd
+++ b/man/res_gamma_pois.Rd
@@ -35,6 +35,7 @@ Other res_dist:
 \code{\link{res_binom}()},
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_gamma_pois_zi.Rd b/man/res_gamma_pois_zi.Rd
index 69aa04c0..21982b30 100644
--- a/man/res_gamma_pois_zi.Rd
+++ b/man/res_gamma_pois_zi.Rd
@@ -44,6 +44,7 @@ Other res_dist:
 \code{\link{res_binom}()},
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_gamma_pois_zt.Rd b/man/res_gamma_pois_zt.Rd
new file mode 100644
index 00000000..201bd1fa
--- /dev/null
+++ b/man/res_gamma_pois_zt.Rd
@@ -0,0 +1,49 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/res.R
+\name{res_gamma_pois_zt}
+\alias{res_gamma_pois_zt}
+\title{Zero-Truncated Gamma-Poisson Residuals}
+\usage{
+res_gamma_pois_zt(x, lambda = 1, theta = 0, type = "dev", simulate = FALSE)
+}
+\arguments{
+\item{x}{A whole numeric vector of values greater than or equal to 1.}
+
+\item{lambda}{A non-negative numeric vector of means.}
+
+\item{theta}{A non-negative numeric vector of the dispersion for the mixture
+models (student, gamma-Poisson and beta-binomial).}
+
+\item{type}{A string of the residual type. 'raw' for raw residuals 'dev' for
+deviance residuals and 'data' for the data.}
+
+\item{simulate}{A flag specifying whether to simulate residuals.}
+}
+\value{
+An numeric vector of the corresponding residuals.
+}
+\description{
+Residuals for the zero-truncated gamma-Poisson (negative binomial)
+distribution. Centring and standardisation use the truncated mean
+\verb{lambda / (1 - P(0; lambda, theta))} and the truncated variance.
+}
+\examples{
+res_gamma_pois_zt(c(1, 2, 5), 2, 1)
+}
+\seealso{
+Other res_dist: 
+\code{\link{res_bern}()},
+\code{\link{res_beta_binom}()},
+\code{\link{res_binom}()},
+\code{\link{res_gamma}()},
+\code{\link{res_gamma_pois}()},
+\code{\link{res_gamma_pois_zi}()},
+\code{\link{res_lnorm}()},
+\code{\link{res_neg_binom}()},
+\code{\link{res_norm}()},
+\code{\link{res_pois}()},
+\code{\link{res_pois_zi}()},
+\code{\link{res_skewnorm}()},
+\code{\link{res_student}()}
+}
+\concept{res_dist}
diff --git a/man/res_lnorm.Rd b/man/res_lnorm.Rd
index 948bbed0..bb4d559b 100644
--- a/man/res_lnorm.Rd
+++ b/man/res_lnorm.Rd
@@ -36,6 +36,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
 \code{\link{res_pois}()},
diff --git a/man/res_neg_binom.Rd b/man/res_neg_binom.Rd
index c8680a53..6c2e8412 100644
--- a/man/res_neg_binom.Rd
+++ b/man/res_neg_binom.Rd
@@ -36,6 +36,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_norm}()},
 \code{\link{res_pois}()},
diff --git a/man/res_norm.Rd b/man/res_norm.Rd
index 0e964dd3..1ac62259 100644
--- a/man/res_norm.Rd
+++ b/man/res_norm.Rd
@@ -35,6 +35,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_pois}()},
diff --git a/man/res_pois.Rd b/man/res_pois.Rd
index 68e0ff4f..7d490785 100644
--- a/man/res_pois.Rd
+++ b/man/res_pois.Rd
@@ -33,6 +33,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_pois_zi.Rd b/man/res_pois_zi.Rd
index 43b3df93..fb42947d 100644
--- a/man/res_pois_zi.Rd
+++ b/man/res_pois_zi.Rd
@@ -35,6 +35,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_skewnorm.Rd b/man/res_skewnorm.Rd
index 0d2c53c5..765e35ec 100644
--- a/man/res_skewnorm.Rd
+++ b/man/res_skewnorm.Rd
@@ -39,6 +39,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/man/res_student.Rd b/man/res_student.Rd
index 52d7d783..d30a4c71 100644
--- a/man/res_student.Rd
+++ b/man/res_student.Rd
@@ -38,6 +38,7 @@ Other res_dist:
 \code{\link{res_gamma}()},
 \code{\link{res_gamma_pois}()},
 \code{\link{res_gamma_pois_zi}()},
+\code{\link{res_gamma_pois_zt}()},
 \code{\link{res_lnorm}()},
 \code{\link{res_neg_binom}()},
 \code{\link{res_norm}()},
diff --git a/tests/testthat/test-dev.R b/tests/testthat/test-dev.R
index 6bd2a5e0..290a0f19 100644
--- a/tests/testthat/test-dev.R
+++ b/tests/testthat/test-dev.R
@@ -316,6 +316,34 @@ test_that("gamma_pois_zi vectorized", {
   expect_equal(dev_gamma_pois_zi(0:3, 3:0, 0:3, seq(0, 1, length.out = 4)), c(6, 1.0464962875291, 2.41568518074605, Inf))
 })
 
+test_that("gamma_pois_zt missing values", {
+  expect_identical(dev_gamma_pois_zt(logical(0), integer(0), numeric(0)), numeric(0))
+  expect_identical(dev_gamma_pois_zt(NA, 1, 1), NA_real_)
+  expect_identical(dev_gamma_pois_zt(1, NA, 1), NA_real_)
+  expect_identical(dev_gamma_pois_zt(1, 1, NA), NA_real_)
+})
+
+test_that("gamma_pois_zt matches NB deviance approximation", {
+  # dev_gamma_pois_zt uses untruncated NB deviance (truncation correction
+  # to the deviance is small at typical x values).
+  expect_identical(dev_gamma_pois_zt(3, 3), dev_gamma_pois(3, 3))
+  expect_equal(dev_gamma_pois_zt(c(1, 3, 4), 3, 2), dev_gamma_pois(c(1, 3, 4), 3, 2))
+  expect_equal(dev_gamma_pois_zt(c(1, 3, 5), 3, 1), dev_gamma_pois(c(1, 3, 5), 3, 1))
+})
+
+test_that("gamma_pois_zt deviance residuals", {
+  # Residuals are signed by sign(x - trunc_mean) and scaled by sqrt(dev_NB).
+  # When x equals lambda, deviance is 0.
+  expect_identical(dev_gamma_pois_zt(3, 3, res = TRUE), 0)
+  expect_equal(dev_gamma_pois_zt(c(1, 2, 5), 2, 1, res = TRUE),
+               c(-0.485351492542021, 0, 0.91931558891504))
+  # Sign of nonzero residuals follows sign(x - trunc_mean). At x = lambda
+  # the untruncated NB deviance is 0 so the residual is exactly 0.
+  trunc_mean_at_2_1 <- 2 / (1 - dnbinom(0, mu = 2, size = 1))  # = 3
+  expect_equal(sign(dev_gamma_pois_zt(c(1, 5), 2, 1, res = TRUE)),
+               sign(c(1, 5) - trunc_mean_at_2_1))
+})
+
 test_that("gamma_pois_zi vectorized missing values", {
   expect_equal(dev_gamma_pois_zi(c(NA, 1), 0:1, 0:1, 0:1), c(NA, Inf))
   expect_equal(dev_gamma_pois_zi(c(0, NA), 0:1, 0:1, 0:1), c(0, NA))
diff --git a/tests/testthat/test-log-lik.R b/tests/testthat/test-log-lik.R
index aca1621c..fd354232 100644
--- a/tests/testthat/test-log-lik.R
+++ b/tests/testthat/test-log-lik.R
@@ -113,6 +113,35 @@ test_that("gamma_pois_zi vectorized", {
   expect_equal(log_lik_gamma_pois_zi(0:3, 3:0, 0:3, seq(0, 1, length.out = 4)), c(-3, -1.90954250488444, -3.43967790223022, -Inf))
 })
 
+test_that("gamma_pois_zt missing values", {
+  expect_identical(log_lik_gamma_pois_zt(numeric(0), numeric(0), numeric(0)), numeric(0))
+  expect_identical(log_lik_gamma_pois_zt(NA, 1, 1), NA_real_)
+  expect_identical(log_lik_gamma_pois_zt(1, NA, 1), NA_real_)
+  expect_identical(log_lik_gamma_pois_zt(1, 1, NA), NA_real_)
+})
+
+test_that("gamma_pois_zt known values", {
+  expect_equal(log_lik_gamma_pois_zt(1, 2), -1.1614393615712)
+  expect_equal(log_lik_gamma_pois_zt(2, 2), -1.1614393615712)
+  expect_equal(log_lik_gamma_pois_zt(1, 2, 0.5), -1.09861228866811)
+  expect_equal(log_lik_gamma_pois_zt(c(1, 3, 4), 3, 1),
+               c(-1.38629436111989, -1.96165850602345, -2.24934057847523))
+})
+
+test_that("gamma_pois_zt vectorized", {
+  expect_equal(log_lik_gamma_pois_zt(1:3, 2, 0),
+               c(-1.1614393615712, -1.1614393615712, -1.56690446967936))
+})
+
+test_that("gamma_pois_zt matches truncated NB conditional", {
+  # log L_zt(x; mu, theta) = log NB(x; mu, theta) - log(1 - P(0; mu, theta))
+  expect_equal(
+    log_lik_gamma_pois_zt(2, 3, 1),
+    dnbinom(2, mu = 3, size = 1, log = TRUE) -
+      log1p(-dnbinom(0, mu = 3, size = 1))
+  )
+})
+
 test_that("gamma missing values", {
   expect_identical(log_lik_gamma(NA), NA_real_)
   expect_identical(log_lik_gamma(1, NA), NA_real_)
diff --git a/tests/testthat/test-ran.R b/tests/testthat/test-ran.R
index 0c375d75..185c0b6e 100644
--- a/tests/testthat/test-ran.R
+++ b/tests/testthat/test-ran.R
@@ -146,6 +146,39 @@ test_that("ran_gamma_pois_zi", {
   expect_identical(ran_gamma_pois_zi(1, c(0.1, 100), theta = 1, prob = 0.5), 0L)
 })
 
+test_that("ran_gamma_pois_zt", {
+  expect_error(ran_gamma_pois_zt(NA_integer_))
+  expect_error(ran_gamma_pois_zt(integer(0)))
+  expect_identical(ran_gamma_pois_zt(0L), integer(0))
+  set.seed(101)
+  expect_identical(ran_gamma_pois_zt(), 1L)
+  expect_identical(ran_gamma_pois_zt(2), c(3L, 2L))
+  expect_identical(ran_gamma_pois_zt(2, 10), c(7L, 7L))
+  expect_identical(ran_gamma_pois_zt(2, c(0.1, 100)), c(1L, 95L))
+  set.seed(101)
+  expect_identical(ran_gamma_pois_zt(2, theta = 10), c(2L, 1L))
+  expect_identical(ran_gamma_pois_zt(2, 10, theta = 10), c(7L, 13L))
+  expect_identical(ran_gamma_pois_zt(2, c(0.1, 100), theta = 1), c(1L, 12L))
+})
+
+test_that("ran_gamma_pois_zt always >= 1", {
+  set.seed(42)
+  # Even with very small lambda where most untruncated draws would be 0,
+  # rejection sampling guarantees all returned values are >= 1.
+  y <- ran_gamma_pois_zt(2000, lambda = 0.5, theta = 1)
+  expect_true(all(y >= 1))
+  expect_length(y, 2000)
+})
+
+test_that("ran_gamma_pois_zt mean approximates truncated mean", {
+  set.seed(42)
+  lambda <- 5
+  theta <- 0.5
+  y <- ran_gamma_pois_zt(20000, lambda = lambda, theta = theta)
+  expected <- lambda / (1 - dnbinom(0, mu = lambda, size = 1 / theta))
+  expect_equal(mean(y), expected, tolerance = 0.05)
+})
+
 test_that("ran_gamma", {
   expect_error(ran_gamma(NA_integer_))
   expect_error(ran_gamma(integer(0)))
diff --git a/tests/testthat/test-res.R b/tests/testthat/test-res.R
index 0bbf882c..bc7bdaae 100644
--- a/tests/testthat/test-res.R
+++ b/tests/testthat/test-res.R
@@ -380,6 +380,50 @@ test_that("res_gamma_pois_zi", {
   expect_equal(sd(res), 0.860406154019735)
 })
 
+test_that("res_gamma_pois_zt", {
+  expect_identical(res_gamma_pois_zt(integer(0), integer(0), integer(0)), numeric(0))
+  expect_identical(res_gamma_pois_zt(NA, 1, 1), NA_real_)
+  expect_identical(res_gamma_pois_zt(1, NA, 1), NA_real_)
+  expect_identical(res_gamma_pois_zt(1, 1, NA), NA_real_)
+  expect_error(res_gamma_pois_zt(1, 1, 1, type = "unknown"))
+
+  # type = "data"
+  expect_identical(res_gamma_pois_zt(1:3, 2, 1, type = "data"), 1:3)
+
+  # type = "raw" — residual is x - truncated mean
+  # truncated mean for lambda = 2, theta = 1 is 2 / (1 - 1/3) = 3
+  expect_equal(res_gamma_pois_zt(c(1, 3, 5), 2, 1, type = "raw"), c(-2, 0, 2))
+  # When theta -> 0, truncated mean approaches lambda for moderate lambda
+  expect_equal(
+    res_gamma_pois_zt(c(1, 10, 20), 10, 0, type = "raw"),
+    c(1, 10, 20) - 10 / (1 - exp(-10))
+  )
+
+  # type = "standardized" — residual is (x - trunc_mean) / sqrt(trunc_var)
+  # For lambda = 3, theta = 1: trunc_mean = 4, trunc_var = 12
+  expect_equal(
+    res_gamma_pois_zt(c(1, 3, 4), 3, 1, type = "standardized"),
+    (c(1, 3, 4) - 4) / sqrt(12)
+  )
+
+  # type = "dev" — uses untruncated NB deviance, signed by sign(x - trunc_mean)
+  expect_equal(
+    res_gamma_pois_zt(c(1, 2, 5), 2, 1),
+    c(-0.485351492542021, 0, 0.91931558891504)
+  )
+
+  # simulate = TRUE replaces x with random draws from ZT distribution
+  set.seed(101)
+  expect_true(all(res_gamma_pois_zt(1:5, 2, 1, simulate = TRUE, type = "data") >= 1))
+
+  # Standardised residuals on simulated data approximate N(0, 1)
+  set.seed(101)
+  res <- res_gamma_pois_zt(rep(2, 10000), 2, 1, simulate = TRUE,
+                           type = "standardized")
+  expect_equal(mean(res), 0, tolerance = 0.05)
+  expect_equal(sd(res), 1, tolerance = 0.05)
+})
+
 test_that("res_gamma", {
   expect_identical(res_gamma(integer(0), integer(0), integer(0)), numeric(0))
   expect_identical(res_gamma(1, 1, 1), 0)

From c861810028b444db793af1cefc5ee02e9e433113 Mon Sep 17 00:00:00 2001
From: "seb.dalgarno@gmail.com" <seb@poissonconsulting.ca>
Date: Fri, 8 May 2026 11:28:52 -0700
Subject: [PATCH 2/3] simplify zero-truncated gamma-Poisson roxygen
 descriptions

Match the brevity of the existing _zi family: keep titles, drop the
multi-paragraph descriptions, keep at most a one-line context note
where it adds something the title does not.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
---
 R/dev.R                      | 9 ++-------
 R/log-lik.R                  | 3 ---
 R/ran.R                      | 4 +---
 R/res.R                      | 4 +---
 man/dev_gamma_pois_zt.Rd     | 9 ++-------
 man/log_lik_gamma_pois_zt.Rd | 3 +--
 man/ran_gamma_pois_zt.Rd     | 4 +---
 man/res_gamma_pois_zt.Rd     | 4 +---
 8 files changed, 9 insertions(+), 31 deletions(-)

diff --git a/R/dev.R b/R/dev.R
index 44bc07a0..43f61b24 100644
--- a/R/dev.R
+++ b/R/dev.R
@@ -246,13 +246,8 @@ dev_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0, res = FALSE) {
 
 #' Zero-Truncated Gamma-Poisson Deviances
 #'
-#' Deviances for the zero-truncated gamma-Poisson (negative binomial)
-#' distribution. The deviance is approximated by the untruncated
-#' gamma-Poisson deviance; deviance residuals are centred on the
-#' truncated mean `lambda / (1 - P(0; lambda, theta))`. The saturated
-#' MLE for the zero-truncated likelihood differs from the observation
-#' and has no closed form; the approximation is accurate for typical
-#' `x` and small zero-mass `P(0; lambda, theta)`.
+#' Approximated by the untruncated gamma-Poisson deviance with deviance
+#' residuals centred on the truncated mean.
 #'
 #' @inheritParams params
 #' @param x A whole numeric vector of values greater than or equal to 1.
diff --git a/R/log-lik.R b/R/log-lik.R
index a3ffbdc6..c9563dfa 100644
--- a/R/log-lik.R
+++ b/R/log-lik.R
@@ -207,9 +207,6 @@ log_lik_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0) {
 
 #' Zero-Truncated Gamma-Poisson Log-Likelihood
 #'
-#' Density of the gamma-Poisson (negative binomial) distribution conditional
-#' on `x >= 1`.
-#'
 #' @inheritParams params
 #' @param x A whole numeric vector of values greater than or equal to 1.
 #'
diff --git a/R/ran.R b/R/ran.R
index aac41de2..d2f9f4b7 100644
--- a/R/ran.R
+++ b/R/ran.R
@@ -99,9 +99,7 @@ ran_gamma_pois_zi <- function(n = 1, lambda = 1, theta = 0, prob = 0) {
 
 #' Zero-Truncated Gamma-Poisson Random Samples
 #'
-#' Draws from the gamma-Poisson (negative binomial) distribution conditional
-#' on the value being greater than or equal to 1. Implemented via rejection
-#' sampling: zero draws are redrawn until accepted.
+#' Drawn via rejection sampling so all returned values are at least 1.
 #'
 #' @inheritParams params
 #' @return A numeric vector of the random samples.
diff --git a/R/res.R b/R/res.R
index aa9361a7..f6119522 100644
--- a/R/res.R
+++ b/R/res.R
@@ -175,9 +175,7 @@ res_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0, type = "dev",
 
 #' Zero-Truncated Gamma-Poisson Residuals
 #'
-#' Residuals for the zero-truncated gamma-Poisson (negative binomial)
-#' distribution. Centring and standardisation use the truncated mean
-#' `lambda / (1 - P(0; lambda, theta))` and the truncated variance.
+#' Centred on the truncated mean and standardised by the truncated variance.
 #'
 #' @inheritParams params
 #' @param x A whole numeric vector of values greater than or equal to 1.
diff --git a/man/dev_gamma_pois_zt.Rd b/man/dev_gamma_pois_zt.Rd
index 0e279209..2dad5ec0 100644
--- a/man/dev_gamma_pois_zt.Rd
+++ b/man/dev_gamma_pois_zt.Rd
@@ -21,13 +21,8 @@ opposed to the deviance.}
 An numeric vector of the corresponding deviances or deviance residuals.
 }
 \description{
-Deviances for the zero-truncated gamma-Poisson (negative binomial)
-distribution. The deviance is approximated by the untruncated
-gamma-Poisson deviance; deviance residuals are centred on the
-truncated mean \verb{lambda / (1 - P(0; lambda, theta))}. The saturated
-MLE for the zero-truncated likelihood differs from the observation
-and has no closed form; the approximation is accurate for typical
-\code{x} and small zero-mass \verb{P(0; lambda, theta)}.
+Approximated by the untruncated gamma-Poisson deviance with deviance
+residuals centred on the truncated mean.
 }
 \examples{
 dev_gamma_pois_zt(c(1, 3, 4), 3, 2)
diff --git a/man/log_lik_gamma_pois_zt.Rd b/man/log_lik_gamma_pois_zt.Rd
index 5e46fd9a..cb9b6c2e 100644
--- a/man/log_lik_gamma_pois_zt.Rd
+++ b/man/log_lik_gamma_pois_zt.Rd
@@ -18,8 +18,7 @@ models (student, gamma-Poisson and beta-binomial).}
 An numeric vector of the corresponding log-likelihoods.
 }
 \description{
-Density of the gamma-Poisson (negative binomial) distribution conditional
-on \code{x >= 1}.
+Zero-Truncated Gamma-Poisson Log-Likelihood
 }
 \examples{
 log_lik_gamma_pois_zt(c(1, 3, 4), 3, 1)
diff --git a/man/ran_gamma_pois_zt.Rd b/man/ran_gamma_pois_zt.Rd
index 7d027d9f..fb2a1064 100644
--- a/man/ran_gamma_pois_zt.Rd
+++ b/man/ran_gamma_pois_zt.Rd
@@ -19,9 +19,7 @@ models (student, gamma-Poisson and beta-binomial).}
 A numeric vector of the random samples.
 }
 \description{
-Draws from the gamma-Poisson (negative binomial) distribution conditional
-on the value being greater than or equal to 1. Implemented via rejection
-sampling: zero draws are redrawn until accepted.
+Drawn via rejection sampling so all returned values are at least 1.
 }
 \examples{
 ran_gamma_pois_zt(10, lambda = 3, theta = 1)
diff --git a/man/res_gamma_pois_zt.Rd b/man/res_gamma_pois_zt.Rd
index 201bd1fa..ec3c4b8b 100644
--- a/man/res_gamma_pois_zt.Rd
+++ b/man/res_gamma_pois_zt.Rd
@@ -23,9 +23,7 @@ deviance residuals and 'data' for the data.}
 An numeric vector of the corresponding residuals.
 }
 \description{
-Residuals for the zero-truncated gamma-Poisson (negative binomial)
-distribution. Centring and standardisation use the truncated mean
-\verb{lambda / (1 - P(0; lambda, theta))} and the truncated variance.
+Centred on the truncated mean and standardised by the truncated variance.
 }
 \examples{
 res_gamma_pois_zt(c(1, 2, 5), 2, 1)

From 11c4aa4c16e834f67b412aa87f5737c795aef111 Mon Sep 17 00:00:00 2001
From: "seb.dalgarno@gmail.com" <seb@poissonconsulting.ca>
Date: Fri, 8 May 2026 13:39:04 -0700
Subject: [PATCH 3/3] compute zero-truncated gamma-Poisson saturated MLE
 numerically

The previous deviance reused the untruncated NB approximation. Replace
with per-observation numerical MLE via parallel_optimize, paralleling
the beta-binomial pattern. Test values updated accordingly.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
---
 R/dev.R                   | 38 +++++++++++++++++++++++++++++++++++---
 man/dev_gamma_pois_zt.Rd  |  4 ++--
 tests/testthat/test-dev.R | 30 +++++++++++++++++-------------
 tests/testthat/test-res.R |  6 ++++--
 4 files changed, 58 insertions(+), 20 deletions(-)

diff --git a/R/dev.R b/R/dev.R
index 43f61b24..1a5f52e5 100644
--- a/R/dev.R
+++ b/R/dev.R
@@ -246,8 +246,8 @@ dev_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0, res = FALSE) {
 
 #' Zero-Truncated Gamma-Poisson Deviances
 #'
-#' Approximated by the untruncated gamma-Poisson deviance with deviance
-#' residuals centred on the truncated mean.
+#' Saturated log-likelihood is found numerically per observation since the
+#' truncated likelihood has no closed-form maximum.
 #'
 #' @inheritParams params
 #' @param x A whole numeric vector of values greater than or equal to 1.
@@ -259,7 +259,28 @@ dev_gamma_pois_zi <- function(x, lambda = 1, theta = 0, prob = 0, res = FALSE) {
 #' @examples
 #' dev_gamma_pois_zt(c(1, 3, 4), 3, 2)
 dev_gamma_pois_zt <- function(x, lambda = 1, theta = 0, res = FALSE) {
-  dev <- dev_gamma_pois(x, lambda = lambda, theta = theta, res = FALSE)
+  force(x)
+  args_not_na <- !is.na(x + lambda + theta)
+  if (length(theta) == 1) {
+    theta_vec <- theta
+  } else {
+    theta_vec <- theta[args_not_na]
+  }
+  opt_lambda <- rep(NA_real_, length(args_not_na))
+  if (any(args_not_na)) {
+    upper <- max(x[args_not_na])
+    if (!is.finite(upper) || upper <= 0) upper <- 1
+    opt_lambda[args_not_na] <- parallel_optimize(
+      f = make_opt_gamma_pois_zt(x[args_not_na], theta_vec),
+      interval = c(.Machine$double.eps, upper),
+      N = sum(args_not_na)
+    )
+  }
+  dev1 <- log_lik_gamma_pois_zt(x = x, lambda = opt_lambda, theta = theta)
+  dev2 <- log_lik_gamma_pois_zt(x = x, lambda = lambda, theta = theta)
+  dev <- dev1 - dev2
+  dev[!is.na(dev) & dev < 0 & dev > -1e-7] <- 0
+  dev <- dev * 2
   if (vld_false(res)) {
     return(dev)
   }
@@ -268,6 +289,17 @@ dev_gamma_pois_zt <- function(x, lambda = 1, theta = 0, res = FALSE) {
   dev_res(x, trunc_mean, dev)
 }
 
+# Objective function for optimization.
+make_opt_gamma_pois_zt <- function(x, theta) {
+  force(x)
+  force(theta)
+  function(lambda) {
+    out <- -log_lik_gamma_pois_zt(x = x, lambda = lambda, theta = theta)
+    out[is.nan(out)] <- Inf
+    out
+  }
+}
+
 #' Log-Normal Deviances
 #'
 #' @inheritParams params
diff --git a/man/dev_gamma_pois_zt.Rd b/man/dev_gamma_pois_zt.Rd
index 2dad5ec0..df87402b 100644
--- a/man/dev_gamma_pois_zt.Rd
+++ b/man/dev_gamma_pois_zt.Rd
@@ -21,8 +21,8 @@ opposed to the deviance.}
 An numeric vector of the corresponding deviances or deviance residuals.
 }
 \description{
-Approximated by the untruncated gamma-Poisson deviance with deviance
-residuals centred on the truncated mean.
+Saturated log-likelihood is found numerically per observation since the
+truncated likelihood has no closed-form maximum.
 }
 \examples{
 dev_gamma_pois_zt(c(1, 3, 4), 3, 2)
diff --git a/tests/testthat/test-dev.R b/tests/testthat/test-dev.R
index 290a0f19..2a620eff 100644
--- a/tests/testthat/test-dev.R
+++ b/tests/testthat/test-dev.R
@@ -323,23 +323,27 @@ test_that("gamma_pois_zt missing values", {
   expect_identical(dev_gamma_pois_zt(1, 1, NA), NA_real_)
 })
 
-test_that("gamma_pois_zt matches NB deviance approximation", {
-  # dev_gamma_pois_zt uses untruncated NB deviance (truncation correction
-  # to the deviance is small at typical x values).
-  expect_identical(dev_gamma_pois_zt(3, 3), dev_gamma_pois(3, 3))
-  expect_equal(dev_gamma_pois_zt(c(1, 3, 4), 3, 2), dev_gamma_pois(c(1, 3, 4), 3, 2))
-  expect_equal(dev_gamma_pois_zt(c(1, 3, 5), 3, 1), dev_gamma_pois(c(1, 3, 5), 3, 1))
+test_that("gamma_pois_zt known values", {
+  expect_equal(dev_gamma_pois_zt(c(1, 3, 4), 3, 2),
+               c(2.69098973, 0.19520565, 0.03152929),
+               tolerance = 1e-6)
+  expect_equal(dev_gamma_pois_zt(3, 3), 0.009522421, tolerance = 1e-6)
+})
+
+test_that("gamma_pois_zt deviance is non-negative", {
+  expect_true(all(dev_gamma_pois_zt(1:20, 5, 1) >= 0))
+  expect_true(all(dev_gamma_pois_zt(1:20, 10, 0.5) >= 0))
+  expect_true(all(dev_gamma_pois_zt(1:20, 1, 2) >= 0))
 })
 
 test_that("gamma_pois_zt deviance residuals", {
-  # Residuals are signed by sign(x - trunc_mean) and scaled by sqrt(dev_NB).
-  # When x equals lambda, deviance is 0.
-  expect_identical(dev_gamma_pois_zt(3, 3, res = TRUE), 0)
+  # Residuals signed by sign(x - trunc_mean) and scaled by sqrt(dev).
   expect_equal(dev_gamma_pois_zt(c(1, 2, 5), 2, 1, res = TRUE),
-               c(-0.485351492542021, 0, 0.91931558891504))
-  # Sign of nonzero residuals follows sign(x - trunc_mean). At x = lambda
-  # the untruncated NB deviance is 0 so the residual is exactly 0.
-  trunc_mean_at_2_1 <- 2 / (1 - dnbinom(0, mu = 2, size = 1))  # = 3
+               c(-1.4823038, -0.4853515, 0.6610002),
+               tolerance = 1e-6)
+  # Sign follows sign(x - trunc_mean). For lambda = 2, theta = 1,
+  # trunc_mean = 2 / (1 - 1/3) = 3.
+  trunc_mean_at_2_1 <- 2 / (1 - dnbinom(0, mu = 2, size = 1))
   expect_equal(sign(dev_gamma_pois_zt(c(1, 5), 2, 1, res = TRUE)),
                sign(c(1, 5) - trunc_mean_at_2_1))
 })
diff --git a/tests/testthat/test-res.R b/tests/testthat/test-res.R
index bc7bdaae..32ad2464 100644
--- a/tests/testthat/test-res.R
+++ b/tests/testthat/test-res.R
@@ -406,10 +406,12 @@ test_that("res_gamma_pois_zt", {
     (c(1, 3, 4) - 4) / sqrt(12)
   )
 
-  # type = "dev" — uses untruncated NB deviance, signed by sign(x - trunc_mean)
+  # type = "dev" — saturated MLE found numerically per observation,
+  # signed by sign(x - trunc_mean)
   expect_equal(
     res_gamma_pois_zt(c(1, 2, 5), 2, 1),
-    c(-0.485351492542021, 0, 0.91931558891504)
+    c(-1.4823038, -0.4853515, 0.6610002),
+    tolerance = 1e-6
   )
 
   # simulate = TRUE replaces x with random draws from ZT distribution