Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
downloads Holly Springs data and plots it, specifically looking at air temperature, soil 10 cm temp, and soil 10 cm H2O
# get hourly Holly Springs and plot it
library("ggplot2")
library("lubridate")
library("dplyr")
library("cowplot")
library("reshape2")
hollyHourly2015 <- read.table("http://www1.ncdc.noaa.gov/pub/data/uscrn/products/hourly02/2015/CRNH0203-2015-MS_Holly_Springs_4_N.txt",
header = FALSE)
# right, got 8760 hours now
# get columns of interest
hollyH <- hollyHourly2015[, c(4, 5, 10, 29, 30, 34, 35)]
colnames(hollyH) <- c("date", "time", "air",
"water5", "water10", "soil5", "soil10")
hollyH$date <- ymd(hollyH$date)
startTime <- ymd_hm(201412311900)
timeVector <- startTime + dhours(0:(length(hollyH$time) - 1))
hollyH$dateTime <- timeVector
hollyH$hour <- hour(hollyH$dateTime)
hollyH$month <- month(hollyH$dateTime)
hollyH <- filter(hollyH, air > -9000)
hollyH <- filter(hollyH, soil10 > -9000)
hollyH <- filter(hollyH, soil5 > -9000)
forPlotholly <- melt(hollyH, id.vars = "dateTime")
forPlotholly <- filter(forPlotholly, variable == "air" |
variable == "soil10")
forPlotholly$variable <- ifelse(forPlotholly$variable == "air", "Air", "Soil 10 cm")
# look at air temp and soil temp at 10 cm by hour
p <- ggplot(data = forPlotholly, aes(x = dateTime, y = value))
holly1 <- p + geom_point(aes(shape = variable, alpha = variable,
colour = variable)) +
scale_shape_manual(values=c(1, 17)) +
scale_alpha_manual(values = c(0.2, 0.4)) +
background_grid(major = "xy") +
scale_x_datetime(date_labels = "%m", date_breaks = "1 month") +
theme(legend.title = element_blank()) +
labs(x = "Month",
y = "Average hourly temperature (°C)",
title = "Holly Springs, Mississippi, 2015")
holly1
save_plot("Desktop/holly_year.png", holly1, base_height = 6, base_aspect_ratio = 1.4)
p <- ggplot(data = hollyH, aes(x = dateTime, y = water10 * 100))
waterHol <- p + geom_point(shape = 1, colour = "blue", alpha = 0.3) +
background_grid(major = "xy") +
scale_y_continuous(limits = c(0, 50), breaks = seq(0, 50, 10)) +
labs(x = "Date",
y = "Average hourly Soil water\nat 10 cm depth (% by volume)",
title = "Holly Springs, Mississippi, 2015")
waterHol
save_plot("Desktop/holly_water_year.png", waterHol, base_height = 6, base_aspect_ratio = 1.4)
# look at hot soil
hotSoil <- filter(hollyH, soil10 >= 23 & month >= 6 & month <= 9)
hotSoil$delta10 <- hotSoil$air - hotSoil$soil10
median(hotSoil$water10)
hotSoil$soilWaterLevel <- ifelse(hotSoil$water10 < .27, "Soil water < 27%", "Soil water >= 27%")
# do boxplot by time of day -- facet it by soil moisture
p <- ggplot(data = hotSoil, aes(x = as.factor(hour), y = delta10))
facetDeltaholly <- p + geom_boxplot(outlier.shape = NA) +
background_grid(major = "xy") +
geom_hline(yintercept = 0, linetype = "dotted") +
geom_jitter(alpha = 0.1, shape = 1) +
facet_wrap(~soilWaterLevel, ncol = 1) +
# scale_y_continuous(limits = c(-8, 6)) +
labs(x = "For the hour ending at this time",
y = "Average air temperature - 10 cm average soil temperature, °C",
title = "Holly Springs, MS, USA, 2015\nall hours from June 1 to Sept 30 when 10 cm soil temperature >= 23°C")
facetDeltaholly
save_plot("Desktop/figure3.png",
facetDeltaholly, base_height = 6, base_aspect_ratio = 1.4)
# try to ascertain something by taking diff of mean Delta air-soil of high
# and low water, done by hour of the day
eto <- hotSoil %>%
group_by(soilWaterLevel, hour) %>%
summarise(meanDelta = mean(delta10))
eto2 <- hotSoil %>%
group_by(soilWaterLevel) %>%
summarise(meanSoil = mean(soil10),
meanAir = mean(air))
dry <- filter(eto, soilWaterLevel == "Soil water < 27%")
wet <- filter(eto, soilWaterLevel == "Soil water >= 27%")
colnames(dry) <- c("level", "hour", "deltaDry")
colnames(wet) <- c("level", "hour", "deltaWet")
forAno <- cbind.data.frame(dry, wet$deltaWet)
colnames(forAno) <- c("level", "hour", "deltaDry", "deltaWet")
p <- ggplot(data = eto, aes(x = hour, y = meanDelta))
meanDiffHolly <- p + geom_point(aes(colour = soilWaterLevel, shape = soilWaterLevel)) +
scale_x_continuous(breaks = 0:23) +
scale_y_continuous(breaks = -5:3) +
background_grid(major = "xy") +
theme(legend.title = element_blank()) +
labs(title = "Holly Springs, MS, USA, 2015\nall hours from June 1 to Sep 30 when 10 cm soil temperature >= 23°C",
x = "Hour ending at",
y = "Avg difference of air temperature - 10 cm soil temperature, °C")
meanDiffHolly
save_plot("Desktop/figure4.png",
meanDiffHolly, base_height = 6, base_aspect_ratio = 1.4)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment