adding visualization code

Code/Visualization/data/adoption_scenarios_5_and_15.csv

@@ -0,0 +1,30 @@
+,NAME,hp_current,hp_current_prop_b,hp_5p,hp_5p_prop_b,hp_15p,hp_15p_prop_b
+0,Aleutians East Borough,0,0,27,0.039764359,339.5,0.5
+1,Aleutians West Census Area,0,0,96,0.072126221,665.5,0.5
+2,Anchorage Municipality,134,0.001129467,217,0.001829063,4195,0.035359069
+3,Bethel Census Area,0,0,25,0.004177807,1247,0.208389037
+4,Bristol Bay Borough,0,0,13,0.015169195,428.5,0.5
+5,Denali Borough,2,0.00124533,19,0.011830635,803,0.5
+6,Dillingham Census Area,0,0,52,0.021621622,1202.5,0.5
+7,Fairbanks North Star Borough,124,0.002906772,147,0.003445932,1248,0.029255257
+8,Haines Borough,59,0.043065693,125,0.091240876,685,0.5
+9,Hoonah-Angoon Census Area,77,0.043331458,192,0.108047271,888.5,0.5
+10,Juneau City and Borough,2000,0.141914426,2192,0.15553821,7046.5,0.5
+11,Kenai Peninsula Borough,66,0.002023175,204,0.006253449,6568,0.201336521
+12,Ketchikan Gateway Borough,600,0.090716662,875,0.132295132,3307,0.5
+13,Kodiak Island Borough,0,0,418,0.071759657,2912.5,0.5
+14,Kusilvak Census Area,0,0,18,0.007685739,887,0.378736123
+15,Lake and Peninsula Borough,0,0,48,0.032675289,734.5,0.5
+16,Matanuska-Susitna Borough,66,0.001292167,196,0.003837344,6148,0.120367289
+17,Nome Census Area,0,0,11,0.002686203,549,0.134065934
+18,North Slope Borough,0,0,1,0.000381825,64,0.024436808
+19,Northwest Arctic Borough,76,0.027972028,79,0.029076187,230,0.08465219
+20,Petersburg Borough,72,0.041909197,197,0.114668219,859,0.5
+21,Prince of Wales-Hyder Census Area,50,0.015494267,347,0.107530214,1613.5,0.5
+22,Sitka City and Borough,400,0.096641701,745,0.179995168,2069.5,0.5
+23,Skagway Municipality,53,0.069736842,78,0.102631579,380,0.5
+24,Southeast Fairbanks Census Area,8,0.002277256,14,0.003985198,320,0.091090236
+25,Valdez-Cordova Census Area,0,0,54,0.008525418,2662,0.42027155
+26,Wrangell City and Borough,54,0.042419482,128,0.100549882,636.5,0.5
+27,Yakutat City and Borough,0,0,6,0.013605442,220.5,0.5
+28,Yukon-Koyukuk Census Area,0,0,0,0,20,0.004982561

Code/Visualization/data/state_estimates.csv

@@ -0,0 +1,19 @@
+,estimate
+total_NPV_current,44083026.58119607
+total_NPV_5p,62133527.08377122
+total_NPV_15p,169797294.20078242
+total_CO2_lbs_current,20009641.53102849
+total_CO2_lbs_5p,30881736.723766647
+total_CO2_lbs_15p,144078754.0879393
+total_CO2_miles_current,22435628.473496106
+total_CO2_miles_5p,34625866.269336246
+total_CO2_miles_15p,161546991.86532924
+total_Heating_Days_Above5_current,1324617.1041045208
+total_Heating_Days_Above5_5p,2223757.7490732973
+total_Heating_Days_Above5_15p,15341798.642153036
+total_Heating_Days_Below5_current,26068.962769193502
+total_Heating_Days_Below5_5p,36617.412259981444
+total_Heating_Days_Below5_15p,460492.9352905859
+total_Heating_Days_Covered_current,0.9806994656948431
+total_Heating_Days_Covered_5p,0.9838002943553924
+total_Heating_Days_Covered_15p,0.9708591040081871

Code/Visualization/modules/vis_borough_data.R

@@ -0,0 +1,210 @@
+require(dplyr)
+require(forcats)
+
+## Load borough-level estimates
+borough_projs <- read.csv("data/borough_estimates_weightedByHuFuel.csv")
+
+## Generate a toy dataframe (simple average within boroughs); replace with real one later
+borough_projs_final <- 
+  borough_projs %>% 
+  rename(name = census_area,
+         NPV = NPV_weighted2,
+         CO2_lbs_saved = CO2_lbs_saved_weighted2,
+         CO2_driving_miles_saved = CO2_driving_miles_saved_weighted2) %>%
+  mutate(name = case_when(name == "Wade Hampton Census Area" ~ "Kusilvak Census Area",
+                          name == "Anchorage municipality" ~ "Anchorage Municipality",
+                          TRUE ~ name),
+         Rebate_dol = case_when(Rebate_dol == 0 ~ "current",
+                                Rebate_dol == 2000 ~ "mid",
+                                Rebate_dol == 4000 ~ "high"),
+         Fuel_Esc_Rate = case_when(Fuel_Esc_Rate == 0.03 ~ "current",
+                                  Fuel_Esc_Rate == 0.06 ~ "mid",
+                                  Fuel_Esc_Rate == 0.09 ~ "high"),
+         Temp_Projection = case_when(Temp_Projection == 0 ~ "current",
+                                     Temp_Projection == 1.5 ~ "mid",
+                                     Temp_Projection == 3 ~ "high")) 
+
+## Write a function to subset data based on scenarios and visualize tilegram
+vis_borough_proj <- function(outcome = "NPV", 
+                             Rebate_dol = "current", 
+                             Fuel_Esc_Rate = "current", 
+                             Temp_Projection = "current") {
+
+  ## Subset data based on combination of 3x3x3 projections 
+  borough_projs_subset <- 
+    borough_projs_final %>%
+    filter(Rebate_dol == !!Rebate_dol &
+             Fuel_Esc_Rate == !!Fuel_Esc_Rate &
+              Temp_Projection == !!Temp_Projection)
+
+  ## Merge adoption data w/ tilegram layout (from 'make_tilegram.R')
+  borough_projs_df <- merge(borough_projs_subset, all_grids, by = 'name')
+
+  ## Rasterize, then convert tile into polygon
+  borough_projs_raster <- 
+    borough_projs_df %>%
+    mutate(name2 = as.integer(factor(name))) %>% # coerce name into integer-factor
+    mutate(row = -row) %>%
+    dplyr::select(col, row, name2) %>%           # select X, Y (for coords), Z (for value)
+    rasterFromXYZ()
+
+  borough_projs_polygon <- 
+    borough_projs_raster %>%
+    rasterToPolygons(dissolve = TRUE) %>%
+    st_as_sf()
+
+  ## Make hover text and title depending on outcome 
+
+  customize_input <- list()
+
+  if (outcome == "NPV") {
+
+      customize_input$hover_text <- ~ paste("Borough:", name, "<br>Net Present Value:", round(NPV))
+      customize_input$legend_title <- "<b>Cost-Effectiveness</b>\n<i>Net Present Value</i>"
+      customize_input$palette <- "RdBu"
+      customize_input$limit <- c(-20000, 20000)
+
+    } else if (outcome == "CO2_lbs_saved") {
+
+      customize_input$hover_text <- ~ paste("Borough:", name, "<br>CO2 Saved:", round(CO2_lbs_saved), "lbs")
+      customize_input$legend_title <- "<b>CO2 Saved<br>(lbs)</b>"
+      customize_input$palette <- "YlGn"
+      customize_input$limit <- NULL
+
+    } else if (outcome == "CO2_driving_miles_saved") {
+
+      customize_input$hover_text <- ~ paste("Borough:", name, "<br>CO2 Saved:", round(CO2_driving_miles_saved), "driving miles")
+      customize_input$legend_title <- "<b>CO2 Saved<br>(in Driving Miles)</b>"
+      customize_input$palette <- "BuGn"
+      customize_input$limit <- NULL
+
+    } else if (outcome == "Heating_Days_Covered") {
+
+      customize_input$hover_text <- ~ paste0("Borough: ", name, "<br>Heating Days Covered: ", round(Heating_Days_Covered, 3)*100, "%")
+      customize_input$legend_title <- '<b>% of Heating Days Covered</b>'
+      customize_input$palette <- "YlOrRd"
+      customize_input$limit <- c(min(borough_projs_final$Heating_Days_Covered), max(borough_projs_final$Heating_Days_Covered))
+
+    }
+
+  ## Make tilegram
+  borough_proj_tilegram <- 
+    # create a plotly object 
+    plot_ly() %>%
+    # add the scatter layer 
+    add_trace(
+      data = borough_projs_df,
+      type = "scatter",
+      mode = 'markers',
+      marker = list(symbol = 'square', size = 9.25),
+      y = ~ (-row),
+      x = ~ col,
+      color = as.formula(paste0("~", outcome)),
+      colors = customize_input$palette,
+      hoverinfo = "text",
+      text = customize_input$hover_text
+    ) %>%
+    # add boundaries 
+    add_sf( 
+      data = borough_projs_polygon,
+      size = I(1),
+      fill = I("transparent"),
+      color = I("black"),
+      hoverinfo = 'skip'
+    ) %>%
+    # customize legend 
+    colorbar(
+      limits = customize_input$limit,
+      title = customize_input$legend_title,
+      orientation = 'h',
+      len = 0.85
+    ) %>%
+    # layout adjustments
+    layout( 
+      xaxis = list(title = ""),
+      yaxis = list(title = "", fixedrange = FALSE, showgrid = FALSE, showline = FALSE, showticklabels = FALSE),
+      showlegend = FALSE,
+      autosize = TRUE,
+      plot_bgcolor = '#D8DEE9',
+      paper_bgcolor = '#D8DEE9',
+      width = 650,
+      height = 650,
+      margin = list(l = 0, r = 0, b = 0, t = 0, pad = 0)
+    )
+
+  return(borough_proj_tilegram)
+
+}
+
+# Write another function to subset data based on scenarios and create bar plots/tables
+vis_borough_barplot <- function(outcome = "NPV",
+                                 Rebate_dol = "current",
+                                 Fuel_Esc_Rate = "current",
+                                 Temp_Projection = "current") {
+
+  ## Subset data based on combination of 3x3x3 projections
+  borough_projs_subset <-
+    borough_projs_final %>%
+    filter(Rebate_dol == !!Rebate_dol &
+             Fuel_Esc_Rate == !!Fuel_Esc_Rate &
+              Temp_Projection == !!Temp_Projection) %>%
+    mutate(name = str_trim(str_remove(name, "Borough|City and Borough|Census Area|Municipality")))
+
+  ## customize data and input 
+  customize_input <- list()
+
+  ## Rank boroughs by outcome
+  if (outcome == "Heating_Days_Covered") {
+
+    ## for Heating_Days_Covered, get bottom 10 boroughs
+    borough_projs_top <- 
+      borough_projs_subset %>%
+      mutate(name = fct_reorder(name, !!sym(outcome))) %>%
+      arrange(desc(name)) %>%
+      slice_tail(n = 10)
+
+    customize_input$xaxis <- list(range = c(0.55, 0.95), title = "% of Heating Days Covered")
+    customize_input$height <- 205
+
+  } else if (outcome == "NPV") {
+
+    borough_projs_top <- 
+      borough_projs_subset %>%
+      mutate(name = fct_reorder(name, !!sym(outcome))) %>%
+      arrange(desc(name)) %>%
+      slice_head(n = 10)
+
+    customize_input$height <- 205
+
+  } else if (outcome %in% c("CO2_lbs_saved", "CO2_driving_miles_saved")) {
+
+    borough_projs_top <- 
+      borough_projs_subset %>%
+      mutate(name = fct_reorder(name, !!sym(outcome))) %>%
+      arrange(desc(name)) 
+
+    customize_input$xaxis <- list(title = ifelse(outcome == "CO2_lbs_saved", "CO2 saved (in lbs)", "CO2 saved (driving miles)"))
+    customize_input$height <- 450
+  }
+
+  ## Plotly bar plot  
+  borough_proj_barplot <- 
+    plot_ly(
+      data = borough_projs_top,
+      x = as.formula(paste0("~", outcome)),
+      y = ~ name
+    ) %>%
+    layout( 
+      xaxis = customize_input$xaxis,
+      yaxis = list(title = ""),
+      height = customize_input$height,
+      plot_bgcolor = '#D8DEE9',
+      paper_bgcolor = '#D8DEE9',
+      margin = list(l = 0, r = 0, b = 0, t = 0, pad = 0)
+    )
+
+  return(borough_proj_barplot)
+
+}
+
+

Code/Visualization/modules/vis_state_data.R

@@ -0,0 +1,108 @@
+require(dplyr)
+require(stringr)
+require(tidyr)
+require(plotly)
+require(raster)
+require(sf)
+
+##### For the three aggregate outcomes #####
+
+## Get state-level aggregates
+state_estimates <- read.csv("data/state_estimates.csv")
+
+## Convert it to wide format
+state_estimates_wide <- 
+  state_estimates %>%
+  filter(str_detect(X, "NPV|CO2|Covered")) %>%
+  separate(X, into = c("total", "outcome", "proj_rate"), sep = "(?<=total)_|_(?=current|5|15)") %>%
+  mutate(proj_rate = paste0("hp_", proj_rate)) %>%
+  pivot_wider(names_from = outcome, values_from = estimate)
+
+##### For tilegram visualization #####
+
+## Read borough 
+adopt_proj <- read.csv("data/adoption_scenarios_5_and_15.csv")
+names(adopt_proj)[2] <- 'name' 
+
+## Merge adoption data w/ tilegram layout (from 'make_tilegram.R')
+adopt_proj_df <- merge(adopt_proj, all_grids, by = 'name')
+
+## Rasterize, then convert tile into polygon
+adopt_proj_raster <- 
+  adopt_proj_df %>%
+  mutate(name2 = as.integer(factor(name))) %>% # coerce name into integer-factor
+  mutate(row = -row) %>%
+  dplyr::select(col, row, name2) %>%           # select X, Y (for coords), Z (for value)
+  rasterFromXYZ()
+
+adopt_proj_polygon <- 
+  adopt_proj_raster %>%
+  rasterToPolygons(dissolve = TRUE) %>%
+  st_as_sf()
+
+## Write a function to visualize adoption rates
+
+vis_adopt_proj <- function(proj_rate, absolute = TRUE) {
+
+  ## Get absolute and relative numbers of HP
+  hp_absolute <- 
+    adopt_proj_df %>% 
+    pull(proj_rate)
+
+  hp_percentage <- 
+    adopt_proj_df %>% 
+    pull(!!sym(paste0(proj_rate, "_prop_b"))) %>% 
+    round(2) * 100
+
+  ## Visualize adoption projection
+  adopt_proj_plotly <- 
+    # create a plotly object 
+    plot_ly() %>%
+    # add the scatter layer 
+    add_trace(
+      data = adopt_proj_df,
+      type = "scatter",
+      mode = 'markers',
+      marker = list(symbol = 'square', size = 9.75),
+      y = ~ (-row),
+      x = ~ col,
+      color = if (absolute == TRUE) as.formula(paste0("~", proj_rate)) else as.formula(paste0("~", paste0(proj_rate, "_prop_b"))),
+      colors = "YlOrRd",
+      hoverinfo = "text",
+      text = ~ paste0("Borough: ", name, "<br>", 
+                      "Number of heat pumps: ", hp_absolute, "<br>", 
+                      "Percentage out of households in borough: ", hp_percentage, "%")
+    ) %>%
+    # add boundaries 
+    add_sf( 
+      data = adopt_proj_polygon,
+      size = I(1),
+      fill = I("transparent"),
+      color = I("black"),
+      hoverinfo = 'skip'
+    ) %>%
+    # customize legend 
+    colorbar(
+      limits = if (absolute == TRUE) c(0,7050) else c(0,0.5), 
+      title = if (absolute == TRUE) "<b>No. of Heat Pumps</b>" else "<b>Prop. w/ Heat Pumps</b>",
+      orientation = 'h',
+      len = 0.85
+    ) %>%
+    # layout adjustments
+    layout( 
+      xaxis = list(title = ""),
+      yaxis = list(title = "", fixedrange = FALSE, showgrid = FALSE, showline = FALSE, showticklabels = FALSE),
+      showlegend = FALSE,
+      autosize = TRUE,
+      plot_bgcolor = '#D8DEE9',
+      paper_bgcolor = '#D8DEE9',
+      width = 650,
+      height = 650,
+      margin = list(l = 0, r = 0, b = 0, t = 0, pad = 0)
+    )
+
+  return(adopt_proj_plotly)
+
+}
+
+vis_adopt_proj(proj_rate = "hp_15p")

Code/Visualization/rsconnect/shinyapps.io/dssg2023heatpump/Code.dcf

@@ -0,0 +1,13 @@
+name: Code
+title: Code
+username: dssg2023heatpump
+account: dssg2023heatpump
+server: shinyapps.io
+hostUrl: https://api.shinyapps.io/v1
+appId: 9602209
+bundleId: 7565395
+url: https://dssg2023heatpump.shinyapps.io/Code/
+version: 1
+asMultiple: FALSE
+asStatic: FALSE
+ignoredFiles: legacy/.Rhistory|legacy/generate_ak_cartogram.R

Code/Visualization/rsconnect/shinyapps.io/dssg2023heatpump/VisApp.dcf

@@ -0,0 +1,12 @@
+name: VisApp
+title: VisApp
+username: dssg2023heatpump
+account: dssg2023heatpump
+server: shinyapps.io
+hostUrl: https://api.shinyapps.io/v1
+appId: 9602329
+bundleId: 7575885
+url: https://dssg2023heatpump.shinyapps.io/VisApp/
+version: 1
+asMultiple: FALSE
+asStatic: FALSE