Trying to map a value for geom_vline, but is not plotting in the correct place on the x axis with ggplot in R

Question

I am currently trying to generate NOAA tide prediction charts (x = datetime, y = water level) with the dawn/sunrise/dusk/sunset times as vertical lines along the x axis timeline.

The rnoaa package calls the data and gives me the prediction date times in POSIXct. The suncalc library provides me a data frame with each date in the range's sunrise, sunset, etc. in POSIXct format as well.

library(rnoaa)
library(tidyverse)
library(ggplot2)
library(suncalc)
march.tides <- as.data.frame(coops_search(station_name = 8551762, 
                                          begin_date = 20200301, end_date = 20200331, 
                                          datum = "mtl", product = "predictions"))
march.tides <- march.tides %>%
  mutate(yyyy.mm.dd = as.Date(predictions.t))
dates <- unique(march.tides$yyyy.mm.dd)
sunlight.times <- getSunlightTimes(date = seq.Date(as.Date("2020/3/1"), as.Date("2020/3/31"), by = 1), 
                lat = 39.5817, lon = -75.5883, tz = "EST")

I then have a loop that spits out separate plots for each calendar date - which works hunky dory. The vertical lines are drawing on the graph without an error, but are definitely in the wrong spot (sunrise is being drawn around 11am when it should be 06:30).

for (i in seq_along(dates)) {
  plot <- ggplot(subset(march.tides, march.tides$yyyy.mm.dd==dates[i])) +
    aes(x = predictions.t, y = predictions.v) +
    geom_line(size = 1L, colour = "#0c4c8a") +
    theme_bw() + 
    geom_vline(xintercept = sunlight.times$sunrise) +
    geom_vline(xintercept = sunlight.times$sunset) +
    geom_vline(xintercept = sunlight.times$dawn, linetype="dotted") +
    geom_vline(xintercept = sunlight.times$dusk, linetype="dotted") +
    ggtitle(dates[i])
  print(plot)
}

I could alternatively facet the separate dates instead of this looping approach. Even when I subset the data to a single date, the vertical lines still did not draw correctly.

I wondered if maybe the issue was a time zone one. If I try to stick a time zone argument onto the tide prediction data call, I get the error:

Error in if (!repeated && grepl("%[[:xdigit:]]{2}", URL, useBytes = TRUE)) return(URL) : 
  missing value where TRUE/FALSE needed

Ben · Accepted Answer · 2020-03-02T15:06:22.940

2

It looks like you want to use EST as your timezone, so you could include in your conversion of predictions.t.

I would be explicit in what you want labeled on your xaxis in ggplot using scale_x_datetime, including the timezone.

library(rnoaa)
library(tidyverse)
library(ggplot2)
library(suncalc)
library(scales)

march.tides <- as.data.frame(coops_search(station_name = 8551762, 
                                          begin_date = 20200301, end_date = 20200331, 
                                          datum = "mtl", product = "predictions"))

march.tides <- march.tides %>%
  mutate(yyyy.mm.dd = as.Date(predictions.t, tz = "EST"))
dates <- unique(march.tides$yyyy.mm.dd)
sunlight.times <- getSunlightTimes(date = seq.Date(as.Date("2020/3/1"), as.Date("2020/3/31"), by = 1), 
                                   lat = 39.5817, lon = -75.5883, tz = "EST")

for (i in seq_along(dates)) {
  plot <- ggplot(subset(march.tides, march.tides$yyyy.mm.dd==dates[i])) +
    aes(x = predictions.t, y = predictions.v) +
    geom_line(size = 1L, colour = "#0c4c8a") +
    theme_bw() + 
    geom_vline(xintercept = sunlight.times$sunrise) +
    geom_vline(xintercept = sunlight.times$sunset) +
    geom_vline(xintercept = sunlight.times$dawn, linetype="dotted") +
    geom_vline(xintercept = sunlight.times$dusk, linetype="dotted") +
    ggtitle(dates[i]) +
    scale_x_datetime(labels = date_format("%b %d %H:%M", tz = "EST"))
  print(plot)
}

Plot

edited Mar 02 '20 at 15:06

answered Feb 29 '20 at 19:39

Ben

28,684
5
23
45

1

Thank you very much - I can't recall how I was attempting the time zone that was causing the error I included in the original question, but this is working perfectly for me now. My only tweak is loading library(scales) for the date_format() function. – M Mason Mar 02 '20 at 14:56
@MMason - glad to hear. I added `library(scales)` to the answer. – Ben Mar 02 '20 at 15:07
I did realize another caveat/remaining issue for this approach since the data call is still in GMT (any other time zone argument in the coops_search() function gives the error), once you convert it you are still working with data that starts and ends 6 hours earlier in EST than desired for a month of data. Their functions are only able to call up to 31 days, so expanding the request dates will not always work. I submitted an issue on GitHub. – M Mason Mar 02 '20 at 16:13
Hmmm...what did you include for `time_zone` in the `coops_search()` call? Was it `lst` or `lst_ldt`? It looks like it refers to the `station_name` - just curious, do you get errors with different stations? – Ben Mar 02 '20 at 17:04
I have tried both, and they both produce the same error. The package documentation says the only options are `gmt`, `lst`, and `lst_ldt`. I already had a nearly identical script that was grabbing verified water levels (post-field work), which only changes the `product` argument and uses `time_zone = "lst_ldt"` without any problems. It may be a data source error with the tide prediction product from CO-OPs, so I have created an issue on the package github in case they might know something. I have not tried different stations however, and will try that next. – M Mason Mar 03 '20 at 19:01

Trying to map a value for geom_vline, but is not plotting in the correct place on the x axis with ggplot in R

1 Answers1