Exploracion de datos

This

library(tidyverse)
## Warning: package 'tidyverse' was built under R version 4.1.3
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.6     v purrr   0.3.4
## v tibble  3.1.7     v dplyr   1.0.9
## v tidyr   1.2.0     v stringr 1.4.0
## v readr   2.1.2     v forcats 0.5.1
## Warning: package 'ggplot2' was built under R version 4.1.3
## Warning: package 'tibble' was built under R version 4.1.3
## Warning: package 'tidyr' was built under R version 4.1.3
## Warning: package 'readr' was built under R version 4.1.3
## Warning: package 'purrr' was built under R version 4.1.3
## Warning: package 'dplyr' was built under R version 4.1.3
## Warning: package 'forcats' was built under R version 4.1.3
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(readxl)
## Warning: package 'readxl' was built under R version 4.1.3
library(sf)
## Warning: package 'sf' was built under R version 4.1.3
## Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1; sf_use_s2() is TRUE
library(mapview)
## Warning: package 'mapview' was built under R version 4.1.3
library(dplyr)
library(lubridate)
## Warning: package 'lubridate' was built under R version 4.1.3
## 
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union
library(reshape)
## Warning: package 'reshape' was built under R version 4.1.3
## 
## Attaching package: 'reshape'
## The following object is masked from 'package:lubridate':
## 
##     stamp
## The following object is masked from 'package:dplyr':
## 
##     rename
## The following objects are masked from 'package:tidyr':
## 
##     expand, smiths

Read data

View map

aves_sf <-  st_as_sf(aves, coords = c("decimalLongitude", "decimalLatitude"), crs = 4326) 

mapview(aves_sf)
sites <- aves %>% dplyr::distinct(decimalLatitude, decimalLongitude) # %>% arrange(habitat) %>% # ordena por sistema luego pone site number
  # mutate(site=c(1:length(aves$county))) # pone site number

aves_sp1 <- aves %>% dplyr::filter (scientificName=="Adelomyia melanogenys") 


sites <- aves %>% distinct (decimalLatitude , decimalLongitude) %>% mutate(site=c(1:19))

Aves_site <- inner_join(aves, sites, by="decimalLatitude")
Aves_site$date <- yday(ymd(Aves_site$eventDate))

spp.melt <- reshape::melt(Aves_site, id.var=c("date" , "eventTime", "site",
                                              "habitat", "decimalLatitude" , "decimalLongitude.x", 
                                              "scientificName", "site"), measure.var="individualCount")

hab <- Aves_site %>% dplyr::select (site, habitat) %>% unique()

############ changed habitat by date
spp1 <- reshape::cast(spp.melt, site ~ date ~ scientificName) # 3D arr)ay 
## Aggregation requires fun.aggregate: length used as default
spnames <- as.vector(dimnames(spp1)$scientificName)

spnames
##  [1] "Adelomyia melanogenys"           "Aglaeactis cupripennis"         
##  [3] "Aglaiocercus kingii"             "Anairetes parulus"              
##  [5] "Anisognathus igniventris"        "Anisognathus somptuosus"        
##  [7] "Atlapetes latinuchus"            "Aulacorhynchus prasinus"        
##  [9] "Buteo platypterus"               "Caracara cheriway"              
## [11] "Cardellina canadensis"           "Catamenia inornata"             
## [13] "Catharus ustulatus"              "Chaetocercus mulsant"           
## [15] "Chamaepetes goudotii"            "Chlorostilbon melanorhynchus"   
## [17] "Cistothorus platensis"           "Cnemoscopus rubrirostris"       
## [19] "Coeligena coeligena"             "Coeligena lutetiae"             
## [21] "Coeligena torquata"              "Colaptes rivolii"               
## [23] "Colibri coruscans"               "Colibri cyanotus"               
## [25] "Coragyps atratus"                "Cyanocorax yncas"               
## [27] "Dendrocolaptes picumnus"         "Diglossa albilatera"            
## [29] "Diglossa caerulescens"           "Diglossa cyanea"                
## [31] "Diglossa humeralis"              "Diglossa lafresnayii"           
## [33] "Elaenia frantzii"                "Elaenia pallatangae"            
## [35] "Ensifera ensifera"               "Falco sparverius"               
## [37] "Grallaria ruficapilla"           "Heliangelus exortis"            
## [39] "Icterus chrysater"               "Lafresnaya lafresnayi"          
## [41] "Lepidocolaptes lacrymiger"       "Leptotila verreauxi"            
## [43] "Lesbia nuna"                     "Margarornis squamiger"          
## [45] "Mecocerculus poecilocercus"      "Melanerpes formicivorus"        
## [47] "Metallura tyrianthina"           "Myadestes ralloides"            
## [49] "Myioborus miniatus"              "Myiodynastes chrysocephalus"    
## [51] "Myiothlypis coronata"            "Myiothlypis nigrocristata"      
## [53] "Ochthoeca cinnamomeiventris"     "Ocreatus underwoodii"           
## [55] "Orochelidon murina"              "Patagioenas fasciata"           
## [57] "Phyllomyias nigrocapillus"       "Piaya cayana"                   
## [59] "Piranga rubra"                   "Piranga rubriceps"              
## [61] "Pygochelidon cyanoleuca"         "Pyrrhomyias cinnamomeus"        
## [63] "Rupornis magnirostris"           "Setophaga fusca"                
## [65] "Synallaxis azarae"               "Tangara labradorides"           
## [67] "Tangara vassorii"                "Tangara xanthocephala"          
## [69] "Troglodytes aedon"               "Troglodytes solstitialis"       
## [71] "Turdus fuscater"                 "Tyrannus melancholicus"         
## [73] "Vanellus chilensis"              "Vireo leucophrys"               
## [75] "Xiphocolaptes promeropirhynchus" "Xiphorhynchus triangularis"     
## [77] "Zenaida auriculata"              "Zonotrichia capensis"

Adelomyia_melanogenys model

##### select 1 species
Adelomyia_melanogenys <- as.data.frame(spp1[, , 1] )
#Setophaga_striata <- as.data.frame(spp1[, , 178] )
#Crotophaga_ani <-  as.data.frame(spp1[, , 55] )


# replace NA by 0
Adelomyia_melanogenys [is.na(Adelomyia_melanogenys )] <- 0
#Setophaga_striata[is.na(Setophaga_striata)] <- 0
#Crotophaga_ani[is.na(Crotophaga_ani)] <- 0


Adelomyia_melanogenys$site <- as.numeric(rownames(Adelomyia_melanogenys))
#Setophaga_striata$site <- as.numeric(rownames(Setophaga_striata))
#Crotophaga_ani$site <- as.numeric(rownames(Crotophaga_ani))


# add cov1
#Adelomyia_melanogenys <- Adelomyia_melanogenys  %>% left_join(hab, by = c("habitat"))

unmarked frame

Abundance Adelomyia_melanogenys

# Load unmarked, format data in unmarked data frame and summarize
library(unmarked)
## Warning: package 'unmarked' was built under R version 4.1.3
umf_ss <- unmarkedFramePCount(
  y=Adelomyia_melanogenys [,1:12],  # remove site                                          # Counts matrix
  # siteCovs= data.frame(hab = Setophaga_striata$Sistema,
   #                    dist_forest = Setophaga_striata$dist_forest) # , # Site covariates
  #obsCovs = list(day = as.data.frame(obs.cov[,2:5]))
  )       # Observation covs
summary(umf_ss)
## unmarkedFrame Object
## 
## 19 sites
## Maximum number of observations per site: 12 
## Mean number of observations per site: 12 
## Sites with at least one detection: 14 
## 
## Tabulation of y observations:
##   0   1   2   3   4 
## 207   9   7   4   1
plot(umf_ss)

fm.Nmix0 <- pcount(~1 ~1, data=umf_ss, control=list(trace=T, REPORT=1), K=364)
## initial  value 268.356718 
## iter   2 value 204.224676
## iter   3 value 172.387371
## iter   4 value 126.175143
## iter   5 value 124.616863
## iter   6 value 123.572625
## iter   7 value 123.416707
## iter   8 value 123.282292
## iter   9 value 123.266077
## iter  10 value 123.242449
## iter  11 value 123.183911
## iter  12 value 123.132823
## iter  13 value 123.122313
## iter  14 value 123.118559
## iter  15 value 123.116397
## iter  16 value 123.112487
## iter  17 value 123.104419
## iter  18 value 123.103888
## iter  19 value 123.100876
## iter  20 value 123.100007
## iter  21 value 123.099977
## iter  22 value 123.099966
## iter  22 value 123.099966
## iter  22 value 123.099966
## final  value 123.099966 
## converged
fm.Nmix0
## 
## Call:
## pcount(formula = ~1 ~ 1, data = umf_ss, K = 364, control = list(trace = T, 
##     REPORT = 1))
## 
## Abundance:
##  Estimate    SE    z  P(>|z|)
##      5.68 0.608 9.33 1.01e-20
## 
## Detection:
##  Estimate    SE     z  P(>|z|)
##     -7.44 0.629 -11.8 2.75e-32
## 
## AIC: 250.1999

model convergence

pb <- parboot(fm.Nmix0, nsim=500, report=10) # goodness of fit 
## t0 = 82.32895
## Running parametric bootstrap in parallel on 3 cores.
## Bootstrapped statistics not reported during parallel processing.
plot (pb) # plot goodness of fit