Índice

1. Los Datasets.

Preparado workspace

# Limpiamos el workspace, por si hubiera algun dataset o informacion cargada
rm(list = ls())

# Limpiamos la consola
cat("\014")

Carga de librerias

# source("loadPackages.R")

packages <- c("ggplot2","ggpubr","readr","plotly","tidyverse","lubridate",
              "magrittr","funModeling","skimr", "rgdal", "zoo", "streamgraph")

new <- packages[!(packages %in% installed.packages()[,"Package"])]
if(length(new)) install.packages(new)
a=lapply(packages, require, character.only=TRUE)

2. Integración y selección de los datos de interés a analizar.

Procedemos a cargar el conjunto de datos

seasons_impute <- read.csv('../Work/seasons_impute.csv', sep = ',', header = TRUE)

waspbusters <- read.csv('../Work/WaspBusters_20210504_LARVAE.csv', sep = ',', header = TRUE)

wb_geo <- read.csv('../Work/WBds01_GEO.csv', sep = ',', header = FALSE,
                   skip =0)
wb_meteo <- read.csv('../Work/WBds02_METEO.csv', sep = ',', header = TRUE)

wb_queen_wasp <- read.csv('../Work/WBds03_all_the_queens_wasps.csv', sep = ',')

estaciones <- read.csv('../data/estaciones.csv', sep = ',')

fruta <- read.csv('../data/fruta.csv', sep = ',')

met <- read.csv('../data/met.csv', sep = ',')

nido <- read.csv('../data/nido.csv', sep = ',')%>%
  dplyr::select(1,4,6,7,8,9,0,11,12,13,14,16,19)

loc_estaciones <- read.csv('../Input_open_data/ds05_LOCALIZACION-ESTACIONES-METEOROLOGICAS.csv', sep = ';')

Local. estaciones metereologicas

loc_estaciones1 <- loc_estaciones[complete.cases(loc_estaciones),]
loc_estaciones2 <- data.frame(x=loc_estaciones1$YUTM,y=loc_estaciones1$XUTM) 
coordinates(loc_estaciones2) <- ~x+y 
class(loc_estaciones2)
## [1] "SpatialPoints"
## attr(,"package")
## [1] "sp"
proj4string(loc_estaciones2) <- CRS("+proj=utm +zone=31+datum=WGS84 +units=m +ellps=WGS84") 
loc_estaciones3 <- spTransform(loc_estaciones2,CRS("+proj=longlat +datum=WGS84"))
df_loc_est <- as.data.frame(loc_estaciones3)

data_loc_est <- cbind(loc_estaciones,df_loc_est)
data_loc_est$latitude <- data_loc_est$x+4.219
data_loc_est$longitude <- data_loc_est$y-6.62
# Seleccionamos los datos que nos interesas
data_loc_est <- data_loc_est %>%
  select(1,2,3,6,9,10)

Visualizaciones

Visualizamos las ubicaciones de las Estaciones Metereologicas´

library(leaflet)
# Create a color palette with handmade bins.
mybins <- seq(4, 6.5, by=0.5)
mypalette <- colorBin( palette="YlOrBr",
                       domain=data_loc_est$cota..m.,
                       na.color="transparent",
                       bins=mybins)
# Prepare the text for the tooltip:
mytext <- paste(   
  "Codigo: ", data_loc_est$CODIGO, "<br/>", 
    "Estacion: ", data_loc_est$ESTACION, "<br/>", 
    "Tipo: ", data_loc_est$TIPO, "<br/>", 
    "Cota: ", data_loc_est$cota..m., "<br/>", 
   "Longitud: ", data_loc_est$longitude, "<br/>", 
   "Latitud: ", data_loc_est$latitude, "<br/>") %>%
  lapply(htmltools::HTML)
# Final Map
leaflet(data_loc_est) %>%
    addTiles() %>%
    addMarkers(popup = paste0(
                 "<b>Estación: </b>"
                 , data_loc_est$ESTACION
                 , "<br>"
                 , "<b>Cota: </b>"
                 , data_loc_est$cota..m.
                 , "<br>"
                 , "<b>Tipo: </b>"
                 , data_loc_est$TIPO
                 , "<br>"
                 , data_loc_est$latitude
                 , data_loc_est$longitude
               ), clusterOptions = markerClusterOptions())%>%
      addProviderTiles(providers$Esri.WorldStreetMap) %>%
      addMiniMap(
        tiles = providers$Esri.WorldStreetMap,
        toggleDisplay = TRUE)

Nido de Avispas

colnames(nido)<-c("id","fecha","usuario", "municipio", "direccion", "especie","altura","diametro","longitude","latitude","estado","agente")
nido$fecha <- as.Date(nido$fecha)
nido$longitude <- gsub(",",".",nido$longitude)
nido$longitude <- as.integer(nido$longitude)
nido$latitude <- gsub(",",".",nido$latitude)
nido$latitude <- as.integer(nido$latitude)

Convertir UTM en

# library(rgdal)
nido <- nido[complete.cases(nido),]
nido1 <- data.frame(x=nido$latitude,y=nido$longitude) 
coordinates(nido1) <- ~x+y 
class(nido1)
## [1] "SpatialPoints"
## attr(,"package")
## [1] "sp"
proj4string(nido1) <- CRS("+proj=utm +zone=31+datum=WGS84 +units=m +ellps=WGS84") 
nido2 <- spTransform(nido1,CRS("+proj=longlat +datum=WGS84"))

Fusionar SpatianPoint a Nido

DF <- as.data.frame(nido2)

data_nido <- cbind(nido,DF)
data_nido$latitude <- data_nido$x+4.219
data_nido$longitude <- data_nido$y-6.62

Obtener Mes y año

data_nido[, "año"] <- as.numeric(format(data_nido[,"fecha"], "%Y"))
data_nido[, "mes"] <- as.numeric(format(data_nido[,"fecha"], "%m"))
mymonths <- c("Jan","Feb","Mar",
              "Apr","May","Jun",
              "Jul","Aug","Sep",
              "Oct","Nov","Dec")
#add abbreviated month name
data_nido <- transform(data_nido, mes = month.abb[mes])
head(data_nido)

Guardamos el fichero

write.csv(data_nido,'data/data_nido.csv', row.names = FALSE)

Visualizaciones

ggplotly(ggplot(data = data_nido, aes(x = mes, y = id, fill = especie)) + 
  geom_bar(stat='identity', position='fill') + 
  facet_wrap(('año')) + 
  xlab("mes") + ylab("% Cantidad de nidos declarados") +
  ggtitle("Especies observadas por mes y año") + 
  scale_fill_brewer(type='seq', palette='YlOrBr'))
ggplotly(ggplot(data = data_nido, aes(x=municipio, y=especie)) + 
  geom_bar(stat='identity') + 
  facet_wrap("año") + 
  coord_flip() + coord_flip()+
  theme(axis.text.x = element_text(angle = 45, hjust = 1, size=7,color="darkred")) + 
  theme(axis.text.y = element_text(hjust = 1, size=3, color="darkblue")) +
  ggtitle("Nidos por municipio y por año") + 
  labs(x = "MUNICIPIO ", y = "ESPACIE  "))
library(ggpmisc)
nidos_fechas <- data_nido%>%
  select(fecha,especie)%>%
  group_by(fecha,especie)%>%
  mutate(num_nidos = n())%>% 
  arrange(fecha)

nidos_fechas <- nidos_fechas[!duplicated(nidos_fechas), ]

ggplotly(ggplot(data = nidos_fechas, aes(x = fecha, y = num_nidos, color = especie)) +
    geom_point() +
  geom_smooth(method=lm) +
    labs(x = "Fecha", y= "Numero de nidos",colour = "Especie",title = "Tendencia nidos/colmenas por expecie") +
  theme_bw())

Vemos ahora la relación entre la altura de los nidos y su diámetro

library(ggvis)
data_nido%>%
  ggvis(~altura, ~diametro, fill = ~mes, opacity := 0.3 ) %>%
  layer_points() %>%
  layer_smooths() %>%
  layer_smooths(span = 0.3, stroke := "red") %>%
 layer_model_predictions(model = "lm", stroke := "green")
nidos_municipio <- data_nido%>%
  select(fecha, municipio)%>%
  group_by(month=floor_date(fecha, "month"),municipio) %>%
  mutate(num_nidos = n())%>% 
  arrange(month)

nidos_municipio <- nidos_municipio[!duplicated(nidos_municipio), ]

ggplotly(ggplot(data = nidos_municipio, aes(x = fecha, y = num_nidos, color = municipio)) +
    geom_line() +
       labs(x = "Fecha", y= "Numero de nidos",colour = "Municipios",title = "Numero de nidos por fecha y población") +
    theme_bw())
library(leaflet)
# Create a color palette with handmade bins.
mybins <- seq(4, 6.5, by=0.5)
mypalette <- colorBin( palette="YlOrBr",
                       domain=data_nido$diametro,
                       na.color="transparent",
                       bins=mybins)
# Prepare the text for the tooltip:
mytext <- paste(   
  "Fecha: ", data_nido$fecha, "<br/>", 
    "Especie: ", data_nido$especie, "<br/>", 
    "Altura: ", data_nido$altura, "<br/>", 
    "Diametro: ", data_nido$diametro, "<br/>", 
   "Estado: ", data_nido$estado, "<br/>", 
   "Agente: ", data_nido$agente, "<br/>", 
   "Municipio: ", data_nido$municipio, sep="") %>%
  lapply(htmltools::HTML)
# Final Map
m <- leaflet(data_nido) %>% 
  addTiles()  %>% 
  setView( lat=43,26271, lng = -2.9334110 , zoom=10) %>%
  addProviderTiles("Esri.WorldImagery") %>%
  addCircleMarkers(~longitude, ~latitude, 
    fillColor = ~mypalette(diametro), fillOpacity = 0.7, color="white", radius=8, stroke=FALSE,
    label = mytext,
    labelOptions = labelOptions( style = list("font-weight" = "normal",
                                              padding = "3px 8px"),
                                 textsize = "13px",
                                 direction = "auto")) %>%
  leaflet::addLegend( pal =mypalette, values =~diametro, 
             opacity=0.9, title = "Diametros", 
             position = "bottomright" )

m

Frutales

by_cod <- fruta %>% group_by(CODIGO.MUNICIPIO)%>% distinct()
head(by_cod)
by_cod_list <- fruta %>% group_by(CODIGO.MUNICIPIO)%>% distinct() %>%
 summarise(alltypes = paste(PRODUCTODES_C, collapse=", ") )
head(by_cod)

Visualizaciones

plot(by_cod)

ggplotly(ggplot(by_cod, aes(PRODUCTODES_C, fill=NOMBRE.LOCALIDAD)) +
  geom_bar(fill = "#0073C2FF") +
  labs(x = "Tipo de Arbol", y= "Numero de municipios",colour = "Municipios",title = "Numero de municipios por tipo de arbol") +
  theme_pubclean())

Tiempo

head(seasons_impute)
seasons_impute$month1 <- str_replace_all(
  seasons_impute$month, # column we want to search
  c("ENE" = "1","FEB" = "2","MAR" = "3","ABR" = "4","MAY" = "5","JUN" = "6",
    "JUL" = "7","AGO" = "8","SEP" = "9","OCT" = "10","NOV" = "11","DIC" = "12") # each string schould be matched with a replacement
)
seasons_impute$month1 <- as.numeric(seasons_impute$month1)
seasons_impute$year <- as.numeric(seasons_impute$year)

seasons_impute$fecha2 <- as.yearmon(paste(seasons_impute$year, seasons_impute$month1), "%Y-%m")
seasons_impute$fecha1<-as.Date(with(seasons_impute,paste(year,month1,sep="-")),"%Y-%m")
weather <- transform(seasons_impute, fecha = as.Date(as.yearmon(paste(year, month1, sep = "-"))))
weather <- weather %>%
  dplyr::select(47,4,5,8:43)
weather[, "año"] <- as.numeric(format(weather[,"fecha"], "%Y"))
weather[, "mes"] <- as.numeric(format(weather[,"fecha"], "%m"))

Fusionamos los datos de las estaciones metereologicas con la localización

weather_df <- merge(weather, data_loc_est, by.x = "codigo",  by.y = "CODIGO")

weather_df <-  weather_df %>%
  select(2,1,40:46,4,6,8,10,12,14,16,18,20,22,24,26,28,30,32,34,36,38)

head(weather_df)

Visualizaciones

data3 <- na.aggregate(weather_df)
data1 <- na.omit(data3)                           # Apply na.omit function

library(dygraphs)
library(xts)
data1[7:27] <- lapply(data1[7:27], as.numeric)
data1$fecha <- as.Date(data1$fecha)
# Then you can create the xts necessary to use dygraph
don <- xts(x = data1[8:27], order.by = data1$fecha)

# Finally the plot
p <- dygraph(don, 
        main = "Variable Metereologicas por fecha", 
        ylab = "Valor", xlab = "Fecha") %>%
  dyOptions(labelsUTC = TRUE, fillGraph=TRUE, fillAlpha=0.1, drawGrid = FALSE, colors="#D8AE5A") %>%
  dyRangeSelector() %>%
  dyCrosshair(direction = "vertical") %>%
  dyHighlight(highlightCircleSize = 5, highlightSeriesBackgroundAlpha = 0.2, hideOnMouseOut = FALSE)  %>%
  dyRoller(rollPeriod = 1)

p

Guardamos el fichero

write.csv(weather_df,'data/weather.csv')