Executive Summary
Predicted customer value over the next 365 days
The short answer
The model projects 17,155 in total customer value over the next 365 days—an average of 14.3 per customer, but a median of only 0.33, with the top 10% capturing 83.5% of the total. The forecast held up well on held-out history: it predicted 733.6 repeat purchases against 648 actually observed over 219 unseen days (a 113% aggregate rate), correlating 0.67 with what each customer really did.
The detail
Total predicted customer lifetime value (365 days): 17,154.77. Average CLV: 14.3 per customer. Median CLV: 0.33 per customer. Top 10% of customers account for 83.5% of predicted value. Expected transactions across the base: 560.5. Holdout validation: the model was refit on the first 510 days of history and asked to predict the following 219 days (time it had never seen). Correlation with actual repeat purchases: 0.668. Aggregate prediction: 733.6 predicted vs 648 actual (113% of actual). Mean absolute error: 0.577 purchases per customer. This back-test rating is good. The forecast assumes past purchasing behaviour continues and excludes costs, margin, and discounting.
What this can't tell you
The average CLV of 14.3 describes almost none of the customer base: only 15.3% of customers are projected above the mean. The median of 0.33 reflects the fact that 88.6% of customers are forecast at effectively zero value because the model judges them very likely to have already churned. These are projections under historical continuity; they cannot account for changes in behaviour, competitive dynamics, or seasonality beyond what the 729-day window shows.
Analysis Overview
BG/NBD + Gamma-Gamma lifetime value across 1,200 customers and 6,259 transactions.
The short answer
The analysis forecasts what 1,200 customers are worth over the next 365 days using two linked models: one that separates quiet customers from churned ones by tracking purchase rates and dropout risk, and another that estimates transaction value. This gives a per-customer projection that assumes purchasing behaviour stays as it has been, and excludes costs, margin, and discounting.
The detail
BG/NBD models the purchase and churn process across 1,200 customers with 6,259 transactions spanning 729 days (2024-01-01 to 2025-12-30). Of these, 1,122 (94%) made more than one purchase—the only evidence the model has for distinguishing a customer who buys rarely from one who has churned. Gamma-Gamma models transaction value by blending each customer's observed average with the population mean. The forecast horizon is 365 days, chosen as the observed history's own span capped at one year. The result is expected customer lifetime value over that next year, assuming past behaviour continues.
What this can't tell you
This is a projection, not a ledger. It cannot account for changes in purchasing behaviour, competitive moves, or seasonality beyond what the historical window reveals. Costs, margin, and discounting are not factored in unless already embedded in 'Transaction Amount'. Customers with very little history (those who joined near the end of the observation window) lean heavily on the population estimate rather than their own signal.
Data Quality
Row accounting, date parsing, and how purchase occasions were formed.
The short answer
Of 6,271 transaction rows loaded, 6,259 were usable: 12 had unparseable dates and were dropped, and 212 same-day purchases by the same customer were merged into single purchase occasions. The data spans 729 days and covers 1,200 distinct customers. Customer age is measured from each customer's first purchase to the final date in the dataset.
The detail
Initial load: 6,271 rows. Rows removed: 12 (unparseable dates in 'Transaction Date'). Rows merged: 212 transactions sharing a day with another purchase by the same customer were consolidated into single purchase occasions. Final usable rows: 6,259. Date range: 2024-01-01 to 2025-12-30 (729 days). Customer age (T) is the span from each customer's first purchase to 2025-12-30, and recency is the age at their most recent purchase—both quantities BG/NBD requires. This measurement means customers who appear only near the end of the observation window have very little history and their individual forecasts rely more on the population average than on their own signal.
What this can't tell you
Merging same-day purchases into single occasions assumes each day is a purchase occasion, not an order occasion. If a customer placed multiple orders on the same day, they are counted as one purchase occasion. Consider requesting a transaction-level export if order-level granularity is needed for operational decisions.
How Predicted Value Is Distributed
Projected 365-day value across every customer.
The short answer
The predicted value distribution is severely right-skewed: the mean (14.3) sits far above the median (0.33), only 15.3% of customers exceed the average, and 88.6% are projected at effectively zero because the model judges them already churned. This concentration means an average CLV figure misrepresents almost all of the customer base.
The detail
Each customer receives a predicted 365-day value. The distribution is strongly right-skewed: mean 14.3 versus median 0.33. Only 15.3% of customers are worth more than the average. 88.6% of customers are forecast at effectively nothing—the model assigns them a high probability of having already churned, given how long they have been quiet relative to their own historical buying rhythm. The practical consequence is that budgets set against the mean CLV will systematically overpay for the long tail (the many customers projected near zero) and underpay for the head (the few customers carrying the forecast).
What this can't tell you
The heavy concentration at zero is a feature of the BG/NBD model: it distinguishes quiet customers from churned ones by measuring recency against each customer's own purchase rate. Customers with long gaps between purchases—relative to their own history—receive high churn probability. This is accurate given the data, but it means the forecast is only as good as the recency signal in the data. If purchase patterns have shifted recently, the model may misclassify quiet-but-active customers as churned. A transaction-level export with more recent purchase dates would help clarify whether recent quietness represents churn or natural buying cycles.
Highest-Value Customers
Top 20 by predicted 365-day value, with P(alive).
| Customer | Predicted Clv | P Alive | Expected Transactions | Expected Avg Value | Past Purchases |
|---|---|---|---|---|---|
| CUST-01071 | 1733 | 0.934 | 27.19 | 63.72 | 59 |
| CUST-00659 | 881.9 | 0.943 | 14.99 | 58.83 | 28 |
| CUST-00258 | 611.3 | 0.85 | 18.82 | 32.48 | 44 |
| CUST-00475 | 534.3 | 0.899 | 9.5 | 56.26 | 23 |
| CUST-00978 | 445.3 | 0.726 | 5.86 | 76.02 | 18 |
| CUST-00637 | 418.6 | 0.931 | 7.28 | 57.48 | 18 |
| CUST-00573 | 410.9 | 0.938 | 13.22 | 31.07 | 32 |
| CUST-00309 | 377 | 0.253 | 6.72 | 56.12 | 48 |
| CUST-00749 | 366.6 | 0.805 | 6.53 | 56.17 | 16 |
| CUST-00758 | 306.1 | 0.933 | 14.05 | 21.79 | 31 |
| CUST-00887 | 304.4 | 0.978 | 22.81 | 13.35 | 60 |
| CUST-00189 | 299.1 | 0.96 | 18.33 | 16.31 | 42 |
| CUST-00151 | 282.4 | 0.757 | 3.84 | 73.55 | 10 |
| CUST-00019 | 249.9 | 0.912 | 9.12 | 27.41 | 20 |
| CUST-00009 | 234.6 | 0.819 | 3.14 | 74.61 | 9 |
| CUST-00043 | 220.9 | 0.926 | 6.72 | 32.89 | 18 |
| CUST-00747 | 219.9 | 0.88 | 6.46 | 34.05 | 16 |
| CUST-00650 | 188.8 | 0.709 | 9.88 | 19.12 | 25 |
| CUST-01173 | 179 | 0.739 | 2.15 | 83.26 | 5 |
| CUST-00328 | 169.9 | 0.941 | 13.06 | 13.01 | 30 |
The short answer
The top 20 customers account for 49.2% of the total 17,155 projected value. Nineteen of the twenty have P(alive) above 0.5, meaning the model expects them to be active. One customer stands out: CUST-00309 has predicted value of 376.97 but P(alive) of only 0.253—valuable if reactivated, but the model doubts they will return on their own.
The detail
Top customer CUST-01071 is projected at 1,732.71 in value (59 past purchases, P(alive) 0.934), driven by 27.19 expected transactions at an average of 63.72 each. CUST-00659 follows at 881.92 (28 past purchases, P(alive) 0.943). Most of the top 20 have high P(alive), meaning the model expects them to remain active. The exception is CUST-00309: 376.97 predicted value, 48 past purchases, but P(alive) only 0.253. This customer is valuable if they return, but the model already doubts they will on their own—a reactivation candidate. Expected transactions and expected average value are the two factors behind each projection; you can see which one is carrying each customer by comparing these two columns.
What this can't tell you
P(alive) is a probability derived from the customer's purchase rate and recency relative to their own buying rhythm. It is not a causal statement about why they have gone quiet. A low P(alive) paired with high value could reflect a customer who bought in bursts and is now in a natural quiet phase, or one who has genuinely churned. Outbound contact or behavioural signals (email opens, site visits) would help distinguish these cases. Consider a transaction-level export with timestamps to see whether recent activity (outside the Transaction Date field) might suggest the customer is still engaged.
Does the Forecast Hold Up?
Predicted vs actual repeat purchases over 219 held-out days.
The short answer
The model was refit on the first 510 days of history and asked to predict repeat purchases over the following 219 unseen days. It predicted 733.6 purchases and 648 actually occurred—a 13% over-prediction in aggregate. The correlation between predicted and actual per-customer purchases is 0.668, indicating a tight upward pattern. This back-test rates good.
The detail
Holdout period: 219 days (2025-05-26 to 2025-12-30). Model refit on first 510 days (2024-01-01 to 2025-05-25). Each of 1,200 customers has one point: predicted repeat purchases on the horizontal axis, actual repeat purchases on the vertical. Correlation: 0.668. Mean absolute error: 0.577 purchases per customer. Aggregate: 733.6 predicted versus 648 actual (model over-predicts by 13%). The diagonal is perfect prediction; points scatter vertically because individual purchase timing is random even when the underlying rate is correct. Vertical scatter is expected and does not indicate model failure. The aggregate total and the correlation are the parts to judge.
What this can't tell you
A 13% over-prediction in aggregate is small and consistent with random purchase timing. However, this holdout period is only 219 days—about 30% of the original 729-day history. A longer held-out window would test whether the model's churn estimates remain stable over time. Additionally, the holdout test uses the same data source and customer base; it does not validate whether the model generalizes to new customers or a different market. To test generalization, consider running the model on a cohort from an earlier period and validating against a later, completely separate cohort.
Fitted Model Parameters
Maximum-likelihood estimates and fit diagnostics for both models.
| Parameter | Value | Interpretation |
|---|---|---|
| r (purchase-rate shape) | 0.7754 | Spread of the underlying purchase rate across customers |
| alpha (purchase-rate scale) | 13.67 | Purchase rate is measured in transactions per day of Transaction Date |
| a (dropout shape) | 1.33 | Shape of the churn-probability distribution |
| b (dropout shape) | 2.441 | Second shape of the churn-probability distribution |
| mean purchases per day | 0.0567 | Model-implied average rate: one purchase every 18 days |
| mean dropout probability | 0.3526 | Model-implied chance a customer stops after any given purchase |
| p (value shape) | 5.195 | Spread of transaction values within a customer |
| q (value shape) | 4.067 | Spread of average transaction value between customers (q>1 required) |
| gamma (value scale) | 17.93 | Scale of the value distribution, in units of Transaction Amount |
| population mean transaction value | 30.38 | What a randomly chosen customer's typical Transaction Amount is worth |
| BG/NBD log-likelihood | -22,497.10 | Maximised log-likelihood of the transaction/churn fit |
| Gamma-Gamma log-likelihood | -4645 | Maximised log-likelihood of the monetary fit |
| frequency vs value correlation | 0.075 | Near zero, as the Gamma-Gamma model assumes |
| customers fitted | 1200 | Distinct values of Customer ID used in the fit |
| repeat customers | 1122 | Customers with more than one purchase — the only ones informing churn |
The short answer
Both models converged. BG/NBD estimates a mean purchase rate of 0.0567 per day (roughly one purchase every 18 days for an active customer) and a 35.3% dropout probability after any purchase. Gamma-Gamma puts the population's typical transaction value at 30.38 (matching the observed average of 30.08 on repeat transactions). Repeat frequency and average transaction value correlate at only 0.075, consistent with the independence assumption the Gamma-Gamma model requires.
The detail
BG/NBD parameters: r = 0.7754, alpha = 13.6709, a = 1.3296, b = 2.4408. Derived quantities: mean purchase rate 0.0567 per day; mean dropout probability 0.3526 (35.26%). Gamma-Gamma parameters: p = 5.1954, q = 4.0674, gamma = 17.9345. Population mean transaction value: 30.38 (observed average on repeat transactions: 30.08). Log-likelihoods: BG/NBD −22497.1; Gamma-Gamma −4645.2. Frequency-value correlation: r = 0.075. Fitted on 1,200 customers; 1,122 repeat customers informed the churn component. The raw shape parameters (r, alpha, a, b) are only weakly identified individually; the ratios in this table (mean rate, mean dropout probability) are the trustworthy summaries.
What this can't tell you
Weak individual identification of r and alpha means small changes in the data can move these parameters, but the derived mean purchase rate (0.0567) is stable. The near-zero frequency-value correlation (0.075) supports the Gamma-Gamma independence assumption, but it does not rule out unmeasured confounding—for example, customer segment or product category might drive both metrics in opposite directions. A segmented fit (e.g., by product line or acquisition cohort) would test whether the population-level parameters hide important heterogeneity.
Predictive Customer Lifetime Value — BG/NBD + Gamma-Gamma
Fits the standard buy-till-you-die model to a raw transaction log and projects each customer forward: how many more purchases they are likely to make, how likely they are to still be active at all, what a typical transaction of theirs is worth, and therefore what they are worth over a stated horizon.
Why This Method?
RFM tells you what a customer HAS done — it is a backward-looking segmentation of observed behaviour. BG/NBD + Gamma-Gamma is the forward-looking counterpart: it estimates the latent purchase rate and the latent churn probability behind that behaviour, so a customer who is quiet because they buy twice a year is separated from one who is quiet because they are gone. The two are complements — segment with RFM, forecast with this.
What This Analysis Covers
- BG/NBD fit by maximum likelihood (transaction + dropout process)
- Gamma-Gamma fit by maximum likelihood (monetary value)
- Per-customer expected transactions, P(alive), expected value, and CLV
- A calibration/holdout back-test: predicted vs actual repeat purchases
- An explicit check of the frequency/monetary independence assumption
Standard Library
Platform standard-library module (LAT-1545): runs on ANY transaction log via the semantic mapping {customer, date, amount}. All narrative is derived from the user's own column names and computed values. The model is fit here by hand (optim on the exact log-likelihoods) — no BTYD package.
suppressPackageStartupMessages(library(DT))
suppressPackageStartupMessages(library(htmlwidgets))
suppressPackageStartupMessages(library(arrow))
suppressPackageStartupMessages(library(knitr))
suppressPackageStartupMessages(library(rmarkdown))
suppressPackageStartupMessages(library(dplyr))
suppressPackageStartupMessages(library(tidyr))
suppressPackageStartupMessages(library(ggplot2))
suppressPackageStartupMessages(library(stringr))
suppressPackageStartupMessages(library(lubridate))
suppressPackageStartupMessages(library(broom))
suppressPackageStartupMessages(library(Matrix))
suppressPackageStartupMessages(library(cluster))
suppressPackageStartupMessages(library(data.table))Core Analysis Pipeline
Multi-format date parser
Tries each candidate format across the whole vector, keeps the format that parses the most values, then back-fills stragglers from the remaining formats. Returns a Date vector (NA where unparseable).
parse_dates_multi <- function(x) {
if (inherits(x, "Date")) return(x)
if (inherits(x, "POSIXct") || inherits(x, "POSIXlt")) return(as.Date(x))
s <- trimws(as.character(x))
s[is.na(s) | s == ""] <- NA_character_
fmts <- c("%Y-%m-%d", "%Y/%m/%d", "%m/%d/%Y", "%d/%m/%Y", "%m-%d-%Y",
"%d-%m-%Y", "%d.%m.%Y", "%b %d, %Y", "%B %d, %Y",
"%d %b %Y", "%d %B %Y")
parsed <- lapply(fmts, function(f) suppressWarnings(as.Date(s, format = f)))
counts <- vapply(parsed, function(p) sum(!is.na(p)), integer(1))
best <- parsed[[which.max(counts)]] # counts is integer, never all-NA
for (p in parsed[order(-counts)]) {
fill <- is.na(best) & !is.na(p)
if (any(fill)) best[fill] <- p[fill]
}
best
}log(exp(u) + exp(v)) without overflow
v is allowed to be -Inf (the x == 0 branch of the BG/NBD likelihood).
log_sum_exp2 <- function(u, v) {
m <- pmax(u, v)
out <- rep(-Inf, length(m))
ok <- is.finite(m)
out[ok] <- m[ok] + log(exp(u[ok] - m[ok]) + exp(v[ok] - m[ok]))
out
}Gaussian hypergeometric 2F1 by its defining power series
2F1(a,b;c;z) = sum_n (a)_n (b)_n / (c)_n * z^n / n!. Convergent for |z| < 1, which always holds here because z = t / (alpha + T + t) with all terms positive. Vectorised over z; iterates until every element stops changing. Returns NA for elements that fail to settle so callers can refuse rather than silently propagate a truncated series.
h2f1_series <- function(a, b, c, z, max_iter = 5000) {
n <- length(z)
term <- rep(1, n)
total <- rep(1, n)
settled <- rep(FALSE, n)
for (j in seq_len(max_iter)) {
prev <- total
term <- term * (a + j - 1) * (b + j - 1) / ((c + j - 1) * j) * z
total <- total + term
settled <- (total == prev)
if (all(settled)) break
}
total[!settled] <- NA_real_
total
}BG/NBD negative log-likelihood (Fader, Hardie & Lee 2005)
Parameters are carried on the log scale so optim searches an unconstrained space and every parameter stays strictly positive. x = number of REPEAT transactions t_x = age at the last transaction (recency, same units as T) T = total observed age (first transaction -> end of observation)
bgnbd_nll <- function(par_log, x, t_x, T_obs) {
pr <- exp(par_log)
if (any(!is.finite(pr)) || any(pr <= 0) || any(pr > 1e6)) return(1e12)
r <- pr[1]; alpha <- pr[2]; a <- pr[3]; b <- pr[4]
ln_a1 <- lgamma(r + x) - lgamma(r) + r * log(alpha)
ln_a2 <- lgamma(a + b) + lgamma(b + x) - lgamma(b) - lgamma(a + b + x)
ln_a3 <- -(r + x) * log(alpha + T_obs)
ln_a4 <- rep(-Inf, length(x))
pos <- x > 0
if (any(pos)) {
ln_a4[pos] <- log(a) - log(b + x[pos] - 1) -
(r + x[pos]) * log(alpha + t_x[pos])
}
ll <- ln_a1 + ln_a2 + log_sum_exp2(ln_a3, ln_a4)
if (any(!is.finite(ll))) return(1e12)
-sum(ll)
}Fit BG/NBD by maximum likelihood
Multi-start Nelder-Mead (the surface is flat and ridged, so a single start is not trustworthy) refined with BFGS. Returns ok = FALSE with a reason whenever the optimiser fails or the estimate lands somewhere the model cannot be believed at — the caller must refuse rather than degrade.
fit_bgnbd <- function(x, t_x, T_obs) {
mean_rate <- mean(x / pmax(T_obs, 1e-6))
alpha_hint <- if (is.finite(mean_rate) && mean_rate > 0) 1 / mean_rate else 10
starts <- list(
log(c(1, 1, 1, 1)),
log(c(0.5, max(alpha_hint, 1e-3), 1.2, 2.5)),
log(c(1, max(alpha_hint / 4, 1e-3), 2, 4)),
log(c(2, 10, 1.5, 1.5)),
log(c(0.25, 50, 0.8, 3))
)
best <- NULL
for (s in starts) {
fit <- tryCatch(
optim(s, bgnbd_nll, x = x, t_x = t_x, T_obs = T_obs,
method = "Nelder-Mead",
control = list(maxit = 20000, reltol = 1e-12)),
error = function(e) NULL)
if (is.null(fit) || !is.finite(fit$value)) next
if (is.null(best) || fit$value < best$value) best <- fit
}
if (is.null(best)) {
return(list(ok = FALSE,
reason = "the optimiser could not evaluate the BG/NBD likelihood at any starting point"))
}
refined <- tryCatch(
optim(best$par, bgnbd_nll, x = x, t_x = t_x, T_obs = T_obs,
method = "BFGS", control = list(maxit = 5000, reltol = 1e-12)),
error = function(e) NULL)
if (!is.null(refined) && is.finite(refined$value) && refined$value <= best$value) {
best <- refined
}
pr <- exp(best$par)
names(pr) <- c("r", "alpha", "a", "b")
if (best$convergence != 0) {
return(list(ok = FALSE, params = pr, convergence = best$convergence,
reason = sprintf("optim did not converge(code %d)", best$convergence)))
}
if (any(!is.finite(pr)) || any(pr < 1e-4) || any(pr > 1e4)) {
return(list(ok = FALSE, params = pr, convergence = best$convergence,
reason = paste0("the fitted parameters ran to an implausible boundary(",
paste(sprintf("%s=%.4g", names(pr), pr), collapse = ", "), ")")))
}
list(ok = TRUE, params = pr, loglik = -best$value,
convergence = best$convergence, reason = "")
}P(customer is still active | x, t_x, T)
bgnbd_palive <- function(pm, x, t_x, T_obs) {
r <- pm[["r"]]; alpha <- pm[["alpha"]]; a <- pm[["a"]]; b <- pm[["b"]]
extra <- rep(0, length(x))
pos <- x > 0
if (any(pos)) {
extra[pos] <- (a / (b + x[pos] - 1)) *
((alpha + T_obs[pos]) / (alpha + t_x[pos]))^(r + x[pos])
}
1 / (1 + extra)
}E[transactions in the next horizon | x, t_x, T] (FHL 2005, eq. 10)
Returns NA where the hypergeometric series or the (a-1) denominator make the closed form untrustworthy; the caller refuses on any NA.
bgnbd_expected_transactions <- function(pm, x, t_x, T_obs, horizon) {
r <- pm[["r"]]; alpha <- pm[["alpha"]]; a <- pm[["a"]]; b <- pm[["b"]]
if (abs(a - 1) < 1e-6) return(rep(NA_real_, length(x)))
cc <- a + b + x - 1
if (any(cc <= 0)) return(rep(NA_real_, length(x)))
z <- horizon / (alpha + T_obs + horizon)
hyp <- h2f1_series(r + x, b + x, cc, z)
num <- (cc / (a - 1)) *
(1 - ((alpha + T_obs) / (alpha + T_obs + horizon))^(r + x) * hyp)
den <- 1 / bgnbd_palive(pm, x, t_x, T_obs) # == 1 + delta * (...)
out <- num / den
out[!is.finite(out)] <- NA_real_
out
}Gamma-Gamma negative log-likelihood (Fader & Hardie 2013)
Fit on customers with at least one repeat transaction, using their mean repeat transaction value mbar and repeat count x.
gg_nll <- function(par_log, x, mbar) {
pr <- exp(par_log)
if (any(!is.finite(pr)) || any(pr <= 0) || any(pr > 1e6)) return(1e12)
p <- pr[1]; q <- pr[2]; g <- pr[3]
ll <- lgamma(p * x + q) - lgamma(p * x) - lgamma(q) +
q * log(g) + (p * x - 1) * log(mbar) + p * x * log(x) -
(p * x + q) * log(g + mbar * x)
if (any(!is.finite(ll))) return(1e12)
-sum(ll)
}Fit Gamma-Gamma by maximum likelihood
q > 1 is required for the expected value to be finite; anything at or below it is an outright refusal, not a caveat.
fit_gg <- function(x, mbar) {
m0 <- mean(mbar)
starts <- list(
log(c(1, 2, max(m0, 1e-3))),
log(c(2, 3, max(2 * m0, 1e-3))),
log(c(6, 4, max(m0 / 2, 1e-3))),
log(c(0.5, 1.5, max(m0, 1e-3)))
)
best <- NULL
for (s in starts) {
fit <- tryCatch(
optim(s, gg_nll, x = x, mbar = mbar, method = "Nelder-Mead",
control = list(maxit = 20000, reltol = 1e-12)),
error = function(e) NULL)
if (is.null(fit) || !is.finite(fit$value)) next
if (is.null(best) || fit$value < best$value) best <- fit
}
if (is.null(best)) {
return(list(ok = FALSE,
reason = "the optimiser could not evaluate the Gamma-Gamma likelihood at any starting point"))
}
refined <- tryCatch(
optim(best$par, gg_nll, x = x, mbar = mbar, method = "BFGS",
control = list(maxit = 5000, reltol = 1e-12)),
error = function(e) NULL)
if (!is.null(refined) && is.finite(refined$value) && refined$value <= best$value) {
best <- refined
}
pr <- exp(best$par)
names(pr) <- c("p", "q", "gamma")
if (best$convergence != 0) {
return(list(ok = FALSE, params = pr, convergence = best$convergence,
reason = sprintf("optim did not converge(code %d)", best$convergence)))
}
if (any(!is.finite(pr)) || any(pr < 1e-4) || any(pr > 1e4)) {
return(list(ok = FALSE, params = pr, convergence = best$convergence,
reason = paste0("the fitted parameters ran to an implausible boundary(",
paste(sprintf("%s=%.4g", names(pr), pr), collapse = ", "), ")")))
}
if (pr[["q"]] <= 1.0001) {
return(list(ok = FALSE, params = pr, convergence = best$convergence,
reason = sprintf(paste0("the monetary shape parameter q estimated at %.3f; ",
"at q <= 1 the model's expected transaction value is ",
"infinite, so no value figure can be believed"),
pr[["q"]])))
}
list(ok = TRUE, params = pr, loglik = -best$value,
convergence = best$convergence, reason = "")
}E[transaction value | p, q, gamma, mbar, x]
A precision-weighted blend of the customer's own observed mean and the population mean; customers with no repeat history fall back entirely to the population mean.
gg_expected_value <- function(pm, x, mbar) {
p <- pm[["p"]]; q <- pm[["q"]]; g <- pm[["gamma"]]
pop_mean <- g * p / (q - 1)
out <- rep(pop_mean, length(x))
pos <- x > 0 & is.finite(mbar) & mbar > 0
if (any(pos)) {
den <- p * x[pos] + q - 1
out[pos] <- (g * p + mbar[pos] * p * x[pos]) / den
}
out[!is.finite(out)] <- NA_real_
out
}Build the BTYD per-customer summary (x, t_x, T, mean repeat value)
Multiple transactions on the SAME DAY are treated as one transaction, the standard BTYD convention: the models describe purchase OCCASIONS.
btyd_summary <- function(tx, obs_end) {
daily <- tx %>%
group_by(customer, txn_date) %>%
summarise(amount = sum(amount), .groups = "drop")
daily %>%
arrange(customer, txn_date) %>%
group_by(customer) %>%
summarise(
first_date = min(txn_date),
last_date = max(txn_date),
n_txn = n(),
x = n() - 1,
t_x = as.numeric(max(txn_date) - min(txn_date)),
T_obs = as.numeric(obs_end - min(txn_date)),
mean_repeat_value = if (n() > 1) {
v <- amount[txn_date > min(txn_date)]
v <- v[is.finite(v) & v > 0]
if (length(v) > 0) mean(v) else NA_real_
} else NA_real_,
total_value = sum(amount),
.groups = "drop"
) %>%
as.data.frame(stringsAsFactors = FALSE)
}
compute_shared <- function(df, params, col_map = list()) {
# === SHARED EXPORTS ===
# initial_rows/final_rows/rows_removed $ row accounting
# cust_name/date_name/amount_name $ humanized user column names
# n_bad_dates/n_missing_id/n_nonpositive/n_same_day_merged $ prep counts
# date_min/date_max/span_days/horizon_days
# n_customers/n_repeaters $ integer
# bg $ list — BG/NBD fit (ok, params r/alpha/a/b, loglik)
# gg $ list — Gamma-Gamma fit (ok, params p/q/gamma, loglik)
# cust_df $ data.frame per customer incl. expected_transactions, p_alive,
# expected_avg_value, predicted_clv
# total_clv/avg_clv/median_clv/top10_share/total_expected_txn
# fm_cor / fm_cor_p / independence_violated
# hold $ list — calibration/holdout back-test (available, cor, mae,
# predicted_total, actual_total, ratio, verdict, n)
# clv_hist_df / top_customers_df / holdout_df / model_params_df
# metrics / json_output
# === /SHARED EXPORTS ===Step 1: Required semantic columns + humanized names
initial_rows <- nrow(df)
cust_name <- humanize_semantic("customer", col_map)[1]
date_name <- humanize_semantic("date", col_map)[1]
amount_name <- humanize_semantic("amount", col_map)[1]
for (need in c("customer", "date", "amount")) {
if (!need %in% names(df)) {
stop(sprintf("Required column '%s' is not mapped.",
humanize_semantic(need, col_map)[1]))
}
}Step 2: Customer id — character, drop blank/missing
cid <- trimws(as.character(df$customer))
keep_id <- !is.na(cid) & cid != ""
n_missing_id <- sum(!keep_id)
df <- df[keep_id, , drop = FALSE]
cid <- cid[keep_id]Step 3: Dates — multi-format parse; refuse if >5% unparseable
dts <- parse_dates_multi(df$date)
n_bad_dates <- sum(is.na(dts))
if (nrow(df) > 0 && n_bad_dates / nrow(df) > 0.05) {
stop(sprintf(
paste0("Could not read %s of %s values in '%s' as dates (%.1f%%). More ",
"than 5%% of the transaction dates are unparseable, so customer ",
"age and recency cannot be measured — please check that column's ",
"date format."),
format(n_bad_dates, big.mark = ","), format(nrow(df), big.mark = ","),
date_name, 100 * n_bad_dates / nrow(df)))
}
keep_dt <- !is.na(dts)
df <- df[keep_dt, , drop = FALSE]; cid <- cid[keep_dt]; dts <- dts[keep_dt]Step 4: Amount — 95% numeric coercion rule; non-positive amounts
(refunds, zero-value rows) cannot be used by the Gamma-Gamma model.
val <- df$amount
if (!is.numeric(val)) {
conv <- suppressWarnings(as.numeric(as.character(val)))
n_orig <- sum(!is.na(val) & trimws(as.character(val)) != "")
if (n_orig == 0 || sum(!is.na(conv)) < 0.95 * n_orig) {
stop(sprintf(paste0("Column '%s' is not numeric enough to use as a ",
"transaction amount — fewer than 95%% of its values ",
"could be read as numbers."), amount_name))
}
val <- conv
}
keep_val <- !is.na(val)
cid <- cid[keep_val]; dts <- dts[keep_val]; val <- val[keep_val]
n_nonpositive <- sum(val <= 0)
if (length(val) > 0 && n_nonpositive / length(val) > 0.5) {
stop(sprintf(
paste0("%s of %s values in '%s' are zero or negative (%.1f%%). The ",
"Gamma-Gamma monetary model is defined only on positive ",
"transaction values, so with the majority non-positive no ",
"lifetime value can be estimated. Supply gross transaction ",
"amounts rather than a net/refund ledger."),
format(n_nonpositive, big.mark = ","), format(length(val), big.mark = ","),
amount_name, 100 * n_nonpositive / max(1, length(val))))
}
final_rows <- length(val)
rows_removed <- initial_rows - final_rows
tx <- data.frame(customer = cid, txn_date = dts, amount = val,
stringsAsFactors = FALSE)
if (final_rows < 1) {
stop(sprintf("No usable transactions remain after cleaning '%s', '%s', and '%s'.",
cust_name, date_name, amount_name))
}Step 5: Observation window + the BTYD per-customer summary
date_min <- min(tx$txn_date); date_max <- max(tx$txn_date)
span_days <- as.numeric(date_max - date_min)
obs_end <- date_max
cust_df <- btyd_summary(tx, obs_end)
n_same_day_merged <- final_rows - sum(cust_df$n_txn)
n_customers <- nrow(cust_df)
n_repeaters <- sum(cust_df$x > 0)Step 6: Identifiability guards — refuse loudly, never degrade
if (n_customers < 50) {
stop(sprintf(
paste0("Only %d distinct customers in '%s'. Predictive lifetime value ",
"needs at least 50 customers before the purchase-rate and churn ",
"distributions mean anything — with fewer, the fitted parameters ",
"are noise."),
n_customers, cust_name))
}
if (n_repeaters < 20) {
stop(sprintf(
paste0("Only %d of %d customers in '%s' made more than one purchase. ",
"The BG/NBD model learns the churn process from REPEAT ",
"behaviour, so with fewer than 20 repeat customers it is ",
"unidentifiable — no honest lifetime-value projection is ",
"possible from this data. This usually means the export covers ",
"too short a window, or '%s' is not a stable customer identifier."),
n_repeaters, n_customers, cust_name, cust_name))
}
if (span_days < 60) {
stop(sprintf(
paste0("The transactions in '%s' span only %.0f days (%s to %s). A ",
"lifetime-value model needs at least 60 days so that part of the ",
"history can be held back to test the forecast; without that ",
"back-test the projection would be untestable."),
date_name, span_days, format(date_min), format(date_max)))
}Step 7: Fit BG/NBD on the full observation window
bg <- fit_bgnbd(cust_df$x, cust_df$t_x, cust_df$T_obs)
if (!isTRUE(bg$ok)) {
stop(sprintf(
paste0("The BG/NBD purchase model could not be fitted to this data: %s. ",
"Rather than present lifetime-value figures the model does not ",
"support, the analysis stops here. This usually means the repeat ",
"purchase pattern in '%s' is too sparse or too irregular for the ",
"buy-till-you-die assumptions."),
bg$reason, cust_name))
}Step 8: Fit Gamma-Gamma on repeat customers with positive value
gg_rows <- cust_df$x > 0 & is.finite(cust_df$mean_repeat_value) &
cust_df$mean_repeat_value > 0
if (sum(gg_rows) < 20) {
stop(sprintf(
paste0("Only %d customers have a positive-value repeat purchase in ",
"'%s'. The Gamma-Gamma monetary model needs at least 20 to ",
"estimate transaction value, so no lifetime value can be ",
"reported."),
sum(gg_rows), amount_name))
}
gg <- fit_gg(cust_df$x[gg_rows], cust_df$mean_repeat_value[gg_rows])
if (!isTRUE(gg$ok)) {
stop(sprintf(
paste0("The Gamma-Gamma monetary model could not be fitted to '%s': %s. ",
"Expected transaction value is the money half of lifetime value, ",
"so the analysis stops rather than report a currency figure it ",
"cannot stand behind."),
amount_name, gg$reason))
}Step 9: Independence check — the Gamma-Gamma assumption, measured
fm_cor <- NA_real_; fm_cor_p <- NA_real_
if (sum(gg_rows) >= 3) {
ct <- tryCatch(cor.test(cust_df$x[gg_rows],
cust_df$mean_repeat_value[gg_rows]),
error = function(e) NULL)
if (!is.null(ct)) { fm_cor <- unname(ct$estimate); fm_cor_p <- ct$p.value }
}
independence_violated <- is.finite(fm_cor) && abs(fm_cor) >= 0.1 &&
is.finite(fm_cor_p) && fm_cor_p < 0.05Step 10: Horizon — the data's own span, capped at one year
horizon_days <- max(30, min(365, round(span_days)))Step 11: Per-customer projection
exp_txn <- bgnbd_expected_transactions(bg$params, cust_df$x, cust_df$t_x,
cust_df$T_obs, horizon_days)
if (any(!is.finite(exp_txn)) || any(exp_txn < 0, na.rm = TRUE)) {
stop(sprintf(
paste0("The fitted BG/NBD model produced impossible transaction ",
"forecasts(non-finite or negative) for %d of %d customers, ",
"which means the parameter estimate(r=%.4g, alpha=%.4g, a=%.4g, ",
"b=%.4g) is not usable for projection. No lifetime-value figures ",
"are reported."),
sum(!is.finite(exp_txn) | exp_txn < 0), n_customers,
bg$params[["r"]], bg$params[["alpha"]], bg$params[["a"]], bg$params[["b"]]))
}
p_alive <- bgnbd_palive(bg$params, cust_df$x, cust_df$t_x, cust_df$T_obs)
exp_val <- gg_expected_value(gg$params, cust_df$x, cust_df$mean_repeat_value)
if (any(!is.finite(exp_val))) {
stop("The fitted Gamma-Gamma model produced non-finite expected transaction values; no lifetime-value figures are reported.")
}
cust_df$expected_transactions <- exp_txn
cust_df$p_alive <- p_alive
cust_df$expected_avg_value <- exp_val
cust_df$predicted_clv <- exp_txn * exp_val
total_clv <- sum(cust_df$predicted_clv)
avg_clv <- mean(cust_df$predicted_clv)
median_clv <- median(cust_df$predicted_clv)
total_expected_txn <- sum(cust_df$expected_transactions)
ord <- order(-cust_df$predicted_clv)
n_top10 <- max(1, round(0.10 * n_customers))
top10_share <- 100 * sum(cust_df$predicted_clv[ord][seq_len(n_top10)]) /
max(total_clv, .Machine$double.eps)
pop_mean_value <- gg$params[["gamma"]] * gg$params[["p"]] /
(gg$params[["q"]] - 1)
observed_total <- sum(tx$amount)Step 12: Calibration / holdout back-test — the honesty backbone
Split the observation window at 70% of its span, refit BG/NBD on the calibration half only, and compare each customer's PREDICTED repeat transactions in the holdout window against what they actually did.
hold <- list(available = FALSE, reason = "")
cal_end <- date_min + as.difftime(round(0.7 * span_days), units = "days")
holdout_days <- as.numeric(date_max - cal_end)
cal_tx <- tx[tx$txn_date <= cal_end, , drop = FALSE]
if (holdout_days < 14) {
hold$reason <- sprintf("the holdout window would be only %.0f days long", holdout_days)
} else if (nrow(cal_tx) < 1) {
hold$reason <- "no transactions fall inside the calibration window"
} else {
cal_sum <- btyd_summary(cal_tx, cal_end)
cal_sum <- cal_sum[cal_sum$T_obs > 0, , drop = FALSE]
if (sum(cal_sum$x > 0) < 20) {
hold$reason <- sprintf(
"only %d customers made a repeat purchase inside the calibration window",
sum(cal_sum$x > 0))
} else {
bg_cal <- fit_bgnbd(cal_sum$x, cal_sum$t_x, cal_sum$T_obs)
if (!isTRUE(bg_cal$ok)) {
hold$reason <- paste0("the calibration-period model did not fit(",
bg_cal$reason, ")")
} else {
pred <- bgnbd_expected_transactions(bg_cal$params, cal_sum$x,
cal_sum$t_x, cal_sum$T_obs,
holdout_days)
hold_tx <- tx[tx$txn_date > cal_end, , drop = FALSE]
actual_tbl <- hold_tx %>%
group_by(customer, txn_date) %>%
summarise(.groups = "drop") %>%
group_by(customer) %>%
summarise(actual = n(), .groups = "drop") %>%
as.data.frame(stringsAsFactors = FALSE)
actual <- actual_tbl$actual[match(cal_sum$customer, actual_tbl$customer)]
actual[is.na(actual)] <- 0
ok <- is.finite(pred) & is.finite(actual)
if (sum(ok) < 20) {
hold$reason <- "the calibration model could not produce forecasts for enough customers"
} else {
pr <- pred[ok]; ac <- actual[ok]
ct <- tryCatch(cor.test(pr, ac), error = function(e) NULL)
hcor <- if (!is.null(ct)) unname(ct$estimate) else NA_real_
mae <- mean(abs(pr - ac))
ptot <- sum(pr); atot <- sum(ac)
ratio <- if (atot > 0) ptot / atot else NA_real_
verdict <- if (!is.finite(hcor) || !is.finite(ratio)) "inconclusive"
else if (hcor >= 0.5 && ratio >= 0.8 && ratio <= 1.25) "good"
else if (hcor >= 0.3 && ratio >= 0.6 && ratio <= 1.6) "fair"
else "poor"
set.seed(42)
sidx <- if (length(pr) > 2000) sample(length(pr), 2000) else seq_along(pr)
hold <- list(
available = TRUE, reason = "",
n = length(pr), cor = hcor, mae = mae,
predicted_total = ptot, actual_total = atot, ratio = ratio,
verdict = verdict, holdout_days = holdout_days,
cal_end = cal_end, cal_params = bg_cal$params,
df = data.frame(predicted_transactions = round(pr[sidx], 4),
actual_transactions = as.numeric(ac[sidx]),
stringsAsFactors = FALSE)
)
}
}
}
}Machine-readable verdict for the Gamma-Gamma independence check above — same fm_cor/fm_cor_p/independence_violated, never a second computation (LAT-1783).
assumption_checks = list(
list(name = "gamma_gamma_independence",
verdict = if (!is.finite(fm_cor) || !is.finite(fm_cor_p)) "warn"
else if (independence_violated) "fail" else "pass",
statistic = if (!is.finite(fm_cor) || !is.finite(fm_cor_p)) "correlation not estimable"
else sprintf("r = %.3g, p = %.3g", fm_cor, fm_cor_p))
)
)
}Your turn
Bring your own data and the question you actually need answered.
CympleData Scientist Send me your data and question, I’ll send you the analytics. ds@mcpanalytics.ai