In patperu/AriaToscana: Data About Air Pollutants from ARPAT Toscana

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.retina = 2,
  fig.height = 8.27, 
  fig.width = 11.69,
  dev=c('png')
)

library('ggplot2')
library('tidyr')
library('dplyr')
library('feather')
library('AriaToscana')

options(stringsAsFactors = FALSE,
        knitr.table.format = 'markdown')

addmargins(table(at_inq$parameter, at_inq$year)) %>% 
  knitr::kable(format.args = list(big.mark = '.', decimal.mark = ','))

addmargins(table(at_pm$parameter, at_pm$year)) %>% 
  knitr::kable(format.args = list(big.mark = '.', decimal.mark = ','))

# If 4 or more measures are 'NA', then 'NA' for the whole time frame (br)
ff <- function(L) ifelse(sum(is.na(L)) >= 4, NA, mean(L, na.rm = TRUE))

x1 <- tbl_df(at_inq) %>%
       filter(year %in% c(2008, 2012, 2014, 2015) & parameter %in% c("O3", "NO2", "CO")) %>%
       mutate(br = factor(cut(hour, breaks = 4), 
                          labels = c("00-06h", "06-12h", 
                                     "12-18h", "18-24h")),
              provincia = substr(stazione, 1, 2)) %>%
       group_by(provincia, stazione, year, month, day, parameter, br) %>%
       summarise(val_na   = sum(is.na(value)),
                 val_mean = ff(value))

p1 <- ggplot(subset(x1, subset = year == 2015 & parameter == "NO2" & provincia == "FI"), 
             aes(as.factor(month), val_mean, group = year, color = year))
p1 <- p1 + geom_point(shape = 1, col ="blue") 
p1 <- p1 + facet_grid(br ~ stazione) 
p1 <- p1 + stat_smooth(aes(group = 1), col = "red")
p1 <- p1 + labs(title = "Biossido di azoto NO2, 2015, Firenze (FI)")
p1 <- p1 + labs(x = "Mese", y = expression(paste(Biossido~di~azoto~(NO2)~µg/m^3)))
p1

p2 <- ggplot(subset(x1, subset = 
                         stazione %in% c("LU-CAPANNORI", "FI-MOSSE", 
                                         "FI-GRAMSCI", "PO-ROMA") & 
                         parameter == "NO2"), aes(as.factor(month), val_mean))
p2 <- p2 + facet_grid(br ~ stazione)
p2 <- p2 + stat_smooth(aes(group = as.factor(year), color = as.factor(year)), se = FALSE)
p2 <- p2 + labs(title = "Biossido di azoto NO2, 2008/2012/2014/2015")
p2 <- p2 + labs(x = "Mese", 
                y = expression(paste(Biossido~di~azoto~(NO2)~µg/m^3)))
p2

PM10

This is table 4.1.2. on page 18 from the last report 2015

m <- tbl_df(at_pm) %>%
       filter(parameter == "PM10" & valid == 1) %>%
       group_by(stazione, year) %>%
       mutate(days_gt50 = cumsum(value > 50)) %>%
       select(stazione, year, days_gt50) %>%
       top_n(1, days_gt50) # %>%
       #distinct(stazione, year) %>%
       #ungroup() %>%
       #arrange(stazione, year)

spread(unique(m), year, days_gt50, fill = ".") %>% knitr::kable()

This is table 4.1.5. on page 24 from the last report 2015

m <- tbl_df(at_pm) %>%
       filter(parameter == "PM10" & valid == 1) %>%
       group_by(stazione, year) %>%
       summarise(val_mean = round(mean(value), 0)) %>%
       select(stazione, year, val_mean) %>%
       ungroup() %>%
       arrange(year, stazione)

spread(m, year, val_mean, fill = ".") %>% knitr::kable()

Percent of missing observations per month

x3 <- tbl_df(x1) %>%
       group_by(provincia, stazione, year, month, parameter) %>%
        summarise(pct_na = round((sum(val_na) / (max(day) * 24)) * 100, 1) )

filter(x3, year == 2015 & 
           stazione == "FI-GRAMSCI" & 
           parameter == "NO2") %>% knitr::kable()

Missing PM10

Percent of missing values of PM10 by year and month

x4 <- tbl_df(at_pm) %>%
       group_by(stazione, year, month, parameter) %>%
       summarise(pct_na = round((sum(is.na(value)) /  (max(day))) * 100, 1) ) %>%
       filter(stazione == "FI-GRAMSCI" & parameter == "PM10") %>%
       ungroup() %>%
       select(year, month, pct_na)

spread(x4, year, pct_na) %>% knitr::kable()

Write/Read 'at_inq' with 'feather'

This was running on a T410 with 8GB memory and a HDD....

str(at_inq)

path <- "at_inq.feather"

system.time(write_feather(at_inq, path))

utils:::format.object_size(file.size(path), units = "Mb")

system.time(x <- read_feather(path))

str(x)

identical(x, at_inq)