Fast G-Computation (Standardization) Point Estimate and Model-Based Inference for a Proportional-Odds Ordinal Model (C++)

Computes the same standardized (G-computation) marginal difference in expected ordinal category score as gcomp_ordinal_proportional_odds_post_fit_cpp — under a fitted proportional-odds (cumulative logit) model \(\mathrm{logit}\,\Pr(Y_i \le k \mid x_i) = \alpha_k - x_i^\top\beta\), \(k = 1, \ldots, K-1\) — but additionally supplies model-based inferential quantities that the other function omits. The fitted linear predictor for each subject is recomputed with the treatment column (j_treat) forced to 1 (\(\eta_{1,i}\)) and to 0 (\(\eta_{0,i}\)), the standardized expected scores mean1, mean0, and their difference md are formed exactly as in gcomp_ordinal_proportional_odds_post_fit_cpp, and then:

• An ordinal-regression log-likelihood Hessian is evaluated at \([\hat\alpha, \hat\beta]\) and inverted (with a symmetrization step and a finiteness check) to give the model-based (non-sandwich) variance-covariance matrix of the full parameter vector; vcov/std_err/z_vals report the \(\hat\beta\)-block of this covariance.

• The delta-method standard error of md, se_md, is obtained by numerically differentiating md (central differences, step \(10^{-5}\)) with respect to every element of \([\alpha,\beta]\) to get a gradient \(g\), then computing \(\sqrt{g^\top V g}\) where \(V\) is the full parameter covariance above; NA if any perturbed md evaluation is non-finite (e.g. a perturbed \(\alpha\) violates monotonicity) or the resulting variance is negative/non-finite.

Unlike the sandwich-based post-fit helpers elsewhere in this package (e.g. gcomp_logistic_post_fit_cpp), the covariance here comes purely from the model's observed information and does not attempt to be robust to misspecification.

Usage

ordinal_gcomp_post_fit_cpp(X_fit, y, coef_hat, alpha_hat, j_treat)

Arguments

X_fit

Numeric matrix of predictors used to fit the model.

y

Numeric vector of the ordinal responses used to fit the model (needed to reconstruct the OrdinalRegression likelihood for the Hessian; only categorization, not numeric coding, matters).

coef_hat

Numeric vector of fitted proportional-odds regression coefficients \(\hat\beta\), same length and column order as X_fit.

alpha_hat

Numeric vector of the \(K-1\) fitted, increasing category thresholds \(\hat\alpha_1, \ldots, \hat\alpha_{K-1}\).

j_treat

1-based column index of the treatment indicator in X_fit.

Value

A list with elements vcov (model-based covariance of \(\hat\beta\)), std_err, z_vals (both length p, one per column of X_fit), mean1, mean0, md (mean1 - mean0), and se_md (delta-method SE of md). Errors (via stop()) if j_treat is out of bounds, the dimensions of y/coef_hat are inconsistent with X_fit, the Hessian is not invertible, or the resulting covariance has any non-finite entry.

See also

gcomp_ordinal_proportional_odds_post_fit_cpp for the point-estimate-only variant (no y, no Hessian, no inference) that this function's point estimates match; fast_ordinal_regression_cpp for the fitting routine that produces coef_hat/alpha_hat.