`statspai.dose_response`¶

dose_response ¶

Continuous Treatment Effects: Dose-Response Estimation.

Estimates the dose-response function E[Y(t)] for continuous treatments using the generalized propensity score (GPS) approach and kernel-based methods.

References

Hirano, K. & Imbens, G. W. (2004). The Propensity Score with Continuous Treatments. Applied Bayesian Modeling and Causal Inference from Incomplete-Data Perspectives, 226164, 73-84. [@hirano2004propensity]

Kennedy, E. H., Ma, Z., McHugh, M. D., & Small, D. S. (2017). Non-parametric methods for doubly robust estimation of continuous treatment effects. JRSS-B, 79(4), 1229-1245. [@kennedy2017parametric]

DoseResponse ¶

Generalized Propensity Score dose-response estimator.

Parameters:

Name	Type	Default
`data`	`DataFrame`	required
`y`	`str`	required
`treat`	`str`	required
`covariates`	`list of str`	required
`n_dose_points`	`int`	`20`
`dose_range`	`tuple`	`None`
`treatment_model`	`sklearn estimator`	`None`
`outcome_model`	`sklearn estimator`	`None`
`n_bootstrap`	`int`	`200`
`alpha`	`float`	`0.05`
`random_state`	`int`	`42`

Examples:

>>> import numpy as np
>>> import pandas as pd
>>> import statspai as sp
>>> rng = np.random.default_rng(0)
>>> n = 150
>>> age = rng.normal(50, 10, n)
>>> weight = rng.normal(70, 12, n)
>>> dosage = 0.5 * age + 0.3 * weight + rng.normal(0, 5, n)
>>> outcome = 0.8 * dosage + 0.2 * age + rng.normal(0, 3, n)
>>> df = pd.DataFrame({"outcome": outcome, "dosage": dosage,
...                    "age": age, "weight": weight})
>>> est = sp.DoseResponse(df, y="outcome", treat="dosage",
...                       covariates=["age", "weight"],
...                       n_dose_points=6, n_bootstrap=8, random_state=0)
>>> res = est.fit()
>>> isinstance(res, sp.CausalResult)
True
>>> list(res.detail.columns)
['dose', 'response', 'se', 'ci_lower', 'ci_upper', 'marginal_effect']

References

[@hirano2004propensity], [@kennedy2017parametric]

fit ¶

fit() -> CausalResult

Estimate the dose-response function.

VCNetResult `dataclass` ¶

Bases: ResultProtocolMixin

Dose-response curve returned by :func:vcnet / :func:scigan.

Attributes:

Name	Type	Description
`t_grid`	`ndarray`	Treatment grid at which the curve is evaluated.
`mu_hat`	`ndarray`	Estimated dose-response :math:`\hat\mu(t)` on `t_grid`.
`se`	`ndarray`	Bootstrap pointwise standard errors.
`ci_lo, ci_hi`	`ndarray`	Pointwise lower / upper confidence band.
`coef_matrix`	`ndarray`	Varying-coefficient matrix of shape `(n_basis, p + 1)`.
`n_obs, n_basis`	`int`
`detail`	`dict`

Examples:

>>> import statspai as sp
>>> import numpy as np
>>> import pandas as pd
>>> df = pd.DataFrame({
...     "y": np.random.randn(200),
...     "dose": np.random.rand(200),
...     "x1": np.random.randn(200),
... })
>>> res = sp.vcnet(df, y="y", treatment="dose",
...                covariates=["x1"])
>>> res.mu_hat[:3]            # dose-response on the t-grid
>>> print(res.summary())     # tidy curve table

dose_response ¶

dose_response(data: DataFrame, y: str, treat: str, covariates: List[str], n_dose_points: int = 20, dose_range: Optional[Tuple[float, float]] = None, treatment_model: Optional[BaseEstimator] = None, outcome_model: Optional[BaseEstimator] = None, n_bootstrap: int = 200, alpha: float = 0.05, random_state: int = 42, method: str = 'gps', bandwidth: Optional[float] = None, n_folds: int = 5) -> CausalResult

Estimate the dose-response function for a continuous treatment.

Parameters:

Name	Type	Description	Default
`data`	`DataFrame`	Input data.	required
`y`	`str`	Outcome variable.	required
`treat`	`str`	Continuous treatment variable.	required
`covariates`	`list of str`	Covariate names.	required
`n_dose_points`	`int`	Number of dose levels to evaluate.	`20`
`dose_range`	`tuple of (float, float)`	(min, max) dose range. If None, uses 5th-95th percentile.	`None`
`treatment_model`	`sklearn estimator`	Model for E[T\|X]. Default: GBM.	`None`
`outcome_model`	`sklearn estimator`	Model for E[Y\|T, GPS]. Default: GBM. The outcome model sees the two columns `[T, GPS]`, so `make_pipeline(PolynomialFeatures(2), LinearRegression())` with `treatment_model=LinearRegression()` is the quadratic specification of Hirano and Imbens (2004). Neither choice is reliable in general. On a linear confounded design (`T = 1 + 0.8x + 0.3x2 + e`, `Y = 0.8T + x + 0.5x2 + u`, n = 1000, 30 draws) the default GBM pair is nearly unbiased for the 25th-75th percentile effect (bias -0.13 against a truth of about 1.9) but noisy (SD 0.43) and its curve was never monotone, while the Hirano-Imbens quadratic is precise (SD 0.07) and badly biased (+1.27, close to the unadjusted slope), because `E[Y \| T, GPS]` is not quadratic there. The default is therefore kept, and a call with both default models warns. For an average slope with valid inference use `sp.dml(model='plr')`.	`None`
`n_bootstrap`	`int`	Bootstrap iterations for confidence bands.	`200`
`alpha`	`float`	Significance level.	`0.05`
`random_state`	`int`		`42`
`method`	`(gps, dr)`	`'gps'`: the generalized-propensity-score regression of Hirano and Imbens (2004) with bootstrap bands (the default, unchanged). `'dr'`: the doubly robust curve of Kennedy, Ma, McHugh and Small (2017): a local linear regression of the pseudo-outcome `(Y - mu(X, T)) / pi(T \| X) * mean_j pi(T \| X_j) + mean_j mu(X_j, T)` on `T`, consistent if either the outcome regression `mu` (`outcome_model`, on `[X, T]`; default gradient boosting) or the normal treatment density `pi` (`treatment_model` for its mean, default linear, and a linear regression of squared residuals for its variance) is right. Nuisances are cross-fitted (`n_folds`); the bandwidth minimises the smoother's leave-one-out risk unless `bandwidth` is given; bands are the robust local-linear sandwich on the pseudo-outcomes (no bootstrap), for the smoothed curve. R `npcausal::ctseff` computes the same estimator with other nuisance and smoothing choices; the curves agree in shape, not to parity.	`'gps'`
`bandwidth`	`float`	`method='dr'`: kernel bandwidth in units of the treatment.	`None`
`n_folds`	`int`	`method='dr'`: cross-fitting folds (1 = no cross-fitting).	`5`

Returns:

Type	Description
`CausalResult`	detail contains the dose-response curve: columns 'dose', 'response', 'se', 'ci_lower', 'ci_upper'. model_info contains 'avg_marginal_effect' (derivative).

Examples:

>>> import statspai as sp, numpy as np, pandas as pd
>>> rng = np.random.default_rng(0)
>>> n = 400
>>> age = rng.normal(40, 10, n)
>>> weight = rng.normal(70, 12, n)
>>> dosage = 0.5 * age + rng.normal(0, 5, n)      # confounded dose
>>> outcome = 0.2 * dosage + 0.1 * weight + rng.normal(0, 1, n)
>>> df = pd.DataFrame({"outcome": outcome, "dosage": dosage,
...                    "age": age, "weight": weight})
>>> result = sp.dose_response(df, y="outcome", treat="dosage",
...                           covariates=["age", "weight"],
...                           n_dose_points=10, n_bootstrap=20)
>>> curve = result.detail            # dose-response curve
>>> curve.shape
(10, 6)
>>> list(curve.columns)
['dose', 'response', 'se', 'ci_lower', 'ci_upper', 'marginal_effect']

scigan ¶

scigan(data: DataFrame, y: str, treatment: str, covariates: Sequence[str], t_grid: Optional[Sequence[float]] = None, propensity_weights: Optional[ndarray] = None, **kwargs: Any) -> VCNetResult

Adversarial dose-response estimator (Bica et al. 2020).

The reference SCIGAN architecture is a GAN — too heavy to ship as a dependency-free algorithm. This entry point exposes a propensity-weighted VCNet fit that captures the essential "balancing" behaviour of SCIGAN: residual imbalance along the treatment axis is re-weighted using a user-provided propensity_weights (g(T | X)^{-1}), which enter the fit as observation weights. With unit weights the result is :func:vcnet's. For the full SCIGAN training loop, plug in your own GAN-generated counterfactual samples and pass the re-weighting through propensity_weights.

Examples:

>>> import statspai as sp
>>> import numpy as np
>>> import pandas as pd
>>> df = pd.DataFrame({
...     "y": np.random.randn(200),
...     "dose": np.random.rand(200),
...     "x1": np.random.randn(200),
... })
>>> w = np.ones(len(df))  # inverse-propensity balancing weights
>>> res = sp.scigan(df, y="y", treatment="dose",
...                 covariates=["x1"],
...                 propensity_weights=w)
>>> res.mu_hat[:3]  # balanced dose-response curve

statspai.dose_response¶

dose_response ¶

DoseResponse ¶

fit ¶

VCNetResult dataclass ¶

dose_response ¶

scigan ¶

`statspai.dose_response`¶

VCNetResult `dataclass` ¶