estimate.default: Influence function based inference
In lava: Latent Variable Models

estimate.default

R Documentation

Influence function based inference

Description

Primary tool for obtaining parameter estimates with robust (sandwich) standard errors, applying the delta method, and testing linear hypotheses. The function returns an object of class estimate which serves as a general container for parameter estimates and their influence functions (IFs). Three calling conventions are supported:

Usage

## Default S3 method:
estimate(
  x = NULL,
  f = NULL,
  ...,
  data,
  id,
  coef,
  IC = TRUE,
  vcov,
  stack = TRUE,
  average = FALSE,
  subset,
  keep,
  use,
  regex = FALSE,
  ignore.case = FALSE,
  print = NULL,
  labels,
  label.width,
  center = FALSE
)

Arguments

`x`	model object (`glm`, `lvmfit`, ...) or an existing `estimate` object. When two model objects are supplied (e.g., `estimate(g, g0)`) a likelihood-ratio test is performed.
`f`	transformation of model parameters. Accepts several input types: A function `f(p)` or `f(p, data)`: applies the delta method. When `f` returns a named list the names are used as parameter labels. A matrix: used as a contrast (linear combination) matrix. A numeric vector of parameter indices: converted to a contrast that selects and differences those parameters. A list of indices: each element selects one parameter. Character expressions: supports wildcards (`"?"`, `""`) and arithmetic on parameter names (e.g., `"z" - "x"`, `2 "z" - 3 * "x"`).
`...`	additional arguments to lower level functions
`data`	an optional data frame, list or environment (or object coercible by `as.data.frame` to a data frame) containing the variables in the model. If not found in `data`, the variables are taken from `environment(formula)`, typically the environment from which `glm` is called.
`id`	(optional) cluster identifier. Can be a vector of cluster IDs, a one-sided formula (evaluated in `data`), a single character column name. When supplied, the IF is aggregated within clusters to produce cluster-robust standard errors. When `average = TRUE`, `id` refers to the rows of `data`, i.e., the values of `f` (see expectation). `id = NULL` removes the id (index) and the row names of the influence function from the returned object.
`coef`	(optional) named parameter vector. Used instead of `coef(x)` when constructing an `estimate` object without a model.
`IC`	if `TRUE` (default) the influence function matrix is estimated and stored in the returned object (extract with the IC method). Can also be a user-supplied IF matrix (one row per observation, one column per parameter), which is used directly instead of estimating it from `x`.
`vcov`	(optional) covariance matrix of parameter estimates, or a logical. If `TRUE`, stats::vcov is used to obtain the (model-based) covariance matrix from `x`, yielding non-robust standard errors. If a matrix is supplied it is used directly. When omitted or `FALSE`, robust standard errors are computed from the influence function.
`stack`	if `TRUE` (default) the influence function contributions are summed within each cluster defined by `id`. Set to `FALSE` to keep the un-stacked (per-observation) decomposition.
`average`	(Deprecated) if `TRUE` the function computes the standardized estimate `\hat\Psi = P_n f(X; \hat\theta)`, i.e., the empirical mean of `f(p, data)`, as defined by the `f` argument, over all rows of `data`. The influence function accounts for both the empirical averaging and the parameter estimation uncertainty (see expectation).
`subset`	an optional vector specifying a subset of observations to be used in the fitting process. (See additional details about how this argument interacts with data-dependent bases in the ‘Details’ below.)
`keep`	(optional) index of parameters to keep from final result. Accepts integer indices, character names, or (with `regex = TRUE`) perl-compatible regular expressions.
`use`	(optional) index of parameters to use in calculations. The selected parameters are first extracted (via `keep`) and then the remaining arguments (`f`, `contrast`, etc.) are applied to this subset.
`regex`	if `TRUE` use perl-compatible regular expressions for `keep` and `use` arguments
`ignore.case`	ignore case in regular expressions
`print`	(optional) custom print function for the resulting `estimate` object
`labels`	(optional) character vector of coefficient names
`label.width`	(optional) max display width of labels
`center`	if TRUE the influence function will automatically be centered (zero mean)

Details

estimate(x, ...) – extract estimates from a model object
estimate(coef=, IC=, ...) – construct from coefficients and IF matrix
estimate(coef=, vcov=, ...) – construct from coefficients and covariance matrix

Value

Object of class estimate with the following elements:

`coef`	Named vector of parameter estimates.
`vcov`	Variance-covariance matrix.
`IC`	Influence function matrix (observations x parameters).
`coefmat`	Formatted coefficient table (estimate, std.err, confidence limits, p-value).
`id`	Cluster/id variable used.
`ncluster`	Number of clusters.
`n`	Number of observations.
`compare`	(When `null` or contrasts are specified) Wald test result.

Influence functions and robust standard errors

An estimator \widehat{\theta} is regular and asymptotically linear (RAL) when it admits the iid decomposition

\sqrt{n}(\widehat{\theta}-\theta) = \frac{1}{\sqrt{n}}\sum_{i=1}^n \mathrm{IC}(Z_i; P) + o_p(1)

where \mathrm{IC} is the unique influence function satisfying E\{\mathrm{IC}(Z; P)\} = 0. By the central limit theorem

\sqrt{n}(\widehat{\theta}-\theta) \overset{d}{\longrightarrow} N(0,\; \mathrm{Var}\{\mathrm{IC}(Z; P)\})

and the asymptotic variance is consistently estimated by the empirical variance of the plugin IF estimate, yielding robust (sandwich) standard errors. The estimated IF can be extracted with the IC method.

Parameter transformations (delta method)

When f is a function \phi: R^p \to R^m, the delta method is applied:

\sqrt{n}\{\phi(\widehat{\theta}) - \phi(\theta)\} = \frac{1}{\sqrt{n}}\sum_{i=1}^n \nabla\phi(\theta)\,\mathrm{IC}(Z_i; P) + o_p(1)

Derivatives are computed numerically via numDeriv::jacobian unless the function returns an attribute "grad" with the analytic Jacobian.

Alternatively, estimate objects support direct arithmetic operations (e.g., a * b, exp(a), a^b) which apply the delta method with exact (analytical) derivatives computed automatically. This influence function calculus allows building complex transformations from simple building blocks without numerical differentiation. See the last example section ("influence function calculus") and vignette("influencefunction", package = "lava") for details.

Standardization (g-computation)

When average = TRUE and f(p, data) depends on covariates, the target parameter is the standardized (g-computation) estimate \Psi = E\{f(X;\theta)\}. The IF for the averaged estimate accounts for both the empirical averaging and parameter estimation uncertainty:

\mathrm{IC}_\Psi(Z; P) = f(X;\theta) - \Psi + [E\nabla_\theta f(X;\theta)]\,\phi(Z; P)

When subset is also specified, the average is conditioned on the subpopulation, yielding a conditional marginalized estimate. The model may be estimated on a different (e.g., smaller) dataset than data. See expectation for details.

Cluster-robust standard errors

When id is supplied, the per-observation IF contributions are summed within clusters (when stack = TRUE), producing the cluster-level IF \widetilde{\mathrm{IC}}(Z_i; P) = \sum_{k=1}^{N_i} \frac{n}{N}\mathrm{IC}(Z_{ik}; P). The resulting variance estimate is equivalent to the GEE working independence sandwich estimator.

For full theoretical background and worked examples see vignette("influencefunction", package = "lava").

Examples


## Simulation from logistic regression model
m <- lvm(y~x+z);
distribution(m,y~x) <- dist_bernoulli("logit")
d <- sim(m,1000)
g <- glm(y~z+x,data=d,family=binomial())
g0 <- glm(y~1,data=d,family=binomial())

## LRT
estimate(g, g0)


estimate(g)

## Testing contrasts
summary(estimate(g), null=0)
estimate(g, rbind(c(1,1,0), c(1,0,2)))
summary(estimate(g, rbind(c(1,1,0), c(1,0,2))), null=c(1,2))
estimate(g, 2:3) ## same as cbind(0,1,-1)
estimate(g, as.list(2:3)) ## same as rbind(c(0,1,0),c(0,0,1))
## Alternative syntax
estimate(g, "z", "z"-"x", 2*"z"-3*"x")
estimate(g, "?")  ## Wildcards
estimate(g, "*Int*", "z")
summary(estimate(g, "1", "2"-"3"), null = c(0,1))
estimate(g, 2, 3)

## Usual (non-robust) confidence intervals
estimate(g, vcov=TRUE)
estimate(g, vcov=vcov(g))

## Transformations
estimate(g, function(p) p[1]+p[2])

## Multiple parameters
e <- estimate(g, function(p) c(p[1]+p[2], p[1]*p[2]))
e
vcov(e)

## Label new parameters
estimate(g, function(p) list("a1"=p[1]+p[2], "b1"=p[1]*p[2]))
#'
## Multiple group
m <- lvm(y~x)
m <- baptize(m)
d2 <- d1 <- sim(m,50,seed=1)
e <- estimate(list(m,m),list(d1,d2))
estimate(e) ## Wrong
ee <- estimate(e, id=rep(seq(nrow(d1)), 2)) ## Clustered
ee
estimate(lm(y~x,d1))

## Standardization (g-computation)
f <- function(p,data)
  list(p0=expit(p["(Intercept)"] + p["z"]*data[,"z"]),
       p1=expit(p["(Intercept)"] + p["x"] + p["z"]*data[,"z"]))
e <- expectation(f, g)
e
estimate(e,diff)
estimate(e,cbind(1,1))
## see also ?expectation

## More examples with clusters:
m <- lvm(c(y1,y2,y3)~u+x)
d <- sim(m,10)
l1 <- glm(y1~x,data=d)
l2 <- glm(y2~x,data=d)
l3 <- glm(y3~x,data=d)

## Some random id-numbers
id1 <- c(1,1,4,1,3,1,2,3,4,5)
id2 <- c(1,2,3,4,5,6,7,8,1,1)
id3 <- seq(10)

## Un-stacked and stacked i.i.d. decomposition
IC(estimate(l1,id=id1,stack=FALSE))
IC(estimate(l1,id=id1))

## Combined i.i.d. decomposition
e1 <- estimate(l1,id=id1)
e2 <- estimate(l2,id=id2)
e3 <- estimate(l3,id=id3)
(a2 <- merge(e1,e2,e3))

## If all models were estimated on the same data we could use the
## syntax:
## Reduce(merge,estimate(list(l1,l2,l3)))

## Same:
IC(a1 <- merge(l1,l2,l3,id=list(id1,id2,id3)))

IC(merge(l1,l2,l3,id=TRUE)) # one-to-one (same clusters)
IC(merge(l1,l2,l3,id=FALSE)) # independence


# ------ influence function calculus -------
ic1 <- scale(rnorm(10), scale=FALSE)
a <- estimate(coef = c("a" = 0.5), IC = ic1, id = 1:10)
b <- estimate(coef = c("b" = 0.8), IC = ic1, id = 1:10)

e <- c(a, b) # merge
merge(a, b)
c(e1=a, b) # naming of par
labels(e, c("p1", "p2")) # renaming parameters
e["a"] # subset
subset(e, "a")

# pipes
# c(a, b) |>
#  transform(function(x) x^2) |>
#  subset("a") |>
#  labels("sq")

# Parameter transformation with automatic calculation of derivatives
a * b
(3 * cos(a) / sqrt(b) + 1) / a
expit(c(a,b))
c(sum=sum(e), sum2=a+b,
  prod=prod(e), prod2=a*b)
e %*% e # inner prod.
c(1, 2) %*% e
c(pow = a^b)
a^c(0.5, 2)
c(b=e["a"] * e["b"] / a, also.b=e["b"])

B <- rbind(c(1,-1), c(1,0), c(0,1))
B %*% e
e == 1 # wald-test, null-hypothesis H0: b=1
e == c(1,2)
B %*% e == 1

lava documentation built on Oct. 23, 2026, 1:09 a.m.