# Calculate change in odds and Pseudo-R-Squared with R

I have a question that touches on both technical solutions in R and statistics. I have a huge dataset with 2,400 respondents in total. I performed a logistic regression in order to analyze views on corruption in local government across different socioeconomic groups. The respondents could either say "not a lot/hardly any corrupt official" or "most/every offical is corrupt".

I am now looking for a way to calculate the change in odds of thinking that local officials are mostly corrupt. So I could e.g. say that there is a x percent decrease in the odds of thinking that corruption is widespread for men as opposed to women.

In addition to that, I would like to calculate Pseudo-R-Squared for each predictive variable, controlling for any other variables. I know how to do that for OLS, but the code does not work if I use on my glm model.

This is my model with "Not a lot/hardly any corrupt official" as the reference category for the dependent variable. Reference category for gender is female, and for education it is basic education.

glm(formula = corruption_local_recoded ~ gender + age + education_cat,
family = binomial(link = "logit"), data = lebanon, subset = (corruption_local_recoded !=
"Don't know" & education_cat != "No formal education"))

Deviance Residuals:
Min      1Q  Median      3Q     Max
-1.671  -1.290   0.896   1.017   1.468

Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept)                       1.274611   0.169750   7.509 5.97e-14 ***
genderMale                        0.169807   0.085740   1.980 0.047650 *
age                              -0.018510   0.002972  -6.228 4.74e-10 ***
education_catSecondary education -0.217526   0.107645  -2.021 0.043302 *
education_catHigher education    -0.402557   0.121817  -3.305 0.000951 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

Null deviance: 3139  on 2327  degrees of freedom
Residual deviance: 3095  on 2323  degrees of freedom
(42 observations deleted due to missingness)
AIC: 3105

Number of Fisher Scoring iterations: 4


This is a sample of the first 250 rows in my dataset:

structure(list(age = c(41L, 36L, 33L, 26L, 28L, 33L, 31L, 45L,
70L, 18L, 23L, 20L, 24L, 44L, 38L, 39L, 23L, 45L, 26L, 54L, 26L,
22L, 33L, 62L, 18L, 67L, 28L, 28L, 26L, 40L, 53L, 36L, 58L, 52L,
43L, 24L, 28L, 29L, 21L, 41L, 33L, 37L, 23L, 21L, 48L, 20L, 65L,
26L, 38L, 24L, 59L, 48L, 26L, 33L, 36L, 39L, 24L, 28L, 75L, 26L,
38L, 32L, 43L, 28L, 63L, 68L, 28L, 32L, 18L, 34L, 20L, 21L, 56L,
31L, 52L, 30L, 26L, 40L, 28L, 38L, 36L, 60L, 56L, 53L, 25L, 66L,
29L, 19L, 33L, 55L, 20L, 40L, 49L, 24L, 47L, 25L, 58L, 31L, 20L,
41L, 71L, 27L, 34L, 19L, 40L, 55L, 36L, 25L, 55L, 38L, 27L, 52L,
21L, 19L, 70L, 38L, 53L, 70L, 22L, 22L, 18L, 18L, 30L, 38L, 45L,
21L, 53L, 48L, 19L, 72L, 35L, 25L, 30L, 58L, 25L, 53L, 47L, 19L,
27L, 28L, 37L, 25L, 48L, 60L, 20L, 21L, 26L, 43L, 38L, 24L, 48L,
26L, 52L, 22L, 21L, 38L, 41L, 30L, 40L, 19L, 55L, 24L, 18L, 18L,
56L, 70L, 43L, 24L, 24L, 18L, 55L, 48L, 36L, 27L, 32L, 28L, 50L,
60L, 27L, 57L, 36L, 31L, 18L, 22L, 45L, 25L, 24L, 29L, 35L, 36L,
48L, 31L, 35L, 30L, 44L, 45L, 37L, 31L, 61L, 58L, 25L, 39L, 18L,
34L, 30L, 36L, 48L, 20L, 21L, 24L, 49L, 61L, 52L, 33L, 45L, 21L,
42L, 28L, 35L, 33L, 25L, 21L, 46L, 52L, 45L, 24L, 34L, 56L, 60L,
36L, 69L, 23L, 63L, 40L, 70L, 70L, 23L, 29L, 29L, 60L, 38L, 65L,
38L, 52L, 28L, 29L, 22L, 26L, 28L, 48L), gender = c("Male", "Female",
"Female", "Male", "Male", "Male", "Male", "Male", "Male", "Male",
"Female", "Male", "Male", "Female", "Male", "Female", "Female",
"Male", "Female", "Female", "Male", "Male", "Male", "Female",
"Male", "Male", "Male", "Male", "Male", "Female", "Male", "Female",
"Male", "Female", "Male", "Male", "Female", "Male", "Male", "Male",
"Male", "Female", "Male", "Male", "Male", "Male", "Female", "Female",
"Female", "Male", "Male", "Female", "Male", "Male", "Female",
"Female", "Male", "Male", "Male", "Female", "Male", "Female",
"Female", "Male", "Male", "Female", "Female", "Male", "Male",
"Female", "Female", "Male", "Male", "Female", "Female", "Male",
"Male", "Male", "Male", "Male", "Female", "Female", "Female",
"Male", "Male", "Male", "Female", "Female", "Male", "Male", "Male",
"Female", "Female", "Female", "Male", "Male", "Female", "Female",
"Female", "Female", "Male", "Male", "Male", "Female", "Female",
"Male", "Female", "Male", "Female", "Male", "Female", "Female",
"Female", "Female", "Male", "Male", "Female", "Female", "Male",
"Female", "Female", "Female", "Male", "Female", "Male", "Female",
"Female", "Female", "Male", "Female", "Female", "Male", "Female",
"Male", "Female", "Female", "Female", "Female", "Female", "Female",
"Male", "Female", "Male", "Female", "Female", "Male", "Male",
"Female", "Female", "Female", "Male", "Male", "Female", "Male",
"Male", "Female", "Female", "Male", "Female", "Female", "Female",
"Female", "Female", "Female", "Female", "Female", "Female", "Male",
"Female", "Male", "Male", "Female", "Male", "Male", "Male", "Female",
"Female", "Male", "Female", "Female", "Female", "Female", "Male",
"Female", "Female", "Male", "Male", "Female", "Male", "Female",
"Female", "Male", "Female", "Female", "Female", "Male", "Female",
"Male", "Female", "Female", "Female", "Female", "Male", "Female",
"Male", "Female", "Female", "Male", "Female", "Male", "Female",
"Male", "Female", "Male", "Female", "Male", "Female", "Female",
"Female", "Male", "Female", "Female", "Male", "Male", "Female",
"Female", "Female", "Male", "Male", "Male", "Female", "Male",
"Female", "Male", "Male", "Female", "Female", "Male", "Male",
"Male", "Female", "Male", "Male", "Female", "Female", "Female",
"Male", "Male", "Male", "Female"), education_cat = structure(c(2L,
2L, 3L, 2L, 3L, 3L, 2L, 2L, 2L, 3L, 4L, 3L, 4L, 2L, 2L, 4L, 3L,
4L, 2L, 3L, 3L, 2L, 2L, 2L, 3L, 3L, 2L, 2L, 2L, 2L, 2L, 2L, 3L,
2L, 2L, 3L, 3L, 3L, 2L, 2L, 2L, 4L, 3L, 4L, 2L, 2L, 1L, 4L, 3L,
3L, 3L, 3L, 2L, 2L, 2L, 2L, 2L, 3L, 2L, 3L, 4L, 3L, 2L, 2L, 2L,
2L, 4L, 2L, 2L, 3L, 3L, 2L, 3L, 3L, 3L, 3L, 3L, 2L, 3L, 3L, 2L,
3L, 2L, 2L, 4L, 2L, 2L, 2L, 3L, 3L, 4L, 2L, 2L, 3L, 3L, 2L, 3L,
4L, 4L, 2L, 2L, 2L, 3L, 4L, 4L, 3L, 3L, 4L, 2L, 3L, 2L, 2L, 2L,
4L, 3L, 3L, 4L, 2L, 3L, 2L, 4L, 2L, 3L, 4L, 2L, 2L, 2L, 3L, 4L,
2L, 3L, 4L, 4L, 2L, 2L, 1L, 2L, 4L, 3L, 2L, 2L, 3L, 2L, 2L, 3L,
2L, 3L, 2L, 2L, 4L, 3L, 3L, 2L, 2L, 3L, 3L, 4L, 2L, 3L, 4L, 3L,
4L, 2L, 2L, 2L, 2L, 2L, 4L, 4L, 3L, 2L, 2L, 4L, 2L, 2L, 3L, 3L,
2L, 4L, 2L, 4L, 3L, 2L, 4L, 2L, 3L, 2L, 4L, 4L, 2L, 1L, 3L, 2L,
2L, 3L, 2L, 3L, 2L, 2L, 2L, 4L, 3L, 2L, 3L, 4L, 3L, 2L, 4L, 4L,
3L, 2L, 4L, 3L, 3L, 2L, 2L, 2L, 3L, 2L, 2L, 4L, 3L, 3L, 2L, 3L,
2L, 2L, 3L, 2L, 3L, 2L, 4L, 2L, 3L, 2L, 2L, 3L, 1L, 4L, 2L, 3L,
2L, 2L, 2L, 4L, 3L, 3L, 3L, 4L, 2L), .Label = c("No formal education",
"Basic education", "Secondary education", "Higher education"), class = "factor"),
corruption_local_recoded = structure(c(1L, 1L, 1L, 1L, NA,
NA, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L,
2L, 2L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L,
2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L,
2L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 2L, 2L, 1L,
1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
1L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, NA, 1L, 2L, 2L, 2L,
2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L,
2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 2L,
2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L,
1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L,
1L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L,
2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 1L,
2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 1L,
2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L,
NA, 2L, 2L, 2L, 2L), .Label = c("Not a lot/hardly any corrupt official",
"Most/every official is corrupt", "Don't know", "Refused to answer"
), class = "factor")), row.names = c(NA, 250L), class = "data.frame")

• Some hints. However I'm afraid that if this is purely a software-related question, it is off-topic. If, on the contrary, you're interested in statistical issues with $R^2$ like measures for GLM, please clarify.
– chl
Oct 13, 2020 at 18:39

I am now looking for a way to calculate the change in odds of thinking that local officials are mostly corrupt. So I could e.g. say that there is a x percent decrease in the odds of thinking that corruption is widespread for men as opposed to women.

Save your regression to a variable called model.

odds = exp(coef(model))


In addition to that, I would like to calculate Pseudo-R-Squared for each predictive variable, controlling for any other variables.

model = glm(formula = corruption_local_recoded ~ gender + age + education_cat,
family = binomial(link = "logit"), data = lebanon, subset = (corruption_local_recoded !=
"Don't know" & education_cat != "No formal education"))


There's step one, you already did it. Now:

nullmodel = glm(formula = corruption_local_recoded ~ 1,
family = binomial(link = "logit"), data = lebanon, subset = (corruption_local_recoded !=
"Don't know" & education_cat != "No formal education"))


Now,

pseudoR = 1 - loglik(model)/loglik(nullmodel)


pseudoR is McFadden's pseudo-R-squared. I'm not sure what you mean when you say you want to calculate it for each individual variable controlling for other variables.

• Thanks! I meant that I want to have the Pseudo-R for each predictive variable, so I can see which variables contribute most to the variances.Both when only one variable is included in the model, and the full model with all variables. Oct 14, 2020 at 7:14
• There are several statistics called "Pseudo-R squared." I'm not familiar with the one you're searching for.
– AJV
Oct 14, 2020 at 12:26
• I don't have any preferences regarding Pseudo-R squared. Oct 14, 2020 at 14:27
• See hbiostat.org/bib/r2.html and also consider the anova.rms function in the R rms package. Aug 20 at 11:58