Influence function based inference

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)

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.

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

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)
#> NULL


estimate(g)
#>              Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept) -0.001888 0.09785 -0.1937 0.1899 9.846e-01
#> z            0.953974 0.08114  0.7949 1.1130 6.481e-32
#> x            1.009058 0.14720  0.7205 1.2976 7.135e-12

## Testing contrasts
summary(estimate(g), null=0)
#> Call: estimate.default(x = x)
#> ────────────────────────────────────────────────────────────
#>              Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept) -0.001888 0.09785 -0.1937 0.1899 9.846e-01
#> z            0.953974 0.08114  0.7949 1.1130 6.481e-32
#> x            1.009058 0.14720  0.7205 1.2976 7.135e-12
#> ────────────────────────────────────────────────────────────
#> Null Hypothesis: 
#>   [(Intercept)] = 0
#>   [z] = 0
#>   [x] = 0 
#>  
#> chisq = 180.732, df = 3, p-value < 2.2e-16
estimate(g, rbind(c(1,1,0), c(1,0,2)))
#>                      Estimate Std.Err   2.5% 97.5%   P-value
#> [(Intercept)] + [z]    0.9521  0.1260 0.7052 1.199 4.075e-14
#> [(Intercept)] + 2[x]   2.0162  0.2404 1.5451 2.487 4.939e-17
summary(estimate(g, rbind(c(1,1,0), c(1,0,2))), null=c(1,2))
#> Call: estimate.default(x = x, f = ..1)
#> ────────────────────────────────────────────────────────────
#>                      Estimate Std.Err   2.5% 97.5% P-value
#> [(Intercept)] + [z]    0.9521  0.1260 0.7052 1.199  0.7037
#> [(Intercept)] + 2[x]   2.0162  0.2404 1.5451 2.487  0.9462
#> ────────────────────────────────────────────────────────────
#> Null Hypothesis: 
#>   [[(Intercept)] + [z]] = 1
#>   [[(Intercept)] + 2[x]] = 2 
#>  
#> chisq = 0.1447, df = 2, p-value = 0.9302
estimate(g, 2:3) ## same as cbind(0,1,-1)
#>           Estimate Std.Err    2.5%  97.5% P-value
#> [z] - [x] -0.05508  0.1537 -0.3563 0.2461    0.72
estimate(g, as.list(2:3)) ## same as rbind(c(0,1,0),c(0,0,1))
#>   Estimate Std.Err   2.5% 97.5%   P-value
#> z    0.954 0.08114 0.7949 1.113 6.481e-32
#> x    1.009 0.14720 0.7205 1.298 7.135e-12
## Alternative syntax
estimate(g, "z", "z"-"x", 2*"z"-3*"x")
#>             Estimate Std.Err    2.5%   97.5%   P-value
#> z            0.95397 0.08114  0.7949  1.1130 6.481e-32
#> [z] - [x]   -0.05508 0.15367 -0.3563  0.2461 7.200e-01
#> 2[z] - 3[x] -1.11922 0.43991 -1.9814 -0.2570 1.095e-02
estimate(g, "?")  ## Wildcards
#>           Estimate Std.Err    2.5%  97.5% P-value
#> [z] - [x] -0.05508  0.1537 -0.3563 0.2461    0.72
estimate(g, "*Int*", "z")
#>              Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept) -0.001888 0.09785 -0.1937 0.1899 9.846e-01
#> z            0.953974 0.08114  0.7949 1.1130 6.481e-32
summary(estimate(g, "1", "2"-"3"), null = c(0,1))
#> Call: estimate.default(x = x, f = "1", ..2)
#> ────────────────────────────────────────────────────────────
#>              Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept) -0.001888 0.09785 -0.1937 0.1899 9.846e-01
#> [z] - [x]   -0.055083 0.15367 -0.3563 0.2461 6.606e-12
#> ────────────────────────────────────────────────────────────
#> Null Hypothesis: 
#>   [(Intercept)] = 0
#>   [[z] - [x]] = 1 
#>  
#> chisq = 77.8686, df = 2, p-value < 2.2e-16
estimate(g, 2, 3)
#>   Estimate Std.Err   2.5% 97.5%   P-value
#> z    0.954 0.08114 0.7949 1.113 6.481e-32
#> x    1.009 0.14720 0.7205 1.298 7.135e-12

## Usual (non-robust) confidence intervals
estimate(g, vcov=TRUE)
#>              Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept) -0.001888 0.09817 -0.1943 0.1905 9.847e-01
#> z            0.953974 0.08318  0.7909 1.1170 1.892e-30
#> x            1.009058 0.14639  0.7221 1.2960 5.469e-12
estimate(g, vcov=vcov(g))
#>              Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept) -0.001888 0.09817 -0.1943 0.1905 9.847e-01
#> z            0.953974 0.08318  0.7909 1.1170 1.892e-30
#> x            1.009058 0.14639  0.7221 1.2960 5.469e-12

## Transformations
estimate(g, function(p) p[1]+p[2])
#>             Estimate Std.Err   2.5% 97.5%   P-value
#> (Intercept)   0.9521   0.126 0.7052 1.199 4.075e-14

## Multiple parameters
e <- estimate(g, function(p) c(p[1]+p[2], p[1]*p[2]))
e
#>                Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept)    0.952086 0.12596  0.7052 1.1990 4.075e-14
#> (Intercept).1 -0.001801 0.09335 -0.1848 0.1812 9.846e-01
vcov(e)
#>             (Intercept) (Intercept)
#> (Intercept)  0.01586612 0.008982930
#> (Intercept)  0.00898293 0.008714966

## Label new parameters
estimate(g, function(p) list("a1"=p[1]+p[2], "b1"=p[1]*p[2]))
#>     Estimate Std.Err    2.5%  97.5%   P-value
#> a1  0.952086 0.12596  0.7052 1.1990 4.075e-14
#> b1 -0.001801 0.09335 -0.1848 0.1812 9.846e-01
#'
## 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
#>        Estimate Std.Err     2.5%  97.5%   P-value
#> y@1      0.1044 0.08278 -0.05785 0.2666 2.073e-01
#> y~x@1    0.9665 0.08727  0.79541 1.1375 1.673e-28
#> y~~y@1   0.6764 0.10626  0.46809 0.8846 1.953e-10
ee <- estimate(e, id=rep(seq(nrow(d1)), 2)) ## Clustered
ee
#>        Estimate Std.Err    2.5%  97.5%   P-value
#> y@1      0.1044  0.1171 -0.1251 0.3338 3.726e-01
#> y~x@1    0.9665  0.1234  0.7246 1.2084 4.853e-15
#> y~~y@1   0.6764  0.1503  0.3818 0.9709 6.771e-06
estimate(lm(y~x,d1))
#>             Estimate Std.Err    2.5%  97.5%   P-value
#> (Intercept)   0.1044  0.1171 -0.1251 0.3338 3.726e-01
#> x             0.9665  0.1234  0.7246 1.2084 4.853e-15

## 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 Std.Err   2.5%  97.5%    P-value
#> p0   0.5010 0.02140 0.4591 0.5429 3.025e-121
#> p1   0.7007 0.01973 0.6620 0.7393 3.516e-276
estimate(e,diff)
#>    Estimate Std.Err  2.5%  97.5%   P-value
#> p1   0.1997  0.0279 0.145 0.2543 8.194e-13
estimate(e,cbind(1,1))
#>             Estimate Std.Err  2.5% 97.5% P-value
#> [p0] + [p1]    1.202 0.03027 1.142 1.261       0
## 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))
#>   (Intercept)           x
#> 1 -0.04904705 -0.14189362
#> 1  0.14744213 -0.15724052
#> 4  1.15168804 -0.54934034
#> 1  0.37775868  0.06124426
#> 3 -0.01647188 -0.03981494
#> 1 -1.75658710  2.38806661
#> 2 -0.48187968 -0.18690182
#> 3 -0.56660619 -0.30858145
#> 4  0.77907752 -0.82047631
#> 5  0.41462553 -0.24506186
#> attr(,"bread")
#>             (Intercept)          x
#> (Intercept)   0.6136517 -0.1000796
#> x            -0.1000796  0.9011828
IC(estimate(l1,id=id1))
#>   (Intercept)           x
#> 1  -0.6402167  1.07508836
#> 4   0.9653828 -0.68490833
#> 3  -0.2915390 -0.17419820
#> 2  -0.2409398 -0.09345091
#> 5   0.2073128 -0.12253093
#> attr(,"bread")
#>             (Intercept)          x
#> (Intercept)   0.6136517 -0.1000796
#> x            -0.1000796  0.9011828
#> attr(,"N")
#> [1] 10

## 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))
#>               Estimate Std.Err     2.5% 97.5%   P-value
#> (Intercept)     0.5418  0.2472  0.05731 1.026 0.0283949
#> x               0.9636  0.2592  0.45568 1.472 0.0002006
#> (Intercept).1   0.2156  0.5105 -0.78505 1.216 0.6728124
#> x.1             1.2369  0.4854  0.28554 2.188 0.0108277
#> (Intercept).2   0.4769  0.2979 -0.10705 1.061 0.1094554
#> x.2             1.1545  0.3714  0.42660 1.882 0.0018795

## 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)))
#>    (Intercept)          x (Intercept).1         x.1 (Intercept).2         x.2
#> 1   -1.2804333  2.1501767    3.53262530 -3.15431642     0.7554987  2.18566521
#> 4    1.9307656 -1.3698167   -2.27428000 -0.36871845    -1.5561326 -0.25228854
#> 3   -0.5830781 -0.3483964    1.08830429 -0.51910711    -0.1215704  0.05798750
#> 2   -0.4818797 -0.1869018    0.06810688 -0.07263298     0.5200220 -0.55458052
#> 5    0.4146255 -0.2450619    0.25612485  0.61909106    -0.4880168 -1.17960761
#> 6    0.0000000  0.0000000   -2.62738927  3.57191544     1.6742259 -2.27609724
#> 7    0.0000000  0.0000000    0.33520853  0.13001396     0.8010558  0.31069745
#> 8    0.0000000  0.0000000   -0.37870058 -0.20624550     0.1460747  0.07955428
#> 9    0.0000000  0.0000000    0.00000000  0.00000000    -1.3102941  1.37992084
#> 10   0.0000000  0.0000000    0.00000000  0.00000000    -0.4208633  0.24874863

IC(merge(l1,l2,l3,id=TRUE)) # one-to-one (same clusters)
#>    (Intercept)           x (Intercept).1         x.1 (Intercept).2         x.2
#> 1  -0.04904705 -0.14189362   0.004626991  0.01338593     0.7554987  2.18566521
#> 2   0.14744213 -0.15724052   0.068106877 -0.07263298     0.5200220 -0.55458052
#> 3   1.15168804 -0.54934034   1.088304293 -0.51910711    -0.1215704  0.05798750
#> 4   0.37775868  0.06124426  -2.274280005 -0.36871845    -1.5561326 -0.25228854
#> 5  -0.01647188 -0.03981494   0.256124851  0.61909106    -0.4880168 -1.17960761
#> 6  -1.75658710  2.38806661  -2.627389265  3.57191544     1.6742259 -2.27609724
#> 7  -0.48187968 -0.18690182   0.335208530  0.13001396     0.8010558  0.31069745
#> 8  -0.56660619 -0.30858145  -0.378700581 -0.20624550     0.1460747  0.07955428
#> 9   0.77907752 -0.82047631   2.342596539 -2.46707796    -1.3102941  1.37992084
#> 10  0.41462553 -0.24506186   1.185401770 -0.70062440    -0.4208633  0.24874863
IC(merge(l1,l2,l3,id=FALSE)) # independence
#>    (Intercept)          x (Intercept).1         x.1 (Intercept).2        x.2
#> 1  -0.14714116 -0.4256809    0.00000000  0.00000000     0.0000000  0.0000000
#> 2   0.44232639 -0.4717216    0.00000000  0.00000000     0.0000000  0.0000000
#> 3   3.45506412 -1.6480210    0.00000000  0.00000000     0.0000000  0.0000000
#> 4   1.13327605  0.1837328    0.00000000  0.00000000     0.0000000  0.0000000
#> 5  -0.04941565 -0.1194448    0.00000000  0.00000000     0.0000000  0.0000000
#> 6  -5.26976131  7.1641998    0.00000000  0.00000000     0.0000000  0.0000000
#> 7  -1.44563904 -0.5607055    0.00000000  0.00000000     0.0000000  0.0000000
#> 8  -1.69981856 -0.9257444    0.00000000  0.00000000     0.0000000  0.0000000
#> 9   2.33723256 -2.4614289    0.00000000  0.00000000     0.0000000  0.0000000
#> 10  1.24387660 -0.7351856    0.00000000  0.00000000     0.0000000  0.0000000
#> 11  0.00000000  0.0000000    0.01388097  0.04015779     0.0000000  0.0000000
#> 12  0.00000000  0.0000000    0.20432063 -0.21789893     0.0000000  0.0000000
#> 13  0.00000000  0.0000000    3.26491288 -1.55732132     0.0000000  0.0000000
#> 14  0.00000000  0.0000000   -6.82284001 -1.10615534     0.0000000  0.0000000
#> 15  0.00000000  0.0000000    0.76837455  1.85727317     0.0000000  0.0000000
#> 16  0.00000000  0.0000000   -7.88216780 10.71574631     0.0000000  0.0000000
#> 17  0.00000000  0.0000000    1.00562559  0.39004188     0.0000000  0.0000000
#> 18  0.00000000  0.0000000   -1.13610174 -0.61873650     0.0000000  0.0000000
#> 19  0.00000000  0.0000000    7.02778962 -7.40123388     0.0000000  0.0000000
#> 20  0.00000000  0.0000000    3.55620531 -2.10187319     0.0000000  0.0000000
#> 21  0.00000000  0.0000000    0.00000000  0.00000000     2.2664960  6.5569956
#> 22  0.00000000  0.0000000    0.00000000  0.00000000     1.5600661 -1.6637416
#> 23  0.00000000  0.0000000    0.00000000  0.00000000    -0.3647111  0.1739625
#> 24  0.00000000  0.0000000    0.00000000  0.00000000    -4.6683977 -0.7568656
#> 25  0.00000000  0.0000000    0.00000000  0.00000000    -1.4640503 -3.5388228
#> 26  0.00000000  0.0000000    0.00000000  0.00000000     5.0226778 -6.8282917
#> 27  0.00000000  0.0000000    0.00000000  0.00000000     2.4031674  0.9320923
#> 28  0.00000000  0.0000000    0.00000000  0.00000000     0.4382241  0.2386628
#> 29  0.00000000  0.0000000    0.00000000  0.00000000    -3.9308824  4.1397625
#> 30  0.00000000  0.0000000    0.00000000  0.00000000    -1.2625898  0.7462459


# ------ 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)
#>   Estimate Std.Err   2.5%  97.5%   P-value
#> a      0.5  0.2346 0.0401 0.9599 0.0331012
#> b      0.8  0.2346 0.3401 1.2599 0.0006511
c(e1=a, b) # naming of par
#>    Estimate Std.Err   2.5%  97.5%   P-value
#> e1      0.5  0.2346 0.0401 0.9599 0.0331012
#> b       0.8  0.2346 0.3401 1.2599 0.0006511
labels(e, c("p1", "p2")) # renaming parameters
#>    Estimate Std.Err   2.5%  97.5%   P-value
#> p1      0.5  0.2346 0.0401 0.9599 0.0331012
#> p2      0.8  0.2346 0.3401 1.2599 0.0006511
e["a"] # subset
#>   Estimate Std.Err   2.5%  97.5% P-value
#> a      0.5  0.2346 0.0401 0.9599  0.0331
subset(e, "a")
#>   Estimate Std.Err   2.5%  97.5% P-value
#> a      0.5  0.2346 0.0401 0.9599  0.0331

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

# Parameter transformation with automatic calculation of derivatives
a * b
#>   Estimate Std.Err    2.5%  97.5% P-value
#> a      0.4   0.305 -0.1979 0.9979  0.1898
(3 * cos(a) / sqrt(b) + 1) / a
#>   Estimate Std.Err   2.5% 97.5% P-value
#> a    7.887   5.319 -2.539 18.31  0.1382
expit(c(a,b))
#>   Estimate Std.Err   2.5%  97.5%   P-value
#> a   0.6225 0.05514 0.5144 0.7305 1.503e-29
#> b   0.6900 0.05019 0.5916 0.7884 5.360e-43
c(sum=sum(e), sum2=a+b,
  prod=prod(e), prod2=a*b)
#>       Estimate Std.Err    2.5%  97.5%  P-value
#> sum        1.3  0.4693  0.3802 2.2198 0.005604
#> sum2       1.3  0.4693  0.3802 2.2198 0.005604
#> prod       0.4  0.3050 -0.1979 0.9979 0.189758
#> prod2      0.4  0.3050 -0.1979 0.9979 0.189758
e %*% e # inner prod.
#>     Estimate Std.Err    2.5% 97.5% P-value
#> sum     0.89  0.6101 -0.3057 2.086  0.1446
c(1, 2) %*% e
#>     Estimate Std.Err   2.5% 97.5%  P-value
#> sum      2.1  0.7039 0.7203  3.48 0.002853
c(pow = a^b)
#>     Estimate Std.Err   2.5%  97.5%   P-value
#> pow   0.5743  0.1222 0.3348 0.8139 2.609e-06
a^c(0.5, 2)
#>    Estimate Std.Err    2.5%  97.5%   P-value
#> p1   0.7071  0.1659  0.3819 1.0323 2.029e-05
#> p2   0.2500  0.2346 -0.2099 0.7099 2.867e-01
c(b=e["a"] * e["b"] / a, also.b=e["b"])
#>        Estimate Std.Err   2.5% 97.5%   P-value
#> b           0.8  0.2346 0.3401  1.26 0.0006511
#> also.b      0.8  0.2346 0.3401  1.26 0.0006511

B <- rbind(c(1,-1), c(1,0), c(0,1))
B %*% e
#>           Estimate Std.Err    2.5%   97.5%   P-value
#> [a] - [b]     -0.3  0.0000 -0.3000 -0.3000 0.0000000
#> a              0.5  0.2346  0.0401  0.9599 0.0331012
#> b              0.8  0.2346  0.3401  1.2599 0.0006511
e == 1 # wald-test, null-hypothesis H0: b=1
#> Call: estimate.default(data = NULL, id = id, coef = coefs, IC = ic0, 
#>     stack = FALSE, keep = keep)
#> ────────────────────────────────────────────────────────────
#>   Estimate Std.Err   2.5%  97.5% P-value
#> a      0.5  0.2346 0.0401 0.9599  0.0331
#> b      0.8  0.2346 0.3401 1.2599  0.3940
#> ────────────────────────────────────────────────────────────
#> Null Hypothesis: 
#>   [a] = 1
#>   [b] = 1 
#>  
#> chisq = 2.2249, df = 1, p-value = 0.1358
e == c(1,2)
#> Call: estimate.default(data = NULL, id = id, coef = coefs, IC = ic0, 
#>     stack = FALSE, keep = keep)
#> ────────────────────────────────────────────────────────────
#>   Estimate Std.Err   2.5%  97.5%   P-value
#> a      0.5  0.2346 0.0401 0.9599 3.310e-02
#> b      0.8  0.2346 0.3401 1.2599 3.153e-07
#> ────────────────────────────────────────────────────────────
#> Null Hypothesis: 
#>   [a] = 1
#>   [b] = 2 
#>  
#> chisq = 13.1221, df = 1, p-value = 0.0002918
B %*% e == 1
#> Call: estimate.default(x = y, f = x)
#> ────────────────────────────────────────────────────────────
#>           Estimate Std.Err    2.5%   97.5% P-value
#> [a] - [b]     -0.3  0.0000 -0.3000 -0.3000  0.0000
#> a              0.5  0.2346  0.0401  0.9599  0.0331
#> b              0.8  0.2346  0.3401  1.2599  0.3940
#> ────────────────────────────────────────────────────────────
#> Null Hypothesis: 
#>   [[a] - [b]] = 1
#>   [a] = 1
#>   [b] = 1 
#>  
#> chisq = 2.2249, df = 1, p-value = 0.1358