I cannot figure out how to extract the standard error "sd(coef)" information from the logistf() regression model. These models are of class logistf and the manual states that data can be extracted this way:
The following generic methods are available for logistf`s output object: print, summary, coef, vcov, confint, anova, extractAIC, add1, drop1, profile, terms, nobs.
However, the standard error is not there.There isn't an object in str(summary(fit)) for the se(coef) and I have had a look at the source code without luck.
Any help would be appreciated!
logistf(formula = newdata[, i] ~ newdata[, j], data = newdata,
firth = TRUE)
Model fitted by Penalized ML
Confidence intervals and p-values by Profile Likelihood
coef se(coef) lower 0.95 upper 0.95 Chisq p
(Intercept) -4.110874e+00 0.8276236 -6.283919 -2.824075 Inf 0
newdata[, j] 1.691332e-08 1.6689839 -4.993849 2.957865 3.552714e-15 1
Likelihood ratio test=3.552714e-15 on 1 df, p=1, n=122
Indeed sd(fit) does not work and I am not sure if that usually works for other model classes.
However, the covariance matrix is available via vcov(fit) assuming the logistf model object is in fit. Then you can get the se(coef) column simply by calculating the square root of the main diagonal of the covariance matrix: sqrt(diag(vcov(fit)))
data(sex2)
fit<-logistf(case ~ age+oc+vic+vicl+vis+dia, data=sex2)
> fit
logistf(formula = case ~ age + oc + vic + vicl + vis + dia, data = sex2)
Model fitted by Penalized ML
Confidence intervals and p-values by Profile Likelihood
coef se(coef) lower 0.95 upper 0.95 Chisq p
(Intercept) 0.12025404 0.4855415 -0.8185574 1.07315110 0.06286298 8.020268e-01
age -1.10598130 0.4236601 -1.9737884 -0.30742658 7.50773092 6.143472e-03
oc -0.06881673 0.4437934 -0.9414360 0.78920202 0.02467044 8.751911e-01
vic 2.26887464 0.5484160 1.2730233 3.43543174 22.93139022 1.678877e-06
vicl -2.11140816 0.5430824 -3.2608608 -1.11773667 19.10407252 1.237805e-05
vis -0.78831694 0.4173676 -1.6080879 0.01518319 3.69740975 5.449701e-02
dia 3.09601263 1.6750089 0.7745730 8.03029352 7.89693139 4.951873e-03
Likelihood ratio test=49.09064 on 6 df, p=7.15089e-09, n=239
> diag(vcov(fit))
[1] 0.2357506 0.1794879 0.1969526 0.3007601 0.2949384 0.1741957 2.8056550
> sqrt(diag(vcov(fit)))
[1] 0.4855415 0.4236601 0.4437934 0.5484160 0.5430824 0.4173676 1.6750089