In MazamaScience/PWFSLSmoke: Utilities for Working with Air Quality Monitoring Data

knitr::opts_chunk$set(fig.width=7, fig.height=5)

Setup

Start by loading the PWFSLSmoke package and loading the data.

suppressPackageStartupMessages({
  library(PWFSLSmoke)
  library(ggplot2)
})

setSpatialDataDir("~/Data/Spatial")

logger.setLevel(ERROR)
wrcc_2015 <- wrcc_load()

Orientation

Create some maps and plots to get an idea about what the data looks like.

monitor_leaflet(wrcc_2015)

days2015 <- sum(!is.na(airsis_2015$data[-1]))/24
days2016 <- sum(!is.na(airsis_2016$data[-1]))/24
days2017 <- sum(!is.na(airsis_2017$data[-1]))/24
deployments2015 <- ncol(airsis_2015$data) - 1
deployments2016 <- ncol(airsis_2016$data) - 1
deployments2017 <- ncol(airsis_2017$data) - 1

Let's see how they were distributed through time:

counts2015 <-  data.frame(time = wrcc_2015$data$datetime, 
                          count = apply(wrcc_2015$data[-1], 1, function(x) sum(!is.na(x))))
ggplot(counts2015, aes(time, count))+
  geom_bar(stat="identity")+
  labs(title='Number of valid data points by hour', subtitle='2015')+
  ylim(0,35)

Now that we've seen when the data was, let's look at where it is. There are several different options for plotting monitors from a ws_monitor object. monitor_esriMap() will plot points for each monitor with AQI colors indicating the maximum PM2.5 value using a map from ESRI, and uses R's base plotting capabilities so it is easier to add points or polygons on top. The other option is monitor_map() which plots the points over polygons instead of a map image.

# Get the bbox for the map
bbox <- RgoogleMaps::qbbox(wrcc_2015$meta$latitude, wrcc_2015$meta$longitude, margin = list(m = c(10,10,10,10), TYPE = "perc"))
basemap <- esriMap_getMap(bboxString = paste0(bbox$lonR[1], ",", bbox$latR[1], ",", bbox$lonR[2], ",", bbox$latR[2]))
layout(matrix(c(1:(4*13)), 13, 4, byrow = TRUE), widths=c(.8,4,4,4), heights = c(1, rep(4, 13)))
par(mar = c(0,0,0,0))
plot.new()
for(year in c("2015", "2016", "2017")){
  plot(0,0,col="transparent", axes = F)
  text(0,0, year, cex = 2)
}

for ( month in 1:12) {
  # First 2015
  tlim2015 <- if (month < 12) {
    paste0('2015', stringr::str_pad(as.character(c(month, month+1)), 2, side = 'left', "0"), '01')
  } else {
    c('20151201', '20160101')
  }
  ws_monitor2015 <- monitor_subset(wrcc_2015, tlim = tlim2015)
  month_str <- strftime(ws_monitor2015$data$datetime[30], format = "%B")
  plot(0,0,col="transparent", axes = F)
  text(0,0,month_str, srt = 90, cex = 1.5)

  monitor_esriMap(ws_monitor2015, mapRaster = basemap)
  box("figure", col = 'gray50')
  plot.new()
  plot.new()
}

Eastern Washington

Eastern Washington tends to be prone to wildfires. Let's take a look at how air quality in Okanogan County in each of the three years to see what stories that might be able to tell us about the fire season in 2015, 2016, and 2017 and how they changed.

sisk2015 <- wrcc_2015

loadSpatialData("USCensusCounties")
monitor_map(sisk2015, cex = 2)
title("2015", cex.main = 3)

Now, let's take a look at the pm2.5 data from each year.

monitor_timeseriesPlot(sisk2015, style='gnats', xlab = "", localTime = FALSE, ylim = c(0,800))
mtext("2015", 4, padj = 1)
title("All pm2.5 readings from monitors in Okanogan County")

Just from looking at these plots, we can say some things about the fire season in Northern California. There is little smoke before July, and PM2.5 levels tend to peak around August or September, and then air quality returns to normal around the beginning of November. PM2.5 levels to not follow a smooth curve. During high-smoke season, PM2.5 levels tend to stay above a baseline, and spike briefly before returning to near the baseline.

We can see that 2015 and 2017 had more smoke, and for a longer time, than 2016. Using this as a proxy for fire intensity and duration, we can say that 2015 and 2017 had longer fire seasons than 2016. It appears that there were two distinct smoke events in 2016, the first much more intense than the second. Upon closer examination, however, we see that there is no data from 2016 between September 23 and October 5.

monitor_timeseriesPlot(sisk2016, tlim = c(20160920, 20161010), type = "l")

Let's look a little closer at the smokiest part of the summer for the three years: July to October.

monitor_dygraph(sisk2015, tlim = c(20150701, 20151001), title = "2015")

Zooming in closer, a periodic pattern emerges between the dramatic smoke events.

monitor_dygraph(title = "2017", monitor_subset(sisk2017), tlim = c(20170708, 20170720))

Let's see how PM2.5 levels change throughout the day for one monitor. This plot shows the hourly PM2.5 values from every day between July 8 - July 20 overlaid on top of each other with the mean in a darker color.

# monitorPlot_timeOfDaySpaghetti() is now defunct, see `?monitorPlot_timeOfDaySpaghetti()`

Air quality in this location tends to peak around 9am, and disappate throughout the day. Smoke will settle in valleys when there is little wind, then disappate when winds pick back up again, so patterns like this are due to the typical weather patterns in the area.

Let's look at 2017 and see how air quality is related to major fires in the region. From the Cal Fire report we can see that there were three fires in Okanogan County which burned over 10,000 acres: The Eclipse Complex, Salmon August Complex, and Orleans Complex. We can find their locations in (InciWeb)[https://inciweb.nwcg.gov/]. Here are their locations, with markers indicating the locations of temporary monitors:

fireNames <- c("Eclipse Complex", "Salmon August Complex", "Orleans Complex")
fireLongitudes <- c(-123.493, -123.099, -123.647)
fireLatitudes <- c(41.792, 41.263, 41.566)
fireStartDate <- strptime(c("20170815", "20170813", "20170725"), format = "%Y%m%d")
fireEndDate <- strptime(c("20171010", "20171020", "20171012"), format = "%Y%m%d")
fireDF <- data.frame(name = fireNames, longitude =fireLongitudes, latitude = fireLatitudes, 
                     startDate = fireStartDate, endDate = fireEndDate)
monitor_map(sisk2017, ylim = c(40.8, 42.15), xlim = c(-124.4, -121), cex = 0)
addMarker(sisk2017$meta$longitude, sisk2017$meta$latitude, color = "blue")
addIcon('orangeFlame', fireDF$longitude, fireDF$latitude, expansion = .001)
text(fireLongitudes, fireLatitudes, fireNames, pos = 1)

Let's see if the timing of the fires lines up with spikes in PM2.5 concentrations.

# Join the monitors into one monitor object with average PM2.5 values
collapsed2017 <- monitor_collapse(sisk2017, monitorID = "Okanogan_2017")
monitor_timeseriesPlot(collapsed2017, tlim = c(20170701, 20171101), type = "l", ylim = c(-78, 400), col = "gray40")
title("Okanogan County, CA 2017 fires and smoke")
abline(0,0,col = 'gray80')

# Add bars below the plot indicating the duration of each fire
rect(fireDF$startDate[1], -33, fireDF$endDate[1], -8, col = 'deepskyblue4', border = NA)
rect(fireDF$startDate[2], -63, fireDF$endDate[2], -38, col = 'darkolivegreen4', border = NA)
rect(fireDF$startDate[3], -93, fireDF$endDate[3], -68, col = 'darkorange3', border = NA)
text(fireDF$startDate, c(-20, -50, -80), fireDF$name, pos = 2, cex = .8)

We start seeing elevated PM2.5 levels shortly after the Orleans Complex is sparked. It's not until the middle of August, however, when two large fires start in the county that we really start seeing high PM2.5 concentrations. The PM2.5 concentrations drop to pre-fire levels around the middle of September and generally stay down, even as the fires continue burning well into October.

MazamaScience/PWFSLSmoke documentation built on July 3, 2023, 11:03 a.m.