diff --git a/habitats/air_tmp_vulnerability.Rmd b/habitats/air_tmp_vulnerability.Rmd new file mode 100644 index 0000000..24d786f --- /dev/null +++ b/habitats/air_tmp_vulnerability.Rmd @@ -0,0 +1,119 @@ +--- +title: "Untitled" +output: html_document +date: "2024-09-12" +--- + +Going to use current pressure data and habitat distributions to determine how vulnerable some key habitats are to the air temperature pressure. + +```{r setup, include=FALSE} +library(tidyverse) +library(terra) +library(here) +library(countrycode) + +``` + +```{r} +heat <- rast("/home/shares/ohi/stressors_2021/combining_pressures/rescaled_pressures/air-heat-index_ssp245_current.tif") +seagrass <- rast(here("habitats/data/seagrass.tif")) + +stack <- c(heat, seagrass) +stack_df <- data.frame(stack) + +seagrass <- stack_df %>% + filter(!is.na(mean)) %>% + mutate(seagrass = ifelse(is.na(seagrass), 0, seagrass)) + +Hmisc::wtd.quantile(seagrass$mean, weights=seagrass$seagrass, probs= c(0.001, 0.01, 0.1, 0.5, 0.90, 0.99, 0.999)) + +tmp <- seagrass %>% + group_by(seagrass) %>% + mutate(mean = ifelse(mean >0.21, 1, 0)) %>% + mutate(seagrass = ifelse(seagrass > 0, 1, 0)) +table(tmp) + +# > 9855/205281 +#[1] 0.04800737 +#> 254*0.048 +#[1] 12.192 +#> (1-(9855/205281))^254 +#[1] 3.740515e-06 + +ggplot(seagrass, aes(x=mean, y=seagrass)) + + geom_jitter(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red") + + ylab("seagrass cover") + + xlab("heat-index pressure") + + + +mangrove <- rast(here("habitats/data/mangroves.tif")) + +stack <- c(heat, mangrove) +stack_df <- data.frame(stack) + +mangrove <- stack_df %>% + filter(!is.na(mean)) %>% + mutate(mangroves = ifelse(is.na(mangroves), 0, mangroves)) + + +ggplot(mangrove, aes(x=mean, y=mangroves)) + + geom_point(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red") + + ylab("mangrove cover") + + xlab("heat-index pressure") + + +ice <- rast(here("habitats/data/ice.tif")) + +stack <- c(heat, ice) +stack_df <- data.frame(stack) + +ice <- stack_df %>% + filter(!is.na(mean)) %>% + mutate(ice = ifelse(is.na(ice), 0, ice)) + +Hmisc::wtd.quantile(ice$mean, weights=ice$ice, probs= c(0.001, 0.01, 0.1, 0.5, 0.90, 0.99, 0.999)) +mean(ice$mean[ice$ice>0]) +ggplot(ice, aes(x=mean, y=ice)) + + geom_point(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red")+ + ylab("seaice cover") + + xlab("heat-index pressure") + + +marsh <- rast(here("habitats/data/salt-marsh.tif")) + +stack <- c(heat, marsh) +stack_df <- data.frame(stack) + +marsh <- stack_df %>% + filter(!is.na(mean)) %>% + mutate(salt_marsh = ifelse(is.na(salt_marsh), 0, salt_marsh)) + +Hmisc::wtd.quantile(marsh$mean, weights=marsh$salt_marsh, probs= c(0.001, 0.01, 0.1, 0.5, 0.90, 0.99, 0.999)) + +tmp <- marsh %>% + group_by(salt_marsh) %>% + mutate(mean = ifelse(mean >0.20, 1, 0)) %>% + mutate(salt_marsh = ifelse(salt_marsh > 0, 1, 0)) +table(tmp) +#> 4997/210048 +#[1] 0.0237898 +#> 345*0.0237898 +#[1] 8.207481 +#> (1 - 4997/210048)^345 +#[1] 0.0002468611 + +ggplot(marsh, aes(x=mean, y=salt_marsh)) + + geom_point(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red") + + ylab("saltmarsh cover") + + xlab("heat-index pressure") + +``` \ No newline at end of file diff --git a/habitats/data/beach.tif b/habitats/data/beach.tif new file mode 100644 index 0000000..f8220b3 Binary files /dev/null and b/habitats/data/beach.tif differ diff --git a/habitats/data/coral-reef.tif b/habitats/data/coral-reef.tif new file mode 100644 index 0000000..78d5e61 Binary files /dev/null and b/habitats/data/coral-reef.tif differ diff --git a/habitats/data/d-h-bottom.tif b/habitats/data/d-h-bottom.tif new file mode 100644 index 0000000..1806563 Binary files /dev/null and b/habitats/data/d-h-bottom.tif differ diff --git a/habitats/data/d-s-benthic.tif b/habitats/data/d-s-benthic.tif new file mode 100644 index 0000000..1743591 Binary files /dev/null and b/habitats/data/d-s-benthic.tif differ diff --git a/habitats/data/deep-waters.tif b/habitats/data/deep-waters.tif new file mode 100644 index 0000000..de74944 Binary files /dev/null and b/habitats/data/deep-waters.tif differ diff --git a/habitats/data/hard-shelf.tif b/habitats/data/hard-shelf.tif new file mode 100644 index 0000000..3eb07b3 Binary files /dev/null and b/habitats/data/hard-shelf.tif differ diff --git a/habitats/data/hard-slope.tif b/habitats/data/hard-slope.tif new file mode 100644 index 0000000..8aa29e6 Binary files /dev/null and b/habitats/data/hard-slope.tif differ diff --git a/habitats/data/ice.tif b/habitats/data/ice.tif new file mode 100644 index 0000000..0a0e48d Binary files /dev/null and b/habitats/data/ice.tif differ diff --git a/habitats/data/inttidalmud.tif b/habitats/data/inttidalmud.tif new file mode 100644 index 0000000..3e3a35b Binary files /dev/null and b/habitats/data/inttidalmud.tif differ diff --git a/habitats/data/kelp.tif b/habitats/data/kelp.tif new file mode 100644 index 0000000..78be28b Binary files /dev/null and b/habitats/data/kelp.tif differ diff --git a/habitats/data/mangroves.tif b/habitats/data/mangroves.tif new file mode 100644 index 0000000..e86bd3b Binary files /dev/null and b/habitats/data/mangroves.tif differ diff --git a/habitats/data/rky-intidal.tif b/habitats/data/rky-intidal.tif new file mode 100644 index 0000000..1fa813e Binary files /dev/null and b/habitats/data/rky-intidal.tif differ diff --git a/habitats/data/rocky-reef.tif b/habitats/data/rocky-reef.tif new file mode 100644 index 0000000..38843f2 Binary files /dev/null and b/habitats/data/rocky-reef.tif differ diff --git a/habitats/data/s-t-s-bottom.tif b/habitats/data/s-t-s-bottom.tif new file mode 100644 index 0000000..c44aecb Binary files /dev/null and b/habitats/data/s-t-s-bottom.tif differ diff --git a/habitats/data/salt-marsh.tif b/habitats/data/salt-marsh.tif new file mode 100644 index 0000000..8a304d8 Binary files /dev/null and b/habitats/data/salt-marsh.tif differ diff --git a/habitats/data/seagrass.tif b/habitats/data/seagrass.tif new file mode 100644 index 0000000..a7dad9a Binary files /dev/null and b/habitats/data/seagrass.tif differ diff --git a/habitats/data/seamounts.tif b/habitats/data/seamounts.tif new file mode 100644 index 0000000..d4def35 Binary files /dev/null and b/habitats/data/seamounts.tif differ diff --git a/habitats/data/soft-shelf.tif b/habitats/data/soft-shelf.tif new file mode 100644 index 0000000..86fb843 Binary files /dev/null and b/habitats/data/soft-shelf.tif differ diff --git a/habitats/data/soft-slope.tif b/habitats/data/soft-slope.tif new file mode 100644 index 0000000..956acab Binary files /dev/null and b/habitats/data/soft-slope.tif differ diff --git a/habitats/data/surface-waters.tif b/habitats/data/surface-waters.tif new file mode 100644 index 0000000..ed7be84 Binary files /dev/null and b/habitats/data/surface-waters.tif differ diff --git a/habitats/data/suspension-reef.tif b/habitats/data/suspension-reef.tif new file mode 100644 index 0000000..bf5dd7a Binary files /dev/null and b/habitats/data/suspension-reef.tif differ diff --git a/habitats/habitat_n.tif b/habitats/habitat_n.tif new file mode 100644 index 0000000..5d89b3c Binary files /dev/null and b/habitats/habitat_n.tif differ diff --git a/habitats/sst_vulnerability.Rmd b/habitats/sst_vulnerability.Rmd new file mode 100644 index 0000000..536c2e8 --- /dev/null +++ b/habitats/sst_vulnerability.Rmd @@ -0,0 +1,96 @@ +--- +title: "Untitled" +output: html_document +date: "2024-09-12" +--- + +Going to use current pressure data and habitat distributions to determine how vulnerable some key habitats are to the air temperature pressure. + +```{r setup, include=FALSE} +library(tidyverse) +library(terra) +library(here) +library(countrycode) + +``` + +```{r} +sst <- rast("/home/shares/ohi/stressors_2021/_dataprep/SST/five_year_average_final/sst_avg_ssp245_2015_2019.tif") +sst <- sst-273.15 + +kelp <- rast(here("habitats/data/kelp.tif")) + +stack <- c(sst, kelp) +stack_df <- data.frame(stack) + +kelp <- stack_df %>% + filter(!is.na(focal_mean)) %>% + mutate(kelp = ifelse(is.na(kelp), 0, kelp)) + +Hmisc::wtd.quantile(kelp$focal_mean, weights=kelp$kelp,probs= c(0.001, 0.01, 0.1, 0.5, 0.90, 0.99, 0.999)) + +ggplot(kelp, aes(x=focal_mean, y=kelp)) + + geom_point(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red") + + ylab("kelp") + + xlab("SST") + + +coral <- rast(here("habitats/data/coral-reef.tif")) + +stack <- c(sst, coral) +stack_df <- data.frame(stack) + +coral <- stack_df %>% + filter(!is.na(focal_mean)) %>% + mutate(coral_reef = ifelse(is.na(coral_reef), 0, coral_reef)) + +Hmisc::wtd.quantile(coral$focal_mean, weights=coral$coral_reef, probs= c(0.001, 0.01, 0.1, 0.5, 0.90, 0.95, 0.99, 0.999)) + +ggplot(coral, aes(x=focal_mean, y=coral_reef)) + + geom_point(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red") + + ylab("coral") + + xlab("SST") + + +marsh <- rast(here("habitats/data/salt-marsh.tif")) + +stack <- c(sst, marsh) + +stack_df <- data.frame(stack) + +marsh <- stack_df %>% + filter(!is.na(focal_mean)) %>% + mutate(salt_marsh = ifelse(is.na(salt_marsh), 0, salt_marsh)) + +Hmisc::wtd.quantile(marsh$focal_mean, weights=marsh$salt_marsh, probs= c(0.001, 0.01, 0.1, 0.5, 0.90, 0.99, 0.999)) + +ggplot(marsh, aes(x=focal_mean, y=salt_marsh)) + + geom_point(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red") + + ylab("marsh") + + xlab("SST") + + +seagrass <- rast(here("habitats/data/seagrass.tif")) + +stack <- c(sst, seagrass) + +stack_df <- data.frame(stack) + +seagrass <- stack_df %>% + filter(!is.na(focal_mean)) %>% + mutate(seagrass = ifelse(is.na(seagrass), 0, seagrass)) + +Hmisc::wtd.quantile(seagrass$focal_mean, weights=seagrass$seagrass, probs= c(0.001, 0.01, 0.1, 0.5, 0.90, 0.99, 0.999)) + +ggplot(seagrass, aes(x=focal_mean, y=seagrass)) + + geom_point(size=0.2, alpha=0.5, height=0.1) + + geom_smooth(method = "gam", se = TRUE) + + geom_hline(yintercept=0, color="red") + + ylab("seagrass") + + xlab("SST")