Commit d950bb7a authored by MathieuGO's avatar MathieuGO Committed by Bakaniko

ACP/CAH SCRIPT intégré au rapport

parent 0e470cea
### ACP réalisée sur un bouquet de variables ###
# Liste des libraries R nécéssaires
library(factoextra)
library(corrplot)
library(readOGR)
library(rCarto)
# Liste des variables retenues pour l'ACP:
# PERTE_INDIGENA_1 [200]-> PERTE_INDIGENA = SI -> Part de la population se considèrant comme indigène
# TENENCIA_1 [161]-> Part des propriétaires
# CLAVIVP_7 [22] -> Local no construido para habitación -> Part de l'habitat précaire ??
# PISOS_2 [27] -> Cemento o firme -> Par des habitations dont le sol est en ciment
# TotPersCuart -> Nombre moyen d'occupant par pièces
# CONACT_10 [2]-> trabajó? -> Part des personnes qui ont un travail
# DESTINO_BASURA_6 [126] -> la tiran en otro lugar? (Calle, baldío, río) -> Part des vivienda dont les déchets sont simplement jetés à l'extérieur de l'habitation
# AGUA_ENTUBADA_3 [89] -> ¿No tienen agua entubada? -> Part des vivienda sans eau courante
# NIVACAD_13 [224] -> Maestría -> Part des personne ayant un niveau de scolarité équivalent au master
# DRENAJE_5 [119] -> ¿No tiene drenaje? -> Part des habitations qui n'ont pas de système d'évacuation des eaux usées
# Import des données brut sur les vivienda
tab_frequence_2015_viv <-read.csv("../data/traitees/tab_frequence_2015_viv.csv",header=TRUE, sep=",")
# Sélection des variables relatives au nombre de personnes par vivienda
Tab_NUMPERS <- tab_frequence_2015_viv[,c(152:178)]
# La formule pour calculer le nombre d'occupant moyen par unité d'habitation:
# somme(nombre de vivienda concerné *nombre de personnes par vivienda)/somme(nombre de vivianda par municipalité)
NbrOccupant<-(Tab_NUMPERS[,1]*1/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,2]*2)/(rowSums(Tab_NUMPERS))+
(Tab_NUMPERS[,3]*3/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,4]*4/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,5]*5/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,6]*6/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,7]*7/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,8]*8/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,9]*9/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,10]*10/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,11]*11/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,12]*12)/(rowSums(Tab_NUMPERS))+
(Tab_NUMPERS[,13]*13/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,14]*14/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,15]*15/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,16]*16/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,17]*17/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,18]*18/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,19]*19/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,20]*20/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,21]*21/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,22]*22)/(rowSums(Tab_NUMPERS))+
(Tab_NUMPERS[,23]*23/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,24]*24/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,25]*25/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,26]*26/(rowSums(Tab_NUMPERS)))+
(Tab_NUMPERS[,27]*30/(rowSums(Tab_NUMPERS)))
# Calcule du nombre moyen de pièce par vivienda
# Sélection des variables relatives au nombre de pièces par vivienda
Tab_TOTCUART <- tab_frequence_2015_viv[,c(16:40)]
# La formule pour calculer le nombre moyen de pièce par unité d'habitation:
# somme(nombre de vivienda concerné *nombre de pièce par vivienda)/somme(nombre de vivianda par municipalité)
NbrPiece<-(Tab_TOTCUART[,1]*1/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,2]*2)/(rowSums(Tab_TOTCUART))+
(Tab_TOTCUART[,3]*3/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,4]*4/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,5]*5/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,6]*6/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,7]*7/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,8]*8/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,9]*9/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,10]*10/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,11]*11/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,12]*12)/(rowSums(Tab_TOTCUART))+
(Tab_TOTCUART[,13]*13/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,14]*14/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,15]*15/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,16]*16/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,17]*17/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,18]*18/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,19]*19/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,20]*20/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,21]*21/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,22]*22)/(rowSums(Tab_TOTCUART))+
(Tab_TOTCUART[,23]*23/(rowSums(Tab_TOTCUART)))+
(Tab_TOTCUART[,24]*24/(rowSums(Tab_TOTCUART)))
# Calcule du nombre moyen d'occupants par pièce et par unité d'habitation
TotPersCuart <- as.data.frame(NbrOccupant/NbrPiece)
# Import des données relatives aux autres variables retenues pour l'ACP
tab_freq_pourc_pers_viv <- readRDS("../data/traitees/tab_freq_pourc_pers_viv.rds")
# Jointure des variables retenues issue de tab_freq_pourc_pers_viv + variable nombre moyen d'occupants par pièce par vivienda
ACP_Data <- merge(TotPersCuart,tab_freq_pourc_pers_viv[,c(1,2,27,22,89,119,126,161,200,224)], by="code_mun")
# Utilisation du code municipios comme nom de ligne
ACP_Data_2 <- data.frame(ACP_Data,row.names=ACP_Data$code_mun)
# Centrage-réduction des données
ACP_Data_3 <- scale(ACP_Data_2[,c(2:11)])
# Utilisation de la fonction princomp pour réaliser l'ACP
PCR<- princomp(ACP_Data_3, cor = TRUE)
# Extraction de la variance et la variance cumulée pour chaque dimensions
eig.val <- get_eigenvalue(PCR)
eig.val
# Représentation graphique de la variance pour chaque dimension sous forme d'histogramme
fviz_eig(PCR, addlabels = TRUE, ylim = c(0, 50))
# Extraction des résultats de l'ACP pour les variables et toutes les dimensions
ResultsPCA <- get_pca_var(PCR)
# Représentation graphique du Cos² obtenu pour chaque variable pour toutes les dimensions
corrplot(ResultsPCA$cos2, is.corr=FALSE)
# Représentation graphique du Cos² obtenu pour chaque variable pour les dimensions 1 et 2
fviz_pca_var(PCR, col.var = "cos2",gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),repel = TRUE
# Représentation graphique des contributions de chaque variable pour toutes les dimensions
corrplot(ResultsPCA$contrib, is.corr=FALSE)
# Représentation graphique des contributions de chaque variable pour la dimension 1
fviz_contrib(PCR, choice = "var", axes = 1, top = 10)
# Représentation graphique des contributions de chaque variable pour la dimension 2
fviz_contrib(PCR, choice = "var", axes = 2, top = 10)
# Représentation des Cos² obtenus par les individus pour la dimension 1 et 2
fviz_pca_ind (PCR, col.ind = "cos2",gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),repel = TRUE
# Représentation graphique des contributions des individus pour les dimensions 1 et 2
fviz_pca_ind (PCR, col.ind = "contrib",gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),repel = TRUE
# Représentation graphique des individus et variables pour les dimensions 1 et 2
fviz_pca_biplot(PCR, repel = TRUE,col.var = "#2E9FDF", col.ind = "#696969"
# Cartographie des résulats de l'ACP
# Import de la géométrie des municipios sous forme de shapefile
Municipios_2015 <- readOGR(dsn= "../data/areas_geoestadisticas_municipales.shp", layer = "areas_geoestadisticas_municipales",encoding= "utf8", stringsAsFactors = FALSE)
# Choix de la projection WGS84
projection <- "+proj=utm +zone=14 +datum=WGS84 +units=m +no_defs + ellps=WGS84"
# Application aux données géométriques de la projection WGS84
Municipios_2015 <- spTransform(Municipios_2015, projection)
# Création d'un identifiant code municipios pour les données géométriques
Municipios_2015$ID <- paste(Municipios_2015$CVE_ENT, Municipios_2015$CVE_MUN, sep="")
# Extraction des résultats de la 1ère dimension de l'ACP
PC1 <- PCR$scores[,1] ; round(PC1, 2)
# Aggrégation des résulats de l'ACP aux données centrées-réduites pour relier code municipios et résultats ACP
PC1tab <- cbind(PC1,ACP_Data)
# Calcule du nombre de classe recommandé
1+(10/3)*log10(76)
# Délimitation des classes qui regroupent les résultats de l'ACP
Discretisation <- c(-7,-3,-2,-1,0,1,2,3,7)
# Discrétisation pour les valeurs négatives
Discretisation2 <- c(-7,-3,-2,-1)
# Discrétisation pour les valeurs positives
Discretisation3 <- c(0,1,2,3,7)
# Détermination du nombre de classes pour les valeurs positives
Nval <- length(Discretisation2)
# Détermination du nombre de classes pour les valeurs négatives
Nval2 <- length(Discretisation3)
choroLayer(spdf = Municipios_2015, spdfid ="ID",
df = PC1tab, dfid = "code_mun",
var = "PC1", breaks = Discretisation,
col= carto.pal(pal1 = "blue.pal",n1 = Nval,pal2 ="red.pal" ,n2 = Nval2, middle = FALSE, transparency = FALSE),
legend.title.txt = "Coordonnées PCA1",legend.pos = "topleft")
layoutLayer(title = "ACP sur les municipios de Mexico: coordonnées de la première composante",
sources = "INEGI - Encuesta Intercensal, 2015",postitle = "center", author = "",scale =0,north = TRUE)
labelLayer(spdf=Municipios_2015, df=PC1tab, spdfid ="ID", dfid ="code_mun", txt="code_mun", col = "black",cex = 0.5)
### CAH réalisée sur le niveau d'équipement et le niveau d'étude des populations des municipios ###
# Liste des variables retenues pour le niveau d'équipement:
#REFRIGERADOR -> Les foyers équipé d'un frigo
#CELULAR -> Les foyers équipé d'un téléphone portable
#TELEVISOR -> Les foyers équipé d'une télévision
#TELEVISOR_PP -> Les foyers équipé d'une télévision à écran plat
#RADIO -> Les foyers équipé d'une radio
#LAVADORA -> Les foyers équipé d'une machine à laver
#HORNO -> Les foyers équipé d'un four
#AUTOPROP -> Les foyers équipé d'une voiture
#INTERNET -> Les foyers équipé d'internet
#COMPUTADORA -> Les foyers équipé d'un ordinateur
#TELEPHONO -> Les foyers équipé d'un téléphone
# Liste des variables retenues pour le niveau d'étude:
#Preescolar -> école maternelle
#Primaria -> école primaire
#Secundaria -> collège
#Bachillerato -> baccalauréat général
#Bachillerato Techno -> baccalauréat technologique
#Licenciatura -> licence
#Especialidad -> dîplome spécialisé
#Maestría -> master
#Doctorado -> doctorat
# Liste des libraries R nécéssaires
library(cluster)
library(cartography)
library(reshape2)
library(ggplot2)
library(rgdal)
library(factoextra)
source("http://www.sthda.com/upload/rquery_cormat.r")
# Création du taux d'équipement des vivienda par municipios
# Import des données brutes sur les vivienda
tab_frequence_2015_viv <- read.csv("../data/traitees/tab_frequence_2015_viv.csv",header=TRUE, sep=",")
# Sélection des variables relatives au taux d'équipement des viviendas dont la réponse est SI
tab_Equipement_SI <- tab_frequence_2015_viv[,c(114,117,120,123,126,129,132,135,138,141,144)]
# Le tableau "tab_Equipement_SI" ne retient que les variables indiquant les ménages équipés.
# Calcule des corrélations entre variables relatives au taux d'équipement
TabCor <- as.data.frame(cor(tab_Equipement_SI))
# Extraction des coefficients de corréalation et des scores de p-value
TabCorAttrib<- rquery.cormat(TabCor, type="flatten", graph=FALSE)
# Représentation graphique des corrélations entre chaque variable étudiée
rquery.cormat(tab_Equipement_SI)
# Sélection de l'ensemble des variables relatives aux équipements
tab_Equipement <- tab_frequence_2015_viv[,c(114:144)]
# Somme de l'ensemble des viviendas par municipios
TotRep <- (rowSums(tab_Equipement))
# Calcul de la part relative des vivienda qui sont concernées par l'équipement étudié
# Calul de de la part des réponse positives (SI) par rapport à l'ensemble des réponses (SI,NO,99) pour tous les équipements
TxEquipement <- as.data.frame((rowSums(tab_Equipement_SI))/TotRep)
# Nom plus explicite donné à la colonne résultante du calcul
colnames(TxEquipement)[1]<-"TxEquipement"
# Import des données persona vivienda sous forme de pourcentage
tab_freq_pourc_pers_viv <- readRDS("../data/traitees/tab_freq_pourc_pers_viv.rds")
# Sélection des données relatives au niveau d'étude
Tab_NIVACAD <- tab_freq_pourc_pers_viv[,c(212:216,222:225)]
# Aggréagation des données d'équipement et des données d'éducation
Data <- cbind(Tab_NIVACAD,TxEquipement)
# Centrage-réduction des données
M <-as.data.frame(scale(Data ))
# Classification des variables par rapport à leur distance euclidenne avec la fonction agnes
ClusterM <- agnes(M,metric = "euclidean", method = "ward")
# Création du dendrogramme
dendroM <- as.dendrogram(ClusterM)
# Représentation graphique du dendogramme
plot(dendroM, type = "triangle", ylab = "Height", xlab = "Municipio")
# Découpage du dendogramme ou choix du nombre de classes issues de la CAH
Tab <- cutree(ClusterM,k=4)
# Aggrégation des données centrées réduites en fonction de leur appartenance à une classe
M$Cluster <- factor(Tab,levels = 1:4,labels = paste("Clus",1:4))
# Création d'un profil moyen pour les données centrées-réduites associées à chaque classe et chaque variable
ClusProfile<- aggregate(M[, 1:10],by=list(M$Cluster),mean)
# Nom plus explicites données aux colonnes
colnames(ClusProfile)[1] <-"CLUSTER"
colnames(ClusProfile)[2] <-"Preescolar"
colnames(ClusProfile)[3] <-"Primaria"
colnames(ClusProfile)[4] <-"Secundaria"
colnames(ClusProfile)[5] <-"Bachillerato"
colnames(ClusProfile)[6] <-"Bachillerato Techno"
colnames(ClusProfile)[7] <-"Licenciatura"
colnames(ClusProfile)[8] <-"Especialidad"
colnames(ClusProfile)[9] <-"Maestría"
colnames(ClusProfile)[10] <-"Doctorado"
# Mise en colonne de la moyenne des valeurs centrées-réduites pour chaque classes
clusLong <-melt(ClusProfile, id.vars = "CLUSTER")
# Représentation graphique des classes en fonction de leurs caractéristiques centrées-réduites
ggplot(clusLong) +
geom_bar(aes(x = variable, y = value, fill = CLUSTER ),stat = "identity") +
scale_fill_manual("Classes :", values = c("dodgerblue", "green3","darkviolet","darkorange")) +
facet_wrap (~ CLUSTER) +
coord_flip() + theme_dark()
# Création d'un tableau qui associe codes municipios et résultats de la CAH
DataMap <- data.frame(tab_freq_pourc_pers_viv$code_mun, row.names(M), M$Cluster)
colnames(DataMap)[1] <-"Code_Mun"
# Import de la géométrie des municipios sous forme de shapefile
Municipios_2015 <- readOGR(dsn= "../data/areas_geoestadisticas_municipales.shp", layer = "areas_geoestadisticas_municipales",encoding= "utf8", stringsAsFactors = FALSE)
# Choix d'une projection
projection <- "+proj=utm +zone=14 +datum=WGS84 +units=m +no_defs + ellps=WGS84"
# Association des données géométriques à la projection WGS84
Municipios_2015 <- spTransform(Municipios_2015, projection)
# Création d'un identifiant code municipios associé au découpage géométrique
Municipios_2015$ID <- paste(Municipios_2015$CVE_ENT, Municipios_2015$CVE_MUN, sep="")
# Cartographie de la CAH
typoLayer(spdf = Municipios_2015,
df = DataMap,
var = "M.Cluster",
spdfid = "ID",
dfid = "Code_Mun",
col = c("dodgerblue", "green3","darkviolet","darkorange"),
legend.title.txt = "Cluster",
legend.frame = FALSE,
legend.pos = "left",
border = "grey20")
layoutLayer(
title = "CAH Municipios de Mexico: Niveau d'étude et taux d'équipement",
postitle = "center",
sources = "INEGI,encuesta intercensal 2015",
author = "",
scale = 0,
north = TRUE)
labelLayer(spdf=Municipios_2015, df=DataMap, spdfid ="ID", dfid ="Code_Mun", txt="Code_Mun", col = "black",cex = 0.5)
# Représentation graphique du dendogramme coloré en fonction des classes choisies
fviz_dend(dendroM, k = 4, # Cut in four grcex = 0.5, # label size
k_colors = c("#0092ff", "#00ff0d","#8c1aff","#F76A2E"),
color_labels_by_k = TRUE, # color labels by groups
rect = TRUE, # Add rectangle around groups
rect_border = c("#0092ff", "#00ff0d","#8c1aff","#F76A2E"),
rect_fill = TRUE,main = "",xlab = "Communes",ylab="Inertie")
This diff is collapsed.
container_name: liftr_container_a11015af_67df_477e_ac5b_a9af83e57f15
image_name: document_test_ekaterina
docker_build_cmd: docker build --no-cache=false --rm=true -t="document_test_ekaterina"
/home/ekaterina/m2_G2M/rapport/docs
docker_run_cmd: docker run --rm --name "liftr_container_a11015af_67df_477e_ac5b_a9af83e57f15"
-u `id -u $USER` -v "/home/ekaterina/m2_G2M/rapport/docs:/liftrroot/" document_test_ekaterina
Rscript -e "library('knitr');library('rmarkdown');library('shiny');setwd('/liftrroot/');render(input
= 'document_test_ekaterina.Rmd')"
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment