Preparación de la base de datos y análisis exploratorio
Grupo 4
2023-03-29
library(rmarkdown) # libreria para mostrar tablas
library(tidyr)
library(dplyr)
library(stringr)
library(writexl)
library(sf)
library(readxl)
library(leaflet)
library(lubridate)Preparación
Carga inicial de los datos
La primera tarea que se realizó después de descargar y agrupar los distintos archivos csv fue la de cargar los datos en un objeto r y modificar el nombre de las columnas. Estos nombre no son los definitivos, puesto que debido a la utilización de varias fuentes de datos, ha sido necesario renombrar algunas columnas.
df = read.csv("bbdd.csv", sep=';', header = FALSE)
colnames(df) <- c('gid', 'name', 'number_', 'address', 'open', 'available',
'free', 'total', 'ticket', 'updated_at',
'{globalid,created_user,created_date,last_edited_user,
last_edited_date}', 'geo_shape', 'geo_point_2d')En una primera instancia, la base de datos tenía la siguiente forma:
Por tanto, para poder trabajar correctamente, eliminamos estas columnas y reorganizamos los nombres de las columnas de la base de datos:
df <- df[,-12:-15]
colnames(df) <- c('gid', 'name', 'number_', 'address', 'open', 'available',
'free', 'total', 'ticket', 'updated_at',
'globalid__created_user__created_date__last_edited_user__last_edited_date',
'geo_shape', 'geo_point_2d')Después de esto, la base de datos resultante es la que sigue:
Como parte de este exploración inicial, pudimos observar como la columna referente al tiempo (“updated_at”) estaba compuesta por la fecha y la hora. Para poder agrupar los datos por horas y unir los diferentes datos de los que disponemos, se decidió separar dicha columna en dos columnas distintas: “fecha” y “hora”.
df[c('fecha', 'hora')] <- str_split_fixed(df$updated_at, ' ', 2)Finalmente, se eliminaron todas las variables de las que no ibamos a hacer uso (updated_at, globalid__created_user__created_date__last_edited_user__last_edited_date, geo_shape) y “geo_point_2d” que, pese a que será muy útil en el estudio y, de hecho, se utilza más adelante en este mismo documento, será introducida de nuevo más adelante.
df = subset(df, select = -c(updated_at, globalid__created_user__created_date__last_edited_user__last_edited_date,geo_shape, geo_point_2d))Agrupación de los datos
Una de las desicisones que se han tomado en este proyecto ha sido la de agrupar los datos, que en un origen estaban cada 15 minutos, por horas:
df$hora_hora = format(strptime(df$hora,"%H:%M:%S"),'%H')
df_agrupado = df %>% group_by(number_, fecha, hora_hora) %>%
reframe(gid = gid, name = name, address = address, open = open,
avg_av=mean(available), avg_free = mean(free),
avg_total = mean(total), ticket = ticket)
df_agrupado = df_agrupado %>% distinct()
df_agrupado$hora_hora = paste(df_agrupado$hora_hora, ":00:00", sep = '') Datos AEMET
Para obtener datos metereológicos se hizo uso de la API de AEMET. El programa Python que realiza las llamadas a la API y crea el csv que posteriormente será cargado es el siguiente:
Una vez generado el csv, se cargo a R:
df_aemet = read.csv("aemet_data.csv")Como se puede ver, en la base de datos de valenbisi (df) las fechas tienen el siguiente formato: “21/01/2023”, mientras que en la base de datos de AEMET (df_aemet) las fechas tienen la siguiente forma: “2023-01-21”. Para poder, crear la base de datos definitiva a partir de concatenaciones, se tuvo que unificar el formato de fecha, eligiendo como formato canónico el de la base de datos de valenbisi. En consecuencia, las fechas del dataframe “df_aement” tuvieron que ser modificadas de la forma que sigue:
df_aemet$fecha = format(as.Date(df_aemet$fecha, '%Y-%m-%d'), "%d/%m/%Y")Además, se eliminaron las variables “indicativo”, “nombre”, “provincia” y “altitud” cuyos datos eran comunes a todas las instancias:
df_aemet = subset(df_aemet, select = -c(indicativo, nombre,provincia, altitud))Una vez depurada al máximo la base de datos metereológicos, se procedió a crear “df_merge”, el dataframe resultado de concatenar “df” y “df_aemet” utilizando como clave común la variable “fecha”:
df_merge = merge(df_agrupado, df_aemet, by="fecha")Datos de festivos
En cuento a los datos de festivos, al tratarse de un cantidad considerablemente reducida de días, fueron generados utilizando simple funciones en Excel. El resultado fue el siguiente:
df_festivos = read.csv("datos_festivos.csv", sep=';')Tal y como se ha hecho con los datos de AEMET, a continuación se concatenan “df_festivos” con, en este caso y a partir de ahora, “df_merge” utilizando como clave la variable “fecha”:
df_merge = merge(df_merge, df_festivos, by="fecha")Tratamiento de los barrios
Por lo que respecta a los datos de los barrios que serán utilizados para realizar análisis como “Comprobar si el comportamiento de los usuarios que utilizan bicicletas varía entre barrios periféricos y barrios del centro”.
Dichos datos fueron obtenidos en el portal de datos abiertos del Ayuntamiento de València en formato geojson puesto que, pese a que estaban disponibles en formato csv, con este formato ha sido más sencillo determinar en qué barrio está cada estación.
Dicha base de datos tiene la siguiente forma:
df_barrios_json = st_read("barris-barrios.geojson")## Reading layer `barris-barrios' from data source
## `C:\Users\Jose\UPV\Javier Luque Saiz - PROYII\rMarkdown\barris-barrios.geojson'
## using driver `GeoJSON'
## Simple feature collection with 88 features and 8 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -0.432535 ymin: 39.27893 xmax: -0.2753685 ymax: 39.56659
## Geodetic CRS: WGS 84Teniendo en cuenta la gran cantidad de filas que tiene “df_merge”, para determinar en qué barrio está cada estación no se hizo uso de dicha base de datos, sino que se utilizó uno de los muchos archivos csv de los que disponiamos donde cada instancia es una estación distinta. Para simplificar las operaciones, el archivo csv original se convirtió a un archivo xlsx y se quitarón todas las variables no necesarias para realizar la concatenación:
df_aux = read_excel("valenbisi.xlsx")Para poder realizar la concatenación, en primer lugar se creó un objeto sf, concretamente un punto (por tanto, a partir de las variables longitud y latitud), para cada una de las estaciones. Para ello, se uso la función st_as_sf de la forma que sigue:
coords_puntos <- st_as_sf(df_aux, coords = c("longitud", "latitud"), crs = 4326)A partir de este nuevo dataframe, que contiene los puntos, y de “df_barrios_json” que contiene los poligonos que delimitan los barrios, utilizando la función “st_join” se puede determinar a qué poligono (barrio) pertenece cada punto (estación):
resultado <- st_join(coords_puntos, df_barrios_json)
df_barrios_json2 = select(df_barrios_json, coddistbar, geometry)
resultado <- merge(resultado %>% as.data.frame(), df_barrios_json2 %>% as.data.frame(), by = "coddistbar")
resultado <- select(resultado, number_, coddistbar, nombre, codbarrio, coddistrit, geometry.x, geometry.y)
names(resultado)[names(resultado) == 'nombre'] <- 'nombre_barrio'
names(resultado)[names(resultado) == 'geometry.x'] <- 'geometry_punto'
names(resultado)[names(resultado) == 'geometry.y'] <- 'geometry_barrio'De la misma forma que se ha hecho anteriormente, a continuación se concatenan “resultado” con “df_merge” utilizando como clave la variable “number_”:
df_merge = merge(df_merge, resultado, by='number_')Para comprobar los resultados de forma gráfica es posible crar un mapa en en que queden indicadas las estaciones dibujadas sobre un mapa con los barrios:
df_aux <- df_aux %>%
mutate(longitud = as.numeric(longitud),
latitud = as.numeric(latitud))
m <- leaflet() %>%
addTiles() %>%
addPolygons(data = df_barrios_json, popup = ~nombre, fillColor = "blue", fillOpacity = 0.5, color = "black", weight = 1) %>%
addCircleMarkers(data = df_aux, lng = ~longitud, lat = ~latitud, label = ~name,
radius = 3, fillColor = "red", fillOpacity = 1, stroke = FALSE)Así, por ejemplo, buscando en el mapa la estación “001_GUILLEN_DE_CASTRO”, que en la base de datos tiene como barrio El Carme, se confirma que, en efecto, las coordenadas de dicha estación se encuentran en el polígono de dicho barrio
Tratamiento de los distritos
Además de por barrios, puede ser interesante realizar análisis por distritos. Para agregar estos datos, en primer lugar se carga el archivo geojson obtenido también del portal de datos abiertos del Ayuntamiento de València:
df_distritos_json = st_read("districtes-distritos.geojson")## Reading layer `districtes-distritos' from data source
## `C:\Users\Jose\UPV\Javier Luque Saiz - PROYII\rMarkdown\districtes-distritos.geojson'
## using driver `GeoJSON'
## Simple feature collection with 22 features and 11 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -0.432535 ymin: 39.27893 xmax: -0.2753685 ymax: 39.56659
## Geodetic CRS: WGS 84A continuación se utilizará el dataframe “df_aux” que ya ha sido cargado y “coords_puntos” que ya había sido generado. Con esto, se procede de forma muy similar que en el caso anterior:
resultado <- st_join(coords_puntos, df_distritos_json)
df_distritos_json2 = select(df_distritos_json, coddistrit, geometry)
resultado <- merge(resultado %>% as.data.frame(), df_distritos_json2 %>% as.data.frame(), by = "coddistrit")
resultado <- select(resultado, number_, nombre, geometry.y)
names(resultado)[names(resultado) == 'nombre'] <- 'nombre_distrito'
names(resultado)[names(resultado) == 'geometry.y'] <- 'geometry_distrito'De la misma forma que se ha hecho anteriormente, a continuación se concatenan “resultado” con “df_merge” utilizando como clave la variable “number_”:
df_merge = merge(df_merge, resultado, by='number_')De la misma forma que en el caso de los barrios, para comprobar los resultados de forma gráfica es posible crar el mismo mapa pero en enste caso indicando los distritos:
df_aux <- df_aux %>%
mutate(longitud = as.numeric(longitud),
latitud = as.numeric(latitud))
m <- leaflet() %>%
addTiles() %>%
addPolygons(data = df_distritos_json, popup = ~nombre, fillColor = "blue", fillOpacity = 0.5, color = "black", weight = 1) %>%
addCircleMarkers(data = df_aux, lng = ~longitud, lat = ~latitud, label = ~name,
radius = 3, fillColor = "red", fillOpacity = 1, stroke = FALSE)Lugares de interés
Por último, falta cargar el archivo en el cual se indica que lugar de interés tiene próximo cada estación:
df_interés = read.csv("lugares_de_interes.csv", sep = ";")Como en todas las ocasiones anteriores, se contatenan “df_interés” y “df_merge” utilizando como clave “number_”:
df_merge = merge(df_merge, df_interés, by='number_')Análisis
Introducción
La base de datos resultante de cruzar los datos de todas las fuentes tiene 545376 observaciones y 35 variables.
Podemos observar que existen dos tipos de variables, las que definen la naturaleza de una estación y las que definen el estado de una estación. Las variables que definen el estado de la estación van cambiando a lo largo del tiempo, como los bornes o bicis disponibles; mientras que las variables que definen la naturaleza de la estación se mantienen constantes, como las coordenadas, el código de barrio, los lugares de interés cercanos, etc.
Resumen de los datos
Empezaremos haciendo un resumen general de todos los datos.
summary(df_merge)## number_ fecha hora_hora gid
## Min. : 1.00 Length:545376 Length:545376 Min. :901581
## 1st Qu.: 69.75 Class :character Class :character 1st Qu.:901650
## Median :138.50 Mode :character Mode :character Median :901719
## Mean :138.50 Mean :901719
## 3rd Qu.:207.25 3rd Qu.:901787
## Max. :276.00 Max. :901856
## name address open avg_av
## Length:545376 Length:545376 Mode :logical Min. : 0.0
## Class :character Class :character FALSE:1976 1st Qu.: 2.5
## Mode :character Mode :character TRUE :543400 Median : 6.5
## Mean : 7.7
## 3rd Qu.:12.0
## Max. :40.0
## avg_free avg_total ticket tmed
## Min. : 0.00 Min. : 9.00 Mode :logical Length:545376
## 1st Qu.: 6.25 1st Qu.:15.00 FALSE:270712 Class :character
## Median :11.75 Median :20.00 TRUE :274664 Mode :character
## Mean :11.83 Mean :19.82
## 3rd Qu.:16.00 3rd Qu.:20.00
## Max. :40.00 Max. :40.00
## prec tmin horatmin tmax
## Length:545376 Length:545376 Length:545376 Length:545376
## Class :character Class :character Class :character Class :character
## Mode :character Mode :character Mode :character Mode :character
##
##
##
## horatmax dir velmedia racha
## Length:545376 Min. : 0.0 Length:545376 Length:545376
## Class :character 1st Qu.:10.0 Class :character Class :character
## Mode :character Median :23.0 Mode :character Mode :character
## Mean :22.9
## 3rd Qu.:26.0
## Max. :99.0
## horaracha festivo vacaciones fin_de_semana
## Length:545376 Min. :0.00000 Min. :0.0000 Min. :0.0000
## Class :character 1st Qu.:0.00000 1st Qu.:0.0000 1st Qu.:0.0000
## Mode :character Median :0.00000 Median :0.0000 Median :0.0000
## Mean :0.03644 Mean :0.2308 Mean :0.2915
## 3rd Qu.:0.00000 3rd Qu.:0.0000 3rd Qu.:1.0000
## Max. :1.00000 Max. :1.0000 Max. :1.0000
## día_semana coddistbar nombre_barrio codbarrio
## Length:545376 Length:545376 Length:545376 Length:545376
## Class :character Class :character Class :character Class :character
## Mode :character Mode :character Mode :character Mode :character
##
##
##
## coddistrit geometry_punto geometry_barrio
## Length:545376 POINT :545376 POLYGON :545376
## Class :character epsg:4326 : 0 epsg:4326 : 0
## Mode :character +proj=long...: 0 +proj=long...: 0
##
##
##
## nombre_distrito geometry_distrito interés metro
## Length:545376 POLYGON :545376 Min. :0.000 Min. :0.0000
## Class :character epsg:4326 : 0 1st Qu.:0.000 1st Qu.:0.0000
## Mode :character +proj=long...: 0 Median :0.000 Median :0.0000
## Mean :1.728 Mean :0.1159
## 3rd Qu.:3.000 3rd Qu.:0.0000
## Max. :8.000 Max. :1.0000Lo primero que se puede observar es que muchas variables no tienen el formato que deberían, por ejemplo, muchas variables categóricas están como numéricas; por lo que pasaremos estas variables a estos formatos:
A categóricas: number, hora_hora (está probablemente podrá utilizarse como numérica en algún momento, pero de momento la pasaremos a categórica), gid, día_semana, codbarrio, coddistrit
A fecha: fecha (mm/dd/yyyy).
A numérica: tmed, prec, tmin, tmax, velmedia, racha.
A lógica: festivo, vacaciones, fin_de_semana, metro.
A hora: horatmin, horatmax
Datos faltantes
na_count = sapply(df_merge, function(y) sum(length(which(is.na(y)))))
na_count = data.frame(na_count)¿Existen aumentos o disminuciones anómalas en los datos?
castro_23_02 = subset(df_merge, df_merge$number_ == 1 & df_merge$fecha == '23/02/2023')
library(ggplot2)
xValue <- sub(':00:00', '' ,castro_23_02$hora_hora)
yValue <- as.numeric(castro_23_02$avg_av)
data <- data.frame(xValue,yValue)
p = ggplot(data, aes(x=xValue, y=yValue, group = 1)) +
geom_line() +
geom_point()
p + (labs(x = "Horas", y = 'Bicis disponibles'))
Como resultado de esta anomalía, se ha creado un script que genera un vector binario que se anexará a la base de datos:
prueba=df_merge[c(1,2,3,8)]
prueba$fecha=as.Date(prueba$fecha, format="%d/%m/%Y")
##ordenacion
prueba=prueba[order(prueba$number_,prueba$fecha,prueba$hora_hora),]
#bucle
longitud=length(prueba[,1])
marcas=c()
for (i in 1:longitud) {
if(as.Date(prueba$fecha[i]) != as.Date('02/12/2022', format = "%d/%m/%Y") &
as.Date(prueba$fecha[i]) != as.Date('23/02/2023', format = "%d/%m/%Y")) {
previo <- i - 1
posterior <- i + 1
if(prueba$number[i] == prueba$number[previo] &
prueba$number[i] == prueba$number[posterior]) {
mediaanterior <- abs(prueba$avg_av[i] - prueba$avg_av[previo])
medioposterior <- abs(prueba$avg_av[i] - prueba$avg_av[posterior])
if(mediaanterior > 8 | medioposterior > 8) {
marcas <- c(marcas, 1)
} else {
marcas <- c(marcas, 0)
next
}
} else if(prueba$number[i] == prueba$number[previo] &
prueba$number[i] != prueba$number[posterior]) {
if(mediaanterior > 8) {
marcas <- c(marcas, 1)
} else {
marcas <- c(marcas, 0)
next
}
} else if(prueba$number[i] != prueba$number[previo] &
prueba$number[i] == prueba$number_[posterior]) {
if(medioposterior > 8) {
marcas <- c(marcas, 1)
} else {
marcas <- c(marcas, 0)
next
}
} else {
marcas <- c(marcas, 0)
}
} else {marcas <- c(marcas, 0)}
}load('df_marcas.Rda')a <- table(df_marcas)
df <- data.frame(check = c("1", "0"),
sum=c(a[names(a)==1], a[names(a)==0]))
g<-ggplot(data=df, aes(x=check, y=sum)) +
geom_bar(stat="identity")
g
Más análisis
estaciones_aux = df_merge %>% select(c('number_', 'name', 'address', 'avg_total','coddistbar', 'nombre_barrio', 'codbarrio', 'coddistrit', 'geometry_punto', 'geometry_barrio', 'nombre_distrito', 'geometry_distrito', 'interés', 'metro'))
estaciones = estaciones_aux %>% distinct(number_, .keep_all = TRUE)Analisis de estaciones por barrio
table(estaciones$nombre_barrio)##
## AIORA ALBORS
## 7 4
## ARRANCAPINS BENICALAP
## 9 11
## BENIFERRI BENIMACLET
## 2 5
## BENIMAMET CABANYAL-CANYAMELAR
## 2 12
## CAMI DE VERA CAMI FONDO
## 2 1
## CAMI REAL CAMPANAR
## 2 8
## CIUTAT DE LES ARTS I DE LES CIENCIES CIUTAT FALLERA
## 7 2
## CIUTAT JARDI CIUTAT UNIVERSITARIA
## 3 5
## EL BOTANIC EL CALVARI
## 2 1
## EL CARME EL GRAU
## 5 3
## EL MERCAT EL PLA DEL REMEI
## 2 6
## ELS ORRIOLS EN CORTS
## 3 2
## EXPOSICIO FAVARA
## 4 2
## JAUME ROIG L'AMISTAT
## 1 2
## L'HORT DE SENABRE L'ILLA PERDUDA
## 2 2
## LA CARRASCA LA CREU COBERTA
## 14 1
## LA CREU DEL GRAU LA FONTETA S.LLUIS
## 2 2
## LA FONTSANTA LA GRAN VIA
## 3 6
## LA LLUM LA MALVA-ROSA
## 1 6
## LA PETXINA LA PUNTA
## 3 1
## LA RAIOSA LA ROQUETA
## 4 2
## LA SEU LA VEGA BAIXA
## 2 2
## LA XEREA LES TENDETES
## 3 1
## MALILLA MARXALENES
## 7 4
## MESTALLA MONTOLIVET
## 6 5
## MORVEDRE NA ROVELLA
## 2 3
## NATZARET NOU MOLES
## 2 8
## PATRAIX PENYA-ROJA
## 3 6
## RUSSAFA SAFRANAR
## 7 3
## SANT ANTONI SANT FRANCESC
## 1 6
## SANT ISIDRE SANT LLORENS
## 4 4
## SANT MARCEL.LI SANT PAU
## 2 11
## SOTERNES TORMOS
## 1 2
## TORREFIEL TRES FORQUES
## 6 3
## TRINITAT VARA DE QUART
## 4 4Analisis de estaciones por distrito
table(estaciones$nombre_distrito)##
## ALGIROS BENICALAP BENIMACLET CAMINS AL GRAU
## 23 13 7 20
## CAMPANAR CIUTAT VELLA EL PLA DEL REAL EXTRAMURS
## 21 18 16 16
## JESUS L'EIXAMPLE L'OLIVERETA LA SAIDIA
## 11 19 16 13
## PATRAIX POBLATS DE L'OEST POBLATS MARITIMS QUATRE CARRERES
## 16 4 23 27
## RASCANYA
## 13Más
etsinf = subset(df_merge, df_merge$number_ ==114)
etsinf_ord=etsinf[order(etsinf$fecha,etsinf$hora_hora),]
dos_semanas = etsinf_ord[1:168, ]
dos_semanas$date = paste(dos_semanas$fecha, dos_semanas$hora_hora, sep=" ")
dos_semanas$date = as.POSIXct(dos_semanas$date, format="%d/%m/%Y %H:%M")
library(ggplot2)
xValue <- dos_semanas$date
yValue <- as.numeric(dos_semanas$avg_av)
data <- data.frame(xValue,yValue)
p = ggplot(data, aes(x=xValue, y=yValue, group = 1)) +
geom_line() +
geom_point()
p + (labs(x = "Horas", y = 'Bicis disponibles'))