Log-Normal Distribution and the Retransformation Problem

In plain language

Healthcare costs, length of hospital stay, and laboratory values like viral loads tend to follow a log-normal distribution — a bell curve when you take the logarithm of the values, but a strongly right-skewed curve on the original dollar or unit scale. When you fit a regression model to log-transformed costs and convert the predicted values back to dollars, you get the geometric mean (the typical value for a median patient), not the arithmetic mean (the true average spend that drives a budget). Correcting for this gap requires a step called Duan smearing: you use the spread of the model's residuals to adjust the back-transformed prediction upward toward the true average — and failing to do so systematically underestimates total population spending.

What the log-normal distribution is and where it comes from

A random variable Y follows a log-normal distribution when its natural logarithm log(Y) is normally distributed. More precisely: if X is normally distributed with mean mu and variance sigma-squared, then Y = exp(X) is log-normal with parameters mu and sigma-squared. The name is counterintuitive — Y itself is not normal, but log(Y) is — and this asymmetry is the source of nearly every practical complication the distribution creates for analysts.

The log-normal arises naturally wherever the data-generating process is multiplicative rather than additive. Many biological and economic processes compound sequentially. A drug's concentration at time t is the dose multiplied by a sequence of absorption, distribution, and elimination factors. A patient's annual healthcare spend is a series of claim amounts that compound through utilization rates, DRG weights, and facility price multipliers. A virus's titre doubles or halves with each replication cycle. When independent multiplicative shocks are logged, the Central Limit Theorem applies to the sum of log-terms, producing a normal distribution on the log scale and a log-normal on the original scale.

In health economics and outcomes research (HEOR), the log-normal is ubiquitous because the cost-generating mechanism is fundamentally multiplicative: a hospitalization multiplies by its DRG weight, which multiplies by facility rates, which compound with complication flags and LOS. The result is a distribution that is right-skewed — a long right tail of catastrophically expensive patients — with a mode well below the mean. Hospital length of stay, pharmacy costs, number of outpatient visits, antibody titers from serology assays, viral loads from PCR, and pharmacokinetic AUC values routinely exhibit this shape.

Arithmetic vs geometric mean: the sigma-squared-over-two term

The most important mathematical property of the log-normal is the relationship between the parameters of the underlying normal (mu, sigma-squared) and the moments on the original scale:

• Arithmetic mean (expected cost): E[Y] = exp(mu + sigma^2 / 2)
• Median (equals geometric mean): exp(mu)
• Variance: Var[Y] = (exp(sigma^2) minus 1) times exp(2*mu + sigma^2)

The sigma^2/2 term in the arithmetic mean formula is the entire story of why log- transformation and back-transformation require careful handling. When sigma^2 is small (distributions that are only mildly skewed on the log scale), the arithmetic mean and the geometric mean are close. As sigma^2 grows — as the distribution becomes more right-skewed — the arithmetic mean increasingly exceeds the geometric mean by the factor exp(sigma^2/2).

For a typical commercial claims cost distribution with sigma = 1.5 on the natural-log scale, the arithmetic mean is exp(sigma^2/2) = exp(1.125) approximately 3.08 times the geometric mean. A naive analysis that reports the geometric mean as "average cost" would understate the true population average spend by more than 200 percent. For budget-impact models and cost-effectiveness analyses, this is not a minor rounding error — it is the difference between a credible submission and a rejected one.

This has a direct consequence for regression. If you fit an ordinary least-squares (OLS) model with log(Y) as the outcome, the predicted value on the log scale, mu-hat, is an estimate of E[log Y] = mu. Back-transforming with exp(mu-hat) gives the geometric mean of Y — the median under log-normality — which is not the same as E[Y], the arithmetic mean. This discrepancy is the retransformation problem.

The retransformation problem and Duan's smearing estimator

The correct estimator for the arithmetic mean under log-OLS is:

E[Y] = exp(mu-hat) times the Duan smearing factor

where the smearing factor is the sample mean of the back-transformed OLS residuals:

Duan factor = (1/n) times the sum of exp(e_i)

and e_i = log(y_i) minus mu-hat_i are the OLS residuals on the log scale. This nonparametric correction — proposed by Duan (1983) — requires no assumption about the residual distribution and is consistent even when the log-scale errors are not exactly normal. Under homoscedastic log-normal errors (equal variance across all covariate values), the smearing factor estimates exp(sigma^2/2), exactly correcting for the sigma^2/2 bias in the arithmetic mean.

Manning and Mullahy (2001) provided the definitive applied assessment of when log-OLS with smearing is adequate and when it fails. The smearing estimator performs well when the log-scale residuals are homoscedastic — that is, when the spread of log-cost around the fitted value is similar across all patient subgroups. It can fail badly under heteroscedasticity: when sicker or more complex patients have more variable costs, sigma^2 varies across the covariate space, the smearing factor should differ by subgroup, and the single pooled smearing factor applied to the full dataset can seriously misestimate group- specific arithmetic means. In the heteroscedastic case, a generalized linear model (GLM) with a log link and gamma or Tweedie variance function is almost universally preferred for primary cost inference: it directly models the arithmetic mean on the original scale without any retransformation step and is robust to variance heterogeneity.

Interpreting the output

Consider a log-OLS regression of total annual costs on a binary treatment indicator. The estimated treatment coefficient is 0.405 on the natural-log scale (95% CI: 0.18 to 0.63).

(1) Formal interpretation. exp(0.405) = 1.50. This is the ratio of geometric means: the median cost in the treated group is approximately 50 percent higher than in the control group under log-normality. Under the additional assumption that the log-scale residuals are homoscedastic (equal variance in both groups), the smearing factors for the two groups are identical and cancel in the ratio, so exp(0.405) also approximates the ratio of arithmetic means — meaning average spending is about 50 percent higher in the treated group. If the residuals are heteroscedastic (unequal variance between groups, which is common when sicker patients are more likely to receive a given treatment), this arithmetic-mean interpretation does not hold. The ratio of arithmetic means must instead be computed as [exp(alpha-hat + beta-hat) times smearing-factor-treated] divided by [exp(alpha-hat) times smearing-factor-control], using group-specific smearing factors.

(2) Practical interpretation. "Patients in the treatment group had typical (median) costs about 50 percent higher than control patients — that is, for a control patient whose costs sit at the median, the corresponding treated patient would be expected to spend 50 percent more. Whether average total spending — the quantity relevant for a budget-impact model — is also 50 percent higher depends on whether cost variability is similar in both groups. A sensitivity analysis using group-specific smearing factors or a gamma GLM is recommended before quoting the arithmetic-mean ratio to a payer or health technology assessment body."

This formal/practical distinction is the core analytical skill for HEOR analysts reporting log-OLS results. The coefficient tells you about the geometric mean. Budget impact requires the arithmetic mean. The smearing correction bridges the two — but only reliably under homoscedasticity.

Geometric means in laboratory and pharmacokinetics reporting

For laboratory assay results — antibody titers from serology, viral loads from PCR, pharmacokinetic AUC or Cmax values, minimum inhibitory concentrations (MICs) from antimicrobial studies — reporting the geometric mean rather than the arithmetic mean is often the scientifically correct choice. These measurements are multiplicative by nature, and the ratio of geometric means is the natural measure of relative magnitude: "the treated group had a titer 4-fold higher." Confidence intervals computed on the log scale and then back- transformed are ratio CIs, directly interpretable as fold-changes with correct coverage properties. Bland and Altman (1996) provide a clear worked demonstration of this approach. For these applications, the geometric mean is the primary target estimand and no smearing correction is needed; the retransformation problem only arises when the target is the arithmetic mean on the original dollar or count scale.

Pros, cons, and trade-offs

Log-OLS (OLS on the log-transformed outcome) with Duan smearing: - Pros: straightforward to fit in any statistical package; residual diagnostics are familiar (Q-Q plots, residual-versus-fitted scatter); interpretable on the log scale as a multiplicative relationship; the coefficient is directly the log of the geometric-mean ratio; smearing correction is simple to compute; widely reported in HEOR literature. - Cons: the naive back-transform gives the geometric mean, not the arithmetic mean; the Duan smearing estimator is biased under heteroscedasticity; zero-cost patients must be excluded (log is undefined at zero) or handled separately by a two-part model; the smearing step is often omitted by analysts who do not know it is required, producing systematically underestimated costs. - When to prefer: when the estimand is the geometric mean or a fold-change ratio (lab/PK data); as a sensitivity analysis alongside a GLM primary analysis; when the log-normality assumption is well-supported and residuals are approximately homoscedastic.

Gamma GLM with log link: - Pros: directly models the arithmetic mean on the original scale; the gamma variance function accommodates heteroscedasticity (variance proportional to mean squared); no back- transformation step or smearing correction is needed; marginally consistent estimator for population mean costs; readily accommodates prediction and covariate adjustment. - Cons: requires specifying the variance function (gamma vs Tweedie vs negative binomial); less familiar to some audiences than log-OLS; the log-link coefficient is a log-mean-ratio, not a log-median-ratio; GLM convergence can be slow for very large datasets. - When to prefer: when the primary estimand is the arithmetic mean cost for budget impact or cost-effectiveness analysis; when cost heteroscedasticity is expected; this is the modern standard for primary cost analyses in HEOR per Manning and Mullahy (2001).

When to use

Use log-transformation and log-scale regression when: - The outcome is continuous, strictly positive, and right-skewed in a multiplicative pattern such that log(Y) is approximately normally distributed. - The primary estimand is the geometric mean or a ratio of geometric means — antibody titers, viral loads, pharmacokinetic parameters, sensitivity analyses on a ratio scale. - The Duan smearing correction is applied and residual heteroscedasticity has been assessed as mild (e.g., by plotting residuals against fitted values or comparing group-specific standard deviations on the log scale). - Exploratory data analysis is needed to understand distributional shape before choosing a primary model specification. - The audience or journal convention strongly favors log-transformed regression over GLMs and Duan smearing is applied correctly.

When NOT to use

Do not use log-transformation when: - The outcome contains zeros. log(0) is undefined, and adding an arbitrary constant (such as log(Y + 1) or log(Y + 0.5)) introduces scale-dependence that biases the back-transformed mean estimate in a direction and magnitude that depend entirely on the chosen constant. For outcomes with a zero spike — as in cost data where a fraction of patients have no claims — use a two-part model (logistic for any use, then log-OLS or gamma GLM conditional on positive use) or a Tweedie GLM designed for semi-continuous data. - The target estimand is the arithmetic mean and the analyst applies only the naive back-transform. This is the single most common error in HEOR cost analyses: reporting exp(mu-hat) as "the mean cost" when it estimates the geometric mean, systematically understating average spend. The arithmetic mean always requires the smearing step. - The outcome is bounded or count-valued. Log-normal is a model for continuous unbounded positive quantities. Bounded scores (0 to 100 quality-of-life instruments), binary outcomes, and count outcomes (number of hospitalizations) have their own distributional families (beta regression, logistic regression, Poisson or negative binomial) that should be used instead. - Residuals on the log scale show strong heteroscedasticity between groups. If sicker or higher-utilizing patients have substantially more variable costs, the pooled smearing factor is biased, the simple ratio exp(beta) is not the arithmetic-mean ratio between arms, and a gamma GLM should be preferred as the primary analysis with log-OLS as a sensitivity check. - The dataset is small (fewer than 30 observations) and log-normality is not supported by prior knowledge. Bootstrap-based mean estimation or a nonparametric approach may be more robust when the distributional assumption cannot be adequately assessed.

Worked example

Scenario

A health economics analyst has three patients from a commercial claims database with total annual costs of $100, $1,000, and $10,000 — values that increase by a factor of ten at each step, a pattern typical of multiplicative cost processes. The analyst wants to understand concretely why fitting a regression on the log-transformed costs and back-transforming the prediction gives the wrong average, and how the Duan smearing factor corrects it to match the true arithmetic mean.

Dataset

Annual total costs (USD) for three patients. The tenfold jumps between consecutive patients illustrate the multiplicative structure common in claims cost distributions, where a small number of high-cost patients dominate average spending.

Steps

• Step 1 — Arithmetic mean (the true average spend a budget analyst needs): (100 + 1000 + 10000) / 3 = 3700. This is the number that, multiplied by population size, gives total budget impact.

• Step 2 — Log-transform each cost using log base 10 (chosen here for clean arithmetic): log10(100) = 2; log10(1000) = 3; log10(10000) = 4.

• Step 3 — Mean of the log10 values, which is what a log-OLS intercept estimates for this one-group dataset: (2 + 3 + 4) / 3 = 3.

• Step 4 — Naive back-transform: 10^3 = 1000. This is the geometric mean — the middle patient on the log scale — NOT the arithmetic mean. Reporting $1,000 as 'average cost' understates true average spend by $2,700.

• Step 5 — Compute OLS residuals on the log10 scale by subtracting the fitted value (3) from each patient's log cost: P1 residual = 2 - 3 = -1; P2 residual = 3 - 3 = 0; P3 residual = 4 - 3 = 1.

• Step 6 — Back-transform each residual from log10 scale: P1 gives 10^(-1) = 0.1; P2 gives 10^0 = 1; P3 gives 10^1 = 10. These are the patient-level smearing terms.

• Step 7 — Duan smearing factor = average of the back-transformed residuals: (0.1 + 1 + 10) / 3 = 3.7. This multiplier captures the asymmetric pull of the upper tail that the geometric mean misses.

• Step 8 — Smeared (corrected) arithmetic mean estimate = geometric mean times smearing factor = 1000 * 3.7 = 3700. The smearing correction exactly recovers the arithmetic mean in this symmetric log-scale example, confirming the estimator is unbiased here.

Result

Arithmetic mean = (100 + 1000 + 10000) / 3 = 3700. Mean log10 = (2 + 3 + 4) / 3 = 3. Geometric mean (naive back-transform) = 10^3 = 1000, which is $2,700 less than the true average. Smearing factor = (0.1 + 1 + 10) / 3 = 3.7. Smeared Duan estimate = 1000 * 3.7 = 3700, which equals the arithmetic mean. Conclusion: the geometric mean is the right answer for "what does a typical (median) patient spend?"; the arithmetic mean is the right answer for "what will this population cost in total?"

Runnable example

import numpy as np
from scipy import stats

# ── Motivating dataset: three patient costs ──────────────────────────────
costs = np.array([100.0, 1000.0, 10000.0])

# 1. Arithmetic mean (the budget target)
arith_mean = costs.mean()
print(f"Arithmetic mean: {arith_mean:.2f}")              # 3700.00

# 2. Geometric mean via log-scale back-transform (natural log throughout)
log_costs = np.log(costs)          # natural log; mean_log = mu-hat from log-OLS
mean_log = log_costs.mean()
geo_mean = np.exp(mean_log)
print(f"Geometric mean (naive back-transform): {geo_mean:.2f}")  # 1000.00

# 3. Duan smearing factor: sample mean of back-transformed residuals
residuals = log_costs - mean_log
smearing_factor = np.exp(residuals).mean()
smeared_mean = geo_mean * smearing_factor
print(f"Smearing factor: {smearing_factor:.4f}")
print(f"Duan smeared estimate of arithmetic mean: {smeared_mean:.2f}")  # 3700.00

# ── Simulated 200-patient treatment comparison ────────────────────────────
rng = np.random.default_rng(42)
n = 200
treat = np.repeat([0, 1], n // 2)
# True model: log(cost) = 6.0 + 0.405 * treat + epsilon, epsilon ~ N(0, 1)
log_y = 6.0 + 0.405 * treat + rng.normal(0, 1.0, n)
y = np.exp(log_y)

# 4. Log-OLS with numpy (design matrix approach)
X = np.column_stack([np.ones(n), treat])
beta_hat = np.linalg.lstsq(X, log_y, rcond=None)[0]
alpha_hat, beta_treat = beta_hat
print(f"\nLog-OLS: alpha = {alpha_hat:.3f}, beta_treat = {beta_treat:.3f}")
print(f"exp(beta_treat) = {np.exp(beta_treat):.3f}  "
      f"<-- ratio of geometric means, NOT necessarily ratio of arithmetic means")

# 5. Pooled smearing factor
fitted = X @ beta_hat
resids = log_y - fitted
pooled_sf = np.exp(resids).mean()
print(f"\nPooled smearing factor: {pooled_sf:.4f}")
for arm, label in [(0, "Control"), (1, "Treated")]:
    mu_arm = alpha_hat + beta_treat * arm
    arith_arm = np.exp(mu_arm) * pooled_sf
    print(f"  {label}: geometric mean = {np.exp(mu_arm):.2f}, "
          f"arithmetic mean = {arith_arm:.2f}")

# 6. Group-specific smearing (correct under heteroscedasticity)
print("\nGroup-specific smearing factors:")
for arm, label in [(0, "Control"), (1, "Treated")]:
    mask = treat == arm
    sf_arm = np.exp(resids[mask]).mean()
    arith_arm = np.exp(alpha_hat + beta_treat * arm) * sf_arm
    print(f"  {label}: smearing factor = {sf_arm:.4f}, "
          f"arithmetic mean estimate = {arith_arm:.2f}")

# 7. Parametric log-normal fit using scipy.stats.lognorm
# scipy uses shape (sigma) and scale (exp(mu)) parameterization
shape, loc, scale = stats.lognorm.fit(costs, floc=0)  # fix loc=0 for strictly positive
mu_fit = np.log(scale)
sigma_fit = shape
print(f"\nLog-normal MLE: mu = {mu_fit:.3f}, sigma = {sigma_fit:.3f}")
print(f"  Estimated arithmetic mean: exp(mu + sigma^2/2) = "
      f"{np.exp(mu_fit + sigma_fit**2 / 2):.2f}")
print(f"  Estimated geometric mean:  exp(mu) = {np.exp(mu_fit):.2f}")

# ── Motivating dataset ───────────────────────────────────────────────────
costs <- c(100, 1000, 10000)

# 1. Arithmetic mean
cat(sprintf("Arithmetic mean: %.2f\n", mean(costs)))          # 3700.00

# 2. Geometric mean via log-scale back-transform (natural log)
log_costs <- log(costs)
mu_hat <- mean(log_costs)
geo_mean <- exp(mu_hat)
cat(sprintf("Geometric mean (naive back-transform): %.2f\n", geo_mean))  # 1000.00

# 3. Duan smearing factor
residuals <- log_costs - mu_hat
smearing_factor <- mean(exp(residuals))
smeared_mean <- geo_mean * smearing_factor
cat(sprintf("Smearing factor: %.4f\n", smearing_factor))
cat(sprintf("Duan smeared estimate: %.2f\n", smeared_mean))   # 3700.00

# Utility: geometric mean function (handles positive values only)
geo_mean_fn <- function(x) exp(mean(log(x[x > 0])))

# ── Simulated 200-patient treatment comparison ────────────────────────────
set.seed(42)
n <- 200
treat <- rep(c(0L, 1L), each = n / 2)
log_y <- 6.0 + 0.405 * treat + rnorm(n, 0, 1.0)
y <- exp(log_y)

# 4. Log-OLS
fit <- lm(log_y ~ treat)
alpha_hat <- coef(fit)[["(Intercept)"]]
beta_hat  <- coef(fit)[["treat"]]
cat(sprintf("\nLog-OLS: alpha = %.3f, beta_treat = %.3f\n",
            alpha_hat, beta_hat))
cat(sprintf("exp(beta_treat) = %.3f  <-- geometric-mean ratio\n",
            exp(beta_hat)))

# 5. Pooled smearing factor
resids <- residuals(fit)
pooled_sf <- mean(exp(resids))
cat(sprintf("Pooled smearing factor: %.4f\n", pooled_sf))
for (arm in c(0, 1)) {
  mu_arm <- alpha_hat + beta_hat * arm
  arith_arm <- exp(mu_arm) * pooled_sf
  cat(sprintf("  Arm %d: geometric mean = %.2f, arithmetic mean = %.2f\n",
              arm, exp(mu_arm), arith_arm))
}

# 6. Group-specific smearing (correct under heteroscedasticity)
cat("\nGroup-specific smearing factors:\n")
for (arm in c(0, 1)) {
  mask <- treat == arm
  sf_arm <- mean(exp(resids[mask]))
  arith_arm <- exp(alpha_hat + beta_hat * arm) * sf_arm
  cat(sprintf("  Arm %d: smearing factor = %.4f, arithmetic mean = %.2f\n",
              arm, sf_arm, arith_arm))
}

# 7. Ratio CI from log-OLS (appropriate for geometric-mean ratio estimands)
ci <- confint(fit)["treat", ]
cat(sprintf("\nGeometric-mean ratio: %.3f (95%% CI: %.3f to %.3f)\n",
            exp(beta_hat), exp(ci[1]), exp(ci[2])))
cat("Note: this CI is for the geometric-mean ratio only.\n")
cat("For the arithmetic-mean ratio, bootstrap the smeared estimates per arm.\n")

# 8. Parametric log-normal fit via MASS::fitdistr
if (requireNamespace("MASS", quietly = TRUE)) {
  fit_ln <- MASS::fitdistr(costs, "lognormal")
  mu_fit    <- fit_ln$estimate["meanlog"]
  sigma_fit <- fit_ln$estimate["sdlog"]
  cat(sprintf("\nLog-normal MLE: mu = %.3f, sigma = %.3f\n", mu_fit, sigma_fit))
  cat(sprintf("  Arithmetic mean: exp(mu + sigma^2/2) = %.2f\n",
              exp(mu_fit + sigma_fit^2 / 2)))
  cat(sprintf("  Geometric mean:  exp(mu) = %.2f\n", exp(mu_fit)))
}