diff --git a/.gitignore b/.gitignore
index ac6cec0ae6148953f7080d9e8bdbc5772b017dda..3acf079a9890636921e00602124682ed4a556f26 100644
--- a/.gitignore
+++ b/.gitignore
@@ -11,5 +11,6 @@ renv/cache/
 
 # Data files
 data/
+plots/
 R/*/
 R/*.html
\ No newline at end of file
diff --git a/R/03_02_absence_preparation.R b/R/03_02_absence_preparation.R
index b72d132a3cd8bb01ef3f856fcf55e22d510b9f3f..5a344071151a53e0b9d599e8e9e0ded56f61737e 100644
--- a/R/03_02_absence_preparation.R
+++ b/R/03_02_absence_preparation.R
@@ -7,9 +7,13 @@ library(furrr)
 library(terra)
 library(CoordinateCleaner)
 library(sf)
-library(MASS)
 library(geos)
 
+# Stat packages
+library(MASS)
+library(dismo)
+library(spatialEco)
+
 # DB packages
 library(Symobio)
 library(DBI)
@@ -49,18 +53,16 @@ range_maps = st_transform(range_maps, proj_string)
 # geographically close (reproduce sampling biases) but environmentally       #
 # dissimilar (avoid false negatives) to the known occurrences                #
 # ---------------------------------------------------------------------------#
-# Geographic distance
 future::plan("multisession", workers = 16)
 pseudo_absences = furrr::future_map(
   .x = target_species, 
   .options = furrr::furrr_options(seed = 42), 
   .f = function(species){
-    # Prepare occurrences
+    # Prepare species occurrence data #####
     occs_spec = dplyr::filter(occs_final, species == !!species)
     
-    occs_poly = occs_spec %>% 
-      summarise(geometry = st_combine(geometry)) %>%
-      st_concave_hull(ratio = 0.5)
+    occs_combined = occs_spec %>% 
+      summarise(geometry = st_combine(geometry))
     
     occs_bff = occs_spec %>% 
       st_buffer(10000) %>% # Buffer by 10 kilometers as exclusion radius for absence sampling
@@ -70,17 +72,47 @@ pseudo_absences = furrr::future_map(
       rename(species = name_matched)
     
     # Define model/sampling region
-    sampling_region = st_union(range_spec, occs_poly) %>% # Union occs poly and expert range
-      st_buffer(dist = 100000) %>% # Buffer by 100 km to extend sampling region into new environments
-      st_intersection(sa_polygon) %>% # Crop to South America 
-      st_difference(occs_bff) # Crop out buffered presences
+    # 1. Union all point coordinates from range polygon and occurrences records
+    # 2. Build concave polygon from unionen points
+    # 3. Buffer by 100 km to extend sampling region into new environments
+    # 4. Crop by South America outline
+    sampling_region = st_cast(range_spec, "MULTIPOINT") %>% 
+      st_union(occs_combined) %>% 
+      geos::geos_concave_hull(ratio = 0.3) %>% # Concave hull in sf requires geos 3.11
+      st_as_sf() %>% 
+      st_set_crs(proj_string) %>% 
+      st_buffer(dist = 100000) %>% 
+      st_intersection(sa_polygon) 
     
-    # Define KDE
+    # Define spatial KDE
     ref = st_bbox(sampling_region)[c(1,3,2,4)]
     min_extent = min(ref[2]-ref[1], ref[4]-ref[3])
     occs_kde = spatialEco::sf.kde(occs_spec, bw = min_extent/2, res = 10000, ref = ref, standardize = T, scale.factor = 1) %>% 
       crop(sampling_region, mask = T, touches = F)
     
+    # Environmental KDE
+    env_sampling_region = raster::raster(terra::crop(raster_data, sampling_region, mask = T))
+    maxent_fit <- maxent(env_sampling_region, st_coordinates(occs_spec))
+
+# Step 3: Predict Across the Raster
+# Predict environmental similarity (continuous output)
+env_pred <- predict(maxent_model, env_sampling_region)
+    
+    
+    svm_model <- e1071::svm(x = env_spec,
+                            nu = 0.01,
+                            gamma = 0.00005,
+                            type = 'one-classification',
+                            kernel = "radial",
+                            scale = FALSE)
+    
+    env_pred <- terra::predict(env_sampling_region, fun = function(x) {
+      dnorm(x, mean = gmm_model$parameters$mean, sd = gmm_model$parameters$variance$sigma)
+    })
+    
+    env_pred = predict(env_sampling_region, svm_model, na.rm = T)
+    plot(env_pred)
+    
     # Sampling of pseudo absence
     abs_spec = st_sf(geometry = st_sfc(), crs = proj_string)
     while(nrow(abs_spec) < nrow(occs_spec)){
@@ -99,10 +131,10 @@ pseudo_absences = furrr::future_map(
     }
     
     plot(sa_polygon)
-    plot(occs_kde, add = T, legend = F)
+    #plot(occs_kde, add = T, legend = F)
     plot(range_spec, add = T, col = "#00000000")
+    plot(occs_spec, add = T, col = "red", pch = 3, cex = 0.5)
     plot(abs_spec, add = T, col = "black", pch = 3, cex = 0.5)
-    plot(occs_spec, add = T, col = "white", pch = 3, cex = 0.5)
     
     
     sample_points = st_as_sf(st_sample(sample_region, nrow(abs_spec))) %>%