On commence par nettoyer un peu les données et rajouter la date précise.

# Read calls
calls <- read_csv("total - V2.csv", col_types = cols()) %>%
  select(-Temps) %>%
  mutate(Id = ifelse(Id == "Chirosp", "ChiroSp", Id)) %>%
  mutate(Id = ifelse(Id == "cris sociaux", "Cris sociaux", Id)) %>%
  mutate(date = make_datetime(Annee, Mois, Jour, Heure, Minute,
                              tz = "Europe/Paris"))
## Warning: 1 parsing failure.
##   row   col   expected           actual             file
## 19971 Temps time like  1,19930555555556 'total - V2.csv'

Observation chronologique

Commençons par observer les données en milieu de la période étudiée sur une durée de 3 jours (à partir du 1er septembre 2018) pour la pipistrelle commune seulement.

calls %>%
  filter(Id == "PippiT") %>%
  filter(date > make_date(2018, 9, 1)) %>%
  filter(date < make_date(2018, 9, 4)) %>%
  ggplot(aes(x = ceiling_date(date, "10 minutes"))) +
  geom_bar()

On constate une activité continue en début de nuit, puis des périodes d’accalmie. Regardons le détail de la nuit du 2 septembre.

calls %>%
  filter(Id == "PippiT") %>%
  filter(date > make_datetime(2018, 9, 1, 12)) %>%
  filter(date < make_datetime(2018, 9, 2, 12)) %>%
  ggplot(aes(x = date)) +
  geom_bar()

La nuit commence entre 20h et 21h, il y a un pic d’activités un peu avant 3h et il peut y avoir jusqu’à 8 cris par minute. S’il s’agit bien d’une activité de chasse, ne considérer que les minutes positives (contact ou non à chaque minute) gomme ce biais.

Début de l’activité

Considérons le début de la nuit pour les pipistrelles sur la période de l’étude.

start_times <- calls %>%
  filter(Id == "PippiT") %>%
  filter(Heure > 12) %>%
  mutate(day = make_date(year = Annee, month = Mois, day = Jour)) %>%
  mutate(hour = make_datetime(hour = Heure, min = Minute)) %>%
  group_by(day) %>%
  summarise(hour = min(hour))
start_times %>%
  ggplot(aes(x = day, y = hour)) +
  geom_point()

L’heure de première sortie avance avec le coucher le soleil d’une façon quasi-parfaite (il faudra le confirmer avec les données de coucher du soleil). Certaine nuit, cependant, les individus retardent leurs sorties, parfois jusqu’à près d’une heure. Il y a des jours sans données (26-27 août et autour du 15 octobre).

Identifions les jours sans données et vérifions que les données sur les nuits impactées sont complètes.

calls %>%
  filter(date < make_date(2018, 8, 27)) %>%
  select(date) %>%
  arrange(desc(date)) %>%
  head()
## # A tibble: 6 x 1
##   date               
##   <dttm>             
## 1 2018-08-26 05:46:00
## 2 2018-08-26 05:46:00
## 3 2018-08-26 04:39:00
## 4 2018-08-26 03:10:00
## 5 2018-08-26 00:02:00
## 6 2018-08-25 22:41:00
calls %>%
  filter(date > make_date(2018, 8, 27)) %>%
  select(date) %>%
  arrange(date) %>%
  head()
## # A tibble: 6 x 1
##   date               
##   <dttm>             
## 1 2018-08-28 20:25:00
## 2 2018-08-28 20:25:00
## 3 2018-08-28 20:25:00
## 4 2018-08-28 20:25:00
## 5 2018-08-28 20:25:00
## 6 2018-08-28 20:26:00
calls %>%
  filter(date < make_date(2018, 10, 15)) %>%
  select(date) %>%
  arrange(desc(date)) %>%
  head()
## # A tibble: 6 x 1
##   date               
##   <dttm>             
## 1 2018-10-14 00:31:00
## 2 2018-10-14 00:31:00
## 3 2018-10-14 00:27:00
## 4 2018-10-14 00:27:00
## 5 2018-10-14 00:27:00
## 6 2018-10-14 00:24:00
calls %>%
  filter(date > make_date(2018, 10, 15)) %>%
  select(date) %>%
  arrange(date) %>%
  head()
## # A tibble: 6 x 1
##   date               
##   <dttm>             
## 1 2018-10-18 19:04:00
## 2 2018-10-18 19:07:00
## 3 2018-10-18 19:23:00
## 4 2018-10-18 19:29:00
## 5 2018-10-18 19:29:00
## 6 2018-10-18 19:31:00
calls %>%
  select(date) %>%
  arrange(desc(date)) %>%
  head()
## # A tibble: 6 x 1
##   date               
##   <dttm>             
## 1 2018-10-30 07:35:00
## 2 2018-10-30 07:25:00
## 3 2018-10-30 07:24:00
## 4 2018-10-30 06:44:00
## 5 2018-10-30 06:40:00
## 6 2018-10-30 06:30:00

La nuit commençant le 25 août semble incomplète (uniquement 5 cris après 22h41) et la nuit commençant le 13 octobre s’arrête à 0h31. Toutes les autres nuits semblent être complètes.

calls_complete_night <- calls %>%
  filter(Date_nuit != "25/08/2018") %>%
  filter(Date_nuit != "13/10/2018")

Revenons aux nuits présentant un retard de l’activité. Sur les 60 nuits, on peut estimer que cela arrive entre 10 et 15 fois. On utilise la régression linéaire robuste pour détecter ces nuits.

start_mod <- MASS::rlm(as.numeric(hour) ~ day, start_times, maxit = 50)
start_times %>%
  add_predictions(start_mod) %>%
  ggplot(aes(x = day, y = hour)) +
  geom_point() +
  geom_line(aes(y = as_datetime(pred)))

start_times %>%
  add_residuals(start_mod) %>%
  arrange(desc(abs(resid))) %>%
  select(day, resid) %>%
  head(5)
## # A tibble: 5 x 2
##   day        resid
##   <date>     <dbl>
## 1 2018-10-21 2779.
## 2 2018-10-22 2766.
## 3 2018-09-21 1370.
## 4 2018-10-11  989.
## 5 2018-10-02  747.

On va analyser les données sur la nuit du 21 septembre car les autres nuits sont très rapprochées et l’anomalie peut se noyer dans les autres données (ça peut être plus facile à trouver dans ce premier cas).

calls %>%
  filter(date > make_datetime(2018, 9, 21, 12)) %>%
  filter(date < make_datetime(2018, 9, 22, 12)) %>%
  ggplot(aes(x = date)) +
  geom_bar()

Il s’agit d’une nuit très calme avec une faible activité en début et fin de nuit seulement (même en considérant toutes les espèces).

Il faudra vérifier les conditions météorologiques pour les dates suivantes: 21/9, 2/10, 11/10, 21/10 et 22/10.

Fin de l’activité

Observons comment se répartie la fin de l’activité pour la pipistrelle. On enlève les deux nuits incomplètes.

end_times <- calls_complete_night %>%
  filter(Id == "PippiT") %>%
  group_by(Date_nuit) %>%
  mutate(hour = max(date)) %>%
  mutate(hour = make_datetime(day = hour(hour) < 12, hour = hour(hour),
                              min = minute(hour))) %>%
  mutate(day = make_date(year = Annee, month = Mois, day = Jour)) %>%
  summarise(day = first(day), hour = first(hour))
end_mod <- MASS::rlm(as.numeric(hour) ~ day, end_times, maxit = 50)
end_times %>%
  add_predictions(end_mod) %>%
  ggplot(aes(x = day, y = hour)) +
  geom_point() +
  geom_line(aes(y = as_datetime(pred)))

L’heure du dernier contact semble moins se décaler pour coller avec l’heure de lever du soleil. Par ailleurs, il y a plus de variations globalement dans la dernière heure d’activité. Sur 8 nuits, l’activité de la pipistrelle s’arrête prématurément.

end_times %>%
  add_residuals(end_mod) %>%
  arrange(desc(abs(resid))) %>%
  select(day, resid) %>%
  head(8)
## # A tibble: 8 x 2
##   day          resid
##   <date>       <dbl>
## 1 2018-10-24 -37894.
## 2 2018-10-02 -34163.
## 3 2018-10-25 -33621.
## 4 2018-10-20 -31107.
## 5 2018-09-29 -30723.
## 6 2018-10-26 -27968.
## 7 2018-10-10 -17918.
## 8 2018-10-18 -16613.

Inspectons la nuit du 29/9.

calls %>%
  filter(Id == "PippiT") %>%
  filter(date > make_datetime(2018, 9, 29, 12)) %>%
  filter(date < make_datetime(2018, 9, 30, 12)) %>%
  ggplot(aes(x = date)) +
  geom_bar()

Il s’agit d’une nuit calme pour la pipistrelle.

Alors que les dernières nuits peuvent correspondre à des baisses d’effectifs liées aux départs pour les gîtes hivernaux, ce que l’on pourra confirmer par une estimation de l’activité, il faudra vérifier les conditions météorologiques pour les nuits suivantes : 29/9, 2/10, 10/10, 18/10 et 20/10.

Quantité d’activités

activity <- calls_complete_night %>%
  mutate(PippiT = Id == "PippiT") %>%
  mutate(Date_nuit = dmy(Date_nuit)) %>%
  group_by(PippiT, Date_nuit, Mois, Jour, Heure) %>%
  summarise(nb_min = n_distinct(Minute), nb_calls = n()) %>%
  group_by(PippiT, Date_nuit) %>%
  summarise(tot_min = sum(nb_min), tot_calls = sum(nb_calls))
activity %>%
  gather(summary, value, tot_min, tot_calls) %>%
  ggplot(aes(x = Date_nuit, y = value, color = PippiT)) +
  geom_point() +
  facet_wrap(~ summary, nrow = 3, scales = "free_y")

activity %>%
  summarise(corr = cor(tot_min, tot_calls))
## # A tibble: 2 x 2
##   PippiT  corr
##   <lgl>  <dbl>
## 1 FALSE  0.797
## 2 TRUE   0.952

Il existe un lien très fort entre les deux indicateurs d’activités pour la pipistrelle (le nombre de cris par nuit et le nombre de minutes positives). Le lien est assez fort pour les autres espèces aussi.

À partir du 15 octobre, l’activité (toutes espèces confondus) risque d’être trop faible pour pouvoir tirer des conclusions. Sans doute le départ en gîte hivernal.

À noter que le calcul des minutes positives réduit davantage le nombre de cris de la pipistrelle que des autres espèces, ce qui semble indiquer que les autres restent moins longtemps à proximité du microphone.

Lien entre nombre de cris et nombre de minutes positives

Réduire les données en considérant les minutes positives permet de gommer l’activité de chasse, productrice de nombreux cris localisés potentiellement dans la zone du micro.

On commence par analyser la relation entre le nombre de cris par heure et le nombre de minutes positives.

# Relation between calls and positive minutes
calls %>%
  group_by(Id, Annee, Mois, Jour, Heure) %>%
  summarise(min_pos = n_distinct(Minute), nb_calls = n()) %>%
  filter(Id == "Nyclei" | Id == "PippiT") %>%
  ggplot(aes(x = min_pos, y = nb_calls, shape = Id, color = Id)) +
  geom_point() +
  stat_smooth()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

On constate que PippiT passe sans doute plus de temps autour du microphone (chasse) et que Nyclei serait plus de passage (souvent un seul cri par minute).

Pour confirmation, étudions le nombre de cris par minute.

calls %>%
  group_by(Id, Annee, Mois, Jour, Heure, Minute) %>%
  summarise(nb_calls = n()) %>%
  filter(Id == "Nyclei" | Id == "PippiT") %>%
  ggplot(aes(x = nb_calls, color = Id)) +
  geom_histogram() +
  facet_wrap(~ Id, scales = "free_y")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

On voit que la présence de nombreux cris dans la même minute est très fréquente pour la pipistrelle. Mais cela pourrait être dû à des périodes “chargées” lors desquelles différentes pipistrelles passent autour du micro sur une longue durée.

Mettons en relation le nombre de cris par minute et les données caractérisant chaque heure (nombres de cris et de minutes positives).

calls %>%
  group_by(Id, Annee, Mois, Jour, Heure) %>%
  mutate(min_pos = n_distinct(Minute), nb_calls = n()) %>%
  group_by(Id, Annee, Mois, Jour, Heure, Minute, min_pos, nb_calls) %>%
  summarise(nb_calls_per_min = n()) %>%
  filter(Id == "Nyclei" | Id == "PippiT") %>%
  ggplot(aes(x = factor(nb_calls_per_min), y = min_pos, color = nb_calls)) +
  geom_jitter(size = 1, alpha = 0.1) +
  scale_color_viridis_c() +
  facet_wrap(~ Id)

Pour Nyclei, l’essentiel des contacts sont des passages (un cri par seconde avec peu d’activités dans l’heure). Pour PippiT, il y a deux situations: peu de contacts par heure et un (ou deux) cri par minute ; plus de trois cris par minutes sur des heures à plus de 200 cris par heure et au moins 50 minutes positives.

En supposant qu’il ne s’agit pas de transit d’essaims de tailles considérables, c’est qu’il doit s’agir des mêmes individus. On pourra utiliser ce critère pour identifier l’activité des pipistrelles.

calls_id <- calls_complete_night %>%
  mutate(rowid = seq_len(nrow(calls_complete_night)))
calls_complete_night <- calls_id %>%
  filter(Id == "PippiT") %>%
  group_by(Annee, Mois, Jour, Heure) %>%
  mutate(min_pos = n_distinct(Minute), nb_calls = n()) %>%
  group_by(Annee, Mois, Jour, Heure, Minute) %>%
  mutate(nb_calls_per_min = n()) %>%
  mutate(transit = nb_calls_per_min == 1 | min_pos < 45) %>%
  ungroup() %>%
  select(rowid, transit) %>%
  right_join(calls_id, by = "rowid") %>%
  select(-rowid)
write_rds(calls_complete_night, "calls.rds", compress = "xz")

Lien entre retard et quantité d’activités

On va vérifier l’hypothèse du lien entre le retard et l’activité sur la nuit. Pour étudier la structure de chaque nuit, nous allons éliminer les nuits incomplètes.

calls_complete_night %>%
  filter(Id == "PippiT") %>%
  mutate(Date_nuit = dmy(Date_nuit)) %>%
  group_by(Date_nuit, Mois, Jour, Heure) %>%
  summarise(nb_min = n_distinct(Minute), nb_calls = n()) %>%
  group_by(Date_nuit) %>%
  summarise(tot_min_inv = 1 / sum(nb_min),
            tot_calls_inv = 1 / sum(nb_calls)) %>%
  left_join(start_times %>% add_residuals(start_mod),
            by = c("Date_nuit" = "day")) %>%
  gather(summary, value, tot_min_inv, tot_calls_inv, resid) %>%
  ggplot(aes(x = Date_nuit, y = value)) +
  geom_point() +
  facet_wrap(~ summary, nrow = 3, scales = "free_y")

On a inversé le nombre total de cris et de minutes positives pour accentuer les nuits calmes. En plus de la nuit du 21/9, il y a effectivement une relation avec le retard du premier contact de chaque nuit car ce lien est encore plus marqué pour les nuits du 21 et 22 octobre qui essuient un retard important. Même chose, en moins marqué, pour d’autres nuits comme les 2, 19 et 20 octobre. Encore moins marqué : le 25 septembre, 3, 18 et 23 octobre. À noter une contre-observation : le 11 octobre, le retard ne s’accompagne pas d’une baisse d’activité.

Il faudra vérifier les conditions météorologiques sur ces nuits.

Conclusion

Il restera à vérifier le lien entre première et dernière activité et heure du coucher et lever du soleil.

summary_activity <- activity %>%
  left_join(start_times %>%
              mutate(start = hour, PippiT = TRUE) %>%
              select(-hour), by = c("Date_nuit" = "day", "PippiT")) %>%
  left_join(end_times %>%
              mutate(end = hour, PippiT = TRUE) %>%
              select(-hour, -Date_nuit), by = c("Date_nuit" = "day", "PippiT"))
write_rds(summary_activity, "summary.rds", compress = "xz")

Il faut aussi inspecter les conditions météorologiques sur les nuits :

À titre de référence, tout le mois de septembre est assez régulier, à l’exception du 21/9 et 29/9.

## R version 3.6.1 (2019-07-05)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.3 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas/libblas.so.3
## LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.2.20.so
## 
## locale:
##  [1] LC_CTYPE=fr_FR.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=fr_FR.UTF-8        LC_COLLATE=fr_FR.UTF-8    
##  [5] LC_MONETARY=fr_FR.UTF-8    LC_MESSAGES=fr_FR.UTF-8   
##  [7] LC_PAPER=fr_FR.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=fr_FR.UTF-8 LC_IDENTIFICATION=C       
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] tidyr_0.8.3     readr_1.3.1     modelr_0.1.4    lubridate_1.7.4
## [5] ggplot2_3.2.0   dplyr_0.8.3    
## 
## loaded via a namespace (and not attached):
##  [1] Rcpp_1.0.1        pillar_1.4.2      compiler_3.6.1   
##  [4] tools_3.6.1       zeallot_0.1.0     digest_0.6.20    
##  [7] viridisLite_0.3.0 evaluate_0.14     tibble_2.1.3     
## [10] gtable_0.3.0      nlme_3.1-141      lattice_0.20-38  
## [13] pkgconfig_2.0.2   rlang_0.4.0       cli_1.1.0        
## [16] yaml_2.2.0        xfun_0.8          withr_2.1.2      
## [19] stringr_1.4.0     knitr_1.23        vctrs_0.2.0      
## [22] hms_0.5.0         generics_0.0.2    grid_3.6.1       
## [25] tidyselect_0.2.5  glue_1.3.1        R6_2.4.0         
## [28] fansi_0.4.0       rmarkdown_1.14    purrr_0.3.2      
## [31] magrittr_1.5      MASS_7.3-51.4     scales_1.0.0     
## [34] backports_1.1.4   htmltools_0.3.6   assertthat_0.2.1 
## [37] colorspace_1.4-1  labeling_0.3      utf8_1.1.4       
## [40] stringi_1.4.3     lazyeval_0.2.2    munsell_0.5.0    
## [43] broom_0.5.2       crayon_1.3.4