Generalized Propensity Scores for Continuous Exposures

A generalized propensity score (GPS) for a continuous exposure is the conditional density of the observed treatment value A given pre-treatment covariates X. It is not the multinomial generalized propensity score used for three or more categorical treatments. It also is not a single adjusted spline coefficient. The causal target must be stated first: for example, the population mean outcome if every eligible person received a specified baseline exposure quantity, compared with another feasible quantity, over a defined baseline exposure window.

Under consistency, conditional exchangeability for the full continuous treatment, and positivity of the relevant exposure density, the GPS can balance measured X at a specified dose. A common implementation models A|X, estimates r(A,X), models E(Y|A,r), and averages predictions over the eligible cohort on a prespecified grid of supported doses. The output is a marginal dose-response function with uncertainty bands. It does not turn an observed dose-response association into a causal effect when severity, contraindications, preferences, or care access remain unmeasured.

Support is the operational constraint

A density estimate exists mathematically even where the data contain almost no comparable patients. Analysts must show the dose distribution overall and by clinically important X, inspect conditional-density tails, and restrict the dose grid to a region with credible support. Restriction changes the target population or intervention domain and must be reported. Do not quietly extrapolate a spline beyond the observed dose range.

When not to use

Do not use this point-treatment card for dose changes driven by interim response, toxicity, or prior treatment: that is treatment-confounder feedback and needs a longitudinal g-method or a clearly defined sustained strategy. Do not use it to claim an effect of impossible dose assignments (for example, a dose forbidden by renal function) or to repair structural nonpositivity. Modified treatment policies and stochastic interventions can define more feasible continuous-treatment questions, but are outside this card's primary GPS workflow and need their own protocol and identification analysis.

Trade-offs

Compared with a descriptive spline, GPS analysis makes an intervention, confounder set, and support requirement explicit. It requires a correctly specified or suitably flexible exposure-density model and more diagnostics. Compared with categorizing dose into low/medium/high groups, it preserves dose information but can become unstable in sparse tails. A binary or multinomial propensity-score approach remains appropriate only after a clinically justified categorical treatment definition.

Data-source operational depth

• Claims: strength and quantity dispensed are proxies for intended supply, not ingestion. Convert them to a prespecified baseline dose using a documented denominator; exclude or separately handle overlapping fills, stockpiling, and inpatient supply that claims do not observe.
• EHR: prescribed, dispensed, administered, and measured concentration are different exposure constructs. Avoid using a dose recorded after response assessment; missed monitoring can make both dose and outcome capture informative.
• Registry: registry dose fields can be clinically granular, but referral severity and protocol-driven titration can create confounding that a sparse covariate set cannot remove.
• Linked data: reconcile order, dispensing, administration, lab, and outcome dates before selecting the baseline exposure window. The linked subset can itself be selected relative to the intended target population.

Decision diagram

%%{init: {'themeVariables': {'fontSize': '16px'}, 'flowchart': {'nodeSpacing': 18, 'rankSpacing': 24}}}%%
flowchart TD
  A["Define feasible baseline<br/>continuous-dose intervention"] --> B["Measure exposure A<br/>and pre-exposure covariates X"]
  B --> C["Estimate conditional density<br/>GPS = f(A given X)"]
  C --> D["Assess support and balance<br/>Restrict grid or population<br/>if needed; report the change"]
  D --> E["Fit outcome surface<br/>Predict at each retained dose<br/>Average across the cohort"]
  E --> F["Dose-response curve<br/>with uncertainty<br/>and stated assumptions"]

Worked example

Scenario

A synthetic linked cohort compares 180-day biomarker change under a continuous baseline exposure quantity, using only exposure values inside a prespecified region of conditional support. Baseline severity, renal function, age, and prior therapy are pre-index covariates.

Dataset

Synthetic baseline rows; exposure is fixed at baseline and outcome change is assessed after follow-up.

Steps

Result

The result is an adjusted 180-day exposure-response function only over the prespecified grid that passes a documented support decision. It is not evidence about exposure changes after follow-up begins.

Trade-offs

Code example

import numpy as np
from scipy.stats import norm
rng = np.random.default_rng(13); n, B = 500, 200
x = rng.normal(size=(n, 2)); mu = 0.35 + .35*x[:, 0] - .20*x[:, 1]
a = np.exp(mu + rng.normal(0, .35, n))
y = 1 + .45*a - .12*a*a + .45*x[:, 0] - .20*x[:, 1] + rng.normal(0, .6, n)
grid = np.quantile(a, [.20, .50, .80])
def curve(a, x, y, grid):
    z = np.column_stack([np.ones(len(a)), x])
    b = np.linalg.lstsq(z, np.log(a), rcond=None)[0]; m = z @ b
    s = np.sqrt(np.sum((np.log(a)-m)**2)/(len(a)-z.shape[1]))
    r = norm.pdf(np.log(a), m, s)/a
    h = np.column_stack([np.ones(len(a)), a, a*a, r, r*r, a*r])
    bh = np.linalg.lstsq(h, y, rcond=None)[0]
    q = []
    for g in grid:
        rg = norm.pdf(np.log(g), m, s)/g
        hg = np.column_stack([np.ones(len(a)), np.repeat(g,len(a)), np.repeat(g*g,len(a)), rg, rg*rg, g*rg])
        q.append([g, np.mean(hg @ bh), np.mean(np.abs((np.log(g)-m)/s)>2)])
    return np.asarray(q)
point = curve(a,x,y,grid)
boot = np.vstack([curve(a[ix],x[ix],y[ix],grid)[:,1] for ix in (rng.integers(0,n,n) for _ in range(B))])
for j, g in enumerate(grid):
    print({"exposure":float(g), "mean_outcome":float(point[j,1]), "bootstrap_95_ci":tuple(np.quantile(boot[:,j],[.025,.975])), "conditional_tail_fraction_abs_z_gt_2":float(point[j,2])})

set.seed(13); n<-500; B<-200; x1<-rnorm(n); x2<-rnorm(n)
a<-exp(.35+.35*x1-.20*x2+rnorm(n,0,.35)); y<-1+.45*a-.12*a^2+.45*x1-.2*x2+rnorm(n,0,.6)
grid<-as.numeric(quantile(a,c(.2,.5,.8)))
f <- function(a, x1, x2, y, g) {
  e <- lm(log(a) ~ x1 + x2)
  m <- fitted(e); s <- summary(e)$sigma
  r <- dlnorm(a, m, s)
  o <- lm(y ~ a + I(a^2) + r + I(r^2) + I(a*r))
  t(sapply(g, function(z) {
    rz <- dlnorm(z, m, s)
    c(mean(predict(o, data.frame(a=z, r=rz))),
      mean(abs((log(z)-m)/s) > 2))
  }))
}
p<-f(a,x1,x2,y,grid); b<-replicate(B,{i<-sample.int(n,n,TRUE);f(a[i],x1[i],x2[i],y[i],grid)[,1]})
print(data.frame(exposure=grid,mean_outcome=p[,1],lcl=apply(b,1,quantile,.025),ucl=apply(b,1,quantile,.975),conditional_tail_fraction_abs_z_gt_2=p[,2]))

data cohort_complete;
  set work.cohort;
  if nmiss(exposure, baseline_severity, egfr, age, outcome)=0 and exposure>0;
  log_exposure=log(exposure);
run;
proc glm data=cohort_complete noprint;
  model log_exposure = baseline_severity egfr age;
  output out=gps_base p=mu r=resid;
run; quit;
proc sql noprint; select sqrt(sum(resid*resid)/(count(resid)-4)) into :sigma from gps_base; quit;
data gps_base; set gps_base;
  gps_obs=pdf('normal', log_exposure, mu, &sigma)/exposure; /* 1/A Jacobian */
run;
proc glm data=gps_base; model outcome = exposure exposure*exposure gps_obs gps_obs*gps_obs exposure*gps_obs; run; quit;

Citations

1. [1]Hirano K, Imbens GW. The propensity score with continuous treatments. In: Applied Bayesian Modeling and Causal Inference from Incomplete-Data Perspectives. 2004;73-84.
2. [2]Imai K, van Dyk DA. Causal inference with general treatment regimes: Generalizing the propensity score. Journal of the American Statistical Association. 2004;99(467):854-866.