Extract variance-covariance matrix from a choicer_fit object
Source:R/methods.R
vcov.choicer_fit.RdWith no arguments, returns the variance-covariance matrix implied by the
fit's own se_method (triggering lazy computation if needed). Passing
type recomputes a different variance estimator post hoc from the
stored data — no refit needed (requires keep_data = TRUE):
"hessian"Inverse of the analytical negated Hessian.
"bhhh"Inverse of the BHHH/OPG information \(\sum_i w_i s_i s_i'\).
"robust"Huber-White sandwich \(A^{-1} (\sum_i w_i^2 s_i s_i') A^{-1}\) — also the valid WESML variance under choice-based weighting.
"cluster"Cluster-robust sandwich \(A^{-1} (\sum_g g_g g_g') A^{-1}\) with \(g_g = \sum_{i \in g} w_i s_i\) the within-cluster sum of weighted scores. Requires
cluster(or a fit made withcluster_col). No small-sample correction is applied.
Here \(i\) indexes choice situations. For repeated choices by the same decision maker (panel data), cluster on the decision maker.
Usage
# S3 method for class 'choicer_fit'
vcov(object, type = NULL, cluster = NULL, ...)Arguments
- object
A choicer_fit object.
- type
NULL(default; return the as-fitted vcov) or one of"hessian","bhhh","robust","cluster".- cluster
Cluster labels for
type = "cluster", one per choice situation. Alignment to the prepared (id-sorted) choice situations is handled as follows:Named (recommended): names are matched against the choice-situation ids, so the vector is safe in any order. Build it by naming your per-situation labels with the id values.
Unnamed: taken to be in the prepared, id-sorted order; a warning flags that assumption. A vector of per-alternative (row-level) length is rejected.
Defaults to the labels stored at fit time via
cluster_col(already aligned). Supplyingclusterwithouttypeimpliestype = "cluster". The safest route is to passcluster_col=at fit time, which sidesteps post-hoc alignment entirely.- ...
Additional arguments (ignored).
Details
Note (mixed logit): clustering repairs the inference, not the
estimand. run_mxlogit() treats each choice situation as an
independent draw from the mixing distribution (a cross-sectional MSL
likelihood, not the panel product form), so on panel data the point
estimates target that cross-sectional model; type = "cluster" makes
their standard errors robust to within-person dependence but does not turn
the fit into a panel mixed logit. For panel random coefficients use
run_hmnlogit (person_col).
Examples
# \donttest{
library(data.table)
set.seed(42)
N <- 50; J <- 3
dt <- data.table(id = rep(1:N, each = J), alt = rep(1:J, N))
dt[, `:=`(x1 = rnorm(.N), x2 = rnorm(.N))]
#> id alt x1 x2
#> <int> <int> <num> <num>
#> 1: 1 1 1.3709584 -0.04069848
#> 2: 1 2 -0.5646982 -1.55154482
#> 3: 1 3 0.3631284 1.16716955
#> 4: 2 1 0.6328626 -0.27364570
#> 5: 2 2 0.4042683 -0.46784532
#> ---
#> 146: 49 2 1.1133860 -0.47733551
#> 147: 49 3 -0.4809928 -0.16626149
#> 148: 50 1 -0.4331690 0.86256338
#> 149: 50 2 0.6968626 0.09734049
#> 150: 50 3 -1.0563684 -1.62561674
dt[, person := rep(1:10, each = 5)[id]]
#> id alt x1 x2 person
#> <int> <int> <num> <num> <int>
#> 1: 1 1 1.3709584 -0.04069848 1
#> 2: 1 2 -0.5646982 -1.55154482 1
#> 3: 1 3 0.3631284 1.16716955 1
#> 4: 2 1 0.6328626 -0.27364570 1
#> 5: 2 2 0.4042683 -0.46784532 1
#> ---
#> 146: 49 2 1.1133860 -0.47733551 10
#> 147: 49 3 -0.4809928 -0.16626149 10
#> 148: 50 1 -0.4331690 0.86256338 10
#> 149: 50 2 0.6968626 0.09734049 10
#> 150: 50 3 -1.0563684 -1.62561674 10
dt[, choice := 0L]
#> id alt x1 x2 person choice
#> <int> <int> <num> <num> <int> <int>
#> 1: 1 1 1.3709584 -0.04069848 1 0
#> 2: 1 2 -0.5646982 -1.55154482 1 0
#> 3: 1 3 0.3631284 1.16716955 1 0
#> 4: 2 1 0.6328626 -0.27364570 1 0
#> 5: 2 2 0.4042683 -0.46784532 1 0
#> ---
#> 146: 49 2 1.1133860 -0.47733551 10 0
#> 147: 49 3 -0.4809928 -0.16626149 10 0
#> 148: 50 1 -0.4331690 0.86256338 10 0
#> 149: 50 2 0.6968626 0.09734049 10 0
#> 150: 50 3 -1.0563684 -1.62561674 10 0
dt[, choice := sample(c(1L, rep(0L, J - 1))), by = id]
#> id alt x1 x2 person choice
#> <int> <int> <num> <num> <int> <int>
#> 1: 1 1 1.3709584 -0.04069848 1 0
#> 2: 1 2 -0.5646982 -1.55154482 1 0
#> 3: 1 3 0.3631284 1.16716955 1 1
#> 4: 2 1 0.6328626 -0.27364570 1 0
#> 5: 2 2 0.4042683 -0.46784532 1 0
#> ---
#> 146: 49 2 1.1133860 -0.47733551 10 0
#> 147: 49 3 -0.4809928 -0.16626149 10 0
#> 148: 50 1 -0.4331690 0.86256338 10 1
#> 149: 50 2 0.6968626 0.09734049 10 0
#> 150: 50 3 -1.0563684 -1.62561674 10 0
fit <- run_mnlogit(dt, "id", "alt", "choice", c("x1", "x2"))
#> Optimization run time 0h:0m:0s
vcov(fit) # as fitted (hessian)
#> x1 x2 ASC_2 ASC_3
#> x1 0.033389428 -0.006452059 -0.003138235 -0.007197134
#> x2 -0.006452059 0.030929336 0.003648743 0.004450948
#> ASC_2 -0.003138235 0.003648743 0.110170532 0.061185377
#> ASC_3 -0.007197134 0.004450948 0.061185377 0.147042595
vcov(fit, type = "robust") # Huber-White, post hoc
#> x1 x2 ASC_2 ASC_3
#> x1 0.0320271769 -0.01290123 -0.0003436706 -0.007503077
#> x2 -0.0129012287 0.03969401 0.0128385440 0.015345820
#> ASC_2 -0.0003436706 0.01283854 0.1116674307 0.056052348
#> ASC_3 -0.0075030769 0.01534582 0.0560523478 0.141521824
# named by situation id -> safe regardless of order
cl <- dt[, person[1L], by = id]
vcov(fit, type = "cluster", cluster = setNames(cl$V1, cl$id))
#> x1 x2 ASC_2 ASC_3
#> x1 0.03443605 -0.01533769 -0.02450414 -0.03281563
#> x2 -0.01533769 0.03462553 0.03170608 0.02136278
#> ASC_2 -0.02450414 0.03170608 0.10860263 0.03377958
#> ASC_3 -0.03281563 0.02136278 0.03377958 0.08791317
# }