Skip to content

Commit

Permalink
Finish static variables and outbreak history
Browse files Browse the repository at this point in the history
  • Loading branch information
@n8layman
n8layman committed Jul 22, 2024
1 parent 442fe88 commit a5a5f8f
Show file tree
Hide file tree
Showing 7 changed files with 497 additions and 217 deletions.
124 changes: 73 additions & 51 deletions R/calc_outbreak_history.R
View file
Original file line number Diff line number Diff line change
@@ -1,44 +1,66 @@

# This is going to be dynamic branching over list of dates. Then a target to convert to raster stacks, one for parquet, and one for animation
get_daily_outbreak_history <- function(dates,
get_daily_outbreak_history <- function(dates_df,
wahis_rvf_outbreaks_preprocessed,
continent_raster_template,
continent_polygon,
country_polygons,
output_dir = "data/outbreak_history_dataset",
output_filename = "outbreak_history.parquet",
beta_dist = .01,
beta_time = 0.5,
beta_cases = 1,
within_km = 500,
max_years = 10,
recent = 1/6) {

daily_outbreak_history <- map_dfr(dates, function(d) calc_outbreak_history(date = d,
wahis_rvf_outbreaks_preprocessed,
continent_raster_template,
continent_polygon,
country_polygons,
beta_dist = .01,
beta_time = 0.5,
beta_cases = 1,
within_km = 500,
max_years = 10,
recent = 1/6))
# Create directory if it does not yet exist
dir.create(output_dir, recursive = TRUE, showWarnings = FALSE)

daily_outbreak_history <- map_dfr(dates_df$date, function(date) calc_outbreak_history(date = date,
wahis_rvf_outbreaks_preprocessed,
continent_raster_template,
continent_polygon,
country_polygons,
beta_dist = .01,
beta_time = 0.5,
beta_cases = 1,
within_km = 500,
max_years = 10,
recent = 1/6))

daily_old_outbreak_history <- terra::rast(daily_outbreak_history$old_outbreaks_rast)
daily_recent_outbreak_history <- terra::rast(daily_outbreak_history$recent_outbreaks_rast)
daily_old_outbreak_history <- terra::rast(daily_outbreak_history$old_outbreaks_rast)

list(daily_old_outbreak_history = daily_old_outbreak_history,
daily_recent_outbreak_history = daily_recent_outbreak_history)
recent_output_filename <- paste0(output_dir, "/", tools::file_path_sans_ext(output_filename), "_recent_", dates_df$year[1], ".", tools::file_ext(output_filename))
old_output_filename <- paste0(output_dir, "/", tools::file_path_sans_ext(output_filename), "_old_", dates_df$year[1], ".", tools::file_ext(output_filename))

if(grepl("\\.parquet", output_filename)) {
# Convert to dataframe
recent <- as.data.frame(daily_recent_outbreak_history, xy = TRUE) |> as_tibble()
arrow::write_parquet(recent, recent_output_filename, compression = "gzip", compression_level = 5)
terra::writeRaster(daily_recent_outbreak_history, filename = paste0(tools::file_path_sans_ext(recent_output_filename), ".tif"), overwrite=T, gdal=c("COMPRESS=LZW"))

old <- as.data.frame(daily_old_outbreak_history, xy = TRUE) |> as_tibble()
arrow::write_parquet(old, old_output_filename, compression = "gzip", compression_level = 5)
terra::writeRaster(daily_old_outbreak_history, filename = paste0(tools::file_path_sans_ext(old_output_filename), ".tif"), overwrite=T, gdal=c("COMPRESS=LZW"))

} else {
terra::writeRaster(daily_recent_outbreak_history, filename = paste0(tools::file_path_sans_ext(recent_output_filename), ".tif"), overwrite=T, gdal=c("COMPRESS=LZW"))
terra::writeRaster(daily_old_outbreak_history, filename = paste0(tools::file_path_sans_ext(old_output_filename), ".tif"), overwrite=T, gdal=c("COMPRESS=LZW"))
}

c(recent_output_filename, old_output_filename)

}

test <- get_daily_outbreak_history(dates = dates,
wahis_rvf_outbreaks_preprocessed = wahis_rvf_outbreaks_preprocessed,
continent_raster_template = continent_raster_template,
continent_polygon = continent_polygon,
country_polygons = country_polygons)

get_outbreak_history_animation(daily_old_outbreak_history)
# test <- get_daily_outbreak_history(dates = dates,
# wahis_rvf_outbreaks_preprocessed = wahis_rvf_outbreaks_preprocessed,
# continent_raster_template = continent_raster_template,
# continent_polygon = continent_polygon,
# country_polygons = country_polygons)
#
# get_outbreak_history_animation(daily_old_outbreak_history)

#' Calculate the proximity in recent history and space of RVF outbreaks
#'
Expand Down Expand Up @@ -73,40 +95,34 @@ calc_outbreak_history <- function(date,
mutate(end_date = pmin(date, coalesce(outbreak_end_date, outbreak_start_date), na.rm = T),
years_since = as.numeric(as.duration(date - end_date), "years"),
weight = ifelse(is.na(cases) | cases == 1, 1, log10(cases + 1))*exp(-beta_time*years_since)) |>
filter(years_since < max_years & years_since > 0)
filter(end_date <= date & years_since < max_years & years_since > 0)

if(!nrow(outbreak_history)) return(NULL)
raster_template <- terra::unwrap(continent_raster_template)
raster_template[] <- 1
raster_template <- terra::crop(raster_template, terra::vect(continent_polygon$geometry), mask = TRUE)
names(raster_template) <- as.Date(date)

outbreak_history <- outbreak_history |> sf::st_as_sf(coords = c("longitude", "latitude"), crs = 4326)

old_outbreaks <- outbreak_history |> filter(years_since >= recent)
recent_outbreaks <- outbreak_history |> filter(years_since < recent)

raster_template <- terra::unwrap(continent_raster_template)
raster_template[] <- 1
raster_template <- terra::crop(raster_template, terra::vect(continent_polygon$geometry), mask = TRUE)
names(raster_template) <- as.Date(date)

recent_outbreak_rast <- get_outbreak_distance_weights(recent_outbreaks, raster_template, within_km)
old_outbreaks_rast <- get_outbreak_distance_weights(old_outbreaks, raster_template, within_km)

# Integrate space and time
# Multiply the time weight by the distance weights
tibble(recent_outbreaks_rast = list(recent_outbreak_rast),
tibble(date = as.Date(date),
recent_outbreaks_rast = list(recent_outbreak_rast),
old_outbreaks_rast = list(old_outbreaks_rast))

# Move to parquets
# 1. Should be a parquet for each data (match weather anomaly)
# of tibble with lat, long, recent_outbreak_weight, old_outbreak_weight
}

# This needs to happen for each outbreak separately otherwise when we calculate weights
# They can blend together. What I mean is that the distance matrix contains distances between
# cells within each circle and every origin not just the origin of that circle.

get_outbreak_distance_weights <- function(outbreaks, raster_template, within_km = 500, beta_dist = 0.01) {

if(!nrow(outbreaks)) {
raster_template[] <- 0
raster_template[!is.na(raster_template)] <- 0
return(raster_template)
}

Expand Down Expand Up @@ -161,10 +177,20 @@ get_outbreak_distance_weights <- function(outbreaks, raster_template, within_km
raster_template
}

get_outbreak_history_animation <- function(outbreak_raster,
get_outbreak_history_animation <- function(input_files,
output_dir = "outputs",
output_filename = "outbreak_history.gif") {

output_filename = "outbreak_history_2007.gif",
layers = NULL,
...) {

# Create directory if it does not yet exist
dir.create(output_dir, recursive = TRUE, showWarnings = FALSE)

# Change resolution to make it easier to make the gif
outbreak_raster <- terra::rast(input_files)

# if(length(layers > 1)) outbreak_raster <- terra::subset(outbreak_raster, layers)

df <- as.data.frame(outbreak_raster, xy=TRUE)

df_long <- df %>%
Expand All @@ -177,28 +203,24 @@ get_outbreak_history_animation <- function(outbreak_raster,
p <- ggplot(df_long, aes(x=x, y=y, fill=value)) +
geom_raster() +
scale_fill_viridis_c(limits=c(min(df_long$value, na.rm=T), max(df_long$value, na.rm=T))) +
labs(title = '{closest_state}', x = "Longitude", y = "Latitude", fill = "Weight") +
labs(title = '{current_frame}', x = "Longitude", y = "Latitude", fill = "Weight") +
theme_minimal() +
theme(text=element_text(size = 14))

# I can't get anim_save to work on my mac. Switching to ImageMagick rather than bother fixing it
# gifs save and render fine but can't be opened once saved. I tried re-installing gifski
# on brew and building the package from source and nothing worked so I gave up.
gganim <- gganimate::animate(p + gganimate::transition_states(states = date,
transition_length = 12,
state_length = 1,
wrap = FALSE),
nframes = 400,
fps = 8.1,
width = 550,
height = 350,
end_pause = 10,
start_pause = 20,
gganim <- gganimate::animate(p + gganimate::transition_manual(frames = date),
fps = 10,
end_pause = 20,
wrap = F,
renderer = gganimate::magick_renderer())

# animate(anim, nframes = 400,fps = 8.1, width = 550, height = 350,
# renderer = gifski_renderer("car_companies_2.gif"), end_pause = 15, start_pause = 25)

magick::image_write(gganim, path=paste(output_dir, output_filename, sep = "/"))

return(output_filename)

}
3 changes: 3 additions & 0 deletions R/get_bioclim_data.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,9 @@ get_bioclim_data <- function(output_dir,
output_filename,
raster_template) {

# Create directory if it does not yet exist
dir.create(output_dir, recursive = TRUE, showWarnings = FALSE)

template <- terra::unwrap(raster_template)
bioclim_data <- geodata::worldclim_global(var = "bio", res = 2.5, path = output_dir)

Expand Down
43 changes: 30 additions & 13 deletions _targets.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -46,19 +46,6 @@ static_targets <- tar_plan(
# modis ndvi = 0.01
tar_target(rsa_polygon, rgeoboundaries::geoboundaries("South Africa", "adm2")),


# OUTBREAK HISTORY -----------------------------------------------------------
tar_target(wahis_outbreak_history, calc_outbreak_history(wahis_rvf_outbreaks_preprocessed,
continent_raster_template,
continent_polygon,
country_polygons)),

tar_target(wahis_daily_outbreak_history, get_daily_outbreak_history(wahis_rvf_outbreaks_preprocessed,
continent_raster_template,
continent_polygon,
country_polygons)),


# SOIL -----------------------------------------------------------
tar_target(soil_directory_raw,
create_data_directory(directory_path = "data/soil")),
Expand Down Expand Up @@ -154,7 +141,37 @@ dynamic_targets <- tar_plan(
tar_target(wahis_rvf_controls_preprocessed,
preprocess_wahis_rvf_controls(wahis_rvf_controls_raw)),

# OUTBREAK HISTORY -----------------------------------------------------------
tar_target(wahis_outbreak_dates, tibble(date = seq(from = min(coalesce(wahis_rvf_outbreaks_preprocessed$outbreak_end_date, wahis_rvf_outbreaks_preprocessed$outbreak_start_date), na.rm = T),
to = max(coalesce(wahis_rvf_outbreaks_preprocessed$outbreak_end_date, wahis_rvf_outbreaks_preprocessed$outbreak_start_date), na.rm = T),
by = "day"),
year = year(date),
month = month(date)) |>
group_by(year) |>
tar_group(),
iteration = "group"),

# Map over year batch over day otherwise too many branches.
tar_target(wahis_outbreak_history, get_daily_outbreak_history(dates_df = wahis_outbreak_dates,
wahis_rvf_outbreaks_preprocessed,
continent_raster_template,
continent_polygon,
country_polygons),
head(wahis_outbreak_dates, 1),
iteration = "vector",
format = "file",
repository = "local"),

tar_target(wahis_outbreak_history_recent_animation, get_outbreak_history_animation(input_files = c("data/outbreak_history_dataset/outbreak_history_recent_2007.tif"),
output_dir = "outputs",
output_filename = "outbreak_history_recent_2007.gif",
wahis_outbreak_history)), # Just included to enforce dependency with wahis_outbreak_history

tar_target(wahis_outbreak_history_old_animation, get_outbreak_history_animation(input_files = c("data/outbreak_history_dataset/outbreak_history_old_2007.tif"),
output_dir = "outputs",
output_filename = "outbreak_history_old_2007.gif",
wahis_outbreak_history)), # Just included to enforce dependency with wahis_outbreak_history

# SENTINEL NDVI -----------------------------------------------------------
# 2018-present
# 10 day period
Expand Down
Loading

0 comments on commit a5a5f8f

Please sign in to comment.