From 7b7ae24f20c54b7325c53f5414bfab442008225d Mon Sep 17 00:00:00 2001
From: aaron <aaron@duckpond.ch>
Date: Wed, 13 Mar 2024 00:02:53 +0100
Subject: [PATCH] working on health plot

---
 ironwood_plots.R | 80 +++++++++++++++++++++++++++++++++++-------------
 1 file changed, 59 insertions(+), 21 deletions(-)

diff --git a/ironwood_plots.R b/ironwood_plots.R
index 1665161..2230549 100644
--- a/ironwood_plots.R
+++ b/ironwood_plots.R
@@ -8,14 +8,15 @@ library(dplyr)
 library(leaflet)
 
 ##
-# parse the input data and declare global values
+# parse the input data, declare global values and auxiliary data
 ##
 
 # read a data frame from the ods document
-df <- read_ods("data/ironwood_data_cleaned.ods", sheet = 1)
-# site base location
-site_lat <- "-33.943917"
-site_lon <- "23.507389"
+df <- read_ods(path = "data/ironwood_data_cleaned.ods", 
+               sheet = 1)
+# site base location (to zero in the map)
+population_lat <- "-33.943917"
+population_lon <- "23.507389"
 # vector of condition names corresponding to the health index numbers
 condition_names <- c("healthy", "light damage", 
                      "medium damage", "severe damage", 
@@ -25,8 +26,8 @@ condition_colors <- c("green", "yellow", "orange", "red", "black")
 
 
 ##
-# 1. asses the tree health
-# - create an overview of the populations health
+# 1.) asses the tree health of the entire population
+#     create an overview of the populations health
 ##
 
 # Calculate the percentage of trees in each health condition
@@ -41,16 +42,41 @@ barplot(percentage,
         col = condition_colors,
         border = "black")
 # Adding a legend
-legend("topright", legend = condition_names, fill = condition_colors)
-# Adding a grid
-#grid(nx = NULL, ny = NULL, col = "lightgray", lty = "dotted")
+legend("topright",
+       legend = condition_names, 
+       fill = condition_colors)
 # Adding a box around the plot
 #box()
 # Add labels with the percentage of trees in each bar
-text(x = barplot(percentage, plot = FALSE), y = percentage, labels = paste0(round(percentage, 1), "%"), pos = 3)
+text(x = barplot(percentage, plot = FALSE),
+     y = percentage,
+     labels = paste0(round(percentage, 1), "%"),
+     pos = 3)
 
 ##
-# 2. perform a shapiro-wilk normality test
+# 2. Create a stacked barchart that represents all site and their health data
+##
+
+# Create data
+set.seed(1124)
+
+# for each site i need a table like this:
+#table(df$tree_health_index[1:14])
+
+colnames(health_data) <- paste("Site", seq(1,20), sep=" ")
+rownames(health_data) <- condition_names
+
+# create color palette:
+library(RColorBrewer)
+coul <- brewer.pal(3, "Pastel2") 
+# Transform this data in %
+data_percentage <- apply(data, 2, function(x){x*100/sum(x,na.rm=T)})
+# Make a stacked barplot--> it will be in %!
+barplot(df$tree_health_index, col=coul , border="white", xlab="group")
+
+
+##
+# 3. perform a shapiro-wilk normality test
 # - the goal is to see if the health is normally distributed
 ##
 
@@ -68,26 +94,38 @@ if (p_value < 0.05) {
 }
 
 ##
-# 3. try to fit health and location data in one plot
+# 4. try to fit health and location data in one plot
 ##
 
+# create a subset of the site locations
+sites <- df[complete.cases(df$site_num),]
+# ensure all coordinates are numeric
+sites$site_lat <- as.numeric(sites$site_lat)
+sites$site_lon <- as.numeric(sites$site_lon)
+
 # create a map from our base location
 map <- leaflet() %>%
-  setView(lng = site_lon,
-          lat = site_lat,
+  setView(lng = population_lon,
+          lat = population_lat,
           zoom = 16) %>%
-  addTiles()
-  #addProviderTiles(providers$CartoDB.Positron)
-
-# Add markers with varying color and size based on population health and number of trees
-map <- map %>%
+  addTiles() %>%
   addCircleMarkers(
     data = df,
     lng = ~tree_lon,
     lat = ~tree_lat,
-    radius = 5,
+    radius = 2,
     color = ~condition_colors[tree_health_index+1],
+    opacity = 1,
     fillOpacity = 1
+  ) %>%
+  addCircleMarkers(
+    data = sites,
+    lng = ~site_lon,
+    lat = ~site_lat,
+    radius = 25,
+    fill = FALSE,
+    color = "green",
+    opacity = 0.1
   )
 
 # show map