edit figure 5

analysis/paper/paper.Rmd

@@ -816,13 +816,52 @@ pdat_final_demand = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
   mutate(total_blah = total_energy_use_tj_new/33417583,
          total_blah = total_blah*1000)
 
-ineq_curr_check = pdat_final_demand$total_blah[10]/pdat_final_demand$total_blah[1]
+ineq_curr = pdat_final_demand$total_blah[10]/pdat_final_demand$total_blah[1]
 
 
 
-df_energy_deciles = read_csv(here("analysis/data/raw/energy_figure_data.csv")) 
+#df_energy_deciles = read_csv(here("analysis/data/raw/energy_figure_data.csv")) 
 
-ineq_curr = df_energy_deciles$pae_energy_use_gj[10]/df_energy_deciles$pae_energy_use_gj[1]
+#ineq_curr = df_energy_deciles$pae_energy_use_gj[10]/df_energy_deciles$pae_energy_use_gj[1]
+
+```
+
+```{r , eval = FALSE}
+# vectorized function that returns scaled quantiles given 
+#quantile index column and column with quantile averages
+#qidx: quantile index
+#qavg_to_scale: column to scale
+#first_target: target value of first quantile
+#mean_target: target mean of scaled quantiles
+scaled_quantiles <- function(.data, 
+                          qidx, 
+                          qavg_to_scale, 
+                          first_target, 
+                          mean_target) {
+  
+  # cumbersomely extract current quantile mean
+  mean_current = .data %>%
+    ungroup() %>%
+    summarise(mean_cur = first(mean({{qavg_to_scale}}))) %>%
+    pull(mean_cur)
+    
+  # cumbersomely extract current first wuantile value
+  first_current = .data %>%
+    ungroup() %>%
+    arrange({{qidx}}) %>%
+    summarise(first_cur = first({{qavg_to_scale}})) %>%
+    pull(first_cur)
+  
+  df_tmp = .data %>%
+    mutate(tmp = {{qavg_to_scale}}*mean_target/mean_current)
+  
+  first_tmp = df_tmp$tmp[1]
+  
+  df_tmp = df_tmp %>%
+    mutate(scaled = mean_target-(mean_target-tmp) * (mean_target-first_target)/(mean_target-first_tmp)) %>%
+    select({{qidx}}, scaled) %>%
+    mutate(v_mean = mean_target, v_first = first_target)
+}
 
 ```
 
@@ -846,6 +885,36 @@ c_mean_mer = round((mer[1]+mer[2])/2)
 
 ```
 
+```{r , eval = FALSE}
+## run once to save file
+df_all = NULL
+for (min_energy in seq(from=mer[1], to=mer[2], by=1)) {
+  for (mean_energy in seq(from=mea[1], to=mea[2], by=1)) {
+    if (min_energy <= mean_energy) {
+      df_all = df_all %>%
+      bind_rows(df_energy_deciles %>%
+                  scaled_quantiles(eu_q_rank, pae_energy_use_gj, min_energy, mean_energy))
+    }
+    
+  }
+}
+
+# thought it was much faster without the bind_rows, turns out that is not the bottleneck
+#n_quantiles = nrow(df_energy_deciles)
+#mea_range = seq(from=mea[1], to=mea[2], by=1)
+#mer_range = seq(from=mer[1], to=mer[2], by=1)
+#df_all2 = data.frame(eu_q_rank = rep(seq(1:n_quantiles),
+#                                     times=length(mer_range)*length(mea_range)),
+#                     v_mean = rep(mea_range, each=n_quantiles, times=length(mer_range)),
+#                     v_first = rep(mer_range, each=length(mea_range)*n_quantiles),
+#                     scaled = 0) %>%
+#  filter(v_first <= v_mean)
+
+
+
+#saveRDS(df_all, here("data/scenarios_fine.rds"))
+```
+
 ```{r , fig.width=7, fig.height=5.5}
 
 round_by = 10
@@ -872,14 +941,14 @@ a = df_all %>%
   geom_tile() +
   geom_hline(yintercept = ineq_curr, alpha=0.8, color="grey", linetype=2) +
   geom_line(data=df_scenario, aes(color=scenario, group=scenario)) +
-  annotate(geom="text", x=max(df_all$v_first)-5.7,y=ineq_curr+0.6,label = "Current 10:10 ratio") +
+  annotate(geom="text", x=max(df_all$v_first)-5.7,y=ineq_curr+0.6,label = "Current (2015) 10:10 ratio") +
   #scale_fill_viridis("Mean energy\navailable") +
   scale_fill_gradient("Mean energy\navailable (GJ/ae)",
                       low=wes_palette("Chevalier1")[3], 
                       high = wes_palette("Rushmore1")[4]) +
   scale_color_manual(values=wes_palette("Darjeeling1")) +
   theme_minimal() +
-  labs(x="Minimum energy requirement (GJ/ae)", y="Maximum inequality (10:10 ratio)", color = "Scenario")+
+  labs(x="Minimum energy requirement (GJ/ae)", y="Maximum energy inequality (10:10 ratio)", color = "Scenario")+
   theme(text=element_text(family="Liberation Sans Narrow"),
         axis.text.x = element_text(size = 13),
         axis.text.y = element_text(size = 13))# +