library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.3.3
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(lubridate)
## Warning: package 'lubridate' was built under R version 4.3.3
##
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
library(plotly)
## Warning: package 'plotly' was built under R version 4.3.3
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
library(readr)
aqidataset <- read_csv("C:/Users/Atul Kumar/Downloads/aqidataset.csv")
## New names:
## • `` -> `...1`
## Rows: 10924 Columns: 15
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (13): Date, Main_Pollutant, Site_Name_(of_Overall_AQI), Site_ID(of_Overa...
## dbl (2): ...1, Overall_AQI_Value
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
View(aqidataset)
data<- aqidataset
#Q1 What is the overall trend in Air Quality Index (AQI) values across all CBSAs in 2023?
# Assuming 'data' is your dataframe and it has columns 'date' and 'AQI'
data$Date <- as.Date(data$Date, format="%m/%d/%Y")
# Plotting the overall AQI trend using the correct column names
ggplot(data, aes(x=Date, y=Overall_AQI_Value)) +
geom_line(color=viridis::viridis(1), size=1) + # You can adjust the color and size of the line
labs(title="Overall AQI Trend in 2023", x="Date", y="Overall AQI Value") +
theme_minimal() +
theme(plot.title = element_text(hjust=0.5), # Center the title
axis.text.x = element_text(angle=45, hjust=1)) # Improve x-axis label readability
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Warning: Removed 4305 rows containing missing values or values outside the scale range
## (`geom_line()`).
#Q2 Which CBSAs consistently reported the worst AQI values in 2023?
data %>%
group_by('Site_Name_of_Overall_AQI' = 'Site_Name_(of_Overall_AQI)') %>%
summarise(Average_AQI = mean(Overall_AQI_Value, na.rm = TRUE)) %>%
top_n(-10, Average_AQI) %>%
ggplot(aes(x=reorder('Site_Name_of_Overall_AQI', Average_AQI), y=Average_AQI)) +
geom_col(fill="red") +
labs(title="Top 10 CBSAs with Worst AQI in 2023", x="CBSA", y="Average AQI") +
theme_minimal() +
coord_flip()
#Q3 How do AQI values compare across different states?
p <- ggplot(data, aes(x=state, y=Overall_AQI_Value, fill=state)) +
geom_boxplot(alpha=0.5) + # Reduced opacity for boxplot to better view jitter
geom_jitter(color="black", size=0.4, alpha=0.5, width=0.15) + # Adding jitter for individual points
labs(title="Comparison of AQI Across Different States in 2023", x="State", y="AQI") +
theme_light() +
theme(axis.text.x = element_text(angle=45, hjust=1))
# Printing the plot
print(p)
#Q4 What are the most common main pollutants identified across CBSAs?
p <- ggplot(data, aes(x=Main_Pollutant, fill=Main_Pollutant)) +
geom_bar() +
labs(title="Most Common Main Pollutants Across CBSAs in 2023", x="Pollutant", y="Count") +
theme_minimal()
# Convert to an interactive plotly object
ggplotly(p)
#Q5 Is there a correlation between Ozone levels and overall AQI values in major CBSAs?
g <- ggplot(data, aes(x=Ozone, y=Overall_AQI_Value)) +
geom_point(alpha=0.5) + # Use semi-transparency for points
geom_smooth(method="lm", color="blue") + # Add a linear model fit line
labs(title="Correlation Between Ozone Levels and AQI", x="Ozone Level", y="AQI") +
theme_light() +
theme(axis.text.x = element_text(angle=45, hjust=1)) # Rotate x-axis labels for better readability
# Convert to an interactive plotly object
ggplotly(g)
## `geom_smooth()` using formula = 'y ~ x'
#Q6 Which CBSAs experienced the most days with AQI categorized as unhealthy or worse?
# Ensure the data is properly filtered and counted
corrected_data <- data %>%
filter(AQI_category %in% c("Unhealthy", "Very Unhealthy", "Hazardous")) %>%
count(`Site_Name_(of_Overall_AQI)`) %>%
arrange(desc(n))
# Ensure we select the top 10 properly
top_sites <- corrected_data %>%
slice_max(order_by = n, n = 10)
# Plot the results
ggplot(top_sites, aes(x=reorder(`Site_Name_(of_Overall_AQI)`, n), y=n)) +
geom_col(fill="darkred") +
labs(title="CBSAs with Most Unhealthy Days in 2023", x="CBSA", y="Number of Days") +
theme_minimal() +
coord_flip()
data$Date <- as.Date(data$Date, format = "%m/%d/%Y") # Adjust the format if your date is in a different format
# Calculate the average PM2.5 by site and month
data_summary <- data %>%
mutate(month = month(Date, label = TRUE, abbr = TRUE)) %>%
group_by(`Site_Name_(of_Overall_AQI)`, month) %>%
summarize(Average_PM25 = mean(PM25, na.rm = TRUE), .groups = 'drop')
## Warning: There were 2305 warnings in `summarize()`.
## The first warning was:
## ℹ In argument: `Average_PM25 = mean(PM25, na.rm = TRUE)`.
## ℹ In group 1: `Site_Name_(of_Overall_AQI) = "A420910013LAT/LON POINT IS OF
## CORNER OF TRAILER"` and `month = Jul`.
## Caused by warning in `mean.default()`:
## ! argument is not numeric or logical: returning NA
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 2304 remaining warnings.
# Ensure that the data_summary is correctly prepared
print(head(data_summary))
## # A tibble: 6 × 3
## `Site_Name_(of_Overall_AQI)` month Average_PM25
## <chr> <ord> <dbl>
## 1 A420910013LAT/LON POINT IS OF CORNER OF TRAILER Jul NA
## 2 A420910013LAT/LON POINT IS OF CORNER OF TRAILER Oct NA
## 3 Air Management Services Laboratory (AMS LAB) Jun NA
## 4 Air Management Services Laboratory (AMS LAB) Nov NA
## 5 Air Management Services Laboratory (AMS LAB) Dec NA
## 6 Air Management Services Laboratory (AMS LAB) <NA> NA
count_data <- data_summary %>%
group_by(`Site_Name_(of_Overall_AQI)`, month) %>%
summarize(Count = n(), .groups = 'drop')
print(count_data)
## # A tibble: 2,305 × 3
## `Site_Name_(of_Overall_AQI)` month Count
## <chr> <ord> <int>
## 1 A420910013LAT/LON POINT IS OF CORNER OF TRAILER Jul 1
## 2 A420910013LAT/LON POINT IS OF CORNER OF TRAILER Oct 1
## 3 Air Management Services Laboratory (AMS LAB) Jun 1
## 4 Air Management Services Laboratory (AMS LAB) Nov 1
## 5 Air Management Services Laboratory (AMS LAB) Dec 1
## 6 Air Management Services Laboratory (AMS LAB) <NA> 1
## 7 Aldino Jan 1
## 8 Aldino Feb 1
## 9 Aldino Mar 1
## 10 Aldino Apr 1
## # ℹ 2,295 more rows
# Check the first few rows of the data summary to understand its structure
print(head(data_summary))
## # A tibble: 6 × 3
## `Site_Name_(of_Overall_AQI)` month Average_PM25
## <chr> <ord> <dbl>
## 1 A420910013LAT/LON POINT IS OF CORNER OF TRAILER Jul NA
## 2 A420910013LAT/LON POINT IS OF CORNER OF TRAILER Oct NA
## 3 Air Management Services Laboratory (AMS LAB) Jun NA
## 4 Air Management Services Laboratory (AMS LAB) Nov NA
## 5 Air Management Services Laboratory (AMS LAB) Dec NA
## 6 Air Management Services Laboratory (AMS LAB) <NA> NA
# Create the ggplot
gg_plot <- ggplot(data_summary, aes(x = month, y = Average_PM25, color = `Site_Name_(of_Overall_AQI)`, group = `Site_Name_(of_Overall_AQI)`)) +
geom_line() + # Connect points by lines
geom_point() + # Show individual data points
labs(title = "Seasonal PM2.5 Levels in Top 5 CBSAs", x = "Month", y = "Average PM2.5") +
theme_minimal()
# Display the ggplot
print(gg_plot)
## Warning: Removed 2305 rows containing missing values or values outside the scale range
## (`geom_line()`).
## Warning: Removed 2305 rows containing missing values or values outside the scale range
## (`geom_point()`).
# Convert ggplot to interactive Plotly plot
interactive_plot <- ggplotly(gg_plot)
# Display the interactive plot
print(interactive_plot)
#Q8 Analyze the variation of NO2 levels in CBSAs with the highest traffic volumes.
ggplot(data, aes(x=NO2, y=Overall_AQI_Value, color=city)) +
geom_point() +
geom_smooth(method="lm") +
labs(title="NO2 Levels in High Traffic CBSAs", x="NO2 Level", y="AQI")
## `geom_smooth()` using formula = 'y ~ x'
#Q9 Which pollutants have shown significant improvement or deterioration during 2023 in specific CBSAs?
# Assuming the data has already been prepared and summarized as data_summary
data_summary <- data %>%
mutate(Date = as.Date(Date, format = "%m/%d/%Y"), # Ensure Date is correctly formatted
Year = year(Date)) %>%
group_by(Year, Main_Pollutant) %>%
summarize(Average_AQI = mean(Overall_AQI_Value, na.rm = TRUE), .groups = 'drop')
# Create the ggplot object
gg_plot <- ggplot(data_summary, aes(x = factor(Year), y = Average_AQI, fill = Main_Pollutant)) +
geom_bar(stat = "identity", position = "dodge") +
labs(title = "Pollutant Trends Over Time in CBSAs", x = "Year", y = "Average AQI") +
theme_minimal()
# Convert to an interactive plotly object
interactive_plot <- ggplotly(gg_plot)
# Print the interactive plot
print(interactive_plot)
#Q10 What is the relationship between carbon monoxide (CO) levels and AQI values across urban centers?
ggplot(data, aes(x=CO, y=Overall_AQI_Value)) +
geom_point(aes(color=city), alpha=0.6) +
geom_smooth(method="lm") +
labs(title="CO Levels vs. AQI in Urban Centers", x="CO Level", y="AQI")
## `geom_smooth()` using formula = 'y ~ x'
#Q11 Compare the monthly average AQI values of CBSAs in coastal versus inland states.
data$month <- format(as.Date(data$Date), "%m")
data %>%
mutate(Type = ifelse(state %in% c("Coastal States Here"), "Coastal", "Inland")) %>%
group_by(month, Type) %>%
summarize(Average_AQI = mean(Overall_AQI_Value, na.rm=TRUE)) %>%
ggplot(aes(x=month, y=Average_AQI, fill=Type)) +
geom_bar(stat="identity", position="dodge") +
labs(title="Monthly Average AQI in Coastal vs. Inland States", x="Month", y="Average AQI")
## `summarise()` has grouped output by 'month'. You can override using the
## `.groups` argument.
#Q12 What are the effects of specific policies or events on the AQI in selected CBSAs?
data$Date <- as.Date(data$Date, format = "%Y-%m-%d") # Adjust the format as necessary
# Categorizing the data into periods based on a specific date
data$Period <- ifelse(data$Date >= as.Date("2023-06-01"), "After Policy", "Before Policy")
# Create the ggplot object
gg_plot <- ggplot(data, aes(x=Period, y=Overall_AQI_Value, fill=Period)) +
geom_boxplot() +
labs(title="Impact of Policy on AQI in Selected CBSAs", x="Period", y="AQI") +
theme_minimal()
# Convert ggplot object to plotly
interactive_plot <- ggplotly(gg_plot)
# Display the interactive plot
print(interactive_plot)
#Q13 How does air quality vary between weekdays and weekends in large CBSAs?
data$Day_of_Week <- weekdays(as.Date(data$Date))
data$Type_Day <- ifelse(data$Day_of_Week %in% c("Saturday", "Sunday"), "Weekend", "Weekday")
ggplot(data, aes(x=Type_Day, y=Overall_AQI_Value, fill=Type_Day)) +
geom_boxplot() +
labs(title="Air Quality: Weekdays vs. Weekends", x="Type of Day", y="AQI")
#Q14 Which CBSAs show the greatest disparity in AQI values from one quarter to the next?
data$Quarter <- quarters(as.Date(data$Date))
data %>%
group_by(Quarter, city) %>%
summarize(Average_AQI = mean(Overall_AQI_Value, na.rm=TRUE)) %>%
ggplot(aes(x=city, y=Average_AQI, fill=Quarter)) +
geom_bar(stat="identity", position="dodge") +
labs(title="Quarterly AQI Disparity in CBSAs", x="City", y="Average AQI")
## `summarise()` has grouped output by 'Quarter'. You can override using the
## `.groups` argument.
#Q15 How does air quality vary by day type (weekday vs. weekend) across different cities?
# Calculate average AQI for each city and day type
aqi_summary <- data %>%
group_by(city, Type_Day) %>%
summarize(Average_AQI = mean(Overall_AQI_Value, na.rm = TRUE), .groups = 'drop')
# Display the first few rows to check the summary
head(aqi_summary)
## # A tibble: 6 × 3
## city Type_Day Average_AQI
## <chr> <chr> <dbl>
## 1 Asheville Weekday 43.5
## 2 Asheville Weekend 42.8
## 3 Athens Weekday 57.9
## 4 Athens Weekend 52.6
## 5 Baltimore Weekday 52.6
## 6 Baltimore Weekend 42.7
# Plotting average AQI by city and day type
gg_plot <- ggplot(aqi_summary, aes(x = city, y = Average_AQI, fill = Type_Day)) +
geom_bar(stat = "identity", position = position_dodge()) +
labs(title = "Average AQI by City and Day Type", x = "City", y = "Average AQI") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Adjust text for better readability
# Print the plot
print(gg_plot)
# Convert ggplot to interactive plotly object
interactive_plot <- ggplotly(gg_plot)
# Display the interactive plot
print(interactive_plot)