Escalamiento Multidimensional
Cargo la data:
folder= "data"
fileName= "idePeru2012.csv"
fileToRead=file.path(folder,fileName)
ide12=read.csv(fileToRead) #recuenden que csv no necesita paquete
head(ide12) #vemos lo que hemos cargado## ubiProv ubiReg DEPARTAMENTO PROVINCIA IDE IDENTIDAD
## 1 10100 10000 AMAZONAS Chachapoyas 0.7737 98.6179
## 2 10200 10000 AMAZONAS Bagua 0.6623 94.6079
## 3 10300 10000 AMAZONAS Bongara 0.6318 97.4681
## 4 10400 10000 AMAZONAS Condorcanqui 0.4598 86.2320
## 5 10500 10000 AMAZONAS Luya 0.6047 96.1927
## 6 10600 10000 AMAZONAS Rodriguez de Mendoza 0.6312 97.3431
## SALUD EDUCACION SANEAMIENTO ELECTRIFICACION POBLACION COSTA CAPITAL
## 1 25.4500 91.4986 70.3454 83.9712 54783 NO SI
## 2 14.6091 79.7902 64.4790 67.9146 77438 NO NO
## 3 9.0102 76.4240 54.8341 72.1693 32317 NO NO
## 4 8.5570 52.2149 37.7145 39.4891 51802 NO NO
## 5 12.4180 74.7260 43.3484 67.3961 52185 NO NO
## 6 14.8787 79.4244 46.5018 67.5461 30236 NO NO
## TAMANO
## 1 2
## 2 2
## 3 1
## 4 2
## 5 2
## 6 1Veo estructura:
str(ide12)## 'data.frame': 195 obs. of 14 variables:
## $ ubiProv : int 10100 10200 10300 10400 10500 10600 10700 20100 20200 20300 ...
## $ ubiReg : int 10000 10000 10000 10000 10000 10000 10000 20000 20000 20000 ...
## $ DEPARTAMENTO : Factor w/ 25 levels "AMAZONAS","ANCASH",..: 1 1 1 1 1 1 1 2 2 2 ...
## $ PROVINCIA : Factor w/ 195 levels "Abancay","Acobamba",..: 46 19 24 59 116 156 186 90 4 11 ...
## $ IDE : num 0.774 0.662 0.632 0.46 0.605 ...
## $ IDENTIDAD : num 98.6 94.6 97.5 86.2 96.2 ...
## $ SALUD : num 25.45 14.61 9.01 8.56 12.42 ...
## $ EDUCACION : num 91.5 79.8 76.4 52.2 74.7 ...
## $ SANEAMIENTO : num 70.3 64.5 54.8 37.7 43.3 ...
## $ ELECTRIFICACION: num 84 67.9 72.2 39.5 67.4 ...
## $ POBLACION : int 54783 77438 32317 51802 52185 30236 118747 161003 7974 16879 ...
## $ COSTA : Factor w/ 2 levels "NO","SI": 1 1 1 1 1 1 1 1 1 1 ...
## $ CAPITAL : Factor w/ 2 levels "NO","SI": 2 1 1 1 1 1 1 2 1 1 ...
## $ TAMANO : int 2 2 1 2 2 1 3 3 1 1 ...Hagmos un truquito para volver a cargar la data:
ide12=read.csv(fileToRead,
colClasses=c(rep("character",4), #las primeras 4 columnas serán chr
rep("numeric",7), # las siguientes 7 columnas serán num
rep("factor",3))) # las siguientes 3 columnas serán factAhora saco estructura nuevamente:
str(ide12)## 'data.frame': 195 obs. of 14 variables:
## $ ubiProv : chr "10100" "10200" "10300" "10400" ...
## $ ubiReg : chr "10000" "10000" "10000" "10000" ...
## $ DEPARTAMENTO : chr "AMAZONAS" "AMAZONAS" "AMAZONAS" "AMAZONAS" ...
## $ PROVINCIA : chr "Chachapoyas" "Bagua" "Bongara" "Condorcanqui" ...
## $ IDE : num 0.774 0.662 0.632 0.46 0.605 ...
## $ IDENTIDAD : num 98.6 94.6 97.5 86.2 96.2 ...
## $ SALUD : num 25.45 14.61 9.01 8.56 12.42 ...
## $ EDUCACION : num 91.5 79.8 76.4 52.2 74.7 ...
## $ SANEAMIENTO : num 70.3 64.5 54.8 37.7 43.3 ...
## $ ELECTRIFICACION: num 84 67.9 72.2 39.5 67.4 ...
## $ POBLACION : num 54783 77438 32317 51802 52185 ...
## $ COSTA : Factor w/ 2 levels "NO","SI": 1 1 1 1 1 1 1 1 1 1 ...
## $ CAPITAL : Factor w/ 2 levels "NO","SI": 2 1 1 1 1 1 1 2 1 1 ...
## $ TAMANO : Factor w/ 5 levels "1","2","3","4",..: 2 2 1 2 2 1 3 3 1 1 ...TAMANO debe ser ordinal:
ide12$TAMANO=as.ordered(ide12$TAMANO)
str(ide12$TAMANO)## Ord.factor w/ 5 levels "1"<"2"<"3"<"4"<..: 2 2 1 2 2 1 3 3 1 1 ...Mapa de similitudes
- Subseteo
- Índice
provMap=ide12[,c(4,6:10)] #ojo: incluyo el nombre
row.names(provMap) = provMap$PROVINCIA ##como sabebemos que el nombre de las provincias no se repite lo podemos usar
provMap = provMap[,-c(1)] #eliminamos la variable del nombre del país, ya no la necesitamos
# Viendo resultado:
head(provMap)## IDENTIDAD SALUD EDUCACION SANEAMIENTO
## Chachapoyas 98.6179 25.4500 91.4986 70.3454
## Bagua 94.6079 14.6091 79.7902 64.4790
## Bongara 97.4681 9.0102 76.4240 54.8341
## Condorcanqui 86.2320 8.5570 52.2149 37.7145
## Luya 96.1927 12.4180 74.7260 43.3484
## Rodriguez de Mendoza 97.3431 14.8787 79.4244 46.5018
## ELECTRIFICACION
## Chachapoyas 83.9712
## Bagua 67.9146
## Bongara 72.1693
## Condorcanqui 39.4891
## Luya 67.3961
## Rodriguez de Mendoza 67.5461- Estandarizo con scale
provMap_s=scale(provMap)
summary(provMap_s)## IDENTIDAD SALUD EDUCACION SANEAMIENTO
## Min. :-6.2335 Min. :-1.2740 Min. :-3.2325 Min. :-2.72782
## 1st Qu.:-0.3082 1st Qu.:-0.7427 1st Qu.:-0.5650 1st Qu.:-0.74504
## Median : 0.2459 Median :-0.2104 Median : 0.1834 Median :-0.09938
## Mean : 0.0000 Mean : 0.0000 Mean : 0.0000 Mean : 0.00000
## 3rd Qu.: 0.6387 3rd Qu.: 0.3862 3rd Qu.: 0.7638 3rd Qu.: 0.78865
## Max. : 0.8962 Max. : 4.4865 Max. : 1.6948 Max. : 2.10498
## ELECTRIFICACION
## Min. :-2.4433
## 1st Qu.:-0.6892
## Median : 0.1391
## Mean : 0.0000
## 3rd Qu.: 0.8236
## Max. : 1.7480- Distancias
provMap_d=dist(provMap_s)- Algoritmo de escalamiento multidimensional con cmdscale:
provMap_r <- cmdscale(provMap_d,eig=TRUE, k=2) # k sugiere dimensiones. El resultado muestra las coordenadas (puntos) de cada provincia. DEJA EL k en 2.
provMap_r$GOF # GOF indica la bondad de ajuste, mientras mas cerca a 1 mejor.## [1] 0.7827708 0.7827708- Mapa de similitudes
titulo="Mapa de Similitudes entre provincias del Peru basado en el IDE"
x <- provMap_r$points[,1]
y <- provMap_r$points[,2]
plot(x, y, main=titulo)
Personalicemos: quiero etiquetas en vez de puntos.
plot(x, y, xlab="Dimension 1", ylab="Dimension 2", main=titulo, type="n") # 'n' evita que se pongan los puntos.
columnForLabels=dimnames(provMap_r[[1]])[[1]] # etiquetas y colores de los puntos
text(x, y,labels = columnForLabels , cex = 0.5) #con cex indicamos el tamaño
Coloreemos las que tienen costa:
plot(x, y, xlab="Dimensión 1", ylab="Dimensión 2", main=titulo, type="n")
# etiquetas y colores de los puntos
columnForLabels=ide12$PROVINCIA
colorForLabels=ide12$COSTA
paleta=c('grey','blue')
text(x, y,labels = columnForLabels, cex = 0.8, col = paleta[colorForLabels])
legend("bottomright", legend = levels(colorForLabels),
fill = paleta,title = "¿Tiene Costa?")
Podemos agregar los datos por región [para que me de la media por departamento]. Con aggregate.
regionIde12=aggregate(cbind(IDENTIDAD,SALUD, EDUCACION,SANEAMIENTO, ELECTRIFICACION) ~ DEPARTAMENTO,data=ide12,FUN=mean)
regionIde12## DEPARTAMENTO IDENTIDAD SALUD EDUCACION SANEAMIENTO ELECTRIFICACION
## 1 AMAZONAS 95.09089 13.576171 75.89234 52.82039 65.94344
## 2 ANCASH 97.20594 11.680745 79.25444 60.00530 75.93114
## 3 APURIMAC 97.02721 8.774300 85.53033 44.30631 72.43510
## 4 AREQUIPA 98.37373 19.257137 85.06124 75.58691 76.48626
## 5 AYACUCHO 99.07324 11.014827 79.90297 44.76749 65.01315
## 6 CAJAMARCA 96.35778 7.549954 77.89242 50.17891 60.18226
## 7 CALLAO 99.49390 26.421800 89.84930 93.04320 99.25290
## 8 CUSCO 97.47522 9.184100 82.17028 51.11976 72.81752
## 9 HUANCAVELICA 98.12110 4.885286 83.45017 29.74453 73.14633
## 10 HUANUCO 95.94294 12.075282 73.67755 40.03158 52.34482
## 11 ICA 99.39004 19.419320 86.02430 84.64438 84.83418
## 12 JUNIN 98.16204 10.321389 87.75080 58.91788 83.14287
## 13 LA LIBERTAD 97.85048 13.855675 73.08089 66.01242 65.89008
## 14 LAMBAYEQUE 97.74597 16.073033 78.45127 78.47643 82.07550
## 15 LIMA 98.96001 20.136650 82.65205 74.67121 83.23084
## 16 LORETO 90.90169 6.885843 62.29814 35.73760 67.34740
## 17 MADRE DE DIOS 97.41320 21.980467 83.19697 59.25963 77.34077
## 18 MOQUEGUA 99.26907 20.488667 90.63790 72.91303 89.43910
## 19 PASCO 97.92057 9.308200 84.91410 47.75177 73.63680
## 20 PIURA 97.35669 9.357750 81.29290 70.67819 72.86267
## 21 PUNO 98.76785 6.629938 86.06802 47.09789 72.34415
## 22 SAN MARTIN 96.09507 10.348420 82.01861 67.37171 72.38856
## 23 TACNA 99.45212 26.117750 91.89978 72.17215 87.01730
## 24 TUMBES 98.59367 10.660467 82.97967 69.49747 83.32957
## 25 UCAYALI 92.78928 9.373175 66.37073 17.39293 62.96545Ahora hago los mismos pasos para departamentos:
regMap=regionIde12 #creo una copia
row.names(regMap)=regionIde12$DEPARTAMENTO #indice
regMap = regMap[,-c(1)]#elimino columna de departamento
regMap_s=scale(regMap) #estandarizamos
regMap_d=dist(regMap_s) #distancias
regMap_r = cmdscale(regMap_d,eig=TRUE, k=2) #calculamos coordenadas (puntos)¿Qué tan bueno fue nuestro resultado?
regMap_r$GOF # mientras mas cerca a 1 mejor.## [1] 0.8677184 0.8677184Mapa:
labelsMap=row.names(regMap)
x <- regMap_r$points[,1]
y <- regMap_r$points[,2]
plot(x, y, main="Mapa IED nivel regiones", type="n")
text(x, y, labels = labelsMap, cex=.5) 
Mapa de modalidades
Creo nuevas columnas donde pondré mis variables categorizadas
ide12$elcOrd=ide12$sanOrd=ide12$eduOrd=ide12$salOrd=ide12$idOrd=NAnames(ide12)## [1] "ubiProv" "ubiReg" "DEPARTAMENTO"
## [4] "PROVINCIA" "IDE" "IDENTIDAD"
## [7] "SALUD" "EDUCACION" "SANEAMIENTO"
## [10] "ELECTRIFICACION" "POBLACION" "COSTA"
## [13] "CAPITAL" "TAMANO" "idOrd"
## [16] "salOrd" "eduOrd" "sanOrd"
## [19] "elcOrd"Doy contenido a mis columnas vacías:
gruposCantidad=3
###
etiquetas=c('bajo','medio','alto')
ide12[,c(15:19)]=lapply(ide12[,c(6:10)],cut,
breaks = gruposCantidad,
labels = etiquetas,
ordered_result = T)str(ide12)## 'data.frame': 195 obs. of 19 variables:
## $ ubiProv : chr "10100" "10200" "10300" "10400" ...
## $ ubiReg : chr "10000" "10000" "10000" "10000" ...
## $ DEPARTAMENTO : chr "AMAZONAS" "AMAZONAS" "AMAZONAS" "AMAZONAS" ...
## $ PROVINCIA : chr "Chachapoyas" "Bagua" "Bongara" "Condorcanqui" ...
## $ IDE : num 0.774 0.662 0.632 0.46 0.605 ...
## $ IDENTIDAD : num 98.6 94.6 97.5 86.2 96.2 ...
## $ SALUD : num 25.45 14.61 9.01 8.56 12.42 ...
## $ EDUCACION : num 91.5 79.8 76.4 52.2 74.7 ...
## $ SANEAMIENTO : num 70.3 64.5 54.8 37.7 43.3 ...
## $ ELECTRIFICACION: num 84 67.9 72.2 39.5 67.4 ...
## $ POBLACION : num 54783 77438 32317 51802 52185 ...
## $ COSTA : Factor w/ 2 levels "NO","SI": 1 1 1 1 1 1 1 1 1 1 ...
## $ CAPITAL : Factor w/ 2 levels "NO","SI": 2 1 1 1 1 1 1 2 1 1 ...
## $ TAMANO : Ord.factor w/ 5 levels "1"<"2"<"3"<"4"<..: 2 2 1 2 2 1 3 3 1 1 ...
## $ idOrd : Ord.factor w/ 3 levels "bajo"<"medio"<..: 3 3 3 1 3 3 3 3 3 3 ...
## $ salOrd : Ord.factor w/ 3 levels "bajo"<"medio"<..: 2 1 1 1 1 1 1 2 1 1 ...
## $ eduOrd : Ord.factor w/ 3 levels "bajo"<"medio"<..: 3 2 2 1 2 2 2 3 3 2 ...
## $ sanOrd : Ord.factor w/ 3 levels "bajo"<"medio"<..: 3 2 2 2 2 2 2 3 2 2 ...
## $ elcOrd : Ord.factor w/ 3 levels "bajo"<"medio"<..: 3 2 2 1 2 2 2 3 3 2 ...Summary de ordinales:
summary(ide12[,c(15:19)])## idOrd salOrd eduOrd sanOrd elcOrd
## bajo : 3 bajo :156 bajo : 18 bajo : 25 bajo :32
## medio: 7 medio: 35 medio: 67 medio:108 medio:89
## alto :185 alto : 4 alto :110 alto : 62 alto :74Sacamos tabla cruzada de tamaño y educación
tablaTE=table(ide12$TAMANO,ide12$eduOrd)
prop.table(tablaTE)##
## bajo medio alto
## 1 0.020512821 0.138461538 0.194871795
## 2 0.056410256 0.128205128 0.143589744
## 3 0.015384615 0.071794872 0.148717949
## 4 0.000000000 0.005128205 0.071794872
## 5 0.000000000 0.000000000 0.005128205Saquemos chi-cuadrado para ver si hay asociación:
chisq.test(tablaTE)## Warning in chisq.test(tablaTE): Chi-squared approximation may be incorrect##
## Pearson's Chi-squared test
##
## data: tablaTE
## X-squared = 18.042, df = 8, p-value = 0.02092Hagamos chi no paramétrico:
chisq.test(tablaTE,simulate.p.value = T)##
## Pearson's Chi-squared test with simulated p-value (based on 2000
## replicates)
##
## data: tablaTE
## X-squared = 18.042, df = NA, p-value = 0.02799Hay relación. Aún no sabemos cómo es.
Hagamos análisis de correspondencias:
Calculemos tabla de correspondencias:
library(ca) #de ser necesario, instalen el paquete!## Warning: package 'ca' was built under R version 3.4.4tablaCA_te=ca(tablaTE)tablaCA_te##
## Principal inertias (eigenvalues):
## 1 2
## Value 0.078998 0.013523
## Percentage 85.38% 14.62%
##
##
## Rows:
## 1 2 3 4 5
## Mass 0.353846 0.328205 0.235897 0.076923 0.005128
## ChiDist 0.140414 0.321618 0.142240 0.746443 0.879049
## Inertia 0.006976 0.033949 0.004773 0.042860 0.003963
## Dim. 1 -0.063547 1.081488 -0.505903 -2.572733 -2.967917
## Dim. 2 1.197668 -0.903582 0.031782 -1.592507 -2.384260
##
##
## Columns:
## bajo medio alto
## Mass 0.092308 0.343590 0.564103
## ChiDist 0.628641 0.266736 0.236666
## Inertia 0.036479 0.024446 0.031596
## Dim. 1 2.001206 0.831911 -0.834179
## Dim. 2 -2.414230 1.103799 -0.277258plot.ca(tablaCA_te, col=c("red","blue"))
Ángulos:
plot.ca(tablaCA_te, col=c("red","blue"), arrows = c(T,T))