diff --git a/1_fetch.R b/1_fetch.R
index 1a8f5af2b117055fdba139fd675e2a898528fd49..d880fe8a7812abf65b9debea0da2456d55c16351 100644
--- a/1_fetch.R
+++ b/1_fetch.R
@@ -39,5 +39,9 @@ p1_targets <- list(
)) |>
filter(group %in% c('CONUS')) |>
rmapshaper::ms_simplify(keep = 0.2) |>
- filter(NAME %in% p1_census_states))
- )
\ No newline at end of file
+ filter(NAME %in% p1_census_states)),
+ # raster data for population density
+ tar_target(p1_pop_density_raster,
+ '1_fetch/in/gpw-v4-population-count-rev11_2020_30_sec_tif/gpw_v4_population_count_rev11_2020_30_sec.tif')
+ )
+
\ No newline at end of file
diff --git a/2_process.R b/2_process.R
index 8dc7aaf94c274f192e9e7d7fd3af5955c13e7c45..0f438ad32c11f5b30e0948d77f1338c2cfc62526 100644
--- a/2_process.R
+++ b/2_process.R
@@ -70,11 +70,12 @@ p2_targets <- list(
# distinct(determinant, .keep_all = TRUE) |>
readr::write_csv('public/determinant_uncertainty.csv')
),
-tar_target(p2_unc_determinant_json,
- read_csv(p2_unc_agg_summary_csv) |>
- toJSON(pretty = TRUE) |>
- write("public/determinant_uncertainty.json")
- ),
+# commented out for now so we don't overwrite spanish names
+#tar_target(p2_unc_determinant_json,
+# read_csv(p2_unc_agg_summary_csv) |>
+# toJSON(pretty = TRUE) |>
+# write("public/determinant_uncertainty.json")
+# ),
tar_target(`p2_unc_agg_summary_ind_csv`,
p2_unc_agg_ind_summary |>
left_join(p2_top_trend_ind_stats) |>
@@ -241,5 +242,19 @@ tar_target(p2_census_acs5sub_disability_data,
survey_var = "acs5",
percent_rename = FALSE),
pattern = map(p2_census_acs5_disability_layers),
- iteration = "list")
+ iteration = "list"),
+# process population density raster data
+tar_target(p2_conus_sf,
+ fetch_conus_sf()),
+tar_target(p2_conus_sf_proj,
+ p2_conus_sf |>
+ st_transform(p1_proj)),
+tar_target(p2_conus_inner,
+ rmapshaper::ms_innerlines(p2_conus_sf_proj)),
+tar_target(p2_pop_density_processed,
+ process_raster(in_raster = p1_pop_density_raster, proj = p1_proj,
+ conus = p2_conus_sf, conus_proj = p2_conus_sf_proj,
+ outfile_path = "2_process/out/pop_density.tif")),
+tar_target(p2_pop_density_filepath,
+ "2_process/out/pop_density.tif")
)
\ No newline at end of file
diff --git a/2_process/src/data_utils.R b/2_process/src/data_utils.R
index dd430c90b2b59cb330dcef6a616222b39cbca603..18bfc86091e2155cb4970f437f010ee53add3041 100644
--- a/2_process/src/data_utils.R
+++ b/2_process/src/data_utils.R
@@ -47,4 +47,47 @@ process_perc <- function(tot_pop, tot_var){
return(joined_perc_df)}
prep_tree_data <- function(data) {
-}
\ No newline at end of file
+}
+
+## data processing for population density raster map
+
+fetch_conus_sf <- function(){
+
+ conus_sf <- tigris::states(cb = TRUE) |>
+ dplyr::filter(NAME %in% c('Washington', 'Oregon', 'California', 'Idaho', 'Nevada',
+ 'Utah', 'Arizona', 'Montana', 'Wyoming', 'Colorado',
+ 'New Mexico', 'North Dakota', 'South Dakota', 'Nebraska', 'Kansas',
+ 'Oklahoma', 'Texas', 'Minnesota', 'Iowa', 'Missouri',
+ 'Arkansas', 'Louisiana')
+ ) |>
+ rmapshaper::ms_simplify(keep = 0.2)
+}
+
+#' @title Process raster data for ploting
+#' @param
+process_raster <- function(in_raster, proj, conus, conus_proj, outfile_path){
+ # in_raster = p1_pop_density_raster, proj = p1_proj, conus = p2_conus_sf, conus_proj = p2_conus_sf_proj
+ raw_data <- rast(in_raster)
+
+ # crop population data to area of interest
+ terra::window(raw_data) <- terra::ext(conus)
+
+ # project population data
+ conus_sf_rast <- rast(conus_proj, resolution = c(1000, 1000))
+ pop_usa_proj <- project(raw_data, conus_sf_rast)
+ # match boundaries of population data to conus data
+ pop_usa_mask <- terra::mask(pop_usa_proj, terra::vect(conus_proj))
+
+ # change to tibble using tidy terra for processing into log scale
+ usa_dat <- as_tibble(pop_usa_mask, xy = TRUE)
+ usa_dat$pop <- ifelse(usa_dat$gpw_v4_population_count_rev11_2020_30_sec < 1
+ & usa_dat$gpw_v4_population_count_rev11_2020_30_sec > 0, 0.1,
+ usa_dat$gpw_v4_population_count_rev11_2020_30_sec)
+ usa_dat$pop_log10 <- log10(usa_dat$pop)
+ usa_dat <- usa_dat |>
+ dplyr::select(-gpw_v4_population_count_rev11_2020_30_sec, -pop)
+ # convert back into raster
+ usa_dat_rast <- as_spatraster(usa_dat)
+
+ writeRaster(usa_dat_rast, filename=outfile_path, overwrite=TRUE)
+}
diff --git a/3_visualize.R b/3_visualize.R
index 6a043398b3d345f56ae414e6377478c9bf658fa5..1bf1be2256f88520d5d1f27ee3e27572000c3d42 100644
--- a/3_visualize.R
+++ b/3_visualize.R
@@ -115,6 +115,15 @@ p3_targets <- list(
),
format = "file"
),
+# population density raster plot
+tar_target(
+ p3_pop_density_plot,
+ plot_raster(in_raster = p2_pop_density_filepath, conus_sf = p2_conus_sf_proj, viz_config_df = p0_viz_config_df,
+ conus_inner = p2_conus_inner, viz_config_pal = p0_viz_config_pal$living_conditions,
+ width = p0_viz_config_df$width_desktop, height = p0_viz_config_df$height_desktop,
+ font_size = p0_viz_config_df$font_size_desktop, barwidth = 20, barheight = 1,
+ outfile_path = "3_visualize/out/pop_density_rast_2020_en.png", low_ramp_col = "#E8F4F4"),
+ format = "file"),
# Spanish version's of maps -----------------------------------------------
tar_target(
p3_med_income_png_es,
diff --git a/3_visualize/src/plot_utils.R b/3_visualize/src/plot_utils.R
index 270e09f2edb5170d5ab2fd97b54c84c0bd6fc328..aaf2e72166279e58ac2d4bd16504b7aa48efe00c 100644
--- a/3_visualize/src/plot_utils.R
+++ b/3_visualize/src/plot_utils.R
@@ -17,7 +17,7 @@
#' @param break_vals, set legend breaks
#' @param low_ramp_col, assign color to low end of color ramp for legend
plot_census_map <- function(census_data, conus_sf, leg_title, outfile_path, var,
- percent_leg, viz_config_df, viz_config_pal, width, height,
+ percent_leg, viz_config_df, viz_config_pal, width, height,
font_size, barwidth, barheight, lim_vals, break_vals,
dollar_leg, low_ramp_col){
@@ -105,4 +105,110 @@ plot_census_map <- function(census_data, conus_sf, leg_title, outfile_path, var,
ggsave(outfile_path, final_map, width = width, height = height, dpi = viz_config_df$dpi, bg = viz_config_df$bg_col, units = "in")
return(outfile_path)
-}
\ No newline at end of file
+}
+
+
+# Map raster data
+#'
+#' @param raster_data .tif file of raster data for specified variable of interest
+#' # notes: update colorpalette according to Elmera's request
+
+plot_raster <- function(in_raster, conus_sf, conus_inner, outfile_path,
+ viz_config_df, viz_config_pal, width, height,
+ font_size, barwidth, barheight,
+ low_ramp_col){
+ # in_raster = p2_pop_density_processed, conus_sf = p2_conus_sf_proj, conus_inner = p2_conus_inner, viz_config_pal = p0_viz_config_pal$living_conditions
+ # width = p0_viz_config_df$width_desktop, height = p0_viz_config_df$height_desktop, font_size = p0_viz_config_df$font_size_desktop, barwidth = 20, barheight = 1,
+ #outfile_path = "3_visualize/out/pop_density_rast_2020_en.png", low_ramp_col = "#eef0ff",
+ raster_data <- rast(in_raster)
+
+ # plot
+ font_legend <- viz_config_df$load_font
+ font_add_google(font_legend)
+ showtext_opts(dpi = 300, regular.wt = 200, bold.wt = 700)
+ showtext_auto(enable = TRUE)
+
+ # doesn't work with cowplot::get_plot_component()
+ legend_label <- expression("Population per 1 km"^2) #bquote('Population per 1'~km^2)
+
+ (pop_map <- ggplot() +
+ #geom_spatraster(data = pop_usa_proj)+
+ geom_sf(data = conus_sf, fill = low_ramp_col, color = low_ramp_col, linewidth = 0.5) +
+ geom_spatraster(data = raster_data)+ #x = x, y = y,
+ geom_sf(data = conus_inner, fill = NA, color = viz_config_df$counties_outline_col, linewidth = 0.5) +
+ scale_fill_gradientn(
+ colors = colorRampPalette(c(low_ramp_col, viz_config_pal))(100),
+ name = legend_label,
+ #limits = lim_vals, # c(0, 100),
+ #breaks = break_vals, # c(0, 25, 50, 75, 100),
+ #labels = if (dollar_leg) scales::label_dollar() else function(x) paste0(x, "%"),
+ labels = c(1, 10, 100, '1k', '10k', '100k'),
+ na.value=NA
+ )+
+ #scale_fill_viridis_c(na.value = NA, option = 'mako', breaks = c(0, 1, 2, 3, 4, 5), direction = -1,
+ # leave a little room for the NAs and Inf (which are 0s)
+ # to be darker than the -1 values
+ # labels = c(1, 10, 100, '1k', '10k', '100k'), begin = 0.05) +
+ labs(fill = legend_label) + #"Population per 1 km^2"
+ theme(plot.background = element_rect(fill = "white", color = "white"),
+ panel.background = element_rect(fill = "white", color = "white"),
+ axis.title = element_blank(),
+ axis.text = element_blank(),
+ axis.ticks = element_blank(),
+ panel.grid = element_blank(),
+ legend.background = element_blank(),
+ legend.direction = "horizontal",
+ legend.position = "bottom") +
+ guides(fill = guide_colorbar(
+ title.position = "top",
+ title.theme = element_text(face = 'bold', family = font_legend, size = font_size),
+ direction = "horizontal",
+ position = "bottom",
+ barwidth = barwidth,
+ barheight = barheight
+ ))
+ # guides(fill = guide_colorbar(
+ # direction = "horizontal",
+ # title = "Population (per km2)",
+ # barwidth = 15,
+ # title.position = "top",
+ # title.theme = element_text(color = "white"),
+ # label.theme = element_text(color = "white")))
+ )
+
+ background_color = "white"
+ plot_margin = 0.025
+
+ # cowplot to get map sizes larger
+ canvas <- grid::rectGrob(
+ x = 0, y = 0,
+ width = 6, height = 6,
+ gp = grid::gpar(fill = background_color, alpha = 1, col = background_color
+ )
+ )
+ pop_map_legend <- get_plot_component(pop_map, 'guide-box', return_all = TRUE)
+ # compose final plot
+ final_map <- ggdraw(ylim = c(0,1),
+ xlim = c(0,1)) +
+ # White background
+ draw_grob(canvas,
+ x = 0, y = 1,
+ height = 6, width = 6,
+ hjust = 0, vjust = 1) +
+ # Add main plot
+ draw_plot(pop_map + theme(legend.position="none"),
+ x = -0.01,
+ y = 0.08,
+ height = 0.98,
+ width = (1-plot_margin)*1.03) +
+ # Add legend
+ draw_plot(pop_map_legend[[3]],
+ x = 0.48,
+ y = 0.02,
+ height = 0.09 ,
+ width = 0.1 - plot_margin)
+
+ ggsave(outfile_path, final_map, width = width, height = height, dpi = viz_config_df$dpi, bg = viz_config_df$bg_col, units = "in")
+
+}
+
\ No newline at end of file
diff --git a/_targets.R b/_targets.R
index 99ea25fbfec38c9ad6cf16b43f0214523ef47039..1f606075799269ffed33bb53de7e5a97fd24875d 100644
--- a/_targets.R
+++ b/_targets.R
@@ -15,7 +15,10 @@ tar_option_set(packages = c('tidyverse',
'readr',
'cowplot',
'rmapshaper',
- 'jsonlite'))
+ 'jsonlite',
+ 'terra',
+ 'tidyterra',
+ 'raster'))
# Phase target makefiles
source("0_config.R")