Snow Water Equivalent Vs ASI

Afghanistan - A closer look

Intro

In code below we try to assess if monthly Snow Water Equivalent (SWE) data has any predictive power on Agriculture Stress Index (ASI) over the March-April-May (MAM) planting season in Afghanistan.

We run the analysis for each province in Afghanistan.

We look at monthly SWE anomaly against MAM ASI: mean, max, max z-score, and max anomaly

Results

Some provinces do show relatively strong correlations between winter SWE (Dec, Jan) and MAM ASI parameters as expected. Some provinces with the stronger observed correlations include Nimroz, Nuristan, Kabul, Kunduz.

However, Faryab, the province of interest does not see a very strong positive correlation

These results concur with the previous analysis of Faryab province that indicated a low correlation/predictive between driest SWE years and least healthy ASI years within

box::use(
  ../R/blob_connect
)

box::use(
  dplyr[...],
  forcats[...],
  ggplot2[...],
  gghdx[...],
  janitor[...],
  lubridate[...],
  purrr[...],
  readr[...],
  scales[...],
  # stats[...],
  stringr[...],
  glue[...],
  tidyr[...]
)

gghdx()


df_snow <- blob_connect$read_blob_file("DF_ADM1_MODIS_SNOW") |> 
  filter(
    # ADM1_NAME == "Faryab",
    parameter == "NDSI_Snow_Cover_mean"
  ) |> 
  clean_names()


df_snow_proc <- df_snow|> 
  # filter(adm1_name == "Faryab") |> 
  group_by(
    adm1_name, adm1_code,
    mo = month(date)
  ) |> 
  mutate(
    swe_anom = (value - mean(value)) ,
    swe_z = swe_anom /sd(value)
  ) |> 
  ungroup()


df_asi <- read_csv("https://www.fao.org/giews/earthobservation/asis/data/country/AFG/MAP_ASI/DATA/ASI_Dekad_Season1_data.csv") |> 
  clean_names()

df_asi_yr <- df_asi |> 
  mutate(
    mo = month(date)
  ) |> 
  filter(
    mo %in% c(3,4,5)
  ) |>
  group_by(
    adm1_code, mo
  ) |> 
  mutate(
    asi_anom = data-mean(data),
    asi_z = asi_anom/sd(data)
  ) |> 
  group_by(
    adm1_code,
    
    yr_date = floor_date(date,"year")
  ) |> 
  summarise(
    asi_max = max(data),
    asi_mean = mean(data),
    asi_max_anom = max(asi_anom),
    asi_max_z = max(asi_z)
    
  ) |> 
  mutate(
    yr = year(yr_date)
  ) |> 
  filter(yr >=2001)


df_snow_long <- df_snow_proc |> 
  select(
    adm1_name,
    adm1_code,
    date, 
    value,
    swe_z,
    swe_anom
  ) |> 
  mutate(
    yr = year(date),
    mo = month(date), 
    yr = ifelse(mo> 5, yr + 1, yr),
    mo = month(mo, abbr = T,label =T)
  ) |> 
    filter(yr>2000)
  
df_snow_monthly_wide <- df_snow_long |>
  pivot_wider(id_cols = c("adm1_name","adm1_code","yr"),
    names_from = mo, values_from = swe_anom
  )


# box::use(corrr[...])
df_swe_asi <- left_join(df_asi_yr,df_snow_monthly_wide, by = c("adm1_code","yr")  )

 df_swe_asi_filt <- df_swe_asi |> 
  filter(yr< 2024) 
  
df_monthly_corrs <- df_swe_asi_filt |> 
  group_by(
    adm1_code, adm1_name
    ) |> 
  summarise(
    across(
      .cols = c("Jun","Jul","Aug","Sep","Oct","Nov","Dec","Jan","Feb","Mar","Apr","May"),
      .fns = ~cor(.,asi_max)
           ),.groups="drop"
  )
asi_stats = list(max = "asi_max",
                 `anomaly` = "asi_max_anom",
                 `z-score`= "asi_max_z",
                 mean = "asi_mean")
ldf_corrs_wide <- asi_stats |> 
  map(
    \(stat_temp){
      df_swe_asi_filt |> 
  group_by(
    adm1_code, adm1_name
    ) |> 
  summarise(
    across(
      .cols = c("Jun","Jul","Aug","Sep","Oct","Nov","Dec","Jan","Feb","Mar","Apr","May"),
      .fns = ~cor(.,!!sym(stat_temp))
           ),.groups="drop"
  )
    }
  )

ldf_corrs_long <- ldf_corrs_wide |> 
  map(
    \(dft){
      dft |> 
      pivot_longer(cols = Jun: May) |> 
        mutate(
          name = fct_relevel(name, "Oct","Nov","Dec","Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep")
        ) 
    }
  )

df_monthly_corrs_long <- df_monthly_corrs |> 
  pivot_longer(cols = Jun: May) |> 
  mutate(
    name = fct_relevel(name, "Oct","Nov","Dec","Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep")
  ) 


lp <- ldf_corrs_long |> 
  imap(
    \(dft,nmt){
      p_title <- glue("Correlation between monthly snow water equivalent (SWE) and ASI ({nmt}) over M-A-M")
                      
      dft |> 
          filter(!name %in% month(6:9, abbr= T, label=T)) |> 
  ggplot(
    aes(x= name,y= adm1_name, fill = value)
  )+
  geom_tile()+
   scale_fill_gradient2(
    low = hdx_hex("tomato-hdx"),      # Color for negative values
    mid = "white",    # Color for zero
    high = hdx_hex("mint-hdx"),   # Color for positive values
    midpoint = 0      # Set midpoint at zero
    
  )+
  geom_tile(
    data= dft |> 
      filter(adm1_name == "Faryab", 
             name %in% c("Oct","Nov","Dec","Jan","Feb","Mar","Apr","May")),
    fill = NA, color ="black", lwd = 1.5
  )+
    
  geom_text(
    aes(label = round(value,2))
  )+
  labs(
    title = p_title,
    subtitle = "Afghanistan by Province",
    y= "Province"
  )+
  theme(
    axis.title.x = element_blank(),
    legend.title= element_blank(),
    plot.title = element_text(size = 12),
    plot.subtitle = element_text(size = 12),
    legend.text = element_text(angle=90)
  )
        
    }
  )

walk(
  lp,
  \(pt){
    plot(pt)
  }
)