Building a new geom_hurricane
I have been meaning to create a new geom for a data set that has been tidied in the following form:
Katrina
# A tibble: 3 x 9
storm_id date latitude longitude wind_speed ne se sw nw
<chr> <dttm> <dbl> <dbl> <fct> <dbl> <dbl> <dbl> <dbl>
1 KATRINA-2005 2005-08-29 12:00:00 29.5 -89.6 34 200 200 150 100
2 KATRINA-2005 2005-08-29 12:00:00 29.5 -89.6 50 120 120 75 75
3 KATRINA-2005 2005-08-29 12:00:00 29.5 -89.6 64 90 90 60 60
I first defined the class and then the actual geom function, however, my output plot turns out to be so miniaturized, So I would appreciate it if you could tell me where I possibly go wrong with the scales.
GeomHurricane <- ggplot2::ggproto("GeomHurricane", Geom,
required_aes = c("x", "y",
"r_ne", "r_se", "r_sw", "r_nw"
),
default_aes = aes(fill = 1, colour = 1,
alpha = 1, scale_radii = 1),
draw_key = draw_key_polygon,
draw_group = function(data, panel_scales, coord) {
coords <- coord$transform(data, panel_scales) %>%
mutate(r_ne = r_ne * 1852 * scale_radii,
r_se = r_se * 1852 * scale_radii,
r_sw = r_sw * 1852 * scale_radii,
r_nw = r_nw * 1852 * scale_radii
)
# Creating quadrants
for(i in 1:nrow(data)) {
# Creating the northeast quadrants
data_ne <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 1:90,
d = data[i,]$r_ne),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
# Creating the southeast quadrants
data_se <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 90:180,
d = data[i,]$r_se),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
# Creating the southwest quadrants
data_sw <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 180:270,
d = data[i,]$r_sw),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
# Creating the northwest quadrants
data_nw <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 270:360,
d = data[i,]$r_nw),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
data_quadrants <- dplyr::bind_rows(list(
data_ne, data_se, data_sw, data_nw
))
data_quadrants <- data_quadrants %>% dplyr::rename(
x = lon,
y = lat
)
data_quadrants$colour <- as.character(data_quadrants$colour)
data_quadrants$fill <- as.character(data_quadrants$fill)
}
coords_data <- coord$transform(data_quadrants, panel_scales)
grid::polygonGrob(
x = coords_data$x,
y = coords_data$y,
default.units = "native",
gp = grid::gpar(
col = coords_data$colour,
fill = coords_data$fill,
alpha = coords_data$alpha
)
)
}
)
and the actual geom function definition:
geom_hurricane <- function(mapping = NULL, data = NULL, stat = "identity",
position = "identity", na.rm = FALSE,
show.legend = NA, inherit.aes = TRUE, ...) {
ggplot2::layer(
geom = GeomHurricane, mapping = mapping,
data = data, stat = stat, position = position,
show.legend = show.legend, inherit.aes = inherit.aes,
params = list(na.rm = na.rm, ...)
)
}
So I went on to plot the following:
ggplot(data = Katrina) +
geom_hurricane(aes(x = longitude, y = latitude,
r_ne = ne, r_se = se, r_sw = sw, r_nw = nw,
fill = wind_speed, colour = wind_speed)) +
scale_colour_manual(name = "Wind speed (kts)",
values = c("red", "orange", "yellow")) +
scale_fill_manual(name = "Wind speed (kts)",
values = c("red", "orange", "yellow"))
The data for this purpose can be found here as Atlantic basin data set 1988 - 2018: http://rammb.cira.colostate.edu/research/tropical_cyclones/tc_extended_best_track_dataset/
For your consideration I used the following codes to tidy the data:
ext_tracks_widths <- c(7, 10, 2, 2, 3, 5, 5, 6, 4, 5, 4, 4, 5, 3, 4, 3, 3, 3,
4, 3, 3, 3, 4, 3, 3, 3, 2, 6, 1)
ext_tracks_colnames <- c("storm_id", "storm_name", "month", "day",
"hour", "year", "latitude", "longitude",
"max_wind", "min_pressure", "rad_max_wind",
"eye_diameter", "pressure_1", "pressure_2",
paste("radius_34", c("ne", "se", "sw", "nw"), sep = "_"),
paste("radius_50", c("ne", "se", "sw", "nw"), sep = "_"),
paste("radius_64", c("ne", "se", "sw", "nw"), sep = "_"),
"storm_type", "distance_to_land", "final")
ext_tracks <- read_fwf("ebtrk_atlc_1988_2015.txt",
fwf_widths(ext_tracks_widths, ext_tracks_colnames),
na = "-99")
storm_observation <- ext_tracks %>%
unite("storm_id", c("storm_name", "year"), sep = "-",
na.rm = TRUE, remove = FALSE) %>%
mutate(longitude = -longitude) %>%
unite(date, year, month, day, hour) %>%
mutate(date = ymd_h(date)) %>%
select(storm_id, date, latitude, longitude, radius_34_ne:radius_64_nw) %>%
pivot_longer(cols = contains("radius"), names_to = "wind_speed",
values_to = "value") %>%
separate(wind_speed, c(NA, "wind_speed", "direction"), sep = "_") %>%
pivot_wider(names_from = "direction", values_from = "value") %>%
mutate(wind_speed = as.factor(wind_speed))
Katrina <- storm_observation %>%
filter(storm_id == "KATRINA-2005", date == ymd_h("2005-08-29-12"))
Solution
Alright, there are 2 problems that I spotted. Problem 1 is that in your draw_group() ggproto method, you convert the radii from nautical miles to meters (I think), but you write this to the coords variable. However, you use the data variable to do the geosphere::destPoint calculation.
Here is a version of that method that I think should work:
draw_group = function(data, panel_scales, coord) {
scale_radii <- if (is.null(data$scale_radii)) 1 else data$scale_radii
data <- data %>%
mutate(r_ne = r_ne * 1852 * scale_radii,
r_se = r_se * 1852 * scale_radii,
r_sw = r_sw * 1852 * scale_radii,
r_nw = r_nw * 1852 * scale_radii
)
# Creating quadrants
for(i in 1:nrow(data)) {
# Creating the northeast quadrants
data_ne <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 1:90, # Should this start at 0?
d = data[i,]$r_ne),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
# Creating the southeast quadrants
data_se <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 90:180,
d = data[i,]$r_se),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
# Creating the southwest quadrants
data_sw <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 180:270,
d = data[i,]$r_sw),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
# Creating the northwest quadrants
data_nw <- data.frame(colour = data[i,]$colour,
fill = data[i,]$fill,
geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
b = 270:360,
d = data[i,]$r_nw),
group = data[i,]$group,
PANEL = data[i,]$PANEL,
alpha = data[i,]$alpha
)
data_quadrants <- dplyr::bind_rows(list(
data_ne, data_se, data_sw, data_nw
))
data_quadrants <- data_quadrants %>% dplyr::rename(
x = lon,
y = lat
)
data_quadrants$colour <- as.character(data_quadrants$colour)
data_quadrants$fill <- as.character(data_quadrants$fill)
}
coords_data <- coord$transform(data_quadrants, panel_scales)
grid::polygonGrob(
x = coords_data$x,
y = coords_data$y,
default.units = "native",
gp = grid::gpar(
col = coords_data$colour,
fill = coords_data$fill,
alpha = coords_data$alpha
)
)
}
The next problem is that you only define 1 x coordinate with the Katrina example. However, the scales don't know about your radius parameters, so they don't adjust the limits to fit your radii in. You can circumvent this by setting xmin, xmax, ymin and ymax bounding box parameters, so that scale_x_continuous() can learn about your radius. (Same thing for the y scale). I'd solve this by using a setup_data method to your ggproto object.
Here is the setup data method that I used to test with, but I'm no spatial genius so you'd have to check if this makes sense.
setup_data = function(data, params) {
maxrad <- max(c(data$r_ne, data$r_se, data$r_sw, data$r_nw))
maxrad <- maxrad * 1852
x_range <- unique(range(data$x))
y_range <- unique(range(data$y))
xy <- as.matrix(expand.grid(x_range, y_range))
extend <- lapply(c(0, 90, 180, 270), function(b) {
geosphere::destPoint(p = xy,
b = b,
d = maxrad)
})
extend <- do.call(rbind, extend)
transform(
data,
xmin = min(extend[, 1]),
xmax = max(extend[, 1]),
ymin = min(extend[, 2]),
ymax = max(extend[, 2])
)
}
After implenting these changes, I get a figure like this:
- Retrieving expected data.frame for testthat expectation
- How to use an anonymous function with pipe operator?
- R passing RSelenium driver environment as function argument
- Keep file names merging list with lapply?
- replace strings that contains one another
- Stop parsing out zeros after decimals in ggplot2's annotate
- Plot data.tree coloring and labelling by level
- Extracting the date and time out of a date and time format
- How to reshape a dataframe into wide format with specified column pairs
- Show partial segments removed after limit scales for `geom_pointrange`
- Control ggplot2 legend look without affecting the plot
- Concatenate strings over all rows in a single column of a data frame
- Is there some way to keep variable names from.SD+.SDcols together with non .SD variable names in data.table?
- How to create mutually dependent checkbox columns in a DT table?
- A problem with an stacked, grouped and faceted graph y gg plot
- Installing R on Linux: configure: error: libcurl >= 7.28.0 library and headers are required with support for https
- Create a column conditional on multiple other character and numerical columns
- data.frame to array (2 columns)
- Changing options() in a function environment without changing options() in global environment in R?
- pivot(), group_by() and summarise() with nested data
- Why legend is not showing up in plot layout using patchwork in R?
- in dplyr::mutate, dplyr::starts_with works for .before but not .after?
- Finding optimal cut-off points for dividing a variable into equal groups considering ties in which quantiles do not work
- How to add a legend to hline?
- R Shiny Suppress slider handle for sliderInput until click
- brmultinom brglm2 vifs and dredge dont work
- R, how to get current time miliseconds, without extra information ? (Sys.time(), format("%X"))
- Installing multiple versions of R
- Command to see 'R' path that RStudio is using
- Fast pairwise longest common substring from start