Load AQI data:

Load the PM2.5 air quality data in Buffalo in 2016: This data comes from EPA Air data that is mentioned in Data above. Some functions in the packages can only recognize “date” as the column name so we first need to change the name.

mytable=read_excel("Data/Buffalo_PM2.5_AQI.xlsx")
names(mytable)[1]<-"date"

According to the method mentioned in reference[1], I need to set some parameters to extract daily and monthly factors to make the heatmap plots below more easier and understandable to be shown:

dat=mytable
dat$month<-as.numeric(as.POSIXlt(dat$date)$mon+1)
dat$monthf<-factor(dat$month,levels=as.character(1:12),labels=c("Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct","Nov","Dec"),ordered=TRUE)
dat$weekday<-as.POSIXlt(dat$date)$wday
dat$weekdayf<-factor(dat$weekday,levels=rev(0:6),labels=rev(c("Sun","Mon","Tue","Wed","Thu","Fri","Sat")),ordered=TRUE)
dat$week <- as.numeric(format(dat$date,"%W"))
dat<-ddply(dat,.(monthf),transform,monthweek=1+week-min(week))

Load Weather data:

Weather data comes from Global historical climatology network (GHCN) and we use three variables as the influential factors, which are average daily wind speed (awnd), daliy maximum temperatures (tmax) and daily minimum temperatures (tmin).

datadir="data"
st = ghcnd_stations()
st=dplyr::filter(st,element%in%c("TMAX","TMIN","AWND"))
coords=as.matrix(geocode("Buffalo, NY"))
dplyr::filter(st,
              grepl("BUFFALO",name)&
              between(latitude,coords[2]-1,coords[2]+1) &
              between(longitude,coords[1]-1,coords[1]+1)&
         element=="TMAX")

## # A tibble: 0 x 11
## # ... with 11 variables: id <chr>, latitude <dbl>, longitude <dbl>,
## #   elevation <dbl>, state <chr>, name <chr>, gsn_flag <chr>,
## #   wmo_id <chr>, element <chr>, first_year <int>, last_year <int>

d=meteo_tidy_ghcnd("USW00014733",
                   var = c("TMAX","TMIN","AWND"),
                   keep_flags=T)

To make results accurate we need to convert correct units of temperatures that is converting it to degrees C.

d_filtered=d%>%
  mutate(tmax=ifelse(qflag_tmax!=" "|tmax==-9999,NA,tmax/10))%>% 
  mutate(tmin=ifelse(qflag_tmin!=" "|tmin==-9999,NA,tmin/10))%>%  
  mutate(awnd=ifelse(qflag_tmin!=" "|awnd==-9999,NA,awnd))%>%  
  arrange(date)

Because we use AQI data in 2016, we need to filter wealth data in 2016. Also, the datatable function is a dynamic visualization method that can show selected data.

d_filtered_2016=filter(d_filtered,date>=as.Date("2016-01-01")&date<=as.Date("2016-12-31"))

datatable(d_filtered_2016, options = list(pageLength=5)) %>% frameWidget(height=400)

Calendar heatmap:

This can take a few seconds to show the plot about Calendar Heatmap of PM2.5 air quality in Buffalo. In this ggplot functions, I use different colours to show different AQI values in which yellow represents low value and red represents high value. Additionally, in Air Quality Index(AQI) leves, it is good day when AQI value is less than 50, and more than 50 but less than 100 means moderate, and more than 100 but less than 150 means unhealthy for sensitive groups. These information can all be shown in this calendar heapmap figure.

ggplot(dat,aes(monthweek, weekdayf,fill=DAILY_AQI_VALUE))+
  
  geom_tile(colour='white') +
  
  facet_wrap(~monthf ,nrow=3) +
  
  scale_fill_gradient(space="Lab",limits=c(0, max(dat$DAILY_AQI_VALUE)),low="yellow", high=" red") +
  
  labs(title="Buffalo PM2.5 AQI Calendar Heatmap in 2016",x="Week of month",y="")

Dygraphics:

Dygraph function provides rich facilities for charting time-series data in R. I use the same AQI data in 2016 in Buffalo by dygraph to compare results with calendat heatmaps.

dt=xts(dat$DAILY_AQI_VALUE,order.by=dat$date)

dygraph(dt, main = "Daily AQI value in Buffalo, NY") %>% 
  dyRangeSelector(dateWindow = c("2016-01-01", "2016-12-31"))%>%
  frameWidget(height =300)

Through the two graphics above, it is clearly to find that calendar heatmap is better for readers to get more inofrmation, such as how many unhealthy days in 2016? or which days are unhealthy days? In dygraphics people can see the trend of the changing about AQI values but it does not provide with daily data visually.