smmargins.Margins

class smmargins.Margins(results, data: DataFrame | None = None, level: float = 0.95, use_t: bool = False, analytic: bool = True, cov_type: str | None = None, vcov: ndarray | None = None, cov_kwds: dict | None = None, weights: ndarray | Sequence[float] | None = None, weight_type: str = 'sampling')

Compute adjusted predictions and marginal effects for a StatsModels fit.

Every statistic is built as a function of \(\beta\) using model.predict(params, exog) (which handles the inverse link). The default delta-method VCE applies J Var(\betâ) J^\top where J is taken analytically when the link supports it (family.link .inverse_deriv or the identity link of OLS/WLS/GLS) and by central finite differences otherwise. Krinsky–Robb and bootstrap VCEs are available via vce= on predict(), dydx(), and did().

Parameters:
  • results (statsmodels results object) – Must expose params, cov_params(), and have model.predict(params, exog) available (OLS, GLM, GLM-like, Logit, Probit, Poisson, NegBin, GEE, mixed, … all qualify).

  • data (pandas.DataFrame, optional) – Original fitting data. If None, we try model.data.frame.

  • level (float) – Confidence level for intervals (default 0.95).

  • use_t (bool) – If True and results.df_resid is available, use t-distribution.

  • analytic (bool) – If True (default), use an analytic outer Jacobian \(\partial g/\partial\beta\) whenever the model exposes a link derivative (any GLM via family.link.inverse_deriv, plus OLS/WLS/GLS via the identity link). Falls back to central finite differences otherwise. Set to False to force FD everywhere.

  • cov_type (str, optional) – Recompute the parameter covariance under this scheme (e.g. "HC0"..``”HC3”, ``"cluster", "HAC"). Passed through to results.get_robustcov_results (or equivalent). If neither cov_type nor vcov is set, results.cov_params() is used as-is — so a results object that was fit with a robust cov_type keeps its sandwich. Setting cov_type here overrides whatever is on results. The same kwargs are also accepted on predict(), dydx(), and did() for per-call overrides.

  • vcov (ndarray, optional) – User-supplied \((k, k)\) parameter covariance, used directly as Var(\betâ). Mutually exclusive with cov_type — passing both raises ValueError.

  • cov_kwds (dict, optional) – Extra keyword arguments forwarded to the covariance computation. Most importantly, {"groups": ...} for cov_type="cluster".

Notes

Formula vs. Raw Exog Mode

If the model was fit using a formula (e.g. smf.ols("y ~ x1 + x2", data)), Margins uses the model’s DesignInfo to rebuild the design matrix. This ensures that interactions (x1:x2) and transformations (log(x1)) are correctly updated when a single variable is perturbed for marginal effects.

If the model was fit using raw matrices (e.g. sm.OLS(y, X).fit()), Margins operates in “raw mode”. In this mode, only the literal column corresponding to the variable is perturbed. Interactions or pre-computed transformations in the design matrix will not be automatically updated. To correctly handle such models, it is recommended to fit using formulas.

Delta method

For a statistic \(g(\beta)\) with Jacobian \(G = \partial g / \partial \beta|_{\hat\beta}\), the delta method approximates

\[\widehat{\mathrm{Var}}[g(\hat\beta)] \approx G \, \widehat V \, G^\top .\]

Inference options

Each of predict(), dydx(), and did() accepts a common block of inference kwargs:

  • vce: "delta" (default), "simulation" (Krinsky–Robb, n_sims draws of \(\beta_s \sim N(\hat\beta, \hat V)\)), or "bootstrap" (refit on resampled rows; supports boot_method="pairs"|"cluster"|"block").

  • cov_type/vcov/cov_kwds: same meaning as on the constructor, but applied per-call and overriding what’s on self.

  • ci_method: "pointwise" (default), "bonferroni", "sidak", or "sup-t". The first three are post-hoc adjustments to the critical value; "sup-t" consumes the simulation/bootstrap draw matrix and so requires vce != "delta".

The point estimate reported is always the analytic \(g(\hat\beta)\) from the original fit, even under simulation or bootstrap — draws contribute only to SEs and intervals.

Examples

A logit with a categorical and an interaction, which would be awkward to differentiate by hand:

import statsmodels.formula.api as smf
from smmargins import Margins

fit = smf.logit(
    "voted ~ age + income + C(educ) + female + age:female",
    data=df,
).fit()
M = Margins(fit)

Adjusted predictions:

M.predict()                                # AAP
M.predict(at="mean")                       # APM
M.predict(atexog={"age": [25, 45, 65]})    # APR

Marginal effects on the response (probability) scale:

M.dydx("age")                              # AME
M.dydx("age", at="mean")                   # MEM
M.dydx("age", atexog={"female": [0, 1]})   # MER, by sex
M.dydx("educ", reference="college")        # discrete contrasts
M.dydx("kids", count=True)                 # x -> x+1 for integers
M.dydx("age", method="eyex")               # full elasticity

Robust SEs and alternative VCEs:

M.dydx("age", cov_type="HC3")              # heteroskedastic-robust
M.dydx("age", cov_type="cluster",
       cov_kwds={"groups": df["clust"]})    # cluster-robust
M.dydx("age", vce="simulation",
       n_sims=2000, sim_seed=0)            # Krinsky–Robb
M.dydx("age", vce="bootstrap",
       n_boot=1000, boot_seed=0)           # pairs bootstrap

Multiple-comparison adjustments for a family of margins:

M.dydx(["age", "income"], ci_method="bonferroni")
M.predict(atexog={"age": [25, 45, 65]},
          vce="simulation", n_sims=2000,
          ci_method="sup-t")

Most calls return a MarginsResult whose __repr__ prints a tidy table of estimates, SEs, z- (or t-) statistics, p-values, and confidence intervals.

A small runnable smoke test:

>>> import numpy as np, pandas as pd, statsmodels.formula.api as smf
>>> from smmargins import Margins
>>> rng = np.random.default_rng(0)
>>> df = pd.DataFrame({
...     "x": rng.standard_normal(200),
...     "g": rng.choice(["A", "B"], 200),
... })
>>> df["y"] = 1.0 + 2.0 * df["x"] + (df["g"] == "B") + rng.standard_normal(200)
>>> fit = smf.ols("y ~ x + C(g)", df).fit()
>>> M = Margins(fit)
>>> aap = M.predict()
>>> ame = M.dydx("x")
>>> aap.estimate.shape, ame.estimate.shape
((1,), (1,))
>>> bool(abs(ame.estimate[0] - 2.0) < 0.2)        # close to truth
True

References

[1] StataCorp. FAQ: How are the standard errors computed with

margins https://www.stata.com/support/faqs/statistics/compute-standard-errors-with-margins/

[2] Williams, R. (2012). Using the margins command to estimate and

interpret adjusted predictions and marginal effects. Stata Journal, 12(2), 308–331 https://www3.nd.edu/~rwilliam/stats/Margins01.pdf

[3] Ai, C., & Norton, E. C. (2003). Interaction terms in logit and

probit models. Economics Letters, 80(1), 123–129 https://doi.org/10.1016/S0165-1765(03)00032-6

See also

MarginsResult

__init__(results, data: DataFrame | None = None, level: float = 0.95, use_t: bool = False, analytic: bool = True, cov_type: str | None = None, vcov: ndarray | None = None, cov_kwds: dict | None = None, weights: ndarray | Sequence[float] | None = None, weight_type: str = 'sampling')

Methods

__init__(results[, data, level, use_t, ...])

contrast([a, b, a_newdata, b_newdata, at, ...])

Joint contrast between two counterfactual specifications.

did(group, condition[, group_levels, ...])

Difference-in-differences on the chosen scale.

dydx(variable[, at, atexog, discrete, ...])

Marginal effect of variable on the chosen scale.

predict([at, atexog, factor_stat, over, ...])

Compute adjusted predictions.

Attributes

n_outcomes

Number of outcome classes (K) for the fitted model.

outcome_labels

Outcome class labels for multi-outcome models, or None.