Skip to content

Instantly share code, notes, and snippets.

@paleolimbot
Created September 18, 2019 13:52
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save paleolimbot/54773c414cdb23aad8669047d2ca84cb to your computer and use it in GitHub Desktop.
Save paleolimbot/54773c414cdb23aad8669047d2ca84cb to your computer and use it in GitHub Desktop.
library(tidyverse)
# best case: your data is already in the right form
perfect_data <- tibble(
day_in_year = c(1, 101, 201, 1, 101, 201, 1, 101, 201),
depth = c(1, 1, 1, 2, 2, 2, 3, 3, 3),
temp = c(15, 16, 17, 14, 15, 16, 13, 14, 15)
)
ggplot(perfect_data, aes(x = day_in_year, y = depth, fill = temp)) +
geom_raster() +
scale_y_reverse()
# imperfect case: your data is only equally spaced in one dimension
# in either the time- or depth- dimension
imerfect_data <- tibble(
day_in_year = c(1, 124, 201, 1, 124, 201, 1, 124, 201),
depth = c(1, 1, 1, 2, 2, 2, 3, 3, 3),
temp = c(15, 16, 17, 14, 15, 16, 13, 14, 15)
)
# create a function that will, given a depth, estimate the temp on any given day
estimate_temp <- function(target_depth, target_day_in_year) {
data_for_depth <- imerfect_data %>% filter(depth == target_depth)
# this is the easiest way to do a linear interpolation
approx(data_for_depth$day_in_year, data_for_depth$temp, xout = target_day_in_year)$y
}
# create a tibble with the time and depth values you want
perfect_data <- crossing(
tibble(day_in_year = seq(1, 201, by = 10)),
# values from the equally spaced dimension must exist in `imperfect_data`
tibble(depth = c(1, 2, 3))
) %>%
# group by the equally spaced dimension
group_by(depth) %>%
# estimate the temperature (our estimator only works at a single depth)
mutate(temp = estimate_temp(depth[1], day_in_year))
ggplot(perfect_data, aes(x = day_in_year, y = depth, fill = temp)) +
geom_raster() +
scale_y_reverse()
# hardest case: unequal in both time and depth (probably what you have)
hard_data <- tibble(
day_in_year = c(1, 124, 201, 1, 124, 201, 1, 124, 201),
depth = c(0.25, 0, 1, 2.5, 2, 2.1, 3, 3, 3),
temp = c(15, 16, 17, 14, 15, 16, 13, 14, 15)
)
# you can always use points to get a feel for what's going on
ggplot(hard_data, aes(day_in_year, depth, color = temp)) +
geom_point(size = 3) +
scale_y_reverse()
# create a function to estimate the temperature at any depth using a smoother
# you will probably have to increase "k" as you get more values
# roughly, k is "wigglyness"...low "k" will give you a smooth surface
estimate_temp <- function(target_depth, target_day_in_year) {
smooth_fit <- mgcv::gam(temp ~ s(day_in_year, depth, k = 4), data = hard_data)
predict(smooth_fit, newdata = tibble(day_in_year = target_day_in_year, depth = target_depth))
}
# create a tibble with the time and depth values you want
perfect_data <- crossing(
tibble(day_in_year = seq(1, 201, by = 10)),
tibble(depth = seq(0, 3, by = 0.2))
) %>%
# estimate the temperature
mutate(temp = estimate_temp(depth, day_in_year))
ggplot(perfect_data, aes(x = day_in_year, y = depth, fill = temp)) +
geom_raster() +
scale_y_reverse()
@paleolimbot
Copy link
Author

(rendered version of the above:)

library(tidyverse)
#> Registered S3 method overwritten by 'rvest':
#>   method            from
#>   read_xml.response xml2

# best case: your data is already in the right form
perfect_data <- tibble(
  day_in_year = c(1, 101, 201, 1, 101, 201, 1, 101, 201),
  depth = c(1, 1, 1, 2, 2, 2, 3, 3, 3),
  temp = c(15, 16, 17, 14, 15, 16, 13, 14, 15)
)

ggplot(perfect_data, aes(x = day_in_year, y = depth, fill = temp)) +
  geom_raster() +
  scale_y_reverse()

# imperfect case: your data is only equally spaced in one dimension
# in either the time- or depth- dimension
imerfect_data <- tibble(
  day_in_year = c(1, 124, 201, 1, 124, 201, 1, 124, 201),
  depth = c(1, 1, 1, 2, 2, 2, 3, 3, 3),
  temp = c(15, 16, 17, 14, 15, 16, 13, 14, 15)
)

# create a function that will, given a depth, estimate the temp on any given day
estimate_temp <- function(target_depth, target_day_in_year) {
  data_for_depth <- imerfect_data %>% filter(depth == target_depth)
  # this is the easiest way to do a linear interpolation
  approx(data_for_depth$day_in_year, data_for_depth$temp, xout = target_day_in_year)$y
}

# create a tibble with the time and depth values you want
perfect_data <- crossing(
  tibble(day_in_year = seq(1, 201, by = 10)),
  # values from the equally spaced dimension must exist in `imperfect_data`
  tibble(depth = c(1, 2, 3))
) %>%
  # group by the equally spaced dimension
  group_by(depth) %>%
  # estimate the temperature (our estimator only works at a single depth)
  mutate(temp = estimate_temp(depth[1], day_in_year))
  
ggplot(perfect_data, aes(x = day_in_year, y = depth, fill = temp)) +
  geom_raster() +
  scale_y_reverse()

# hardest case: unequal in both time and depth (probably what you have)
hard_data <- tibble(
  day_in_year = c(1, 124, 201, 1, 124, 201, 1, 124, 201),
  depth = c(0.25, 0, 1, 2.5, 2, 2.1, 3, 3, 3),
  temp = c(15, 16, 17, 14, 15, 16, 13, 14, 15)
)

# you can always use points to get a feel for what's going on
ggplot(hard_data, aes(day_in_year, depth, color = temp)) +
  geom_point(size = 3) +
  scale_y_reverse()

# create a function to estimate the temperature at any depth using a smoother
# you will probably have to increase "k" as you get more values
# roughly, k is "wigglyness"...low "k" will give you a smooth surface
estimate_temp <- function(target_depth, target_day_in_year) {
  smooth_fit <- mgcv::gam(temp ~ s(day_in_year, depth, k = 4), data = hard_data)
  predict(smooth_fit, newdata = tibble(day_in_year = target_day_in_year, depth = target_depth))
}

# create a tibble with the time and depth values you want
perfect_data <- crossing(
  tibble(day_in_year = seq(1, 201, by = 10)),
  tibble(depth = seq(0, 3, by = 0.2))
) %>%
  # estimate the temperature
  mutate(temp = estimate_temp(depth, day_in_year))

ggplot(perfect_data, aes(x = day_in_year, y = depth, fill = temp)) +
  geom_raster() +
  scale_y_reverse()

Created on 2019-09-18 by the reprex package (v0.2.1)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment