Environment prep

library(tidyverse)      
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ ggplot2 3.4.0     ✔ purrr   1.0.1
## ✔ tibble  3.1.8     ✔ dplyr   1.1.0
## ✔ tidyr   1.3.0     ✔ stringr 1.5.0
## ✔ readr   2.1.3     ✔ forcats 1.0.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
library(spdep)        # for spatial autocorrelation tests
## Loading required package: spData
## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`
## Loading required package: sf
## Linking to GEOS 3.11.0, GDAL 3.5.3, PROJ 9.1.0; sf_use_s2() is TRUE
library(sp)
load('DB_RSdata.Rda')
data_M <- DB_RSdata %>% filter(sex=="M")

Variable notes

  • territorial definition: ifelse(recap>3 & call>1, 1, 0)
  • home.X and home.Y: average of all capture locations
  • call.X and call.Y: average of all calling locations
  • color: B, IB, I, IR, R
  • color2: B, I, R ([IB, I, IR] -> I)
  • color_num: 0,1,2,3,4 based on color
  • color2_num: 0, 1, 2 based on color2

Plots

All males

data_sp <- data_M %>% select(home.X, home.Y, color_num) %>%
  na.omit() %>% # remove those with NA coordinates
  mutate(home.X = jitter(home.X), home.Y = jitter(home.Y))  #add some noise to avoid overlapping points


ggplot(data_sp, aes(x = home.X, y = home.Y, color = color_num)) +
  geom_point() +
  scale_color_gradient(low = "cornflowerblue", high = "indianred") +
  labs(title = "Male locations", x = "X Coordinate", y = "Y Coordinate", color = "Color score") +
  theme_bw(base_size = 15) +
  coord_fixed()

Territorial males

data_sp <- data_M %>% filter(territorial == 1) %>%
  select(call.X, call.Y, color_num) %>%
  na.omit() %>% # remove those with NA coordinates
  mutate(call.X = jitter(call.X), call.Y = jitter(call.Y))  #add some noise to avoid overlapping points

ggplot(data_sp, 
       aes(x = call.X, y = call.Y, color = color_num)) +
  geom_point() +
  scale_color_gradient(low = "cornflowerblue", high = "indianred") +
  labs(title = "Territorial male locations", x = "X Coordinate", y = "Y Coordinate", color = "Color score") +
  theme_bw(base_size = 15) +
  coord_fixed()

Closest male color

# Compute pairwise distances using the Euclidean distance metric
distances <- data_M %>% select(home.X, home.Y, color2) %>%
  dist(method = "euclidean") %>% as.matrix
## Warning in dist(., method = "euclidean"): NAs introduced by coercion
# Set the diagonal of the distance matrix to a large value to exclude self-distance
diag(distances) <- Inf

# Find the index of the nearest neighbor for each data point
nearest_neighbor_indices <- apply(distances, 1, which.min)

# Extract the colors of the nearest neighbors
nearest_neighbor_colors <- data_M$color2[nearest_neighbor_indices]

# Create a data frame with the original data, the nearest neighbor coordinates, and colors
data_NN <- data_M %>%
  mutate(NN.X = data_M$home.X[nearest_neighbor_indices],
         NN.Y = data_M$home.Y[nearest_neighbor_indices],
         NN.color2 = nearest_neighbor_colors)

not really different in terms of proportion

data_NN %>% select(color2) %>% table()%>% prop.table 
## color2
##         B         I         R 
## 0.1184211 0.5394737 0.3421053
data_NN %>% select(color2, NN.color2) %>% table() %>% prop.table(1)
##       NN.color2
## color2          B          I          R
##      B 0.16666667 0.55555556 0.27777778
##      I 0.09756098 0.56097561 0.34146341
##      R 0.13461538 0.57692308 0.28846154

Closest territorial male color

# Compute pairwise distances using the Euclidean distance metric
distances <- data_M %>% filter(territorial==1) %>%
  select(call.X, call.Y, color2) %>%
  dist(method = "euclidean") %>% as.matrix
## Warning in dist(., method = "euclidean"): NAs introduced by coercion
# Set the diagonal of the distance matrix to a large value to exclude self-distance
diag(distances) <- Inf

# Find the index of the nearest neighbor for each data point
nearest_neighbor_indices <- apply(distances, 1, which.min)

# Extract the colors of the nearest neighbors
nearest_neighbor_colors <- data_M$color2[nearest_neighbor_indices]

# Create a data frame with the original data, the nearest neighbor coordinates, and colors
data_NN <- data_M %>% filter(territorial==1)%>%
  mutate(NN.X = data_M$call.X[nearest_neighbor_indices],
         NN.Y = data_M$call.Y[nearest_neighbor_indices],
         NN.color2 = nearest_neighbor_colors)

sample size a bit too small

data_NN %>% select(color2) %>% table()%>% prop.table 
## color2
##    B    I    R 
## 0.15 0.55 0.30
data_NN %>% select(color2, NN.color2) %>% table() %>% prop.table(1)
##       NN.color2
## color2         B         I         R
##      B 0.0000000 0.5000000 0.5000000
##      I 0.1818182 0.5909091 0.2272727
##      R 0.1666667 0.4166667 0.4166667

Spatial autocorrelation for all males

not different from random distribution

data_sp <- data_M %>% select(home.X, home.Y, color_num) %>%
  na.omit() %>% # remove those with NA coordinates
  mutate(home.X = jitter(home.X), home.Y = jitter(home.Y))  #add some noise to avoid overlapping points

# Create a spatial object with X and Y coordinates
spatial_data <- SpatialPoints(data_sp[, c("home.X", "home.Y")])

# Create a distance-based spatial weights matrix (use knn, 3 nearest neighbors)
W <- knn2nb(knearneigh(spatial_data, k = 3))

# Convert the neighbor list (W) to a listw object
W_listw <- nb2listw(W)

# Calculate Moran's I for the "color" variable
moran.test(data_sp$color_num, W_listw)
## 
##  Moran I test under randomisation
## 
## data:  data_sp$color_num  
## weights: W_listw    
## 
## Moran I statistic standard deviate = 0.19218, p-value = 0.4238
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.005071310      -0.006711409       0.003758949

Spatial autocorrelation for territorial males using call locations

color as numeric scores (B = 0, R = 4, and everything in between) not different from random distribution

data_sp <- data_M %>% filter(territorial == 1) %>%
  select(call.X, call.Y, color_num) %>%
  na.omit() %>% # remove those with NA coordinates
  mutate(call.X = jitter(call.X), call.Y = jitter(call.Y))  #add some noise to avoid overlapping points

# Create a spatial object with X and Y coordinates
spatial_data <- SpatialPoints(data_sp[, c("call.X", "call.Y")])

# Create a distance-based spatial weights matrix (use knn, 3 nearest neighbors)
W <- knn2nb(knearneigh(spatial_data, k = 3))

# Convert the neighbor list (W) to a listw object
W_listw <- nb2listw(W)

# Calculate Moran's I for the "color" variable
moran.test(data_sp$color_num, W_listw)
## 
##  Moran I test under randomisation
## 
## data:  data_sp$color_num  
## weights: W_listw    
## 
## Moran I statistic standard deviate = -0.15084, p-value = 0.5599
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.04275336       -0.02564103        0.01287056

Spatial autocorrelation for territorial males using home locations

not different from random distribution

data_sp <- data_M %>% filter(territorial == 1) %>%
  select(home.X, home.Y, color_num) %>%
  na.omit() %>% # remove those with NA coordinates
  mutate(home.X = jitter(home.X), home.Y = jitter(home.Y))  #add some noise to avoid overlapping points

# Create a spatial object with X and Y coordinates
spatial_data <- SpatialPoints(data_sp[, c("home.X", "home.Y")])

# Create a distance-based spatial weights matrix (use knn, 3 nearest neighbors)
W <- knn2nb(knearneigh(spatial_data, k = 3))

# Convert the neighbor list (W) to a listw object
W_listw <- nb2listw(W)

# Calculate Moran's I for the "color" variable
moran.test(data_sp$color_num, W_listw)
## 
##  Moran I test under randomisation
## 
## data:  data_sp$color_num  
## weights: W_listw    
## 
## Moran I statistic standard deviate = -0.64614, p-value = 0.7409
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       -0.10015482       -0.02564103        0.01329916