Pour continuer sur les manipulation de données publiques, je voulais m’inspirer du projet de Cédric, de la formation Data Science pour l’Actuariat sur les résultats au baccalauréat. Les données nécessaires à cette étude sont disponibles sur plusieurs sites,
- Ministère de l’éducation nationale
- Bases INSEE
- NosDonnees.fr
Il ne s’agit aucunement d’une analyse poussée des résultats, juste un peu de visualisation, sans aucune autre prétention ! Ah oui, même si on ne va pas faire de carte (je les trouve peu lisibles) on va quand même utiliser les données spatiales : les établissements scolaires sont géolocalises, et on peut obtenir des informations locales, sur le taux de chômage, ou le revenu médian. Et faire des graphiques.
Ce préambule passé, on peut commencer.
library(dplyr) library(readxl) library(sp) library(ggmap) library(raster) library(leaflet) library(DT) library(cowplot) library(gstat) library(tmap) |
On va commencer par récupérer par établissement, les résultats au bac.
url_resultat_etab = "https://data.education.gouv.fr/explore/dataset/fr-en-indicateurs-de-resultat-des-lycees-denseignement-general-et-technologique/download/?format=csv&timezone=Europe/Berlin&use_labels_for_header=true" download.file(url_resultat_etab,destfile = paste0(librairie,"import_resultat_etab.csv"), method="curl") df_resultat_etab = read.csv("import_resultat_etab.csv",header=TRUE, sep= ";", encoding="UTF-8") |
Comme bien souvent avec les données des administrations françaises, on a souvent des soucis de typographie. Pour simplifier, on va supprimer les accents, et uniformiser un peu les noms
MiseEnForme_Colonnes = function(text) { text <- gsub("è", "e", text) text <- gsub("é", "e", text) text <- gsub("_", ".", text) text <- gsub("serie.", "", text) text <- gsub("Effectif.Presents.", "Effectif.", text) text <- gsub("Taux.","Tx.",text) text <- gsub("Brut.de.Reussite.", "Admis.Etab.", text) text <- gsub("Reussite.Attendu.", "Admis.", text) text <- gsub("brut", "Etab", text) text <- gsub("attendu", "Academie", text) text <- gsub("toutes.", "TOTAL", text) text <- gsub("Total.", "TOTAL", text) text <- gsub("..Etablissement", ".Etab", text) text <- gsub("Pourcentage", "Tx", text) return(text) } for(i in 1:ncol(df_resultat_etab)){ colnames(df_resultat_etab)[i] <- MiseEnForme_Colonnes(names(df_resultat_etab)[i]) } |
On va ensuite supprimer les départements et régions d’outre-mer,
df_resultat_etab = df_resultat_etab[-which(toupper(df_resultat_etab$Departement) %in% c("GUADELOUPE","MAYOTTE","MARTINIQUE","REUNION","GUYANE")),] |
récupérer les noms des colonnes
Colonnes = colnames(df_resultat_etab) |
et comme on s’intéresse aux premières variables
Colonnes_Generiques = Colonnes[1:8] |
on les recupere, pour construire quelques statistiques pour colonnes relatives aux séries L, ES et S
Colonnes_Series = Colonnes[grepl("([a-zA-Z]*?.)*\\.S$|([a-zA-Z]*?.)*\\.ES$|([a-zA-Z]*?.)*\\.L$|([a-zA-Z]*?.)*\\.TOTAL$",Colonnes)] |
Et on finit avec les autres
Colonnes_Autres = Colonnes[grepl("(Tx.Bacheliers.*)|(Tx.acces.*)|(Effectif.de.*)|(libelle.region)|(code.region)|(element)",Colonnes)] df_resultat_etab = cbind(df_resultat_etab[Colonnes_Generiques],df_resultat_etab[Colonnes_Series],df_resultat_etab[Colonnes_Autres]) |
On peut aussi localiser les établissements
url_carto_etab <- "https://www.data.gouv.fr/s/resources/adresse-et-geolocalisation-des-etablissements-denseignement-du-premier-et-second-degres/20160526-143453/DEPP-etab-1D2D.csv" download.file(url_carto_etab,destfile=paste0(librairie,"import_carto_etab.csv")) df_carto_etab = read.csv2("import_carto_etab.csv",header=TRUE |
On récupère ici la géolocalisation de 66556 établissements nationaux ! On peut croiser avec des données socio-économiques des communes
nom_base_emploi = "base-cc-emploi-pop-act-2014" url_baseemploi_popactive = paste0("https://www.insee.fr/fr/statistiques/fichier/2862207/",nom_base_emploi,".zip") download.file(url_baseemploi_popactive,destfile=paste0(librairie,nom_base_emploi,".zip")) unzip(paste0(nom_base_emploi,".zip"),overwrite = TRUE) df_base_emploi_source = read_excel(paste0(nom_base_emploi,".xls"),sheet="COM_2014",skip=5) |
On va exclure les territoires d’outre-mer ici
df_base_emploi_source <- df_base_emploi_source[-which(df_base_emploi_source$DEP %in% c("971","972","973","974","975")),] df_base_emploi_colonnes = c("CODGEO","P14_POP1564","P14_H1564","P14_F1564","P14_ACT1564","P14_ACTOCC1564","P14_CHOM1564","P14_INACT1564", "P14_ETUD1564", "P14_RETR1564", "P14_AINACT1564", "P14_HCHOM1524", "P14_FCHOM1524", "C14_ACT1564","C14_ACT1564_CS1","C14_ACT1564_CS2","C14_ACT1564_CS3","C14_ACT1564_CS4","C14_ACT1564_CS5","C14_ACT1564_CS6","P14_POP15P") df_base_emploi = df_base_emploi_source[,names(df_base_emploi_source) %in% df_base_emploi_colonnes] |
et corriger les soucis classiques de la Corse,
MiseEnForme_CodeGeo = function(text) { text <- gsub("2A", "20", text) text <- gsub("2B", "20", text) } df_base_emploi$CODGEO = MiseEnForme_CodeGeo(df_base_emploi$CODGEO) |
On peut aussi utiliser des données de revenus, par communes
nom_base_revenus = "indic-struct-distrib-revenu-2014-COMMUNES" url_baserevenus = paste0("https://www.insee.fr/fr/statistiques/fichier/3126151/",nom_base_revenus,".zip") download.file(url_baserevenus,destfile=paste0(librairie,nom_base_revenus,".zip")) unzip(paste0(nom_base_revenus,".zip"),overwrite = TRUE) df_base_revenus = read_excel("FILO_DISP_COM.xls",sheet="ENSEMBLE",skip=5)[,c(1,4,7)] df_base_revenus$CODGEO = MiseEnForme_CodeGeo(df_base_revenus$CODGEO) |
On recupere des donnees spatiales relatives aux communes
url_geoloc_communes = "http://www.nosdonnees.fr/wiki/images/b/b5/EUCircos_Regions_departements_circonscriptions_communes_gps.csv.gz" download.file(url_geoloc_communes,destfile=paste0(librairie,"geoloc_communes.csv.gz")) df_geoloc_communes = read.csv2(gzfile("geoloc_communes.csv.gz"),header=TRUE, stringsAsFactors = FALSE,encoding="UTF-8") |
et comme toujours, un peu de corrections s’imposent
df_geoloc_communes = df_geoloc_communes[-which(df_geoloc_communes$numéro_département %in% c("971","972","973","974","975")),] df_geoloc_communes = df_geoloc_communes[,names(df_geoloc_communes) %in% c("code_insee","latitude","longitude","codes_postaux")] df_geoloc_communes_nb <- nrow(df_geoloc_communes) |
On va ensuite creer une fonction de remplacement des valeurs manquantes, et de correction des séparateurs décimaux
MiseEnForme_CoordonneesGeo = function(valeur){ pretraitement = ifelse(as.character(valeur)=="-","0",as.character(valeur)) traitement = as.numeric(ifelse(pretraitement==".","0",gsub(pattern=",",replacement=".",pretraitement))) return(traitement) } df_geoloc_communes$latitude = MiseEnForme_CoordonneesGeo(df_geoloc_communes$latitude) df_geoloc_communes$longitude = MiseEnForme_CoordonneesGeo(df_geoloc_communes$longitude) |
On passe ensuite a l’élimination des lignes en double
df_geoloc_communes = unique(df_geoloc_communes) |
On va ensuite changer les noms des colonnes pour harmoniser avec les autres bases
names(df_geoloc_communes) = c("Codes_Postaux","CODGEO","coordonnee_y","coordonnee_x") |
On peut ensuite rechercher les lignes en double sur les codes insee
liste_CODGEO2 = aggregate(x=df_geoloc_communes$Codes_Postaux,by=list(df_geoloc_communes$CODGEO),FUN="length") list_geoloc_communes_CODGEO2 = liste_CODGEO2[liste_CODGEO2$x>1,1] df_geoloc_communes_CODGEO2 = df_geoloc_communes[df_geoloc_communes$CODGEO %in% list_geoloc_communes_CODGEO2,1:2] |
Ici, un correction manuelle s’impose pour 4 configurations : les données propres à Lyon, Paris et Marseille ne sont pas géolocalisées, les données propres à la ville de Laguépie sont géolocalisées en doubles
df_geoloc_communes_propre = df_geoloc_communes[!df_geoloc_communes$CODGEO %in% list_geoloc_communes_CODGEO2,] df_geoloc_communes_corrige = data.frame(Codes_Postaux=c("13001","69001","75001","82250"), CODGEO=c("13055","69123","75056","82088"), coordonnee_y=c(43.3,45.75,48.85,44.15), coordonnee_x=c(5.4,4.85,2.31,1.97)) df_geoloc_communes = rbind(df_geoloc_communes_propre,df_geoloc_communes_corrige) |
On peut enfin fusionner les bases
df_etab = merge(df_resultat_etab,df_carto_etab,by="Cod.Etab") |
Certains établissements ne peuvent être géolocalisés pour certaines années
df_etab_total_nongeolocalises <- df_resultat_etab[!df_resultat_etab$Cod.Etab %in% df_carto_etab$Cod.Etab,] |
Comme l’étude ne porte que sur les seuls lycées d’enseignement général et technologique, le dataframe est réduit aux observations relatives d’une part aux lycées, d’autre part aux établissements d’enseignement polyvalent, général ou général et technologique.
df_etab = df_etab[grep("LYCÉE",toupper(df_etab$nature_uai_libe)),] df_etab = df_etab[grep("GÉNÉRAL|POLYVALENT",toupper(df_etab$nature_uai_libe)),] df_etab_nongeolocalises = df_etab[df_etab$Cod.Etab %in% df_etab_total_nongeolocalises$Cod.Etab,] df_etab_geolocalise = df_etab[!is.na(df_etab$coordonnee_x),] df_etab_geolocalise = df_etab_geolocalise[!is.na(df_etab_geolocalise$coordonnee_y),] |
Enfin, on va convertir les code géographiques (ici reconnus comme facteur) en chaines de caractères (de 5 caractères) pour pouvoir fusionner les tables
ConvertCODGEO = function(code) { if(is.character(code)) { code_character = ifelse(nchar(code)<5, paste0("0",code), code) return(code_character) } else if(is.factor(code)){ code_character = ifelse(code<10000, paste0("0",as.numeric(as.character(code))), as.numeric(as.character(code))) return(code_character) } else if(is.numeric(code)){ code_character = ifelse(code<10000, paste0("0",code), as.character(code)) return(code_character) } } df_etab_geolocalise$Code.commune = ConvertCODGEO(df_etab_geolocalise$Code.commune) df_etab_geolocalise$Secteur.Public.Prive = sapply(df_etab_geolocalise$Secteur.Public.Prive,function(nature) {ifelse(nature=="PU","Lycées Publics","Lycées Privés")}) |
On conservation alors les établissements dont la commune n’est pas manquante
df_etab_geolocalise = df_etab_geolocalise[!is.na(df_etab_geolocalise$Code.commune),] |
Pour finir, on va creer une base, pour ensuite faire une graphique
tbl_etab_nature_res_source = df_etab_geolocalise[,c(3,8,9,10,11,13,14,15)] for(i in c(6,7,8)){ temp = tbl_etab_nature_res_source[!is.na(tbl_etab_nature_res_source[i]),c(1,2,i-3,i)] temp$Serie = ifelse(i==6,"L",ifelse(i==7,"ES","S")) names(temp)[2:4] = c("Nature","Effectif","Tx.Admis") if(i==6){ tbl_etab_nature_result = temp } else{ tbl_etab_nature_result = rbind(tbl_etab_nature_result,temp) } } graph = ggplot(tbl_etab_nature_result,aes(x=Effectif,y=Tx.Admis,colour=factor(Annee))) graph = graph + geom_point(alpha=0.45) graph = graph + facet_grid(Serie~Nature) graph = graph + xlab("Effectifs de l'établissement en terminale (par série)") + ylab("Taux d'admission (%)") graph = graph + scale_color_discrete(name="Année des\nrésultats") graph = graph + theme(legend.title = element_text(size=9,face="bold"), legend.text = element_text(size=9), strip.background = element_rect(colour="black", fill="gray95"), panel.border = element_rect(linetype = "solid"), panel.grid.major = element_line(colour = "gray75",linetype = "dashed"), panel.grid.minor = element_line(colour = "gray95",linetype = "dashed"), axis.title.x = element_text(size=9, face="bold"), axis.text.x = element_text(size=8), axis.title.y = element_text(size=9, face="bold"), axis.text.y = element_text(size=8)) graph |
On a ici l’evolution des resultats en fonction de la taille des etablissements.
df_communes_CorrNaN = df_communes_Corr[which(!df_communes_Corr$TxChomage == "NaN" & !df_communes_Corr$TxCadres == "NaN" & !df_communes_Corr$TxOuvriers == "NaN" & !df_communes_Corr$NbPopulation == "NaN" & !df_communes_Corr$TxSenior == "NaN" & !df_communes_Corr$RevenusMedians == "NaN"),] df_communes_sp = SpatialPointsDataFrame(coords = df_communes_CorrNaN[, c("coordonnee_x", "coordonnee_y")], data = df_communes_CorrNaN) Grille = as.data.frame(makegrid(df_communes_sp, nsig=2, cellsize = 0.1)) names(Grille) = c("X", "Y") coordinates(Grille) = c("X", "Y") gridded(Grille) = TRUE fullgrid(Grille) = TRUE proj4string(Grille) = proj4string(df_communes_sp) |
On peut ensuite faire du krigeage, histoire de lisser un peu nos donnees de chomage et de revenu
df_communes_sp.TxChomage = krige(TxChomage ~ 1, df_communes_sp, Grille, nmax=1) df_communes_sp.RevenusMedians = krige(RevenusMedians ~ 1, df_communes_sp, Grille, nmax=1) sp_lycee_WGS84@data$TxChomage = extract(R.TxChomage,sp_lycee_WGS84) sp_lycee_WGS84@data$TxCadres = extract(R.TxCadres,sp_lycee_WGS84) sp_lycee_WGS84@data$TxOuvriers = extract(R.TxOuvriers,sp_lycee_WGS84) sp_lycee_WGS84@data$NbPopulation = extract(R.NbPopulation,sp_lycee_WGS84) sp_lycee_WGS84@data$TxSenior = extract(R.TxSenior,sp_lycee_WGS84) sp_lycee_WGS84@data$RevenusMedians = extract(R.RevenusMedians,sp_lycee_WGS84) |
On peut enfin conclure, en faisant une fonction generique de visualisation
Creation_Graphique = function(df, AnneeObs_Ouv, AnneeObs_Clo, Effectifs, Abscisses, Ordonnees, TitreAbs, TitreOrd, CouleurGraph, CouleurLiss, Serie) { df_temp = df[which(df$Annee>=AnneeObs_Ouv & df$Annee<=AnneeObs_Clo),] df_temp = df_temp[which(!is.na(df_temp[,Effectifs]) & !is.na(df_temp[,Abscisses]) & !is.na(df_temp[,Ordonnees])),] df_temp = df_temp[!df_temp[,Effectifs]==0,] df_temp = df_temp[,c(Abscisses,Ordonnees)] graphique = ggplot(df_temp,aes(x = df_temp[,Abscisses],y = df_temp[,Ordonnees])) graphique = graphique + geom_point(data = df_temp, aes(x = df_temp[,Abscisses],y = df_temp[,Ordonnees]),size=1, color=CouleurGraph,alpha=0.25) graphique = graphique + geom_density2d(aes(colour=..level..),show.legend=F) + scale_colour_gradient(low="gray55",high="gray25") graphique = graphique + scale_y_continuous(breaks= seq(80,100,by=2), limits = c(80,100)) graphique = graphique + xlab(TitreAbs) + ylab(TitreOrd) graphique = graphique + ggtitle(Serie) graphique = graphique + theme(plot.title = element_text(size=13,color=CouleurLiss, face="bold", hjust=0), axis.title.x = element_text(size=8, face="bold"), axis.text.x = element_text(size=8), axis.title.y = element_text(size=8, face="bold"), axis.text.y = element_text(size=8), panel.border = element_rect(linetype = "solid"), panel.grid.major = element_line(colour = "gray55",linetype = "dashed"), panel.grid.minor = element_line(colour = "gray75",linetype = "dashed")) graphique = graphique + stat_smooth(method = "loess",fill=CouleurLiss,color=CouleurLiss) return(graphique) } Production_Graphique_VI_1 = function(df, Titre_General, Axe_Abscisses, Titre_Abscisses, Annee_Observee_Ouv, Annee_Observee_Clo){ Graph_S = Creation_Graphique(df, Annee_Observee_Ouv, Annee_Observee_Clo, "Effectif.S", Axe_Abscisses, "Tx.Admis.Etab.S", Titre_Abscisses, "Taux d'admission (%)", "dodgerblue3","dodgerblue4","Série S") Graph_ES = Creation_Graphique(df, Annee_Observee_Ouv, Annee_Observee_Clo, "Effectif.ES", Axe_Abscisses, "Tx.Admis.Etab.ES", Titre_Abscisses, "Taux d'admission (%)","darkorange2","darkorange3","Série ES") Graph_L = Creation_Graphique(df, Annee_Observee_Ouv, Annee_Observee_Clo, "Effectif.L", Axe_Abscisses, "Tx.Admis.Etab.L", Titre_Abscisses, "Taux d'admission (%)", "chartreuse4","darkgreen","Série L") Graph_TS = Creation_Graphique(df, Annee_Observee_Ouv, Annee_Observee_Clo, "Effectif.Etab", Axe_Abscisses, "Tx.Admis.Etab", Titre_Abscisses, "Taux d'admission (%)", "indianred1","red4","Toutes séries") p = plot_grid(Graph_S, Graph_ES, Graph_L, Graph_TS, ncol = 2, nrow = 2,align = 'hv', scale = c(0.95, 0.95, 0.95, 0.95),vjust = 0.9, hjust=-0.5) titre <- ggdraw() + draw_label(Titre_General,fontface="bold", size=10) plot_grid(titre, p, ncol = 1, rel_heights=c(.25,5)) } |
On note ici
df_lycee <- sp_lycee_WGS84@data |
et on peut faire un premier graphique, avec le taux de chomage
Production_Graphique_VI_1(df = df_lycee, Titre_General = "Taux d'admission par série en fonction du taux de chômage \n dans la population active - Tous lycées confondus", Axe_Abscisses = "TxChomage", Titre_Abscisses = "Taux de chômage dans la population active (%)", Annee_Observee_Ouv = "2013", Annee_Observee_Clo = "2015") |
On ne va pas enfoncer les portes ouvertes de l’inference ecologique, en affirmant des choses aussi stupides que “on a moins de chances d’avoir le bac quand on est au chômage”. Mais on peut noter que dans les zones avec un fort taux de chômage, les résultats au bac sont moins bons.
On peut ensuite regarder en fonction du revenu de la commune du lycée
Production_Graphique_VI_1(df = df_lycee, Titre_General = "Taux d'admission par série en fonction du niveau des revenus disponibles médians - Tous lycées confondus", Axe_Abscisses = "RevenusMedians", Titre_Abscisses = "Quantile du niveau des revenus disponibles médians (%)", Annee_Observee_Ouv = "2013", Annee_Observee_Clo = "2015") |
Fascinant, non ? Mais c’est clairement juste une première approche… il faudrait aller plus loin ensuite !
