1
$\begingroup$

Suppose I run a binomial GLM (in R) with response variable [0,1] and 2 predictor variables that are both categorical. Let's call them a and b where a has 3 factor levels (a1,a2,a3) and b has 2 factor levels (b1,b2).

Therefore:

mod <- glm(y~a+b, family=binomial(logit),data=pretend)

Let's say the summary output of the coefficients looks like this:

Coefficients:
                  Estimate Std. Error z value Pr(>|z|)    
(Intercept)        -0.62049    0.06988  -8.879  < 2e-16 ***
a a2                -2.24304    0.16111 -13.923  < 2e-16 ***
a.a3                -1.76965    0.14147 -12.509  < 2e-16 ***
s.b2                -0.79545    0.07918 -10.046  < 2e-16 ***

So I understand a1 and b1 are constrained in the intercept and the estimate values presented are log-odds ratios (I am familiar with interpreting them).

Let's say now I want to determine the estimate when a2 is true or = 1 and b1 is =1 (true), then I would just calculate it as B0 (intercept) + B1 (coeff of a2) = -0.62049 + -2.24304

Or if a3 is true and b2 then I would calculate it as B0 (intercept) + B3 (coeff of a3) + B4(coeff of b2). = -0.62049 + -1.76965 + -0.79545

But now how do I calculate the Standard errors and/or 95% Confidence intervals for the last two scenarios? Can I just add the standard errors in the same way? Or I would think I leave out the SE's of B0 (the intercept)?

So that in my first scenario where a2 is true and b1 is true, then I would just calculate it as SE of B1 (SE of a2) = 0.16111

And in my second scenario where a3 is true and b2 then I would calculate it as SE of B3 (SE of a3) + SE of B4(SE of b2). = 0.14147 + 0.07918.

Am I correct in my assumptions of interpreting the model output or am I completely wrong here?

Improve this question
$\endgroup$
  • $\begingroup$ you're completely wrong (but this is a fine question!) You can't add standard errors, you have to add variances (and you have to take the covariances into account as well). @Jarle's answer is correct; you can also compute this directly along the lines of stats.stackexchange.com/questions/105519/… $\endgroup$ – Ben Bolker Jun 6 '16 at 20:28
2
$\begingroup$

You just do something like

predict(mod,newdata=data.frame(a="a2",b="b1"),se.fit=TRUE)

See ?predict.glm for further details.

Improve this answer
$\endgroup$
  • $\begingroup$ Thank you. I did not want to use the 'predict' function because my response variable is a proportion value between 0 and 1 (fractional logit model) and I am calculating robust standard errors as suggested here: [link]stats.stackexchange.com/questions/117052/… Therefore, I am calculating my SE's after the fact, separate from the model using function 'coeftest' in package 'sandwich'. I am not sure how to bring those SE's produced together with the predict function or using the formulas supplied by @BenBolker's post. $\endgroup$ – MiMi Jun 7 '16 at 3:40
  • $\begingroup$ If you want the standard error of $p=logit^{-1}(\eta)$ instead of at the scale of the linear predictor $\eta$, just add the argument type="response" to predict. $\endgroup$ – Jarle Tufto Jun 7 '16 at 9:20
  • $\begingroup$ No I am calculating different types of SE's - heteroskedastic robust SE's - which cannot be solved by adding 'type=response'. Thank you though. $\endgroup$ – MiMi Jun 7 '16 at 9:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.