Anova para Datos de ARBIMON

Cargar Paquetes

library(readr) # lee tabla
library(tidyverse) # maneja datos
library(lubridate) # fechas
library(hrbrthemes) # temas adicionales 
library(viridis) 
library(xtable)
library(gt) #great tables
library(rstatix) #great tables
library(report) # easy reporting result

Leer Datos

Desde un archivo tipo CSV

recordings <- read_csv("D:/BoxFiles/Box Sync/CodigoR/Acustica_Aves_Andrea/data/recordings-export_v2.csv") # lee datos
xtable(recordings)

Los datos tienen una estructura desordenada.

Poner los datos en formato Tidy

Una sola observacion por fila y columna con variables.

recordings_ext <- recordings  %>% # fix messy data
  pivot_longer(
    cols = Antropofonia__Avion:Geofonia__Viento,# starts_with(c("Anthropophony")),
    names_to = "Source",
    values_to = c("set"),
    values_drop_na = FALSE
  ) # %>%  add_column(Type = "Anthropophony")

# new tidy data object
recordings_tidy <- recordings_ext  %>% separate(., Source, into = c("Fuente", "Tipo"), sep = "__", remove = TRUE) %>% drop_na()
# recordings_tidy <- recordings_ext  %>% separate(., site, into = c("Sistema", "Finca"), sep = "_", remove = TRUE) %>% drop_na()
recordings_tidy$dia <- yday(recordings_tidy$time)  # extract day

# Grabaciones por sitio, fuente, tipo y dia
recordings_count <-  recordings_tidy  %>% 
        count (site, Fuente, Tipo, dia)

recordings_tidy2 <- recordings_count  %>% separate(., site, into = c("Sistema", "Finca"), sep = "_", remove = TRUE) %>% drop_na()

xtable(recordings_tidy2)

Grafica de fuente por sistema

# Boxplot basic
recordings_tidy2 %>%
  ggplot( aes(x=Fuente, y=n)) +
  geom_boxplot() +
  scale_fill_viridis(discrete = TRUE, alpha=0.6, option="A") +
  labs(y = "Número de eventos          ") +
  theme_ipsum() + # theme_bw() + 
  theme(
    #legend.position="none",
    plot.title = element_text(size=11), axis.text.x=element_text (angle=90, hjust=1)
  ) + facet_wrap(~Sistema) +
    # ggtitle("Basic boxplot") +
  xlab("") 

#Guardar esta figura en buena resolución .tiff a 300 dpi
#ggsave("fig/PA_PT.tiff", units="in", width=8, height=6, dpi=300, compression = 'lzw')


# Boxplot basic
recordings_tidy2 %>%
  ggplot( aes(x=Fuente, y=n, fill=Sistema)) +
  geom_boxplot() +
  scale_fill_viridis(discrete = TRUE, alpha=0.6, option="A") +
  labs(y = "Número de eventos          ") +
  theme_bw() + #theme_ipsum() + # theme_bw() + 
  theme(
    #legend.position="none",
    plot.title = element_text(size=11) #, axis.text.x=element_text (angle=90, hjust=1)
  ) + #facet_wrap(~Sistema) +
  # ggtitle("Basic boxplot") +
  xlab("") 

#Guardar esta figura en buena resolución .tiff a 300 dpi
# ggsave("fig/PA_PT.tiff", units="in", width=8, height=6, dpi=300, compression = 'lzw')

Grafica en barras de la media

# the names of the new data frame and the data frame to be summarised
aves_summary <- recordings_tidy2 %>% 
  group_by(Sistema, Fuente) %>%   # the grouping variable
  summarise(mean_PL = mean(n),  # calculates the mean of each group
            sd_PL = sd(n), # calculates the standard deviation of each group
            n_PL = n(),  # calculates the sample size per group
            SE_PL = sd(n)/sqrt(n())) # calculates the standard error of each group



avesPlot <- ggplot(aves_summary, aes(Fuente, mean_PL)) + 
                   geom_col( position = "dodge") +  
                   geom_errorbar(aes(ymin = mean_PL - SE_PL, ymax = mean_PL + SE_PL), width=0.2) + facet_wrap(~Sistema)

avesPlot + labs(y="Eventos ± s.e.     ", x="") + theme_ipsum()

Grafica detallada

# Boxplot basic
recordings_tidy2 %>%
  ggplot( aes(x=Tipo, y=n, color=Sistema  )) +
  geom_boxplot() +
  scale_fill_viridis(discrete = TRUE, alpha=0.6, option="A") +
  labs(y = "Número de eventos") +
  theme_bw() +
  # theme_ipsum() + #theme_bw() + #
  facet_wrap( ~ Fuente, dir="v" ) +
  theme(
     axis.text.x=element_text (angle=90, hjust=1)
    # plot.title = element_text(size=11)
  ) +  
  # ggtitle("Basic boxplot") +
  xlab("")

#Guardar esta figura en buena resolución .tiff a 300 dpi
#ggsave("fig/by_Tipo.tiff", units="in", width=9, height=8, dpi=300, compression = 'lzw')



# Boxplot basic2
recordings_tidy2 %>% #reorder(Fuente, Tipo, Sistema) +
  ggplot( aes(x=Tipo, y=n, color=Sistema  )) +
  geom_boxplot() +
  #scale_x_discrete(limits= recordings_tidy2$Tipo) +
  labs(y = "Número de eventos          ") +
  theme_ipsum() +
  theme(
     axis.text.x=element_text (angle=90, hjust=1, vjust=0)
    # plot.title = element_text(size=11)
  ) + facet_grid( ~ Fuente ) +
  # ggtitle("Basic boxplot") +
  xlab("")

#Guardar esta figura en buena resolución .tiff a 300 dpi
#ggsave("fig/by_Tipo.tiff", units="in", width=9, height=8, dpi=300, compression = 'lzw')

ANOVA de factores

Los datos son conteos. Asi que usamos Poisson.

###  ANOVA ###
# Son datos de conteos <-  distribucion poisson
fm =glm(n ~    Sistema * Fuente , data =recordings_tidy2, family = poisson )
anova1 <- aov( fm) 
anova1 %>% anova_summary %>%  gt()  # Hay diferencias entre sistemas y tipos

Effect	DFn	DFd	F	p	p<.05	ges
Sistema	1	507	16.959	4.46e-05	*	0.032
Fuente	2	507	9.922	5.93e-05	*	0.038
Sistema:Fuente	2	507	0.326	7.22e-01		0.001

# see interpretation in easy language.
report (anova1)

## The ANOVA suggests that:
## 
##   - The main effect of Sistema is significant (F(1, 507) = 16.96, p < .001) and can be considered as small (partial omega squared = 0.03).
##   - The main effect of Fuente is significant (F(2, 507) = 9.92, p < .001) and can be considered as small (partial omega squared = 0.03).
##   - The interaction between Sistema and Fuente is not significant (F(2, 507) = 0.33, p = 0.722) and can be considered as very small (partial omega squared = 0.00).

# cuáles son los grupos que generan estas diferencias?
intervals = TukeyHSD(anova1)
par(mar=c(6,18,3,1)) # adjust as needed
plot(intervals, las=2)

# dev.off()