#Load Data

# 1. Load world map
world <- ne_countries(scale = "medium", returnclass = "sf")

# 2. Load temperature raster data (BIO01: Annual Mean Temperature)
raster_file <- "wc2/wc2.1_2.5m_bio_1.tif"
if (!file.exists(raster_file)) stop("Raster file does not exist in the specified path.")

temperature_raster <- rast(raster_file)
temperature_df <- as.data.frame(temperature_raster, xy = TRUE)
colnames(temperature_df) <- c("Longitude", "Latitude", "Temperature")

setwd("/Users/katiaavinapadilla/Desktop/PAPA")

# 3. Load Phytophthora infestans occurrence data
df <- read.csv("phytophtoradata.csv")
if (!all(c("decimalLongitude", "decimalLatitude") %in% colnames(df))) {
  stop("The file 'phytophtoradata.csv' must contain the columns 'decimalLongitude' and 'decimalLatitude'.")
}

# 4. Load wild potato species occurrence data
all_occurrences <- read.csv("wild_potato_species_gbif.csv")
if (!all(c("decimalLongitude", "decimalLatitude", "species") %in% colnames(all_occurrences))) {
  stop("The file 'wild_potato_species_gbif.csv' must contain the columns 'decimalLongitude', 'decimalLatitude', and 'species'.")
}

#Map 1: Temperature and Phytophthora infestans Distribution

map_temp_phytophthora <- ggplot() +
  # Temperature layer
  geom_tile(data = temperature_df, aes(x = Longitude, y = Latitude, fill = Temperature), alpha = 0.7) +
  scale_fill_viridis(option = "A", name = "Temperature (°C)") +
  
  # Base map
  geom_sf(data = world, fill = "transparent", color = "black", size = 0.3) +
  
  # Phytophthora infestans points
  geom_point(data = df, aes(x = decimalLongitude, y = decimalLatitude), 
             color = "red", size = 1, alpha = 0.8) +
  
  # Labels and theme
  labs(
    title = "Global Temperature and Phytophthora infestans Distribution",
    subtitle = "Annual Mean Temperature (BIO01) overlayed with pathogen occurrence",
    x = "Longitude",
    y = "Latitude",
    caption = "Data sourced from WorldClim and GBIF"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, face = "bold", hjust = 0.5),
    plot.subtitle = element_text(size = 12, hjust = 0.5),
    legend.position = "right",
    legend.key.size = unit(0.5, "cm")
  )

print(map_temp_phytophthora)

#Map 2: Wild Potato Species and Phytophthora infestans with Andean Region Highlighted

# Define Andean region coordinates for the rectangle
andean_region <- data.frame(
  xmin = -80, xmax = -60,
  ymin = -20, ymax = 0
)

map_potato_phytophthora <- ggplot(data = world) +
  # Base map
  geom_sf(fill = "gray95", color = "gray85") +
  
  # Wild potato species points
  geom_point(data = all_occurrences, 
             aes(x = decimalLongitude, y = decimalLatitude, color = species, shape = species), 
             size = 1, alpha = 0.8) +
  
  # Phytophthora infestans points
  geom_point(data = df, aes(x = decimalLongitude, y = decimalLatitude), 
             color = "red", size = 1, alpha = 0.8, inherit.aes = FALSE) +
  
  # Andean region highlighted
  annotate("rect", xmin = andean_region$xmin, xmax = andean_region$xmax,
           ymin = andean_region$ymin, ymax = andean_region$ymax,
           fill = NA, color = "blue", size = 1, linetype = "dashed") +
  
  # Scales for species
  scale_color_manual(values = rep(viridis::viridis(12), length.out = length(unique(all_occurrences$species)))) +
  scale_shape_manual(values = rep(c(16, 17, 18, 19), length.out = length(unique(all_occurrences$species)))) +
  
  # Labels and theme
  labs(
    title = "Distribution of Wild Potato Species and Phytophthora infestans",
    subtitle = "Andean region highlighted with pathogen occurrence overlayed",
    x = "Longitude",
    y = "Latitude",
    caption = "Data sourced from GBIF"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 18, face = "bold", hjust = 0.5),
    plot.subtitle = element_text(size = 12, hjust = 0.5),
    legend.position = "bottom",
    legend.key.size = unit(0.5, "cm")
  )
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
print(map_potato_phytophthora)

#Map 3: Wild Potato Species and climate data

library(ggplot2)
library(viridis)        # Improved color palette
library(ggnewscale)     # Multiple color scales
library(rnaturalearth)  # World map
library(rnaturalearthdata)
library(terra)          # Raster handling
library(sf)             # Spatial data
library(grid)           # For annotation of the rectangle
# Revisar los valores mínimos y máximos de temperatura
summary(temperature_df$Temperature)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -54.7592 -23.5127  -0.7347  -4.3388  18.8702  31.1667
# Definir paleta de colores para especies (mover aquí)
species_colors <- c(
  "Oxalis tuberosa" = "#1b9e77",    # Verde azulado
  "Solanum commersonii" = "#d95f02", # Naranja
  "Solanum stenotomum" = "#7570b3",  # Azul
  "Solanum brevicaule" = "#e7298a",  # Rosa
  "Solanum phureja" = "#66a61e",     # Verde
  "Ullucus tuberosus" = "#e6ab02"    # Amarillo
)

range_temp <- range(temperature_df$Temperature, na.rm = TRUE)
print(range_temp)  # Esto mostrará los valores mínimo y máximo
## [1] -54.75917  31.16667
map_potato_species <- ggplot() +

  geom_raster(data = temperature_df, aes(x = Longitude, y = Latitude, fill = Temperature)) +
  scale_fill_gradientn(
    name = "Temperature (°C)", 
    colours = rev(terrain.colors(10)),  
    limits = range_temp  
  ) +
  

  geom_point(data = all_occurrences, 
             aes(x = decimalLongitude, y = decimalLatitude, color = species), 
             size = 1, alpha = 0.7) +
  
 
  geom_sf(data = world, fill = "transparent", color = "gray80", size = 0.2) +
  

  scale_color_manual(
    name = "Species", 
    values = species_colors
  ) +
  

  labs(
    title = "Wild Potato Species and Climate Data",
    subtitle = "Species occurrences overlaid with climate data",
    x = "Longitude",
    y = "Latitude"
  ) +
  
  theme_minimal() +
  theme(
    legend.position = "top",
    legend.title = element_text(size = 10, face = "bold"),
    legend.text = element_text(size = 8),
    plot.title = element_text(size = 18, face = "bold", hjust = 0.5),
    plot.subtitle = element_text(size = 12, hjust = 0.5)
  )

# Mostrar el mapa
print(map_potato_species)
## Warning: Raster pixels are placed at uneven horizontal intervals and will be shifted
## ℹ Consider using `geom_tile()` instead.

Mapa 4

Save Maps

ggsave("map_temperature_phytophthora.png", map_temp_phytophthora, width = 10, height = 6, dpi = 300)
ggsave("map_potato_phytophthora.pdf", map_potato_phytophthora, width = 10, height = 6, dpi = 300)
ggsave("map_potato_species.pdf", map_potato_species, width = 10, height = 6, dpi = 300)
## Warning: Raster pixels are placed at uneven horizontal intervals and will be shifted
## ℹ Consider using `geom_tile()` instead.