knitr::opts_chunk$set(fig.width = 7, fig.height = 5)
knitr::opts_knit$set(root.dir = "C:/Users/astri/Mirror/Mazamascience/Projects/AirSensor")

Introduction

This tutorial demonstrates how to create files containing 'pas' and 'pat' data for a particular community and how to save and access them in a local directory. Target audiences include grad students, researchers and any member of the public concerned about air quality and comfortable working with R and RStudio.

Tutorials in this series include:

PAS and PAT Objects

pas objects in the AirSensor package are described in Purple Air Synoptic Data. They contain per instrument metadata for all Purple Air sensors operating at a particular moment in time. This will include both "spatial metadata" like longitude, latitude, timezone, etc. as well as per-instrument keys allowing us to obtain timeseries data from web services.

pat objects are described in Purple Air Timeseries Data. These contain both temporally invariant spatial metadata for a single sensor as well as the actual time series measurements made by that sensor.

In order to successfully work with PurpleAir sensor data, we will need to create both a 'pas' object and a collection of 'pat' objects for the sensors we are interested in.

Goal

Our goal in this tutorial is to create 'pas' and 'pat' data for a single month for the Methow Valley -- a community in north-central Washington State. Clean Air Methow operates as a project of the Methow Valley Citizens Council and began deploying Purple Air Sensors in 2018:

In the summer of 2018, Clean Air Methow launched the Clean Air Ambassador Program, an exciting citizen science project, and one of the largest, rural networks of low-cost sensors in the world!

This tutorial will demonstrate how to create 'pas' and 'pat data for this collection of sensors for September, 2020 when the Methow Valley experienced poor air quality due to Wildfire smoke.

Find Sensors

The first task is to identify the senors we wish to use. This can be done by loading the default 'pas' object from the Mazama Science maintained data archvie, zooming in on the Methow Valley and clicking on some of the sensors. It should be quickly apparent that many of those sensors have a label containing "MV Clean Air Ambassador".

# AirSensor package
library(AirSensor)

# Set the archiveBaseUrl so we can get a 'pas' object
setArchiveBaseUrl("https://airsensor.aqmd.gov/PurpleAir/v1")

# Load the default 'pas' (today for the entire US)
pas <- pas_load()

# Subset by state
wa <- pas_filter(pas, stateCode == "WA")

# Look at it
pas_leaflet(wa)

With this information, we can now create a script to create 'pat' objects for each of the Methow Valley Clean Air Ambassador (MVCAA) sensors.

NOTE: We could use the pas_filterArea() function to define a bounding box and get all sensors in the area. But for this tutorial we will limit ourselves to those labeled with "MV Clean Air Ambassador".

R Script

The following R script will take several minutes to run and will create 'pas' and 'pat' data files on your computer.

After running the script, a final section will demonstrate how to load and work with these local data files.

This R script can be used as a starting point for those wishing to create small collections of data for other communities and other dates.

# Methow Valley local data archive

# ----- Setup ------------------------------------------------------------------
# Create an archive directory underneath ~/Data
wd <- getwd()
archiveDir <- "Data/MVCAA"
dir.create(file.path(wd,archiveDir), recursive = TRUE)
archiveDir <- paste0(wd, "/", archiveDir)

# AirSensor package
library(AirSensor)

# Set the archiveBaseUrl so we can get a 'pas' object
setArchiveBaseUrl("https://airsensor.aqmd.gov/PurpleAir/v1")

# ----- Create PAS object ------------------------------------------------------

# Create a 'pas' object limited to MVCAA sensors
#   - load most recent 'pas' for the entire country
#   - subset to include sensors labeled MVCAA
mvcaa <-
  pas_load() %>%
  pas_filter(stringr::str_detect(label, "MV Clean Air Ambassador"))
  # NOTE: Could have filtered by area with:
  #pas_filterArea(w = -120.5, e = -120.0, s = 48.0, n = 49.0)

# Look at it
pas_leaflet(mvcaa)

# Save it in our archive directory
save(mvcaa, file = paste0(archiveDir, "/mvcaa.rda"))

# ----- Create PAT objects -----------------------------------------------------

# Get all the deviceDeploymentIDs
mvcaa_ids <- pas_getDeviceDeploymentIDs(mvcaa)

# Specify times
startdate <- "2020-09-01"
enddate <- "2020-10-01"
timezone <- "America/Los_Angeles"

# Create an empty List to store things
patList <- list()

# Loop over ids and get data (might take a while). Start with just a few to test.
mvcaa_ids_test <- c("ab5dca99422f2c0d_13669", "f6c44edd41c941c7_10182")

#---- Initialize counters first!
idCount<- length(mvcaa_ids_test)
count <- 0
successCount <- 0

for (id in mvcaa_ids_test[1:idCount]) {

  count <- count + 1
  print(sprintf("Working on %s (%d/%d) ...", id, count, idCount))

  # Use a try-block in case you get "no data" errors
  result <- try({

    # It's nice to copy-paste the full function signature so you can see all possible arguments
    patList[[id]] <- pat_createNew(
      id = id,
      label = NULL,        # not needed if you have the id
      pas = mvcaa,
      startdate = startdate,
      enddate = enddate,
      timezone = timezone,
      baseUrl = "https://api.thingspeak.com/channels/",
      verbose = FALSE
    )
    successCount <- successCount + 1

  }, silent = FALSE)

  if ( "try-error" %in% class(result) ) {
    print(geterrormessage())
  }

}

# How many did we get?
print(sprintf("Successfully created %d/%d pat objects.", successCount, idCount))


# Loop over all ids and get data (might take a while).
idCount <- length(mvcaa_ids)
count <- 0 
successCount <- 0

for (id in mvcaa_ids[1:idCount]) {

  count <- count + 1
  print(sprintf("Working on %s (%d/%d) ...", id, count, idCount))

  # Use a try-block in case you get "no data" errors
  result <- try({

    # It's nice to copy-paste the full function signature so you can see all possible arguments
    patList[[id]] <- pat_createNew(
      id = id,
      label = NULL,        # not needed if you have the id
      pas = mvcaa,
      startdate = startdate,
      enddate = enddate,
      timezone = timezone,
      baseUrl = "https://api.thingspeak.com/channels/",
      verbose = FALSE
    )
    successCount <- successCount + 1

  }, silent = FALSE)

  if ( "try-error" %in% class(result) ) {
    print(geterrormessage())
  }

}


# How many did we get?
print(sprintf("Successfully created %d/%d pat objects.", successCount, idCount))

# explore patList: each sensor has metadata and data that can be explored simply by saving them as objects. Have a look to the tables in the R environment. 

meta <- patList$ab5dca99422f2c0d_13669$meta
data <- patList$ab5dca99422f2c0d_13669$data

# Save it in our archive directory
save(patList, file = file.path(archiveDir, "patList.rda"))

# ----- Evaluate patList -------------------------------------------------------

# We can use sapply() to apply a function to each element of the list
sapply(patList, function(x) { return(x$meta$label) })

# How big is patList in memory?
print(object.size(patList), units="MB")
# 37.9 Mb

# How big patList.rda on disk (as compressed binary) 
fileSize <- file.size(file.path(archiveDir, "patList.rda"))
sprintf("%.1f Mb", fileSize/1e6)
# 7.3 MB

Loading Local Data

Now we have a local copy of all of the 'pas' and 'pat' data we wanted, precompiled into 'pas' and 'pat' objects.

# Empty current environment to ensure we're using our local archive
rm(list = ls())

# AirSensor package
library(AirSensor)

# Examine archive directory
list.files("~/Data/mvcaa")

# Load files
mvcaa <- get(load("~/Data/mvcaa/mvcaa.rda"))
patList <- get(load("~/Data/mvcaa/patList.rda"))

# Interactive map
pas_leaflet(mvcaa)

# Print site names and associated ids
sapply(patList, function(x) { return(x$meta$label) })

# Pull out Balky Hill as a separate 'pat' object
Balky_Hill <- patList[["ab5dca99422f2c0d_13669"]]

# Basic plot
pat_multiplot(Balky_Hill)

Best of luck assessing air quality in your community!



MazamaScience/AirSensor documentation built on April 28, 2023, 11:16 a.m.