I have a dataframe in R where each row is an sf-formatted geospatial grid. For each grid, I have the number of venomous snakes and the number of non-venomous snakes seen in that area. I also have data on the diets of birds in each grid area. Specifically, I know how many non-venomous snakes birds ate, and how many venomous snakes birds ate. My goal is to evaluate whether birds are eating fewer venomous snakes than would be expected given their abundance in a given area. That is, are birds avoiding potentially dangerous prey.
See below for an example dataframe (though without the column describing each spatial polygon; I'm not sure how to recreate that for this example).
df <- data.frame(grid_id = c(1,2,3,4,5),
n_venomous = c(50,50,50,50,50),
n_non_venomous = c(50,50,50,50,50),
n_eat_venomous = c(10, 10, 10, 10, 10),
n_eat_non_venomous = c(40, 40, 40, 40, 40))
> df
grid_id n_venomous n_non_venomous n_eat_venomous n_eat_non_venomous
1 1 50 50 10 40
2 2 50 50 10 40
3 3 50 50 10 40
4 4 50 50 10 40
5 5 50 50 10 40
This seems like a relatively simply case of calculating differences between observed and expected values (i.e., if birds weren't choosy, we would expect them to eat 25 venomous and 25 non-venomous snakes in each area). But I imagine I shouldn't just run dozens of individual chi-square tests (I have 50+ different spatial grids). There is also the issue of spatial autocorrelation.
Given my research goal, how should I approach this analysis in R?