Représenter des données
Marie Havret, Ugo Parment, Olivier Gillet
1 Gérer son espace de travail
# Reinitialiser (nettoyer) son espace de travail
# Lister les objets definis dans l'espace de travail actuel puis supprimer
rm(list=ls())
# Dans le répertoire "FORMATION_R", enregistrer un nouveau script R
# Fixer le repertoire courant (attention, indiquer un nom de repertoire valide!)
setwd('.')
# Afficher le repertoire de travail actuel
getwd()... [1] "C:/GIT/R_FORMATION_DEBUTANT"# Packages necessaires pour cette formation :
# install.packages('ggplot2')
library(ggplot2)
# install.packages('plotly')
library(plotly)
# install.packages("reshape")
library(reshape)
# install.packages("tidyr")
library(tidyr)
# install.packages("sf")
library(sf)
# install.packages("leaflet")
library(leaflet)
# install.packages("osmdata")
library(osmdata)2 Ouvrir les données
# Import des donnees precedemment creees, fonction 'load'
load('NYPD_Arrests_Data_2017.Rdata')3 Résumé de la session précédente
- Requêter et parcourir un tableau de données
- Manipuler et créer un tableau de données
- Calculer de nouvelles variables
# Compter le nombre d'infractions par jour, fonction 'table', et convertir en tableau de données (dataframe)
# Creation d'une variable temporaire
tmp <- as.data.frame(table(dfFilter$ARREST_DATE))
# Changer le nom des colonnes de la table
colnames(tmp) <- c('Date', 'nOfArrests')
# de la colonne 'Date'
tmp$Date <- as.POSIXct(as.character(tmp$Date),format="%Y-%m-%d")
# Representation graphique en utilisant les fonctions basiques de R
plot(tmp$Date, tmp$nOfArrests, type = "l", frame = TRUE, pch = NA,
col = "red", xlab = "Date", ylab = "Number of Arrests")
# Representation graphique en utilisant les fonctions de package ggplot2
ggplot(data=tmp,aes(x=Date,y=nOfArrests))+
geom_line()+
labs(x = "Date", y = "Number of Arrests")+
theme_bw() + theme(panel.grid.minor = element_blank())
4 Statistiques descriptives (2017)
# Utilisation de la fonction 'aggregate' pour generer quelques statistiques descriptives (minimum, moyenne, maximum ...)
# La fonction retourne les jours avec le plus et le moins d'infractions, la moyenne quotidienne d'infraction et le total mensuel d'infractions.
# Creer une nouvelle colonne contenant le numéro du mois
tmp[['Month']] <- format(tmp$Date, "%m")
# Realiser plusieurs 'aggregate', plusieurs fonctions (minimum, maximum, moyenne)
minArrests <- aggregate(nOfArrests ~ Month, data = tmp, FUN = min)
maxArrests <- aggregate(nOfArrests ~ Month, data = tmp, FUN = max)
sumArrests <- aggregate(nOfArrests ~ Month, data = tmp, FUN = sum)
meanArrests <- aggregate(nOfArrests ~ Month, data = tmp, FUN = mean)
# Creer un tableau de donnees avec toutes ces statistiques
statsByMonths <- cbind(minArrests, meanArrests[,2],maxArrests[,2], sumArrests[,2])
# Convertir les donnees du tableau
for(i in 2:ncol(statsByMonths)){
statsByMonths[,i] <- as.numeric(statsByMonths[,i])
}
# Renommer les colonnes du tableau
colnames(statsByMonths) <- c('Month', 'Min', 'Mean','Max', 'Sum')
# Autre manière de faire grâce a la fonction cast()
# Elle permet de combiner les fonctions statistiques precedentes
statsByMonths <- cast(Month ~ ., value = 'nOfArrests', data = tmp, fun = c(min, mean, max,sum))
statsByMonths <- as.data.frame(statsByMonths)
colnames(statsByMonths) <- c('Month', 'Min', 'Mean','Max', 'Sum')
# Supprimer les variables inutiles pour la suite du script
rm(tmp, minArrests, meanArrests, maxArrests, sumArrests)
# Mettre en forme la variable qui definit notre abscisse
# Concatener la variable 'Month' avec l'année et un jour du mois
statsByMonths$Month <- paste0('2017-',statsByMonths$Month,'-01')
# Convertir cette chaine de caractere en date
statsByMonths$Month <- as.POSIXct(as.character(statsByMonths$Month),format="%Y-%m-%d")
# Representer graphiquement les résultats
plot(statsByMonths$Month, statsByMonths$Min, type = "b", frame = TRUE,pch = 19, col = "blue",
xlab = "Date", ylab = "Number of Arrests",ylim=c(0,max(statsByMonths$Max)))
lines(statsByMonths$Month, statsByMonths$Mean, pch = 19, col = "grey", type = "b", lty = 2)
lines(statsByMonths$Month, statsByMonths$Max, pch = 19, col = "red", type = "b", lty = 2)
legend("bottomright", legend=c("Min", "Mean","Max"),col=c("blue", "grey", "red"), lty = 1:3, cex=0.5)
Le package “tidyr” a été développé pour manipuler et mettre en forme des tableaux de données. La fonction “gather” issue de ce package permet créer des modalités d’une nouvelle variable à partir de colonnes.
# Autre manière de représentation en utilisant ggplot
tmp <- statsByMonths[,c("Month","Min","Mean","Max")]
# Appel de la foncion 'gather' pour mettre en forme les donnees
tmp <- gather(tmp, key = "Variable", value = "Value", -Month)
ggplot(tmp, aes(x = Month, y = Value)) +
geom_line(aes(color = Variable)) +
scale_color_manual(values = c("steelblue", "grey", "darkred"))+
labs(x = "Date", y = "Number of Arrests")+
theme_bw() + theme(panel.grid.minor = element_blank())
# Supprimer les variables inutiles pour la suite du script
rm(statsByMonths)5 Exercice
Objectif: Calculer un résumé statistique pour le mois de juillet 2017.
- Sélectionner un mois particulier : ici le mois de juillet
- Compter le nombre d’infractions selon les genres et par quartier pour juillet
- Compter le nombre d’infractions selon les classes d’âge et par quartier
- Calculer le pourcentage d’homme et de femme ainsi que le total de crime
- Joindre les deux tables
- Identifier le crime majoritaire pour chaque quartier
dfjuly... Quartiers Femme perFemme Homme perHomme <18 18-24 25-44 45-64 65+ Total
... 1 Bronx 949 19.1 4011 80.9 248 1315 2491 883 23 4960
... 2 Brooklyn 1138 17.1 5531 82.9 370 1600 3309 1323 67 6669
... 3 Manhattan 1097 17.8 5049 82.2 287 1442 3010 1331 76 6146
... 4 Queens 810 17.7 3769 82.3 267 1151 2313 794 54 4579
... 5 Staten Island 215 20.8 817 79.2 86 226 524 188 8 1032
... Crime
... 1 DANGEROUS DRUGS
... 2 DANGEROUS DRUGS
... 3 DANGEROUS DRUGS
... 4 ASSAULT 3 & RELATED OFFENSES
... 5 DANGEROUS DRUGS6 Représentations graphiques
#Creation d'une colonne "WEEK" pour connaitre le numéro de la semaine
dfFilterjuly$WEEK <- format(dfFilterjuly$ARREST_DATE, "%V")
# Compter le nombre d'infractions par jour
nbInfPerDay <- aggregate(INFRACTION ~ ARREST_DATE, data = dfFilterjuly, FUN = length)
# Créer un graphique avec ggplot
ggplot(nbInfPerDay,aes(x = ARREST_DATE,y = INFRACTION))+
geom_line()+
geom_point()+
scale_x_datetime(name = "",date_breaks = "3 day", date_labels = "%a %d")+
ggtitle("Daily number of arrests (july 2017)")+
labs(x = "Date", y = "Number of Arrests")+
theme_bw() + theme(panel.grid.minor = element_blank())
# Stocker cet objet dans une nouvelle variable
c1 <- ggplot(nbInfPerDay,aes(x = ARREST_DATE,y = INFRACTION))+
geom_line()+
geom_point()+
scale_x_datetime(name = "",date_breaks = "3 day", date_labels = "%a %d")+
ggtitle("Daily number of arrests (july 2017)")+
labs(x = "Date", y = "Number of Arrests")+
theme_bw() + theme(panel.grid.minor = element_blank())
# Afficher ce graphique/Objet dans un viewer interactif (plotly package)
ggplotly(c1)# Compter le nombre d'infractions par jour et par genre
nbInfPerDay_Sex <- aggregate(INFRACTION ~ ARREST_DATE + SEX, data = dfFilterjuly, FUN = length)
# Stocker ce graphique dans une nouvelle variable
c2 <- ggplot(nbInfPerDay_Sex,aes(x = ARREST_DATE, y = INFRACTION, colour = SEX))+
geom_line()+
geom_point()+
scale_x_datetime(name = "",date_breaks = "3 day", date_labels = "%a %d")+
ggtitle("Daily number of arrests (july 2017)") +
labs(x = "Date", y = "Number of Arrests")+
theme_bw() + theme(panel.grid.minor = element_blank())
# Afficher ce graphique/Objet dans un viewer interactif (plotly package)
ggplotly(c2)7 Exercice
- Créer le graphique du nombre d’infractions commis la première semaine de juillet 2017
- Représenter graphiquement les infractions commises par tranche d’âge la première semaine de juillet 2017
- Représenter graphiquement les infractions des semaines 27 à 30 (x=jours / y= nb infractions)
8 Quelques améliorations graphiques
c1 <- ggplot(nbInfPerDay,aes(x = ARREST_DATE,y = INFRACTION))+
geom_line()+
geom_point()+
scale_x_datetime(name = "",date_breaks = "2 day",date_labels = "%A %d %B")+
ggtitle("Number of arrests by day, on July 2017")+
theme_bw()+
theme(axis.text.x = element_text(family = "cursive",
colour = "#A9A9A9",
angle = 75,
size=9,
face="bold"),
axis.text.y = element_text(colour = "#A9A9A9",
size=9,
face="italic"),
axis.title = element_text(colour = "black",
size=13,
face=3),
plot.title = element_text(colour = "black",
size=15,
face="bold",
vjust = 0.5),
panel.grid.minor = element_blank())
ggplotly(c1)c1 <- ggplot(nbInfPerDay_Sex,aes(x = ARREST_DATE, y = INFRACTION, colour = SEX))+
geom_line()+
geom_point()+
scale_x_datetime(name = "",date_breaks = "2 day",date_labels = "%A %d %B")+
ggtitle("Number of arrests by gender and by day, on July 2017")+
theme_bw()+
theme(axis.text.x = element_text(family = "cursive",
colour = "#A9A9A9",
angle = 75,
size=9,
face="bold"),
axis.text.y = element_text(colour = "#A9A9A9",
size=9,
face="italic"),
axis.title = element_text(colour = "black",
size=13,
face=3),
plot.title = element_text(colour = "black",
size=15,
face="bold",
vjust = 0.5),
panel.grid.minor = element_blank())
ggplotly(c1)9 Les données spatialisées
Les données spatialisées du package ‘sf’ correspondent à une table dont une colonne est dédiée à la géométrie. Chaque ligne correspond à un individu (ou entités), qui est associé à des variables (i.e un identifiant, nom du quartier, etc.) et à une forme. Ce dernier peut être un polygone, une ligne ou un point. Le package ‘sf’ met à disposition une large gamme de fonctions permettant les traitements spatiaux. Elles sont toutes de la forme ‘st_X’ où X correspond au traitement que l’on souhaite réaliser, par exemple une intersection ou un buffer. Dans notre cas, nous avons à disposition plusieurs table. Par exemple, la table ‘shpQuartier’, où chaque ligne correspond à un quartier associé à une géométrie.
# Reinitialiser (nettoyer) son espace de travail
# Lister les objets definis dans l'espace de travail actuel puis supprimer
rm(list=ls())
# Ouvrir le fichier avec les statistiques -
load('NYPD_Arrests_Statistics_2017.Rdata')
# Ouvrir le shapefile
shpQuartier <- read_sf("Quartiers.shp", as_tibble = F)# Afficher le systeme de projection du shapefile
st_crs(shpQuartier)... Coordinate Reference System:
... User input: 4326
... wkt:
... GEOGCS["WGS84(DD)",
... DATUM["WGS84",
... SPHEROID["WGS84",6378137.0,298.257223563]],
... PRIMEM["Greenwich",0.0],
... UNIT["degree",0.017453292519943295],
... AXIS["Geodetic longitude",EAST],
... AXIS["Geodetic latitude",NORTH],
... AUTHORITY["EPSG","4326"]]# Modifier le systeme de projection
shpQuartier <- st_transform(shpQuartier,4326)
# Afficher l'emprise spatiale du jeu de donnees
st_bbox(shpQuartier$geometry)... xmin ymin xmax ymax
... -74.25559 40.49612 -73.70001 40.91553# Afficher le shapefile - toutes les colonnes de la table attributaire
plot(shpQuartier, col="white")
# Afficher le shapefile - seulement la deuxieme colonne
plot(shpQuartier[2])
# Afficher seulement la geometrie du shapefile
plot(st_geometry(shpQuartier), col="white")
# Afficher les centroïdes des zones
plot(st_centroid(shpQuartier$geometry), pch = 20, col = "lightblue", add = T)
# Joindre les statistiques au shapefile
shpQuartier <- merge(shpQuartier, df, by.x="boro_name", by.y='Quartiers')# Afficher le debut de la table
head(shpQuartier[,c('boro_name','Crime','geometry')])... Simple feature collection with 5 features and 2 fields
... geometry type: MULTIPOLYGON
... dimension: XY
... bbox: xmin: -74.25559 ymin: 40.49612 xmax: -73.70001 ymax: 40.91553
... CRS: EPSG:4326
... boro_name Crime geometry
... 1 Bronx DANGEROUS DRUGS MULTIPOLYGON (((-73.89681 4...
... 2 Brooklyn DANGEROUS DRUGS MULTIPOLYGON (((-73.86706 4...
... 3 Manhattan DANGEROUS DRUGS MULTIPOLYGON (((-74.01093 4...
... 4 Queens ASSAULT 3 & RELATED OFFENSES MULTIPOLYGON (((-73.83668 4...
... 5 Staten Island DANGEROUS DRUGS MULTIPOLYGON (((-74.05051 4...# Afficher le type de l'objet
class(shpQuartier)... [1] "sf" "data.frame"# Afficher le nombre total de crimes pour 2017
plot(shpQuartier[14])
9.1 Collecter des données
# Construire une requête OSM
q <- getbb("New York")%>%
opq()%>%
add_osm_feature("amenity", "police")
# Lancer la requête, fonction 'osmdata_sf'
osmPolice<- osmdata_sf(q)
# Obtenir un objet SF contenant les points
shpPolice <- osmPolice$osm_points
# Afficher le type de l'objet
class(shpPolice)
# Supprimer les variables inutiles pour la suite du script
rm(q, osmPolice)
# Enregistrer le fichier
save(shpPolice, file='Commissariat.Rdata')
# /!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\
# /!\ Charger les données si absence de connexion internet
# load('DONNEES//Commissariat.Rdata')
# /!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\/!\
# Afficher seulement la geometrie du shapefile
plot(st_geometry(shpQuartier), col="white")
# Afficher le shapefile des commissariats
plot(st_geometry(shpPolice), pch = 20, add=T, col='blue')
# Reprojeter les donnees dans le même système de projection
shpPolice <- st_transform(shpPolice, st_crs(shpQuartier))
# Compter le nombre de commissariats pour chaque quartier
shpQuartier$Commissariat <- lengths(st_intersects(shpQuartier, shpPolice))9.2 Requêtes spatiales
# Ouvrir le shapefile des districts policiers
shpDistrict <- read_sf("Precincts_Boundaries.shp", as_tibble = F)# Reprojeter les donnees dans le meme systeme de projection
shpDistrict <- st_transform(shpDistrict, st_crs(shpQuartier))
# Compter le nombre de commissariat pour chaque district
shpDistrict$Commissariat <- lengths(st_intersects(shpDistrict, shpPolice))# Import de donnees, table 'df' et 'dfFilter', fonction 'load'
load('NYPD_Arrests_Data_2017.Rdata')# Creer des donnees spatialisées à partir d'un tableau de donnees
# La fonction st_as_sf() permet de transformer un tableau au format sf, en objet spatialise.
# Il convient egalement de définir le systeme de coordonnees
shpCrime <- st_as_sf(dfFilter,coords = c("LONGITUDE","LATITUDE"),crs=4326)
# Reprojeter les donnees dans le meme systeme de projection
shpCrime <- st_transform(shpCrime, st_crs(shpQuartier))
# Compter le nombre d'infractions pour chaque district
shpDistrict$Crime <- lengths(st_intersects(shpDistrict, shpCrime))9.3 Leaflet, viewer interactif
Le package Leaflet a été developpé pour manipuler des données spatialisées et créer des cartes interactives (viewer)
# Preparation de la carte leaflet
# Utilisation du pipe %>% pour chainer les fonctions du package leaflet
leaflet()%>%
addTiles() #Fond de carte OpenStreetMap par défaut# Centrer la carte sur New York
coord <- st_bbox(shpQuartier$geometry)# Permet de connaitre la couverture spatiale des quartiers de NY
lng_NY <- as.numeric((coord["xmin"] + coord['xmax']) / 2) #(Longitude min + Longitude max) / 2
lat_NY <- as.numeric((coord["ymin"] + coord['ymax']) / 2) #(Latitude min + Latitude max) / 2
#Creer une liste reprenant les centroides
shpQuartier_Centroids <- st_centroid(shpQuartier)
# Changement du fond de carte
leaflet()%>%
setView(lng = lng_NY, lat =lat_NY, zoom = 10)%>%
addProviderTiles("Stamen.Toner")# Modification du fond cartographique
leaflet()%>%
setView(lng = lng_NY, lat = lat_NY, zoom = 10)%>%
addProviderTiles("Esri.NatGeoWorldMap")# Creer un objet permettant d'afficher les quartiers de NY en rouge
map_NY <- leaflet() %>%
addTiles() %>%
addPolygons(data = shpQuartier,
fillColor = "red", # Couleur du polygone
weight = 2, # Epaisseur des bordures
opacity = 1, # Parametre obligatoire - Opacite des bordures
color = "black", # Couleur des bordures
dashArray = "3", # Type de bordure (1: continu, 2 & 3:pointillés)
fillOpacity = 0.4, # Opacite du fond
popup = shpQuartier$boro_name) # Popup avec nom des quartiers
#Afficher la carte
map_NY# Completer la carte en affichant le centroide de chaque quartier
map_NY <- addCircles(map = map_NY,
data = shpQuartier_Centroids,
color = "yellow",
fill = "yellow",
opacity = 1)
#Afficher la carte
map_NY# Ajout de fenetre interactive à la carte (pop-up)
# Creation du texte a afficher lors d'un clic sur l'objet cible
popup_info <- paste0("Quartier : ", shpQuartier$boro_name,"<br>", #Nom du quartier
"Infractions : ",shpQuartier$Total) # Nombre d'infractions
map_NY <- leaflet() %>%
addTiles() %>%
addPolygons(data = shpQuartier,
fillColor = "red",
weight = 2,
opacity = 1,
color = "black",
dashArray = "3",
fillOpacity = 0.4,
popup = popup_info) # Appel de popup-info créé au-dessus
#Afficher la carte
map_NY# Carte des densites d'infractions par quartiers
# Utilisation d'une palette graduee ("YlOrRd")
# Voir les palettes disponibles : http://www.sthda.com/french/wiki/couleurs-dans-r
pal <- colorBin("YlOrRd", # Palette choisie
domain = shpQuartier$Total, # Variable ciblee par la palette
bins = 4) # Nombre de classes
map_NY <- leaflet() %>%
addTiles() %>%
addPolygons(data = shpQuartier,
fillColor = ~pal(Total),
weight = 2,
opacity = 1,
color = "black",
dashArray = "3",
fillOpacity = 0.4,
popup = popup_info) %>%
#Ajouter la legende
addLegend(pal = pal,
values = shpQuartier$Total,
position = "bottomright")
# Afficher la carte
map_NY9.4 Exercice
- Réaliser une carte des densités d’infractions par district