13 Dummy Variables Part 2
There are several places in this file where you need to do things (e.g., write R code, write Latex). Some are labeled “YOUR CODE GOES HERE.” Some places tell you to “UN-COMMENT-OUT” code after you’ve done something else (un-comment-out means remove the hash tags from code in code chunks so that it becomes code instead of comments). Some places say “ReplaceWithCoefficientFromBetaModel”. Others are labeled “PLACEHOLDER”. Make sure to go through the file and find all the places you need to do something.
This chapter goes with LN4.2. In LN4.2 we examined two models of the relationship between education and wages. We examined two different ways we could define sets of dummy variables to represent people with a high school degree, a BA, a Master’s, or a PhD as their highest level of education. We saw how both models could be used to answer questions such as “what is the average wage for someone with a Master’s degree?” or “What is the expected difference in average wage for a person with a Master’s degree compared to a BA?” In this BP chapter you will explore these models empirically using micro-level data (i.e., data on individual people).
Specifically, the data we’ll use comes from the 2021 ACS 1-year Public Use Microdata Sample (PUMS) (the 2022 ACS was only released January 25, 2024 and the API doesn’t seem to have been fully updated yet, so 2021 is the most recently available data). This data has the responses of individuals that the Census Bureau uses to calculate the county- and state-level statistics they publish as part of the ACS (e.g., what you’re using for the CP). Formally to analyze this data properly we need to use survey weights to make calculations done using the individual responses representative of the population as a whole. We’re not going to worry about that. However, note that our calculations will be off somewhat from the true average wages for the groups we explore in this chapter.
We’ll obtain our data using the tidycensus package, as described here. If you don’t already have a Census API key, get one and use it, as described here. For our purposes, it’s fine to leave census_api_key("YOUR API KEY GOES HERE") in your code.
library(tidyverse)
library(tidycensus)
library(stargazer)
library(pander)
## This turns off scientific notation when there are fewer then 4 digits. Otherwise pander displays very small p values with way too many decimal places
options(scipen = 4)13.1 Obtain and prepare data
To work with the PUMS data via tidycensus, you can pull the table of variables to find what variables are included. I’ve set this code chunk to not be evaluated (eval = FALSE) when you build because while doing this is helpful while you’re working with this data, it doesn’t belong in your HTML output
## To work with the PUMS data via tidycensus, you can pull the table of variables to find what variables are included
## I've set this code chunk to not be evaluated (eval = FALSE) when you build. You might want to do this on your own, but it doesn't belong in your HTML output
pums_vars_2021 <- pums_variables %>%
filter(year == 2021, survey == "acs1")Let’s pull the following variables from the 2021 ACS-1 year Public Use Microdata (i.e., data on individual people) for the state of Wisconsin (we’re limiting our sample to Wisconsin so that it doesn’t take too long to load or work with).
| Variable | Description |
|---|---|
PUMA |
Public Use Microdata Areas |
WAGP |
WAGe of Person |
AGEP |
AGE of Person |
SCHL |
educational attainment of person (i.e., highest level obtained) |
sex |
SEX of person |
The results='hide' code chunk option keeps it from displaying status updates while it downloads (which is sometimes about 100 lines of output, so you definitely don’t want it in your HTML output).
## Download 2021 acs1 for Wisconsin.
wiPUMS <- get_pums(
variables = c("PUMA","WAGP", "AGEP", "SCHL", "SEX"),
state = "wi",
survey = "acs1",
year = 2021,
recode = TRUE
)Let’s rename WAGP as wage so that we remember what it is.
Let’s see what levels of education exist in our data, and how many observations we have for each
| SCHL | SCHL_label | n |
|---|---|---|
| 01 | No schooling completed | 1555 |
| 02 | Nursery school, preschool | 654 |
| 03 | Kindergarten | 708 |
| 04 | Grade 1 | 671 |
| 05 | Grade 2 | 611 |
| 06 | Grade 3 | 764 |
| 07 | Grade 4 | 653 |
| 08 | Grade 5 | 765 |
| 09 | Grade 6 | 780 |
| 10 | Grade 7 | 789 |
| 11 | Grade 8 | 1377 |
| 12 | Grade 9 | 966 |
| 13 | Grade 10 | 1143 |
| 14 | Grade 11 | 1329 |
| 15 | 12th grade - no diploma | 859 |
| 16 | Regular high school diploma | 14503 |
| 17 | GED or alternative credential | 1618 |
| 18 | Some college, but less than 1 year | 3734 |
| 19 | 1 or more years of college credit, no degree | 6587 |
| 20 | Associate’s degree | 5398 |
| 21 | Bachelor’s degree | 8796 |
| 22 | Master’s degree | 3209 |
| 23 | Professional degree beyond a bachelor’s degree | 672 |
| 24 | Doctorate degree | 483 |
| bb | N/A (less than 3 years old) | 1657 |
We’re interested in examining how wages differ by educational attainment. We’d like to limit our data to people who are old enough to have completed high school (or equivalent) and college, and then be working. So, let’s limit our sample to people age 25 or older. Let’s also only examine people who have a wage.
We should look at what education levels (SCHL) remain after filtering
| SCHL | SCHL_label | n |
|---|---|---|
| 01 | No schooling completed | 158 |
| 02 | Nursery school, preschool | 4 |
| 03 | Kindergarten | 6 |
| 05 | Grade 2 | 2 |
| 06 | Grade 3 | 4 |
| 07 | Grade 4 | 6 |
| 08 | Grade 5 | 6 |
| 09 | Grade 6 | 38 |
| 10 | Grade 7 | 8 |
| 11 | Grade 8 | 108 |
| 12 | Grade 9 | 79 |
| 13 | Grade 10 | 126 |
| 14 | Grade 11 | 195 |
| 15 | 12th grade - no diploma | 341 |
| 16 | Regular high school diploma | 6340 |
| 17 | GED or alternative credential | 823 |
| 18 | Some college, but less than 1 year | 1860 |
| 19 | 1 or more years of college credit, no degree | 3273 |
| 20 | Associate’s degree | 3627 |
| 21 | Bachelor’s degree | 5882 |
| 22 | Master’s degree | 2130 |
| 23 | Professional degree beyond a bachelor’s degree | 455 |
| 24 | Doctorate degree | 328 |
After filtering, the value “bb” that indicates a person who is under 3 years old and thus can’t have any schooling is now gone. Now all of the levels of SCHL are numeric, which allows us to convert them from a character to a numeric data type. This will then allow us to use inequalities to define dummy variables (e.g., SCHL <16).
We’re comparing wages of people who have a high school degree (or equivalent) or higher, so we also need to drop everyone who has less than a high school degree. Because SCHL is now numeric, we can do this using an inequality (rather than listing all of the options separately).
In order to match with LN4.2, we also need to drop people with a professional degree. In practice we would probably want to examine them too, but for the sake of matching with what is in LN4.2 we’re going to drop them.
Whenever you do something like filtering data, it’s a very good idea to look at the data and make sure it worked. In previous years of 380, students have wasted many hours trying to get models to work, only to finally go back and look at the data to find that they actually messed up an earlier step that seems easy. So even if it’s something easy you think you know how to do, look at the data. You can display summary measures as we’ll do here, but it’s also a good idea to click on it in the Environment tab and actually scroll through it quickly. Typically you don’t display the results in a paper, but for the purposes of the BP, I want to demonstrate what you might do (e.g., get a count by education levels). Here, we could just re-run this:
| SCHL | SCHL_label | n |
|---|---|---|
| 16 | Regular high school diploma | 6340 |
| 17 | GED or alternative credential | 823 |
| 18 | Some college, but less than 1 year | 1860 |
| 19 | 1 or more years of college credit, no degree | 3273 |
| 20 | Associate’s degree | 3627 |
| 21 | Bachelor’s degree | 5882 |
| 22 | Master’s degree | 2130 |
| 24 | Doctorate degree | 328 |
Note that this variable is lists the highest level of education attained. That means it generally builds. We can assume someone who reached a 21 has a 16, and a person who reached a 22 has a 21 (and thus also a 16), and a person who reached 24 has a 22 (and thus also a 21 and 16).
13.2 Define dummy variables
In LN4.2 we had an “alpha model” and “beta model” which defined education in different ways.
| Variable | Alpha Model Definition | Beta Model Definition |
|---|---|---|
BA |
=1 if have B.A. degree, =0 if don’t |
=1 if highest degree is B.A. degree, =0 otherwise |
Masters |
=1 if have Master’s degree, =0 if don’t |
=1 if highest degree is Master’s, =0 otherwise |
PhD |
=1 if have PhD, =0 if don’t |
=1 if highest degree is PhD, =0 otherwise |
13.2.1 Alpha Model
In LN4.2, we defined the “Alpha Model” as
\[ wage = \alpha_0 + \alpha_1 BA + \alpha_2 Masters + \alpha_3 PhD + u \]
We need to create the dummy variables used for the “Alpha Model”. We’ll prefix the variables with “a” for “Alpha”. Later you’ll define “Beta Model” variables and prefix them with “b”.
In our data, some people have some college or an Associates Degree. While we might be interested in differences between people with some college and a only a high school degree, for the purposes of what we’re doing here (learning about dummy variables), let’s classify anyone without a BA as having high school as their highest degree (even if they have some college or an Associates Degree).
# Create dummy variables for "alpha model" (starting with a)
wiPUMS <- wiPUMS %>%
mutate(aBA = ifelse(SCHL >= 21, 1, 0),
aMasters = ifelse(SCHL >= 22, 1, 0),
aPhD = ifelse(SCHL == 24, 1, 0)
)Now let’s estimate the regression shown as “Alpha Model” in LN4
| Estimate | Std. Error | t value | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 44840 | 400.7 | 111.9 | 0 |
| aBA | 22914 | 771.5 | 29.7 | 1.873e-190 |
| aMasters | 11608 | 1279 | 9.079 | 1.179e-19 |
| aPhD | 26332 | 2999 | 8.78 | 1.742e-18 |
| Observations | Residual Std. Error | \(R^2\) | Adjusted \(R^2\) |
|---|---|---|---|
| 24263 | 50562 | 0.07036 | 0.07025 |
13.2.2 Beta Model
In LN4.2, we defined the “Beta Model” (version 1…we won’t use the “version 2” at all) as:
\[ wage = \beta_0 + \beta_1 BA + \beta_2 Masters + \beta_3 PhD + u \]
We need to create the dummy variables for the “Beta Model” from LN4.2. Use the same categories as you did for the Alpha Model (so we’re considering anyone with some college or an associates degree as having high school as their highest education). Prefix these variables with a “b” for “Beta”.
YOUR CODE GOES HERE: add variables bBA, bMasters, and bPhD to wiPUMS
# Create dummy variables for "Beta Model" (prefix with b)
# YOUR CODE GOES HERE
wiPUMS <- wiPUMS %>%
mutate(bBA = ifelse(SCHL == 21, 1, 0),
bMasters = ifelse(SCHL == 22, 1, 0),
bPhD = ifelse(SCHL == 24, 1, 0)
)Now let’s estimate the regression shown as “Beta Model” in LN4.2
UN-COMMENT-OUT the following code after you estimate the Beta Model
| Estimate | Std. Error | t value | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 44840 | 400.7 | 111.9 | 0 |
| bBA | 22914 | 771.5 | 29.7 | 1.873e-190 |
| bMasters | 34523 | 1167 | 29.59 | 4.048e-189 |
| bPhD | 60854 | 2820 | 21.58 | 2.747e-102 |
| Observations | Residual Std. Error | \(R^2\) | Adjusted \(R^2\) |
|---|---|---|---|
| 24263 | 50562 | 0.07036 | 0.07025 |
13.3 Compare the regressions side-by-side
In LN4.2 we talk about how different models can lead to equivalent results theoretically. Here we want to examine that theoretical equivalence numerically. Specifically, we want to write out various conditional expectations that should be equivalent theoretically and show that the estimated values are indeed equivalent.
Let’s start by displaying the two models on the same table. You have to be very careful doing this (so you don’t accidentally mis-label variables), but we can re-name the variable labels for each model so that they are the same (e.g., rename aBA and bBA to both be BA) and thus display on the same row of the stargazer table. We’ll make a copy of the model results that we use specifically for this purpose. Then we’ll rename the coefficients to have the same names. Then we’ll display them using stargazer.
UN-COMMENT-OUT the following code after you estimate the Beta Model
alphaModelStargazer <- alphaModel
betaModelStargazer <- betaModel
names(alphaModelStargazer$coefficients) <- c("(Intercept)", "BA", "Masters", "PhD")
names(betaModelStargazer$coefficients) <- c("(Intercept)", "BA", "Masters", "PhD")
stargazer(alphaModelStargazer, betaModelStargazer,
type = "html",
report=('vc*p'),
notes = "<em>*p<0.1;**p<0.05;***p<0.01</em>",
notes.append = FALSE,
model.numbers = FALSE,
column.labels = c("Alpha","Beta"))| Dependent variable: | ||
| wage | ||
| Alpha | Beta | |
| BA | 22,914.380*** | 22,914.380*** |
| p = 0.000 | p = 0.000 | |
| Masters | 11,608.400*** | 34,522.790*** |
| p = 0.000 | p = 0.000 | |
| PhD | 26,331.610*** | 60,854.400*** |
| p = 0.000 | p = 0.000 | |
| Constant | 44,840.080*** | 44,840.080*** |
| p = 0.000 | p = 0.000 | |
| Observations | 24,263 | 24,263 |
| R2 | 0.070 | 0.070 |
| Adjusted R2 | 0.070 | 0.070 |
| Residual Std. Error (df = 24259) | 50,562.000 | 50,562.000 |
| F Statistic (df = 3; 24259) | 612.030*** | 612.030*** |
| Note: | *p<0.1;**p<0.05;***p<0.01 | |
13.4 Compare the predictions of each model
In LN4.2 we claimed that the Alpha Model and Beta Model were theoretically equivalent. Here we are demonstrating that they are also empirically equivalent. The first thing you notice is that the \(R^2\) and adjusted \(R^2\) are identical. We haven’t talked about the other statistics (Residual Std. Error and F Statistic), but they’re identical too. That’s what we would expect if the models are indeed equivalent. Recall the definition of \(R^2\). You’d expect equivalent models to explain the same fraction of the variation in wages.
Now let’s check if real-world situations that are supposed to be equivalent in theory are indeed equivalent numerically.
First, let’s store our coefficients in variables so they’re easier to work with:
## Alpha Model
a0 <- coef(alphaModel)["(Intercept)"]
a1 <- coef(alphaModel)["aBA"]
a2 <- coef(alphaModel)["aMasters"]
a3 <- coef(alphaModel)["aPhD"]
## Beta Model
## FILL THESE IN WITH COEFFICIENTS FROM YOUR BETA MODEL
b0 <- coef(betaModel)["(Intercept)"]
b1 <- coef(betaModel)["bBA"]
b2 <- coef(betaModel)["bMasters"]
b3 <- coef(betaModel)["bPhD"]
## I'll use this to let you know where you need to fill in answers
PLACEHOLDER <- "PLACEHOLDER"We’ll check a variety of situations explored by these models. For each situation, we’ll write out the following 6 steps (if the steps aren’t clear, also look at the examples below):
First, we’ll write out the situation in English (e.g., “Expected/average wage for a person whose highest degree is high school”)
Next, we’ll write out the conditional expectation corresponding to that situation (e.g., \(E(wage|HS)\)) on the left hand side of the equation. To te right of the
|in theE(...|...), write it out as a short description of the situation in the world (e.g.,HS) instead of writing out all the variable values (e.g.,BA=0, Masters=0, PhD=0).Then, on the right hand side of the equation, we’ll write out the conditional expectation in terms of model parameters using parameters from the Alpha Model (e.g., \(\alpha_0\))…
…and in terms of model parameters using parameters from the Beta Model (e.g., \(\beta_0\)).
Last, for the third and fourth rows on the right hand side of the equation, we’ll calculate the value estimated using the Alpha Model (e.g., 44840.08)…
…and the value estimated using the Beta Model (e.g., 44840.08).
Note that we’re using format(...,nsmall=2) instead of round(...,2) to display the coefficients to two decimal places.
I’ve filled in a few parts for you below to show you what you’re supposed to do. Anyplace you see “PLACEHOLDER” you need to fill in whatever should go there. (Note that I filled in the coefficient b0 for you, but until you do the parts above, it will display as ReplaceWithCoefficientFromBetaModel.)
Use LaTex for parameters (e.g., \(\alpha_0\), \(\beta_0\)). For values, each value should be displayed to 2 decimal places using format(...,nsmall=2). Note that you can add together multiple coefficients inside the format function (e.g. format(a0 + a1,nsmall=2)). This structure is already set up for you below (i.e., you just need to replace everywhere that says PLACEHOLDER).
When you start filling it in, make sure to build frequently so you don’t mess it up and then are unable to figure out what isn’t working. For the first few, I’d build each time you make a change (remember to set it in _bookdown.yml to only build this chapter while you’re working on it, and then just build the entire book when you’re finished before committing/pushing to GitHub).
For each scenrio, your parameters in rows one and two on the right hand side should match the slides, and the values in rows three and four should be identical. If the values aren’t identical, you did something wrong.
Here are the scenarios:
Expected/average wage for a person whose highest degree is high school (HS)
\[ \begin{aligned} E(wage|HS) &= \alpha_0 \\&= \beta_0 \\&= 44840.08 \\&= 44840.08 \end{aligned} \]
Expected/average wage for a person whose highest degree is BA
\[ \begin{aligned} E(wage|BA) &= \alpha_0 + \alpha_1 \\&= \beta_0 + \beta_1 \\&= 67754.47 \\&= 67754.47 \end{aligned} \]
Expected/average wage for a person whose highest degree is Master’s
\[ \begin{aligned} E(wage|Masters) &= \alpha_0 + \alpha_1 + \alpha_2 \\&= \beta_0 + \beta_2 \\&= 79362.87 \\&= 79362.87 \end{aligned} \]
Expected/average wage for a person whose highest degree is PhD
\[ \begin{aligned} E(wage|Masters) &= \alpha_0 + \alpha_1 + \alpha_2 + \alpha_3 \\&= \beta_0 + \beta_3 \\&= 105694.48 \\&= 105694.48 \end{aligned} \]
How much higher do you expect the average wage to be for someone whose highest degree is a BA compared to someone whose highest degree is high school (HS)
\[ \begin{aligned} E(wage|BA) - E(wage|HS) &= \alpha_0 + \alpha_1 - (\alpha_0) = \ alpha_1 \\&= \beta_0 + \beta_1 - (\beta_0) = \beta_1 \\&= 22914.38 \\&= 22914.38 \end{aligned} \]
How much higher do you expect the average wage to be for someone whose highest degree is a Master’s compared to someone whose highest degree is high school (HS)
\[ \begin{aligned} E(wage|Masters) - E(wage|HS) &= \alpha_0 + \alpha_1 + \alpha_2 - (\alpha_0) = \alpha_1 + \alpha_2 \\&= \beta_0 + \beta_2 - (\beta_0) = \beta_2 \\&= 34522.79 \\&= 34522.79 \end{aligned} \]
How much higher do you expect the average wage to be for someone whose highest degree is a PhD compared to someone whose highest degree is high school (HS)
\[ \begin{aligned} E(wage|PhD) - E(wage|HS) &= \alpha_0 + \alpha_1 + \alpha_2 + \alpha_3 - (\alpha_0) = \alpha_1 + \alpha_2 + \alpha_3 \\&= \beta_0 + \beta_3 - (\beta_0) = \beta_3 \\&= 60854.40 \\&= 60854.40 \end{aligned} \]
How much higher do you expect the average wage to be for someone whose highest degree is a Master’s compared to someone whose highest degree is a BA
\[ \begin{aligned} E(wage|Masters) - E(wage|BA) &= \alpha_0 + \alpha_1 + \alpha_2 + - (\alpha_0 + \alpha_1) = \alpha_2 \\&= \beta_0 + \beta_2 - (\beta_0 + \beta_1) = \beta_2 - \beta_1 \\&= 11608.40 \\&= 11608.40 \end{aligned} \]
How much higher do you expect the average wage to be for someone whose highest degree is a PhD compared to someone whose highest degree is a BA
\[ \begin{aligned} E(wage|PhD) - E(wage|BA) &= \alpha_0 + \alpha_1 + \alpha_2 + \alpha_3 - (\alpha_0 + \alpha_1) = \alpha_2 + \alpha_3 \\&= \beta_0 = \beta_3 - (\beta_0 + \beta_1) = \beta_3 - \beta_1 \\&= 37940.02 \\&= 37940.02 \end{aligned} \]
How much higher do you expect the average wage to be for someone whose highest degree is a PhD compared to someone whose highest degree is a Master’s
\[ \begin{aligned} E(wage|PhD) - E(wage|Masters) &= \alpha_0 + \alpha_1 + \alpha_2 + \alpha_3 - (\alpha_0 + \alpha_1 + \alpha_2) = \alpha_3 \\&= \beta_0 = \beta_3 - (\beta_0 + \beta_2) = \beta_3 - \beta_2 \\&= 26331.61 \\&= 26331.61 \end{aligned} \]
13.5 Other things you should think about (and possibly do)
Doing the things below (and doing them correctly) is required to get a check plus. If you don’t do them, you can still get a check.
Remember that a check plus is a 100 and a check is a 93, so the difference is small. Only do the following if you actually have time (e.g., the priority is doing well on the group projects). But if you have time, do them. It will help you better understand.
13.5.1 Group averages
Linear regression when all explanatory variables are dummy variables produces the same estimates of average values for different groups as if you simply took the average of the \(y\) value for each group. Try confirming this using group_by() and summarize() to get average wage for each education level (make sure not to save wiPUMS as grouped it still works as expected below), and then compare them with what we found above in the regressions (i.e., make sure they are the same).
Why then do we use regression? One reason is because we easily get standard errors, confidence intervals, and p-values for hypothesis tests of the group averages and differences between groups. Another reason is because we can also control for other variables (such as age…see below).
13.5.2 Causal estimates?
Make sure you can explain why we should not interpret our results as causal. For example, above we estimated that the average wage is 22914.38 higher for someone whose highest degree is a BA compared to someone whose highest degree is high school (HS). Why can we not claim that getting a BA causes this increase in wage? In other words, why is \(\hat{\alpha}_1\) (or \(\hat{\beta}_1\)) not an unbiased estimate of the true effect of getting a BA on wages?
13.5.3 What about age?
How do the results change if you control for age (i.e., add age as a variable)? You can use either the alpha or beta model and display the results side-by-side (i.e., the original model and the new model including age) using stargazer.
13.5.4 What about sex?
How might you add sex to the model?
In the data, male is coded as sex=1 and female is coded as sex=2 (see the PUMS data dictionary). Create a female dummy variable (i.e., female=1 indicates female, female=0 indicates male because that’s the only other option in the data).
Estimate two models (explained below) and display the results side-by-side using stargazer, along with the original model for comparison (either alpha or beta).
For the first model, simply add female as a variable.
For the second model, allow how education affects wages to be different for males versus females.
What does each model allow you to say about differences in wages by sex? In your answer, include specific examples using your coefficients.