Exploring Mammal Diversity in Peru with perumammals

order_summary <- pm_list_orders() |> 
  dplyr::arrange(dplyr::desc(n_species)) |> 
  dplyr::mutate(percentage = round(n_species / sum(n_species) * 100, 1))

order_summary
#> # A tibble: 13 × 6
#>    order            n_families n_genera n_species n_endemic percentage
#>    <chr>                 <int>    <int>     <int>     <int>      <dbl>
#>  1 Rodentia                 11       67       194        56       33.9
#>  2 Chiroptera                8       64       189         9       33  
#>  3 Didelphimorphia           1       15        47        10        8.2
#>  4 Artiodactyla             10       31        46         0        8  
#>  5 Primates                  3       12        42         9        7.3
#>  6 Carnivora                 7       20        33         0        5.8
#>  7 Pilosa                    4        5         7         0        1.2
#>  8 Cingulata                 2        3         5         1        0.9
#>  9 Eulipotyphla              1        1         3         2        0.5
#> 10 Lagomorpha                1        1         2         0        0.3
#> 11 Paucituberculata          1        2         2         0        0.3
#> 12 Perissodactyla            1        1         2         0        0.3
#> 13 Sirenia                   1        1         1         0        0.2

ggplot(order_summary, 
       aes(x = reorder(order, n_species),
           y = n_species)) +
  geom_col(fill = "steelblue") +
  geom_text(aes(label = n_species),
            hjust = -0.2, size = 3) +
  coord_flip() +
  labs(
    title = "Mammalian Diversity in Peru by Order",
    subtitle = "Based on Pacheco et al. (2021)",
    x = "Order",
    y = "Number of Species"
  ) +
  theme_minimal()

# Obtener datos de ecoregiones
ecoregion_diversity <- pm_list_ecoregions(include_endemic = TRUE)
#> ── Peruvian Mammal Ecoregions (Brack-Egg, 1986) ────────────────────────────────
#> ℹ Number of ecoregions: 10
#> ℹ Total mammal species in Peru: 573
#> 
#> Ecoregions by species richness:
#> 
#> SB - Selva Baja: 320 species (55.8%), 18 endemic (5.6%)
#> YUN - Yungas: 256 species (44.7%), 48 endemic (18.8%)
#> SP - Sabana de Palmera: 83 species (14.5%), 0 endemic (0%)
#> BSE - Bosque Seco Ecuatorial: 81 species (14.1%), 4 endemic (4.9%)
#> VOC - Vertiente Occidental: 72 species (12.6%), 15 endemic (20.8%)
#> PUN - Puna: 71 species (12.4%), 14 endemic (19.7%)
#> BPP - Bosque Pluvial del Pacífico: 69 species (12%), 0 endemic (0%)
#> COS - Costa: 66 species (11.5%), 16 endemic (24.2%)
#> OCE - Oceánica: 30 species (5.2%), 0 endemic (0%)
#> PAR - Páramo: 26 species (4.5%), 4 endemic (15.4%)
#> 
#> Use pm_by_ecoregion() to filter species by ecoregion
#> ────────────────────────────────────────────────────────────────────────────────

# Preparar datos para el gráfico tipo waffle
waffle_data <- ecoregion_diversity |> 
  select(ecoregion_label, pct_endemic) |> 
  mutate(
    total_pct = paste0(round(pct_endemic, 1), "% endemic"),
    ecoregion_full = paste0(ecoregion_label, "\n", total_pct)
  ) |> 
  rowwise() |> 
  mutate(
    grid_data = list(
      expand.grid(x = 1:10, y = 1:10) |> 
        mutate(
          id = row_number(),
          type = ifelse(id <= round(pct_endemic), "Endemic", "Common")
        )
    )
  ) |> 
  unnest(grid_data)

# Crear el gráfico
ggplot(waffle_data, aes(x = x, y = y, fill = type)) +
  geom_tile(color = "white", linewidth = 0.4) +
  facet_wrap(~ecoregion_full, ncol = 2) +
  scale_fill_manual(
    name = NULL,
    values = c("Common" = "#CBD5E1", "Endemic" = "#E07A5F"),
    labels = c("Common species", "Endemic species")
  ) +
  coord_equal() +
  labs(
    title = "Mammal Endemism by Ecoregion in Peru",
    subtitle = "Each square represents 1% of species. The percentage indicates the proportion of endemic species"
  ) +
  theme_void() +
  theme(
    plot.title = element_text(
      face = "bold", 
      size = 16, 
      hjust = 0.5, 
      margin = margin(b = 5)
    ),
    plot.subtitle = element_text(
      hjust = 0.5, 
      size = 10, 
      color = "gray30", 
      margin = margin(b = 20)
    ),
    legend.position = "bottom",
    legend.text = element_text(size = 10),
    legend.key.size = unit(0.7, "cm"),
    legend.margin = margin(t = 15),
    strip.text = element_text(
      size = 9,
      face = "bold", 
      lineheight = 1.1,
      margin = margin(t = 5, b = 5)  # Aumentar márgenes
    ),
    strip.clip = "off",  # Permitir que el texto se extienda fuera del panel
    plot.background = element_rect(fill = "white", color = NA),
    panel.spacing.x = unit(2, "lines"),  # Aumentar espacio horizontal
    panel.spacing.y = unit(2.5, "lines"),  # Aumentar espacio vertical
    plot.margin = margin(15, 15, 15, 15)
  )

endemic_primates <- peru_mammals %>%
  filter(order == "Primates", endemic == TRUE) %>%
  select(genus, species, common_name)

endemic_primates
#> # A tibble: 9 × 3
#>   genus      species       common_name                                
#>   <chr>      <chr>         <chr>                                      
#> 1 Lagothrix  flavicauda    Mono choro de cola amarilla                
#> 2 Aotus      miconax       Mono nocturno peruano, mono nocturno andino
#> 3 Saguinus   illigeri      Pichico pardo de Illiger                   
#> 4 Saguinus   leucogenys    Pichico andino                             
#> 5 Saguinus   nigrifrons    Pichico de frente negra                    
#> 6 Callicebus oenanthe      Tocón del río Mayo                         
#> 7 Callicebus urubambensis  Tocón del Urubamba                         
#> 8 Pithecia   aequatorialis Huapo ecuatorial                           
#> 9 Pithecia   isabela       Huapo de Isabel

field_data <- data.frame(
  site = rep(c("Site_A", "Site_B", "Site_C"), each = 3),
  species = c(
    "Tremarctos ornatos", "Mazama rufina", "Odocoileus virginianus",
    "Leopardus pardalis", "Puma concolor", "Leopardus jacobita",
    "Lagothrix flavicauda", "Ateles belzebuth", "Alouatta seniculus"
  ),
  count = c(2, 1, 3, 1, 2, 1, 5, 3, 2)
)

cleaned_data <- field_data %>%
  left_join(
    validate_peru_mammals(field_data$species),
    by = c("species" = "Orig.Name")
  ) %>%
  select(site, original = species, validated = Matched.Name, count, Match.Level)

cleaned_data
#>     site               original            validated count
#> 1 Site_A     Tremarctos ornatos   Tremarctos ornatus     2
#> 2 Site_A          Mazama rufina        Mazama rufina     1
#> 3 Site_A Odocoileus virginianus           Odocoileus     3
#> 4 Site_B     Leopardus pardalis   Leopardus pardalis     1
#> 5 Site_B          Puma concolor        Puma concolor     2
#> 6 Site_B     Leopardus jacobita   Leopardus jacobita     1
#> 7 Site_C   Lagothrix flavicauda Lagothrix flavicauda     5
#> 8 Site_C       Ateles belzebuth     Ateles belzebuth     3
#> 9 Site_C     Alouatta seniculus   Alouatta seniculus     2
#>              Match.Level
#> 1             Exact rank
#> 2             Exact rank
#> 3 Matched at higher rank
#> 4             Exact rank
#> 5             Exact rank
#> 6             Exact rank
#> 7             Exact rank
#> 8             Exact rank
#> 9             Exact rank