You can download this .qmd file from here. Just hit the Download Raw File button.
Determinants of COVID vaccination rates
First, a little detour to describe several alternatives for reading in data:
If you navigate to my Github account, and find the 264_fall_2026 repo, there is a Data folder inside. You can then click on vacc_Mar21.csv to see the data we want to download. This link should also get you there, but it’s good to be able to navigate there yourself.
Approach 1: create a Data folder in the same location where this .qmd file resides, and then store vaccinations_2021.csv in that Data folder
2
Approach 2: give R the complete path to the location of vaccinations_2021.csv, starting with Home (~)
3
Approach 3: link to our course webpage, and then know we have a Data folder containing all our csvs
4
Approach 4: navigate to the data in GitHub, hit the Raw button, and copy that link
A recent Stat 272 project examined determinants of covid vaccination rates at the county level. Our data set contains 3053 rows (1 for each county in the US) and 14 columns; here is a quick description of the variables we’ll be using:
state = state the county is located in
county = name of the county
region = region the state is located in
metro_status = Is the county considered “Metro” or “Non-metro”?
rural_urban_code = from 1 (most urban) to 9 (most rural)
perc_complete_vac = percent of county completely vaccinated as of 11/9/21
tot_pop = total population in the county
votes_Trump = number of votes for Trump in the county in 2020
votes_Biden = number of votes for Biden in the county in 2020
perc_Biden = percent of votes for Biden in the county in 2020
ed_somecol_perc = percent with some education beyond high school (but not a Bachelor’s degree)
ed_bachormore_perc = percent with a Bachelor’s degree or more
unemployment_rate_2020 = county unemployment rate in 2020
median_HHincome_2019 = county’s median household income in 2019
Consider only Minnesota and its surrounding states (Iowa, Wisconsin, North Dakota, and South Dakota). We want to examine the relationship between the percentage who voted for Biden and the percentage of complete vaccinations by state. Generate two plots to examine this relationship:
A scatterplot with points and smoothers colored by state. Make sure the legend is ordered in a meaningful way, and include good labels on your axes and your legend. Also leave off the error bars from your smoothers.
`geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Describe which plot you prefer and why. What can you learn from your preferred plot?
I prefer the second plot because it allows us to look at each state separately. This makes it easier to see the relationship between the percentage voting for Biden and the percentage completely vaccinated within each state. The first plot makes it easier to compare the states directly, but the second plot makes the patterns within each state easier to see.
From the second plot, I can see that counties with higher percentages of votes for Biden generally tend to have higher vaccination rates. However, the relationship is not exactly the same in every state. Looking at the states separately makes these differences easier to see.
Produce 3 different plots for illustrating the relationship between the rural_urban_code and percent vaccinated. Hint: you can sometimes turn numeric variables into categorical variables for plotting purposes (e.g. as.factor(), ifelse()).
`geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
Warning: Failed to fit group -1.
Caused by error in `smooth.construct.cr.smooth.spec()`:
! x has insufficient unique values to support 10 knots: reduce k.
State your favorite plot, why you like it better than the other two, and what you can learn from your favorite plot. Create an alt text description of your favorite plot, using the Four Ingredient Model. See this link for reminders and references about alt text.
y favorite plot is the boxplot because it makes it easier to compare vaccination rates across the different rural-urban codes. The scatterplot shows all of the individual counties, but some of the points overlap. The average vaccination plot shows the average for each category, but it does not show how much the counties vary within each category.
The boxplot shows that vaccination rates are generally different across the rural-urban codes. More urban counties tend to have higher vaccination rates, while more rural counties tend to have lower vaccination rates. However, there is still variation within each rural-urban category.
Alt text: The boxplot compares the percent of counties that are completely vaccinated across rural-urban codes from 1 to 9. The x-axis shows rural-urban code, with lower numbers representing more urban counties and higher numbers representing more rural counties. The y-axis shows the percent completely vaccinated. Each box shows the distribution of vaccination rates for counties in that rural-urban category. Overall, the boxes show that vaccination rates tend to be higher for the more urban categories and lower for the more rural categories.
BEFORE running the code below, sketch the plot that will be produced by R. AFTER running the code, describe what conclusion(s) can we draw from this plot?
vaccine_data |>filter(!is.na(perc_Biden)) |>mutate(big_states =fct_lump(state, n =10)) |>group_by(big_states) |>summarize(IQR_Biden =IQR(perc_Biden)) |>mutate(big_states =fct_reorder(big_states, IQR_Biden)) |>ggplot() +geom_point(aes(x = IQR_Biden, y = big_states))
Sketch: There will be the 10 largest states plus an “Other” group because fct_lump(state, n = 10) keeps the 10 most common states and combines the rest into “Other”.
The plot shows that the variability in the percentage of votes for Biden differs among the states. The IQR represents the spread of Biden’s vote percentages within each state, so states with a larger IQR have more variation in Biden vote percentages across their counties, while states with a smaller IQR have less variation. The “Other” group represents all states that were not among the 10 most common states in the data.
In this question we will focus only on the 12 states in the Midwest (i.e. where region == “Midwest”).
Create a tibble with the following information for each state. Order states from least to greatest state population.
number of different rural_urban_codes represented among the state’s counties (there are 9 possible)
Use your tibble in (a) to produce a plot of the relationship between proportion of Metro counties and median unemployment rate. Points should be colored by the number of different rural_urban_codes in a state, but a single linear trend should be fit to all points. What can you conclude from the plot?
midwest_states |>ggplot(aes(x = proportion_metro,y = median_unemployment,color = rural_codes )) +geom_point() +geom_smooth(aes(group =1), se =FALSE) +labs(x ="Proportion of Metro counties",y ="Median unemployment rate",color ="Number of Rural-Urban Codes" )
`geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Warning: The following aesthetics were dropped during statistical transformation:
colour.
ℹ This can happen when ggplot fails to infer the correct grouping structure in
the data.
ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
variable into a factor?
The plot shows a general negative relationship between the proportion of Metro counties and median unemployment rate. As the proportion of Metro counties increases, the median unemployment rate tends to decrease. However, there is still some variation among the states.
Describe the tibble temp created above. What would be the dimensions? What do rows and columns represent?
The tibble temp has 3 rows and 4 columns. The three rows represent the three groups created from people_per_MD. The columns are MD_group, n, mean_perc_vac, and mean_white. The n column gives the number of counties in each group, mean_perc_vac gives the average vaccination percentage, and mean_white gives the average percentage of residents who identify as white.
What would happen if we replaced new_perc_vac = ifelse(perc_complete_vac > 95, NA, perc_complete_vac) with new_perc_vac = ifelse(perc_complete_vac > 95, perc_complete_vac, NA)?
The new code would do the opposite of the original code. It would keep vaccination percentages greater than 95 and change all vaccination percentages that are 95 or lower to NA.
What would happen if we replaced mean_white = mean(perc_white, na.rm = TRUE) with mean_white = mean(perc_white)?
If we removed na.rm = TRUE, R would not ignore missing values. If a group contains an NA value for perc_white, the mean for that group would be NA instead of being calculated from the available values.
What would happen if we removed group_by(MD_group)?
If we removed group_by(MD_group), the summary calculations would be done for the entire data set instead of separately for each MD group. The resulting tibble would have one row instead of three rows.
Hypothetical R chunk #2:
# Hypothetical R chunk 2ggplot(data = vaccine_data) +geom_point(mapping =aes(x = perc_over_65, y = perc_complete_vac, color = HR_party)) +geom_smooth()temp <- vaccine_data |>group_by(HR_party) |>summarise(var1 =n()) |>arrange(desc(var1)) |>slice_head(n =3)vaccine_data |>ggplot(mapping =aes(x =fct_reorder(HR_party, perc_over_65, .fun = median), y = perc_over_65)) +geom_boxplot()
Why would the first plot produce an error?
The first plot produces an error because the x and y variables are only included in the mapping for geom_point(). The geom_smooth() layer does not know which variables to use. The mapping should be placed inside ggplot() so that both geom_point() and geom_smooth() can use the same x and y variables.
Describe the tibble temp created above. What would be the dimensions? What do rows and columns represent?
The tibble temp has 3 rows and 2 columns. The three rows represent the three parties with the largest numbers of counties. The HR_party column gives the party, and var1 gives the number of counties for that party.
What would happen if we replaced fct_reorder(HR_party, perc_over_65, .fun = median) with HR_party?
If we replaced fct_reorder(HR_party, perc_over_65, .fun = median) with HR_party, the parties would no longer be ordered according to their median percentage of residents over 65. They would instead use their normal factor order.
Hypothetical R chunk #3:
# Hypothetical R chunk 3vaccine_data |>filter(!is.na(people_per_MD)) |>mutate(state_lump =fct_lump(state, n =4)) |>group_by(state_lump, rural_urban_code) |>summarise(mean_people_per_MD =mean(people_per_MD)) |>ggplot(mapping =aes(x = rural_urban_code, y = mean_people_per_MD, colour =fct_reorder2(state_lump, rural_urban_code, mean_people_per_MD))) +geom_line()
Describe the tibble piped into the ggplot above. What would be the dimensions? What do rows and columns represent?
The tibble contains one row for each combination of state_lump and rural_urban_code that occurs in the data. There are five state groups because the four largest states are kept and the remaining states are combined into “Other.” Since there are nine possible rural-urban codes, there can be up to 45 rows.
The tibble has three columns: state_lump, rural_urban_code, and mean_people_per_MD. The last column gives the average number of people per doctor for each state group and rural-urban code.
Carefully describe the plot created above.
The plot shows the relationship between rural-urban code and the mean number of people per doctor. The x-axis shows the rural-urban code, from more urban counties to more rural counties. The y-axis shows the mean number of people per doctor. Each line represents one of the state groups created by fct_lump().
The plot allows us to compare how the number of people per doctor changes as the rural-urban code increases and to compare the patterns between the different state groups.
What would happen if we removed filter(!is.na(people_per_MD))?
If we removed filter(!is.na(people_per_MD)), counties with missing values for people_per_MD would remain in the data. Because the mean function does not use na.rm = TRUE, a group containing missing values could have a mean of NA. This could cause some parts of the plot to be missing.
What would happen if we replaced fct_reorder2(state_lump, rural_urban_code, mean_people_per_MD) with state_lump?
If we replaced fct_reorder2(state_lump, rural_urban_code, mean_people_per_MD) with state_lump, the lines would still be colored according to the state groups. However, the state groups would no longer be reordered based on the rural-urban code and mean number of people per doctor. They would use their normal factor order instead.
