1

I want to compute a logit regression for rare events. I decided to use the Zelig package (relogit function) to do so.

Usually, I use stargazer to extract and save regression results. However, there seem to be compatibility issues with these two packages (Using stargazer with Zelig).

I now want to extract the following information from the Zelig relogit output:

Coefficients, z values, p values, number of observations, log likelihood, AIC

I have managed to extract the p-values and coefficients. However, I failed at the rest. But I am sure these values must be accessible somehow, because they are reported in the summary() output (however, I did not manage to store the summary output as an R object). The summary cannot be processed in the same way as a regular glm summary (https://stats.stackexchange.com/questions/176821/relogit-model-from-zelig-package-in-r-how-to-get-the-estimated-coefficients)

A reproducible example:

##Initiate package, model and data

require(Zelig)

data(mid)

z.out1 <- zelig(conflict ~ major + contig + power + maxdem + mindem + years,
                data = mid, model = "relogit")

##Call summary on output (reports in console most of the needed information)

summary(z.out1)

##Storing the summary fails and only produces a useless object

summary(z.out1) -> z.out1.sum

##Some of the output I can access as follows

z.out1$get_coef() -> z.out1.coeff
z.out1$get_pvalue() -> z.out1.p
z.out1$get_se() -> z.out1.se

However, I did not find similar commands for other elements, such as z values, AIC etc. However, as they are shown in the summary() call, they should be accessible somehow.

The summary call result:

Model: 

Call:
z5$zelig(formula = conflict ~ major + contig + power + maxdem + 
    mindem + years, data = mid)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-3.0742  -0.4444  -0.2772   0.3295   3.1556  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.535496   0.179685 -14.111  < 2e-16
major        2.432525   0.157561  15.439  < 2e-16
contig       4.121869   0.157650  26.146  < 2e-16
power        1.053351   0.217243   4.849 1.24e-06
maxdem       0.048164   0.010065   4.785 1.71e-06
mindem      -0.064825   0.012802  -5.064 4.11e-07
years       -0.063197   0.005705 -11.078  < 2e-16

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 3979.5  on 3125  degrees of freedom
Residual deviance: 1868.5  on 3119  degrees of freedom
AIC: 1882.5

Number of Fisher Scoring iterations: 6

Next step: Use 'setx' method

1 Answer 1

2

Use from_zelig_model for deviance, AIC.

m <- from_zelig_model(z.out1)
m$aic
...

Z-values are coefficient / sd.

z.out1$get_coef()[[1]]/z.out1$get_se()[[1]]
5
  • Thank you, that works! Just two follow-ups: 1) Do you think z.out1.dat.loglik <- (z.out1.dat$deviance)/-2 is correct for converting to log likelihood? 2) Do you know how I can extract the number of valid observations, excluding NA (n)?
    – deca
    Feb 7, 2018 at 16:22
  • 1
    1) No, I'm not sure how to get the likelihood here unfortunately. Deviance is calculated using loglik of a saturated and your model, but the one for the saturated model is not given here.
    – erocoar
    Feb 7, 2018 at 19:57
  • 1
    2) Do you mean the number of rows used for the regression? If so, length(complete.cases(m$model)) should be a safe way to do that
    – erocoar
    Feb 7, 2018 at 19:59
  • Great, thanks! Your solution works perfectly, so I accept your answer.
    – deca
    Feb 7, 2018 at 20:17
  • Glad it helped!
    – erocoar
    Feb 7, 2018 at 22:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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