Night Sky Weather
Ralph Rogge (RRO)
11/6/2016
Output directory name
output.dir <- "Output/"Additional R packages
library(data.table)
library(dplyr)
library(ggplot2)
library(lubridate)
library(scales)## Warning: package 'scales' was built under R version 3.3.2Define weather station parameters
station <- list(name="02712")Data pool used for cloudiness, precipitaion, and humidity.
data = list()Cloudiness
Read data from file
data$cloudiness <- read.csv("Analytical Data/Cloudiness 02712.csv", as.is = TRUE)
nrow(data$cloudiness)## [1] 379883Reorder day part
data$cloudiness$daypart <- factor(data$cloudiness$daypart, levels = c("Daytime", "Twilight", "Nighttime"))Add season
data$cloudiness$season <- season(data$cloudiness$month)Calculate cloudiness by season
df <- data$cloudiness %>%
filter(daypart != "Twilight") %>%
group_by(season, daypart) %>%
summarize(value = 12.5* mean(cloudiness, na.rm = TRUE))Differences of diurnal cloudiness
data.table(diff.diurnal(df))## season daypart1 value1 daypart2 value2 delta
## 1: Winter Daytime 80.94934 Nighttime 77.79297 3.156370
## 2: Spring Daytime 68.24326 Nighttime 61.67766 6.565601
## 3: Summer Daytime 60.44897 Nighttime 54.84742 5.601547
## 4: Autumn Daytime 73.37257 Nighttime 68.40463 4.967943Differences of seasonal cloudiness
data.table(diff.seasonal(df))## daypart season1 value1 season2 value2 delta
## 1: Daytime Winter 80.94934 Spring 68.24326 12.706078
## 2: Nighttime Winter 77.79297 Spring 61.67766 16.115309
## 3: Daytime Spring 68.24326 Summer 60.44897 7.794295
## 4: Nighttime Spring 61.67766 Summer 54.84742 6.830241
## 5: Daytime Summer 60.44897 Autumn 73.37257 -12.923604
## 6: Nighttime Summer 54.84742 Autumn 68.40463 -13.557208
## 7: Daytime Autumn 73.37257 Winter 80.94934 -7.576769
## 8: Nighttime Autumn 68.40463 Winter 77.79297 -9.388342Calculate average cloud cover by day
df <- data$cloudiness %>%
# Group by month, and day.
group_by(month, day, daypart) %>%
# Average cloudiness in %.
summarize(cloudiness.mean = 12.5* mean(cloudiness, na.rm = TRUE),
cloudiness.sd = sd(cloudiness, na.rm = TRUE)) %>%
# Add date column (year value is a dummy).
mutate(date=as.Date(paste(2016, month, day, sep = "-"))) %>%
# Select columns.
ungroup() %>%
select(date, cloudiness.mean, cloudiness.sd, daypart)Plot average cloud cover by day
p <- ggplot(df) +
geom_point(aes(x = date, y = cloudiness.mean, color = daypart), alpha = 0.5) +
geom_errorbar(aes(x = date,
ymin = cloudiness.mean - cloudiness.sd,
ymax = cloudiness.mean + cloudiness.sd, color = daypart),
alpha = 0.3) +
geom_smooth(aes(x = date, y = cloudiness.mean, color = daypart),
span = 0.3,
alpha = 0.3,
method = "loess",
se = FALSE) +
scale_color_discrete(name = "") +
scale_x_date(breaks = date_breaks("3 months"),
labels = date_format("%b"),
minor_breaks = date_breaks("1 month")) +
labs(title = paste("Average Cloud Cover / Weather Station", station$name,
paste0("(",
min(data$cloudiness$year, na.rm = TRUE),
"-",
max(data$cloudiness$year, na.rm = TRUE),
")")),
x = "Month",
y = "Cloud Cover [%]") +
theme_bw()
p
Write result to files
if (!dir.exists(output.dir)) dir.create(output.dir)
write.csv(df, paste0(output.dir, "Cloudiness 02712.csv"), row.names = FALSE)
write.plot(p, paste0(output.dir, "Cloudiness 02712.png"))## quartz_off_screen
## 2Clean up
rm(df, p)Precipitation
Read data from file
data$precipitation <- read.csv("Analytical Data/Precipitation 02712.csv", as.is = TRUE)
nrow(data$precipitation)## [1] 177712Reorder day part
data$precipitation$daypart <- factor(data$precipitation$daypart, levels = c("Daytime", "Twilight", "Nighttime"))Add season
data$precipitation$season <- season(data$precipitation$month)Calculate sum of precipitation height by season
df <- data$precipitation %>%
filter(daypart != "Twilight") %>%
group_by(season, daypart) %>%
summarize(value = sum(height, na.rm = TRUE))Differences of diurnal precipitation
data.table(diff.diurnal(df))## season daypart1 value1 daypart2 value2 delta
## 1: Winter Daytime 962.4 Nighttime 1510.5 -548.1
## 2: Spring Daytime 1985.8 Nighttime 1059.8 926.0
## 3: Summer Daytime 3259.0 Nighttime 968.7 2290.3
## 4: Autumn Daytime 1611.0 Nighttime 1777.9 -166.9Differences of seasonal precipitation
data.table(diff.seasonal(df))## daypart season1 value1 season2 value2 delta
## 1: Daytime Winter 962.4 Spring 1985.8 -1023.4
## 2: Nighttime Winter 1510.5 Spring 1059.8 450.7
## 3: Daytime Spring 1985.8 Summer 3259.0 -1273.2
## 4: Nighttime Spring 1059.8 Summer 968.7 91.1
## 5: Daytime Summer 3259.0 Autumn 1611.0 1648.0
## 6: Nighttime Summer 968.7 Autumn 1777.9 -809.2
## 7: Daytime Autumn 1611.0 Winter 962.4 648.6
## 8: Nighttime Autumn 1777.9 Winter 1510.5 267.4Calculate average precipitation height by day
df <- data$precipitation %>%
# Sum precipitation height per day.
group_by(year, month, day, daypart) %>%
summarize(height = sum(height, na.rm = TRUE)) %>%
# Average precipitation height by day over all years.
group_by(month, day, daypart) %>%
summarize(height.mean = mean(height, na.rm = TRUE),
height.sd = sd(height, na.rm = TRUE)) %>%
# Add date column (year value is a dummy).
mutate(date = as.Date(paste(2016, month, day, sep = "-"))) %>%
# Select columns.
ungroup() %>%
select(date, height.mean, height.sd, daypart)Plot average precipitation height by day
p <- ggplot(df) +
geom_point(aes(x = date, y = height.mean, color = daypart), alpha = 0.5) +
geom_errorbar(aes(x = date,
ymin = ifelse(height.mean - height.sd >= 0, height.mean - height.sd, 0),
ymax = height.mean + height.sd,
color = daypart),
alpha = 0.3) +
geom_smooth(aes(x = date, y = height.mean, color = daypart),
span = 0.3,
alpha = 0.3,
method = "loess",
se = FALSE) +
facet_grid(. ~ daypart) +
scale_color_discrete(guide = FALSE) +
scale_x_date(breaks = date_breaks("3 months"),
labels = date_format("%b"),
minor_breaks = date_breaks("1 month")) +
labs(title = paste("Average Precipitation / Weather Station", station$name,
paste0("(",
min(data$precipitation$year, na.rm = TRUE),
"-",
max(data$precipitation$year, na.rm = TRUE),
")")),
x = "Month",
y = "Height [mm]") +
theme_bw()
p
Write result to files
if (!dir.exists(output.dir)) dir.create(output.dir)
write.csv(df, paste0(output.dir, "Precipitation 02712.csv"), row.names = FALSE)
write.plot(p, paste0(output.dir, "Precipitation 02712.png"))## quartz_off_screen
## 2Temperature
Read data from file
data$temperature <- read.csv("Analytical Data/Temperature 02712.csv", as.is = TRUE)
nrow(data$temperature)## [1] 394190Reorder day part
data$temperature$daypart <- factor(data$temperature$daypart, levels = c("Daytime", "Twilight", "Nighttime"))Add season
data$temperature$season <- season(data$temperature$month)Calculate average temperature by season
df <- data$temperature %>%
filter(daypart != "Twilight") %>%
group_by(season, daypart) %>%
summarize(value = mean(temperature, na.rm = TRUE))Differences of diurnal temperature
data.table(diff.diurnal(df))## season daypart1 value1 daypart2 value2 delta
## 1: Winter Daytime 1.874211 Nighttime 0.7895658 1.084645
## 2: Spring Daytime 10.950334 Nighttime 6.8197910 4.130543
## 3: Summer Daytime 19.379825 Nighttime 16.3336091 3.046216
## 4: Autumn Daytime 11.266041 Nighttime 8.2169834 3.049057Differences of seasonal temperature
data.table(diff.seasonal(df))## daypart season1 value1 season2 value2 delta
## 1: Daytime Winter 1.8742110 Spring 10.9503344 -9.076123
## 2: Nighttime Winter 0.7895658 Spring 6.8197910 -6.030225
## 3: Daytime Spring 10.9503344 Summer 19.3798252 -8.429491
## 4: Nighttime Spring 6.8197910 Summer 16.3336091 -9.513818
## 5: Daytime Summer 19.3798252 Autumn 11.2660405 8.113785
## 6: Nighttime Summer 16.3336091 Autumn 8.2169834 8.116626
## 7: Daytime Autumn 11.2660405 Winter 1.8742110 9.391829
## 8: Nighttime Autumn 8.2169834 Winter 0.7895658 7.427418Calculate average temperature by day
df <- data$temperature %>%
# Group by month, and day.
group_by(month, day, daypart) %>%
# Average temperature.
summarize(temperature.mean = mean(temperature, na.rm = TRUE),
temperature.sd = sd(temperature, na.rm = TRUE)) %>%
# Add date column (year value is a dummy).
mutate(date=as.Date(paste(2016, month, day, sep = "-"))) %>%
# Select columns.
ungroup() %>%
select(date, temperature.mean, temperature.sd, daypart)Plot average temperature by day
p <- ggplot(df) +
geom_point(aes(x = date, y = temperature.mean, color = daypart), alpha = 0.5) +
geom_errorbar(aes(x= date,
ymin = temperature.mean - temperature.sd,
ymax = temperature.mean + temperature.sd,
color = daypart),
alpha = 0.3) +
geom_smooth(aes(x = date, y = temperature.mean, color = daypart),
span = 0.3,
alpha = 0.3,
method = "loess",
se = FALSE) +
facet_grid(. ~ daypart) +
scale_color_discrete(guide = FALSE) +
scale_x_date(breaks = date_breaks("3 months"),
labels = date_format("%b"),
minor_breaks = date_breaks("1 month")) +
labs(title = paste("Average Temperature / Weather Station", station$name,
paste0("(",
min(data$temperature$year, na.rm = TRUE),
"-",
max(data$temperature$year, na.rm = TRUE),
")")),
x = "Month",
y = "Temperature [°C]") +
theme_bw()
p
Write result to files
if (!dir.exists(output.dir)) dir.create(output.dir)
write.csv(df, paste0(output.dir, "Temperature 02712.csv"), row.names = FALSE)
write.plot(p, paste0(output.dir, "Temperature 02712.png"))## quartz_off_screen
## 2Humidity
Read data from file
data$humidity <- read.csv("Analytical Data/Humidity 02712.csv", as.is = TRUE)
nrow(data$humidity)## [1] 394190Reorder day part
data$humidity$daypart <- factor(data$humidity$daypart, levels = c("Daytime", "Twilight", "Nighttime"))Add season
data$humidity$season <- season(data$humidity$month)Calculate average humidity by season
df <- data$humidity %>%
filter(daypart != "Twilight") %>%
group_by(season, daypart) %>%
summarize(value = mean(humidity, na.rm = TRUE))Differences of diurnal humidity
data.table(diff.diurnal(df))## season daypart1 value1 daypart2 value2 delta
## 1: Winter Daytime 82.02312 Nighttime 87.88503 -5.861908
## 2: Spring Daytime 67.03168 Nighttime 79.35380 -12.322118
## 3: Summer Daytime 67.51415 Nighttime 82.62562 -15.111467
## 4: Autumn Daytime 79.27744 Nighttime 89.31034 -10.032907Differences of seasonal humidity
data.table(diff.seasonal(df))## daypart season1 value1 season2 value2 delta
## 1: Daytime Winter 82.02312 Spring 67.03168 14.9914363
## 2: Nighttime Winter 87.88503 Spring 79.35380 8.5312272
## 3: Daytime Spring 67.03168 Summer 67.51415 -0.4824685
## 4: Nighttime Spring 79.35380 Summer 82.62562 -3.2718175
## 5: Daytime Summer 67.51415 Autumn 79.27744 -11.7632886
## 6: Nighttime Summer 82.62562 Autumn 89.31034 -6.6847294
## 7: Daytime Autumn 79.27744 Winter 82.02312 -2.7456793
## 8: Nighttime Autumn 89.31034 Winter 87.88503 1.4253197Calculate aerage humidity by day
df <- data$humidity %>%
# Group by month, and day.
group_by(month, day, daypart) %>%
# Average humidity.
summarize(humidity.mean = mean(humidity, na.rm = TRUE),
humidity.sd = sd(humidity, na.rm = TRUE)) %>%
# Add date column (year value is a dummy).
mutate(date=as.Date(paste(2016, month, day, sep = "-"))) %>%
# Select columns.
ungroup() %>%
select(date, humidity.mean, humidity.sd, daypart)Plot average humidity by day of year
p <- ggplot(df) +
geom_point(aes(x = date, y = humidity.mean, color = daypart), alpha = 0.5) +
geom_errorbar(aes(x = date,
ymin = humidity.mean - humidity.sd,
ymax = humidity.mean + humidity.sd,
color = daypart),
alpha = 0.3) +
geom_smooth(aes(x = date, y = humidity.mean, color = daypart),
span = 0.3,
alpha = 0.3,
method = "loess",
se = FALSE) +
facet_grid(. ~ daypart) +
scale_color_discrete(guide = FALSE) +
scale_x_date(breaks = date_breaks("3 months"),
labels = date_format("%b"),
minor_breaks = date_breaks("1 month")) +
labs(title = paste("Average Humidity / Weather Station", station$name,
paste0("(",
min(data$humidity$year, na.rm = TRUE),
"-",
max(data$humidity$year, na.rm = TRUE),
")")),
x = "Month",
y = "Humidity [%]") +
theme_bw()
p
Write result to files
if (!dir.exists(output.dir)) dir.create(output.dir)
write.csv(df, paste0(output.dir, "Humidity 02712.csv"), row.names = FALSE)
write.plot(p, paste0(output.dir, "Humidity 02712.png"))## quartz_off_screen
## 2Wind
Read data from file
data$wind <- read.csv("Analytical Data/Wind 02712.csv", as.is = TRUE)
nrow(data$wind)## [1] 494772Reorder day part
data$wind$daypart <- factor(data$wind$daypart, levels = c("Daytime", "Twilight", "Nighttime"))Add season
data$wind$season <- season(data$wind$month)Calculate average wind by season
df <- data$wind %>%
filter(daypart != "Twilight") %>%
group_by(season, daypart) %>%
summarize(value = mean(speed, na.rm = TRUE))Differences of diurnal wind
data.table(diff.diurnal(df))## season daypart1 value1 daypart2 value2 delta
## 1: Winter Daytime 2.361793 Nighttime 2.152585 0.2092075
## 2: Spring Daytime 2.373692 Nighttime 1.910112 0.4635804
## 3: Summer Daytime 2.031986 Nighttime 1.497790 0.5341961
## 4: Autumn Daytime 1.919832 Nighttime 1.682902 0.2369304Differences of seasonal wind
data.table(diff.seasonal(df))## daypart season1 value1 season2 value2 delta
## 1: Daytime Winter 2.361793 Spring 2.373692 -0.01189947
## 2: Nighttime Winter 2.152585 Spring 1.910112 0.24247341
## 3: Daytime Spring 2.373692 Summer 2.031986 0.34170611
## 4: Nighttime Spring 1.910112 Summer 1.497790 0.41232183
## 5: Daytime Summer 2.031986 Autumn 1.919832 0.11215355
## 6: Nighttime Summer 1.497790 Autumn 1.682902 -0.18511222
## 7: Daytime Autumn 1.919832 Winter 2.361793 -0.44196019
## 8: Nighttime Autumn 1.682902 Winter 2.152585 -0.46968302Calculate average wind speed by day
df <- data$wind %>%
# Group by month, and day.
group_by(month, day, daypart) %>%
# Average wind
summarize(speed.mean = mean(speed, na.rm = TRUE),
speed.sd = sd(speed, na.rm = TRUE)) %>%
# Add date column (year value is a dummy).
mutate(date=as.Date(paste(2016, month, day, sep = "-"))) %>%
# Select columns.
ungroup() %>%
select(date, speed.mean, speed.sd, daypart)Plot average wind speed by day
p <- ggplot(df) +
geom_point(aes(x = date, y = speed.mean, color = daypart), alpha = 0.5) +
geom_errorbar(aes(x = date,
ymin = speed.mean - speed.sd,
ymax = speed.mean + speed.sd,
color = daypart),
alpha = 0.3) +
geom_smooth(aes(x = date, y = speed.mean, color = daypart),
span = 0.3,
alpha = 0.3,
method = "loess",
se = FALSE) +
facet_grid(. ~ daypart) +
scale_color_discrete(guide = FALSE) +
scale_x_date(breaks = date_breaks("3 months"),
labels = date_format("%b"),
minor_breaks = date_breaks("1 month")) +
labs(title = paste("Average Wind Speed / Weather Station", station$name,
paste0("(",
min(data$wind$year, na.rm = TRUE),
"-",
max(data$wind$year, na.rm = TRUE),
")")),
x = "Month",
y = "Wind Speed [m/s]") +
theme_bw()
p
Write result to files
if (!dir.exists(output.dir)) dir.create(output.dir)
write.csv(df, paste0(output.dir, "Wind 02712.csv"), row.names = FALSE)
write.plot(p, paste0(output.dir, "Wind 02712.png"))## quartz_off_screen
## 2Day Part Length
Read data from file
data$daypartlength <- read.csv("Analytical Data/Day Part Length 02712.csv", as.is = TRUE)
nrow(data$daypartlength)## [1] 8760Reorder day part
data$daypartlength$daypart <- factor(data$daypartlength$daypart, levels = c("Daytime", "Twilight", "Nighttime"))Add season
data$daypartlength$season <- season(data$daypartlength$month)Calculate average day part length by season
df <- data$daypartlength %>%
filter(daypart != "Twilight") %>%
group_by(season, month, day, daypart) %>%
summarize(length = n()) %>%
group_by(season, daypart) %>%
summarize(value = mean(length))Differences of diurnal day part length
data.table(diff.diurnal(df))## season daypart1 value1 daypart2 value2 delta
## 1: Winter Daytime 8.833333 Nighttime 11.333333 -2.500000
## 2: Spring Daytime 13.434783 Nighttime 6.402174 7.032609
## 3: Summer Daytime 15.228261 Nighttime 3.532609 11.695652
## 4: Autumn Daytime 10.725275 Nighttime 9.560440 1.164835Differences of seasonal day part length
data.table(diff.seasonal(df))## daypart season1 value1 season2 value2 delta
## 1: Daytime Winter 8.833333 Spring 13.434783 -4.601449
## 2: Nighttime Winter 11.333333 Spring 6.402174 4.931159
## 3: Daytime Spring 13.434783 Summer 15.228261 -1.793478
## 4: Nighttime Spring 6.402174 Summer 3.532609 2.869565
## 5: Daytime Summer 15.228261 Autumn 10.725275 4.502986
## 6: Nighttime Summer 3.532609 Autumn 9.560440 -6.027831
## 7: Daytime Autumn 10.725275 Winter 8.833333 1.891941
## 8: Nighttime Autumn 9.560440 Winter 11.333333 -1.772894Calculate average day part length by month
df <- data$daypartlength %>%
# Group by month, and day.
group_by(month, day, daypart) %>%
# Calculate day part length by day.
summarize(length = n()) %>%
# Group by month.
group_by(month, daypart) %>%
# Calcluate average day part length by month.
summarize(length = mean(length)) %>%
# Add date column (year value is a dummy).
mutate(date=as.Date(paste(2016, month, 1, sep = "-"))) %>%
# Select columns.
ungroup() %>%
select(date, length, daypart)Plot average day part length by month
p <- ggplot(df) +
geom_bar(aes(x = date, y = length, fill = daypart), alpha = 0.5, stat = "identity", position="dodge") +
facet_grid(. ~ daypart) +
scale_fill_discrete(name = "") +
scale_x_date(breaks = date_breaks("3 months"),
labels = date_format("%b"),
minor_breaks = date_breaks("1 month")) +
labs(title = paste("Average Day Part Length / Weather Station",
station$name,
"(2015)"),
x = "Month",
y = "Day Part Length [h]") +
theme_bw()
p
Write result to files
if (!dir.exists(output.dir)) dir.create(output.dir)
write.csv(df, paste0(output.dir, "Day Part Length 02712.csv"), row.names = FALSE)
write.plot(p, paste0(output.dir, "Day Part Length 02712.png"))## quartz_off_screen
## 2Clean up
rm(df, p)Clear Nights 2015
clear.sky = list()Cloudiness
cloudiness.threshold <- 2 # 2 * 12.5% = 25%Nightly hours with cloud cover of 25% or less
clear.sky$cloudiness <- data$cloudiness %>%
filter(year == 2015) %>%
filter(daypart == "Nighttime") %>%
filter(cloudiness <= cloudiness.threshold) %>%
select(-quality, -daypart)
nrow(clear.sky$cloudiness)## [1] 964Precipitation
Nightly hours without percipitation
clear.sky$precipitation <- data$precipitation %>%
filter(year == 2015) %>%
filter(daypart == "Nighttime") %>%
filter(indicator == FALSE) %>%
select(-height, -quality, -daypart)
nrow(clear.sky$precipitation)## [1] 2175Humidity
humidity.threshold <- 85Nightly hours with humidity of 85% or less
clear.sky$humidity <- data$humidity %>%
filter(year == 2015) %>%
filter(daypart == "Nighttime") %>%
filter(humidity <= humidity.threshold) %>%
select(-quality, -daypart)
nrow(clear.sky$humidity)## [1] 1292Clear Sky
I want at least 3 hours of clear sky per night
length.threshold <- 3Nights without precipitations and a cloud cover of 25% or less and a humidity of 85% or less for at least 3 hours
df <- merge(merge(clear.sky$cloudiness, clear.sky$precipitation), clear.sky$humidity) %>%
group_by(year, month, day) %>%
summarize(length = n()) %>%
filter(length >= length.threshold)
nrow(df)## [1] 83Calculate hours of clear sky by month
df<- df %>% group_by(year, month) %>%
summarize(length = sum(length)) %>%
mutate(date = as.Date(paste(year, month, 1, sep = "-"))) %>%
ungroup %>%
select(date, length)
months <- seq(1,12)
missing <- months[!(months %in% (as.POSIXlt(df$date)$mon+1))]
if (length(missing > 0)) {
df <- rbind(df, data.frame(date=as.Date(paste(2015, missing, 1, sep = "-")), length=c(0)))
df <- df[order(df$date),]
}
data.table(df %>% mutate(month=strftime(date,"%b")))## date length month
## 1: 2015-01-01 11 Jan
## 2: 2015-02-01 31 Feb
## 3: 2015-03-01 98 Mar
## 4: 2015-04-01 60 Apr
## 5: 2015-05-01 7 May
## 6: 2015-06-01 0 Jun
## 7: 2015-07-01 17 Jul
## 8: 2015-08-01 69 Aug
## 9: 2015-09-01 66 Sep
## 10: 2015-10-01 18 Oct
## 11: 2015-11-01 42 Nov
## 12: 2015-12-01 8 DecPlot hours of clear sky by month
p <- ggplot(df) +
geom_bar(aes(x = date, y = length), stat = "identity") +
scale_x_date(breaks = date_breaks("3 months"),
labels = date_format("%b"),
minor_breaks = date_breaks("1 month")) +
labs(title = paste("Clear Sky / Weather Station", station$name, "(2015)"),
x = "Month",
y = "Clear Sky [h]") +
theme_bw()
p
Write result to files
if (!dir.exists(output.dir)) dir.create(output.dir)
write.csv(df, paste0(output.dir, "Clear Sky 02712.csv"), row.names = FALSE)
write.plot(p, paste0(output.dir, "Clear Sky 02712.png"))## quartz_off_screen
## 2Appendix
diff.diurnal <- function(data) {
dayparts <- c("Daytime", "Nighttime")
seasons <- c("Winter", "Spring", "Summer", "Autumn")
diffs <- data.frame(character(), character(), numeric(), numeric(), character(), numeric())
for (i in seq(1, length(dayparts)-1)) {
for (j in seq(1, length(seasons))) {
x <- data %>% filter(daypart == dayparts[i]) %>% filter(season == seasons[j])
y <- data %>% filter(daypart == dayparts[i+1]) %>% filter(season == seasons[j])
diffs <- rbind(diffs,data.frame(season = x$season,
daypart1 = x$daypart,
value1 = x$value,
daypart2 = y$daypart,
value2 = y$value,
delta = x$value - y$value))
}
}
diffs
}diff.seasonal <- function(data) {
diffs <- data.frame(character(), character(), numeric(), numeric(), character(), numeric())
seasons <- c("Winter", "Spring", "Summer", "Autumn", "Winter")
dayparts <- c("Daytime", "Nighttime")
for (i in seq(1, length(seasons)-1)) {
for (j in seq(1, length(dayparts))) {
x <- data %>% filter(season == seasons[i]) %>% filter(daypart == dayparts[j])
y <- data %>% filter(season == seasons[i+1]) %>% filter(daypart == dayparts[j])
diffs <- rbind(diffs,data.frame(daypart = x$daypart,
season1 = x$season,
value1 = x$value,
season2 = y$season,
value2 = y$value,
delta = x$value - y$value))
}
}
diffs
}season <- function(x) {
y <- cut(x, breaks=c(-1,2,5,8,11,12), labels=c("A","B","C","D","E"))
levels(y) <- list("Winter"=c("A", "E"), "Spring"="B", "Summer"="C", "Autumn"="D")
y
}write.plot <- function(p, filename, width = 1024, height = 576, font.size = 16) {
png(filename,width, height)
print(p + theme_bw(base_size=font.size))
dev.off()
}