Airquality Assignment
Raunak Singh Chhabra
DATA 110
Airquality Tutorial and Homework Assignment
Load in the library
Load library tidyverse in order to access dplyr and ggplot2
library(tidyverse)The source for this dataset is the New York State Department of Conservation and the National Weather Service of 1973 for five months from May to September recorded daily.
Load the dataset into your global environment
Because airquality is a pre-built dataset, we can write it to our data directory to store it for later use.
data("airquality")Look at the structure of the data
str(airquality)## 'data.frame': 153 obs. of 6 variables:
## $ Ozone : int 41 36 12 18 NA 28 23 19 8 NA ...
## $ Solar.R: int 190 118 149 313 NA NA 299 99 19 194 ...
## $ Wind : num 7.4 8 12.6 11.5 14.3 14.9 8.6 13.8 20.1 8.6 ...
## $ Temp : int 67 72 74 62 56 66 65 59 61 69 ...
## $ Month : int 5 5 5 5 5 5 5 5 5 5 ...
## $ Day : int 1 2 3 4 5 6 7 8 9 10 ...View the data using the “head” function
The function, head, will only display the first 6 rows of the dataset. There are 153 observations (rows) and 6 variables.
head(airquality)## Ozone Solar.R Wind Temp Month Day
## 1 41 190 7.4 67 5 1
## 2 36 118 8.0 72 5 2
## 3 12 149 12.6 74 5 3
## 4 18 313 11.5 62 5 4
## 5 NA NA 14.3 56 5 5
## 6 28 NA 14.9 66 5 6Notice that all the variables are classified as either integers
Calculate Summary Statistics
If you want to look at specific statistics, here are some variations on coding. Here are 2 different ways to calculate “mean.”
mean(airquality$Temp)## [1] 77.88235mean(airquality[,4])## [1] 77.88235For the second way to calculate the mean, the matrix \[row,column\] is looking for column #4, which is the Temp column and we use all rows
Calculate Median, Standard Deviation, and Variance
median(airquality$Temp)## [1] 79sd(airquality$Wind)## [1] 3.523001var(airquality$Wind)## [1] 12.41154Rename the Months from number to names
Sometimes we prefer the months to be numerical, but here, we need them as the month names. There are MANY ways to do this. Here is one way to convert numbers 5 - 9 to May through September
airquality1 <- airquality |>
mutate(month_name = case_when(Month == 4 ~ "April",
Month == 5 ~ "May",
Month == 6 ~ "June",
Month == 7 ~ "July",
Month == 8 ~ "August",
Month == 9 ~ "September"))Now look at the summary statistics of the dataset
See how Month has changed to have characters instead of numbers (it is now classified as “character” rather than “integer”) in the new airquality1 dataframe.
summary(airquality1$month_name)## Length N.unique N.blank Min.nchar Max.nchar
## 153 5 0 3 9Month is a categorical variable with different levels, called factors.
This is one way to reorder the Months so they do not default to alphabetical (you will see another way to reorder DIRECTLY in the chunk that creates the plot below in Plot #1
airquality1$Month <- factor(airquality1$month_name,
levels=c("May", "June","July", "August",
"September"))Plot 1: Create a histogram categorized by Month
Here is a first attempt at viewing a histogram of temperature by the months May through September. We will see that temperatures increase over these months. The median temperature appears to be about 75 degrees.
fill = Month colors the histogram by months between May - Sept.
scale_fill_discrete(name = “Month”…) provides the month names on the right side as a legend in chronological order. This is a different way to order than what was shown above.
labs allows us to add a title, axes labels, and a caption for the data source
Plot 1 Code
p1 <- airquality1 |>
ggplot(aes(x=Temp, fill=month_name)) +
geom_histogram(position="identity")+
scale_fill_discrete(name = "Month",
labels = c("May", "June","July", "August", "September")) +
labs(x = "Monthly Temperatures from May - Sept",
y = "Frequency of Temps",
title = "Histogram of Monthly Temperatures from May - Sept, 1973",
caption = "New York State Department of Conservation and the National Weather Service")Plot 1 Output
p1
Is this plot useful in answering questions about monthly temperature
values? Not very. With position = "identity" the five
months are drawn directly on top of one another with no transparency, so
the bars in front hide the bars behind them and it is impossible to read
any single month’s distribution.
Plot 2: Improve the histogram of Average Temperature by Month
Outline the bars in white using the color = “white” command
Use alpha to add some transparency (values between 0 and 1)
Change the binwidth
Add some transparency and white borders around the histogram bars.
Plot 2 Code
p2 <- airquality1 |>
ggplot(aes(x=Temp, fill=month_name)) +
geom_histogram(position="identity", alpha=0.5, binwidth = 5, color = "white")+
scale_fill_discrete(name = "Month", labels = c("May", "June","July", "August", "September")) +
labs(x = "Monthly Temperatures from May - Sept",
y = "Frequency of Temps",
title = "Histogram of Monthly Temperatures from May - Sept, 1973",
caption = "New York State Department of Conservation and the National Weather Service")Plot 2 Output
p2
Here July stands out for having high frequency of 85 degree temperatures. The dark purple color indicates overlaps of months due to the transparency.
Did this improve the readability of the plot? Yes — the transparency and the wider bins make the overlapping months much easier to separate than in Plot 1.
Plot 3: Create side-by-side boxplots categorized by Month
We can see that August has the highest temperatures based on the boxplot distribution.
p3 <- airquality1 |>
ggplot(aes(Month, Temp, fill = month_name)) +
labs(x = "Months from May through September", y = "Temperatures",
title = "Side-by-Side Boxplot of Monthly Temperatures",
caption = "New York State Department of Conservation and the National Weather Service") +
geom_boxplot() +
scale_fill_discrete(name = "Month", labels = c("May", "June","July", "August", "September"))Plot 3 Output
p3
Notice that the points above and below the boxplots in June and July are outliers.
Plot 4: Side by Side Boxplots in Gray Scale
Make the same side-by-side boxplots, but in grey-scale
Use the scale_fill_grey command for the grey-scale legend, and again, use fill=Month in the aesthetics.
Plot 4 Code
Here we just changed the color palette to gray scale using scale_fill_grey
p4 <- airquality1 |>
ggplot(aes(Month, Temp, fill = month_name)) +
labs(x = "Monthly Temperatures", y = "Temperatures",
title = "Side-by-Side Boxplot of Monthly Temperatures",
caption = "New York State Department of Conservation and the National Weather Service") +
geom_boxplot()+
scale_fill_grey(name = "Month", labels = c("May", "June","July", "August", "September"))Plot 4 Output
p4
A Meaningful Change: Fixing the Mislabeled Legend
While reproducing the plots above I found a bug that runs through all four of them, and fixing it is the meaningful change I am submitting.
Every plot maps fill = month_name, which is a
character variable. ggplot2 sorts a character variable
alphabetically, so the fill groups are ordered:
levels(factor(airquality1$month_name))## [1] "August" "July" "June" "May" "September"That is August, July, June, May, September — not chronological. The
tutorial code then overrides the legend text with
labels = c("May", "June", "July", "August", "September"),
which does not reorder anything. It simply renames the
alphabetical groups in that order, so the August group gets labeled
“May,” the July group gets labeled “June,” and so on. Every legend in
Plots 1 through 4 is wrong, and the colors do not match the months they
claim to.
The fix is to map fill to the Month factor
we already built with the correct chronological levels, and then delete
the manual labels argument so ggplot2 uses the real level
names:
p2_fixed <- airquality1 |>
ggplot(aes(x = Temp, fill = Month)) +
geom_histogram(position = "identity", alpha = 0.5, binwidth = 5, color = "white") +
scale_fill_brewer(name = "Month", palette = "YlOrRd") +
labs(x = "Temperature (degrees Fahrenheit)",
y = "Number of days",
title = "Daily Temperatures by Month, May - September 1973",
subtitle = "Legend now in chronological order, with a palette that runs cool to warm",
caption = "New York State Department of Conservation and the National Weather Service")
p2_fixed
Two further improvements are included here. The legend is now in
calendar order, which lets a reader track the seasonal progression down
the legend instead of hunting through an alphabetical list. And
scale_fill_brewer(palette = "YlOrRd") replaces the default
rainbow hues with a sequential yellow-to-red ramp, so warmer months are
literally warmer colors. Month is ordered, so a sequential palette
matches the data better than unordered categorical colors.
The same correction applied to the boxplots:
p3_fixed <- airquality1 |>
ggplot(aes(x = Month, y = Temp, fill = Month)) +
geom_boxplot() +
scale_fill_brewer(name = "Month", palette = "YlOrRd") +
labs(x = "Month", y = "Temperature (degrees Fahrenheit)",
title = "Side-by-Side Boxplots of Monthly Temperatures, 1973",
subtitle = "July and August are the warmest months; May is both coolest and most variable",
caption = "New York State Department of Conservation and the National Weather Service")
p3_fixed
Now the box colors and the legend agree with the x-axis.
Plot 5: Ozone vs. Temperature
For my own plot I moved away from temperature-by-month and looked
instead at the relationship between two numerical variables we had not
yet explored: Ozone (parts per billion) and
Temp, with Wind speed added as a third
dimension.
Ozone has missing values, so I filtered them out first
and reported how many were dropped.
sum(is.na(airquality$Ozone))## [1] 37p5 <- airquality1 |>
filter(!is.na(Ozone)) |>
ggplot(aes(x = Temp, y = Ozone)) +
geom_point(aes(color = Month, size = Wind), alpha = 0.75) +
geom_smooth(method = "loess", se = TRUE, color = "grey25", linewidth = 0.8) +
scale_color_brewer(name = "Month", palette = "YlOrRd") +
scale_size_continuous(name = "Wind (mph)", range = c(1.5, 6)) +
labs(x = "Temperature (degrees Fahrenheit)",
y = "Ozone (parts per billion)",
title = "Ozone Rises Sharply with Temperature in New York, 1973",
subtitle = "Each point is one day; larger points are windier days. 37 days with missing ozone readings excluded.",
caption = "New York State Department of Conservation and the National Weather Service")
p5## `geom_smooth()` using formula = 'y ~ x'
Brief Essay
The plot type. Plot 5 is a scatterplot of two continuous numerical variables, ozone concentration against temperature, with a LOESS smoothing curve laid over it. It is meaningfully different from Plots 1 through 4: those were univariate displays of a single variable (temperature) broken out by a categorical variable (month), using histograms and boxplots. This plot shows the relationship between two numerical variables, which neither a histogram nor a boxplot can do. It also encodes four variables at once — temperature on the x-axis, ozone on the y-axis, month as color, and wind speed as point size — which follows the “show multivariate data” principle from the Unit 1 notes.
Insights. The relationship is clearly positive but not linear. Below roughly 80 degrees, ozone stays low and fairly flat, mostly under 50 ppb. Above 80 degrees the curve bends sharply upward and the spread widens dramatically, with the hottest days reaching over 150 ppb. That accelerating pattern is why the LOESS curve is more honest here than a straight regression line would be: a line would understate ozone on the hottest days and overstate it on cool ones.
The color and size channels add a second finding. The high-ozone points are almost entirely July and August, which is what the earlier boxplots would predict since those were the warmest months. The size channel shows the opposite relationship for wind: the large points, meaning windy days, cluster along the bottom of the plot at low ozone, while the extreme ozone days are small points, meaning nearly still air. That makes physical sense, since wind disperses pollutants while stagnant hot air lets ozone accumulate. Heat alone does not produce a bad air day; heat plus calm does.
One honest limitation: 37 of the 153 days have no ozone reading and are dropped here. If those readings failed disproportionately on particular kinds of days, the pattern could be biased in a way this plot cannot reveal.
Special code. Several things beyond the tutorial’s toolkit:
filter(!is.na(Ozone))removes the 37 missing ozone readings up front, so the count is reported explicitly rather than appearing as a silent ggplot2 warning.geom_smooth(method = "loess", se = TRUE)fits a locally weighted curve with a shaded confidence band. Unlikegeom_histogramorgeom_boxplot, this geom computes a model from the data and plots its predictions.- Mapping
size = Windinsideaes()adds a continuous third variable, andscale_size_continuous(range = c(1.5, 6))controls the point-size span so the smallest points stay visible while the largest do not swamp the plot. alpha = 0.75keeps overlapping points readable in the dense low-temperature cluster.scale_color_brewer(palette = "YlOrRd")reuses the cool-to-warm palette from my corrected plots so that color means the same thing throughout the document.