Fast Negative Binomial Regression with Variance Calculation (R Wrapper)

This function provides a fast implementation of negative binomial regression, wrapping a C++ backend (fast_neg_bin_with_var_cpp; see that page, and fast_dnbinom_mu_vec_cpp, for the full mean/dispersion negative-binomial model). Columns 1 and 2 of X (conventionally the intercept and treatment indicator) are always kept; this wrapper adds automatic, silent handling of rank-deficient or numerically unstable covariate sets beyond those first two columns: (1) upfront, any covariate columns (3 onward) found rank-deficient by qr are dropped before the first fit attempt; (2) if the C++ fit still fails (e.g. an L-BFGS line-search failure), covariates are dropped one at a time, in reverse QR-pivot order (most redundant first), retrying after each drop, until the fit succeeds or only the intercept and treatment columns remain — at which point, if it still fails, this function stop()s with an explicit error rather than returning a corrupted fit.

Usage

fast_negbin_regression_with_var(X, y, j = 2, optimization_alg = "lbfgs")

Arguments

X

A numeric matrix of predictor variables. It is assumed that an intercept column (e.g., a column of ones) is already included in X if desired.

y

A numeric vector of the response variable, representing count data.

j

The index of the coefficient to compute the variance for. Defaults to 2.

optimization_alg

Optimization algorithm: "lbfgs" (default) or "newton_raphson".

Value

A list containing the following components:

b

A numeric vector of the estimated negative binomial regression coefficients \(\hat\beta\), from whichever (possibly covariate-reduced) fit in the retry sequence ultimately succeeded.

ssq_b_j

The variance of the j-th estimated coefficient, computed by inverting the C++ backend's returned Hessian/Fisher-information matrix directly (via solve, not the backend's own vcov); NA if that inversion fails or j exceeds the number of columns remaining after covariate-dropping.

ssq_b_2

The variance of the second estimated coefficient specifically (typically the treatment effect), computed the same way, regardless of what j is.

See also

fast_neg_bin_with_var_cpp for the underlying backend (no automatic rank-deficiency retry) and full model documentation; fast_negbin_regression for the estimate-only, similarly retry-hardened wrapper.

Examples

X = matrix(rnorm(100), 10, 10)
y = rpois(10, 2)
fast_negbin_regression_with_var(X, y)
#> $b
#>  [1] -0.1074797  0.8617488 -1.1355726 -1.5249041  1.3885218 -0.9843362
#>  [7]  0.6125664 -1.3329996  2.2358959 -0.3133249
#> 
#> $ssq_b_j
#> [1] 56.39152
#> 
#> $ssq_b_2
#> [1] 56.39152
#>