pm_spectrometer_tests_20160223.R

# PM2.5 spectrometer tests
library('ggplot2')
library('openair')

## Loading required package: lazyeval
## Loading required package: 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
## 
## Loading required package: maps

library('lubridate')
require('corrplot')

## Loading required package: corrplot

# PMS1003 data
pms1003.data <- read.delim("/home/gustavo/data/Dust_tests/minispectrometer/WithArduino/20160223_pms1003.txt")
pms1003.data$date <- as.POSIXct(paste(pms1003.data$Day,pms1003.data$Time),tz='Pacific/Auckland')
pms1003.data$Day<-NULL
pms1003.data$Time<-NULL
# Number concentrations to pt/cc (raw comes as pt/100cc)
pms1003.data$N300 <-pms1003.data$N300 / 100
pms1003.data$N500 <-pms1003.data$N500 / 100
pms1003.data$N1000 <-pms1003.data$N1000 / 100
pms1003.data$N2500 <-pms1003.data$N2500 / 100
pms1003.data$N5000 <-pms1003.data$N5000 / 100
pms1003.data$N10000 <-pms1003.data$N10000 / 100

# GRIMM data
grimm <- read.csv("/home/gustavo/data/Dust_tests/minispectrometer/WithArduino/grimm_20160223.parsed", sep="")
grimm$date <- ISOdatetime(year = grimm$Year+2000,
                          month = grimm$Month, 
                          day = grimm$Day,
                          hour = grimm$Hour,
                          min = grimm$Minute, 
                          sec = grimm$Second, tz = "GMT")
# Convert counts data into pt/cc (raw comes as pt/1000cc)
grimm[,17:47] <- grimm[,17:47] / 1000
merged_data <- merge(pms1003.data,grimm,by = 'date', all = TRUE)



# Correlation calculations to find Shinyei's size information
for_correl <- timeAverage(merged_data,avg.time = '30 min')

## Warning in checkPrep(mydata, vars, type = "default", remove.calm = FALSE, :
## Detected data with Daylight Saving Time, converting to UTC/GMT

find_channels <- cor(for_correl[,c(2:10,29:59)],use = 'pairwise')

corrplot(find_channels,method = 'circle',diag=FALSE)

# Plotting

ggplot(pms1003.data,aes(x=date))+
  geom_line(aes(y = PM1),colour = 'red')+
  geom_line(aes(y = PM2.5),colour = 'blue')+
  geom_line(aes(y = PM10),colour = 'black')

ggplot(pms1003.data,aes(x=date))+
  geom_line(aes(y = N300),colour = 'red')+
  geom_line(aes(y = N500),colour = 'blue')+
  geom_line(aes(y = N1000),colour = 'black')+
  geom_line(aes(y = N2500),colour = 'green')

timePlot(merged_data,pollutant = c('N300','n265','n290','n324','N1000'),avg.time = '60 min')

## Warning in checkPrep(mydata, vars, type, remove.calm = FALSE): Detected
## data with Daylight Saving Time, converting to UTC/GMT

scatterPlot(merged_data,x="n265",y="N300",avg.time = '60 min')

## Warning in checkPrep(mydata, vars, type = "default", remove.calm = FALSE, :
## Detected data with Daylight Saving Time, converting to UTC/GMT

# Concentrations seem to be off by a factor of 10 so multiplying N300 by 10 we obtain
merged_data$N300_scale <- merged_data$N300*10
timePlot(merged_data,pollutant = c('N300_scale','n265'),avg.time = '60 min',group = TRUE)

## Warning in checkPrep(mydata, vars, type, remove.calm = FALSE): Detected
## data with Daylight Saving Time, converting to UTC/GMT

# Attempting a linear fit of N300 based on the 60min data
summary(n300_linear <- lm(data = timeAverage(merged_data,avg.time = '60 min'),n265~N300))

## Warning in checkPrep(mydata, vars, type = "default", remove.calm = FALSE, :
## Detected data with Daylight Saving Time, converting to UTC/GMT

## 
## Call:
## lm(formula = n265 ~ N300, data = timeAverage(merged_data, avg.time = "60 min"))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.1279 -0.8516 -0.1143  0.7120  6.1791 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   0.1498     0.3495   0.429    0.669    
## N300         10.0751     0.2499  40.316   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.7 on 92 degrees of freedom
## Multiple R-squared:  0.9464, Adjusted R-squared:  0.9458 
## F-statistic:  1625 on 1 and 92 DF,  p-value: < 2.2e-16

merged_data$N300_corr <- predict(n300_linear,newdata = merged_data)
timePlot(merged_data,pollutant = c('N300_corr','n265'),avg.time = '60 min',group = TRUE)

## Warning in checkPrep(mydata, vars, type, remove.calm = FALSE): Detected
## data with Daylight Saving Time, converting to UTC/GMT