#SOLUCION CLAVE A

Desarrolle el siguiente ejercicio: Se necesita construir un indicador multivariado sintético, que mida la “Seguridad Municipal” Para ello se dispone de la siguiente información: Correlación con la variable compleja

##1. (25%) A través del análisis de componentes principales, identifique para un modelo de 3 factores: a. Los ponderadores normalizados para cada factor.##**

library(magrittr)
data_parcial_2 <-
load("C:/Users/MINEDUCYT/Desktop/MAE2/data_parcial_2_A_rev.RData")

head(datos_parcial_2)
##   ID  Municipio X1 X2       X3       X4       X5       X6 X7       X8
## 1  1 ATIQUIZAYA  9  2 20.00000 20.00000  0.00000 0.000000  2  56.4000
## 2  2  EL CARMEN 10  6 62.50000 50.00000 37.50000 3.947368 11 147.3750
## 3  3    ALEGRIA 10 20 50.00000 50.00000 50.00000 2.564103 16 135.0000
## 4  4 SAN JULIAN  8  3 42.85714 42.85714 14.28571 1.351351 35 121.1429
## 5  5    TEJUTLA  7  7 75.00000 75.00000 75.00000 9.090909  8 202.5000
## 6  6  PASAQUINA  6 13 30.00000 30.00000 30.00000 8.108108 25  81.0000
library(dplyr)
library(magrittr)
library(dplyr)

norm_directa <- function(x) {
(x - min(x)) / (max(x) - min(x))
}
norm_inversa <- function(x) {
(max(x) - x) / (max(x) - min(x))
}

#Seleccionando las variables con correlación positiva

datos_parcial_2 %>%
select(X1, X2, X3, X5, X7, X8) %>%
apply(MARGIN = 2, FUN = norm_directa) %>% as.data.frame() -> variables_corr_positiva

#Seleccionando las variables con correlación negativa

datos_parcial_2 %>%
select(X4, X6) %>%
apply(MARGIN = 2, FUN = norm_inversa) %>% as.data.frame() -> variables_corr_negativa

#Juntando y reordenando las variables
variables_corr_positiva %>%
bind_cols(variables_corr_negativa) %>%
select(X1, X2, X3, X4, X5, X6, X7, X8) -> datos_financieros_normalizados

head(datos_financieros_normalizados)
##           X1          X2        X3        X4        X5        X6         X7
## 1 0.19354839 0.000000000 0.0400000 0.8000000 0.0000000 1.0000000 0.00000000
## 2 0.22580645 0.017167382 0.5500000 0.5000000 0.4285714 0.7844130 0.09890110
## 3 0.22580645 0.077253219 0.4000000 0.5000000 0.5714286 0.8599606 0.15384615
## 4 0.16129032 0.004291845 0.3142857 0.5714286 0.1632653 0.9261954 0.36263736
## 5 0.12903226 0.021459227 0.7000000 0.2500000 0.8571429 0.5034965 0.06593407
## 6 0.09677419 0.047210300 0.1600000 0.7000000 0.3428571 0.5571726 0.25274725
##          X8
## 1 0.1582266
## 2 0.5167488
## 3 0.4679803
## 4 0.4133709
## 5 0.7339901
## 6 0.2551724

##Ahora se verificarán las pruebas

Matriz de Correlación & Pruebas de Barlett y KMO para saber si se puede usar análisis factorial

#Matriz de correlación
library(PerformanceAnalytics)
chart.Correlation(as.matrix(datos_financieros_normalizados),histogram = TRUE,pch=12)

#KMO
library(rela)
KMO<-paf(as.matrix(datos_financieros_normalizados))$KMO
print(KMO)
## [1] 0.67931

Se puede realizar el analisis factorial por ser KMO> 0.5

#Prueba de Barlett
library(psych)
options(scipen = 99999)
Barlett<-cortest.bartlett(datos_financieros_normalizados)
print(Barlett)
## $chisq
## [1] 1025.9
## 
## $p.value
## [1] 0.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000046951
## 
## $df
## [1] 28

El p-value es menor por lo tanto existe correlacion. Podemos iniciar el analisis factorial

##Analisi factorial

library(FactoMineR)
library(factoextra)
library(kableExtra)

Rx<-cor(datos_financieros_normalizados)
PC<-princomp(x = datos_financieros_normalizados,cor = TRUE,fix_sign = FALSE)
variables_pca<-get_pca_var(PC)
factoextra::get_eig(PC) %>% kable(caption="Resumen de PCA",
        align = "c",
        digits = 2) %>% 
  kable_material(html_font = "sans-serif") %>% 
  kable_styling(bootstrap_options = c("hover"))
Resumen de PCA
eigenvalue variance.percent cumulative.variance.percent
Dim.1 3.90 48.72 48.72
Dim.2 1.96 24.55 73.27
Dim.3 0.84 10.52 83.78
Dim.4 0.50 6.24 90.03
Dim.5 0.45 5.68 95.70
Dim.6 0.28 3.45 99.16
Dim.7 0.07 0.82 99.98
Dim.8 0.00 0.02 100.00 
fviz_eig(PC,
         choice = "eigenvalue",
         barcolor = "red",
         barfill = "red",
         addlabels = TRUE, 
       )+labs(title = "Gráfico de Sedimentación",subtitle = "Usando princomp, con Autovalores")+
  xlab(label = "Componentes")+
  ylab(label = "Autovalores")+geom_hline(yintercept = 1)

Con la solucion mostrada se retienen 3 factores

library(corrplot)
#Modelo de 3 Factores (Rotada)
numero_de_factores<-3
modelo_factores<-principal(r = Rx,
                             nfactors = numero_de_factores,
                             covar = FALSE,
                             rotate = "varimax")
modelo_factores
## Principal Components Analysis
## Call: principal(r = Rx, nfactors = numero_de_factores, rotate = "varimax", 
##     covar = FALSE)
## Standardized loadings (pattern matrix) based upon correlation matrix
##      RC1   RC2   RC3   h2     u2 com
## X1 -0.16  0.80 -0.03 0.67 0.3316 1.1
## X2  0.08  0.84 -0.03 0.71 0.2879 1.0
## X3  0.93 -0.09  0.28 0.95 0.0493 1.2
## X4 -0.95  0.05 -0.26 0.98 0.0208 1.2
## X5  0.43 -0.06  0.80 0.83 0.1742 1.5
## X6 -0.25  0.03 -0.91 0.89 0.1142 1.2
## X7 -0.07  0.83 -0.04 0.69 0.3107 1.0
## X8  0.96 -0.06  0.27 0.99 0.0087 1.2
## 
##                        RC1  RC2  RC3
## SS loadings           2.97 2.05 1.68
## Proportion Var        0.37 0.26 0.21
## Cumulative Var        0.37 0.63 0.84
## Proportion Explained  0.44 0.31 0.25
## Cumulative Proportion 0.44 0.75 1.00
## 
## Mean item complexity =  1.2
## Test of the hypothesis that 3 components are sufficient.
## 
## The root mean square of the residuals (RMSR) is  0.06 
## 
## Fit based upon off diagonal values = 0.98

Se verifican las correlaciones para observar las dimensiones en cada factor

correlaciones_modelo<-variables_pca$coord
rotacion<-varimax(correlaciones_modelo[,1:numero_de_factores])
correlaciones_modelo_rotada<-rotacion$loadings

corrplot(correlaciones_modelo_rotada[,1:numero_de_factores],
         is.corr = FALSE,
         method = "square",
         addCoef.col="grey",
         number.cex = 0.75)

Respuesta

library(kableExtra)
cargas<-rotacion$loadings[1:6,1:numero_de_factores]
ponderadores<-prop.table(apply(cargas^2,MARGIN = 2,sum))
t(ponderadores) %>% kable(caption="Ponderadores de los Factores Extraídos",
        align = "c",
        digits = 2) %>% 
  kable_material(html_font = "sans-serif") %>% 
  kable_styling(bootstrap_options = c("striped", "hover"))
Ponderadores de los Factores Extraídos
Dim.1 Dim.2 Dim.3
0.41 0.27 0.32

Para verificar la contribucion de cada variable a cada factor

contribuciones<-apply(cargas^2,MARGIN = 2,prop.table)
contribuciones %>% kable(caption="Contribución de las variables en los Factores",
        align = "c",
        digits = 2) %>% 
  kable_material(html_font = "sans-serif") %>% 
  kable_styling(bootstrap_options = c("striped", "hover"))
Contribución de las variables en los Factores
Dim.1 Dim.2 Dim.3
X1 0.01 0.47 0.00
X2 0.00 0.52 0.00
X3 0.42 0.01 0.05
X4 0.44 0.00 0.04
X5 0.09 0.00 0.40
X6 0.03 0.00 0.51

##b. Las variables incluidas en cada factor. Luego de realizar el analisis se puede concluir que las variables quedarian de la siguiente forma

COMPONENTE 1: X3, X4, X8 COMPONENTE 2: X1, X2, X7 COMPONENTE 3: X5, X6

##2. (25%) Para el factor 1, utilice el método CRITIC para obtener los ponderadores normalizados para cada variable.##

COMPONENTE 1: X3, X4, X8

# Normalización de los datos
library(dplyr)

datos_parcial_2 %>%
  select(X3, X8) %>%
  apply(MARGIN = 2, FUN = norm_directa) %>% as.data.frame() -> data_factor

datos_parcial_2 %>%
  select(X4) %>%
  apply(MARGIN = 2, FUN = norm_inversa) %>% as.data.frame() -> data_factorr

data_factor %>% #peguele el otro dataframe
  bind_cols(data_factorr) %>%
  select(X3, X4, X8) -> data_factor_1

head(data_factor_1)
##        X3      X4      X8
## 1 0.04000 0.80000 0.15823
## 2 0.55000 0.50000 0.51675
## 3 0.40000 0.50000 0.46798
## 4 0.31429 0.57143 0.41337
## 5 0.70000 0.25000 0.73399
## 6 0.16000 0.70000 0.25517
#Cálculo de las desviaciones estándar de cada variable

data_factor_1 %>% dplyr::summarise(S3=sd(X3),S4=sd(X4),S8=sd(X8))-> sd_vector
print(sd_vector)
##        S3      S4      S8
## 1 0.24628 0.20111 0.20874
#Cálculo de la matriz de correlación

cor(data_factor_1)->mat_R_F1
print(mat_R_F1)
##          X3       X4       X8
## X3  1.00000 -0.93872  0.95904
## X4 -0.93872  1.00000 -0.99585
## X8  0.95904 -0.99585  1.00000
#Cálculo de los ponderadores brutos
1-mat_R_F1->sum_data
colSums(sum_data)->sum_vector
sd_vector*sum_vector->vj
print(vj)
##        S3      S4      S8
## 1 0.48755 0.79129 0.42517
#Cálculo de los ponderadores netos

vj/sum(vj)->wj
print(wj)
##        S3      S4      S8
## 1 0.28612 0.46437 0.24951
#Ponderadores:
print(round(wj*100,2))
##      S3    S4    S8
## 1 28.61 46.44 24.95

##3. (25%) Para el factor 2, utilice el método de Entropía para obtener los ponderadores normalizados para cada variable.

COMPONENTE 2: X1, X2, X7

datos_parcial_2 %>% dplyr::select(X1,X2,X7)->data_norm
apply(data_norm,2,prop.table)->data_norm
head(data_norm)
##             X1         X2        X7
## [1,] 0.0078125 0.00070734 0.0013986
## [2,] 0.0086806 0.00212202 0.0076923
## [3,] 0.0086806 0.00707339 0.0111888
## [4,] 0.0069444 0.00106101 0.0244755
## [5,] 0.0060764 0.00247569 0.0055944
## [6,] 0.0052083 0.00459770 0.0174825
#Fórmula de entropía
entropy<-function(x){
  return(x*log(x)) 
}
apply(data_norm,2,entropy)->data_norm_2
head(data_norm_2)
##             X1         X2        X7
## [1,] -0.037906 -0.0051310 -0.009192
## [2,] -0.041204 -0.0130618 -0.037443
## [3,] -0.041204 -0.0350233 -0.050270
## [4,] -0.034513 -0.0072664 -0.090806
## [5,] -0.031010 -0.0148572 -0.029013
## [6,] -0.027383 -0.0247457 -0.070744
#Número de variables en el factor:
ncol(data_norm)->m
#Constante de entropía:
-1/log(m)->K
print(K)
## [1] -0.91024
#Cálculo de las entropías
K*colSums(data_norm_2)->Ej
print(Ej)
##     X1     X2     X7 
## 4.1805 3.2029 3.7019
#Cálculo de las especificidades:
1-Ej->vj
print(vj)
##      X1      X2      X7 
## -3.1805 -2.2029 -2.7019
#Cálculo de los ponderadores:
prop.table(vj)->wj #es igual a usar vj/sum(vj)
print(wj)
##      X1      X2      X7 
## 0.39337 0.27245 0.33417

##4. (25%) Para el factor 3, utilice el método de Ranking para obtener los ponderadores normalizados para cada variable (utilice la numeración de las variables para establecer la jerarquía).

COMPONENTE 3: En primer lugar estaria X5 y el segundo lugar seria para la variable X6

##Rango de suma

library(magrittr)
#Vector de Jerarquías
rj<-c(1:2)
names(rj) <- c ("X5","X6")

#Función para generar los pesos
ponderadores_subjetivos_rank_suma<-function(vector_jerarquias){
  n<-length(vector_jerarquias)
  vector_pesos<-n-vector_jerarquias+1
  list(w_brutos=vector_pesos,w_normalizados=vector_pesos/sum(vector_pesos))
}
#Aplicando la función:
pesos_ranking_suma<-ponderadores_subjetivos_rank_suma(rj)

#Pesos brutos
pesos_ranking_suma$w_brutos
## X5 X6 
##  2  1
#Pesos normalizados
pesos_ranking_suma$w_normalizados %>% round(digits = 3)
##    X5    X6 
## 0.667 0.333
#Gráfico de los pesos normalizados
barplot(pesos_ranking_suma$w_normalizados,
        main = "Ponderadores Ranking de Suma",
        ylim = c(0,0.8),col = "orange")

##Jerarquia reciproca##

library(magrittr)
#Vector de Jerarquías
rj<-c(1:2)
names(rj)<-c("X5", "X6")

#Función para generar los pesos
ponderadores_subjetivos_rank_reciproco<-function(vector_jerarquias){
  vector_pesos<-1/vector_jerarquias
  list(w_brutos=vector_pesos,w_normalizados=vector_pesos/sum(vector_pesos))
}
#Aplicando la función:
pesos_ranking_reciproco<-ponderadores_subjetivos_rank_reciproco(rj)

#Pesos brutos
pesos_ranking_reciproco$w_brutos
##  X5  X6 
## 1.0 0.5
#Pesos normalizados
pesos_ranking_reciproco$w_normalizados %>% round(digits = 3)
##    X5    X6 
## 0.667 0.333
#Gráfico de los pesos normalizados
barplot(pesos_ranking_reciproco$w_normalizados,
        main = "Ponderadores Ranking Recíproco",
        ylim = c(0,0.8),col = "orange")

##Jerarquia exponencial

library(magrittr)
#Vector de Jerarquías
rj<-c(1:2)
names(rj)<-c("X5","X6")

#Función para generar los pesos
ponderadores_subjetivos_rank_exponencial<-function(vector_jerarquias,p=2){
  n<-length(vector_jerarquias)
  vector_pesos<-(n-vector_jerarquias+1)^p
  list(w_brutos=vector_pesos,w_normalizados=vector_pesos/sum(vector_pesos))
}
#Aplicando la función:
pesos_ranking_exponencial<-ponderadores_subjetivos_rank_exponencial(rj)

#Pesos brutos
pesos_ranking_exponencial$w_brutos
## X5 X6 
##  4  1
#Pesos normalizados
pesos_ranking_exponencial$w_normalizados %>% round(digits = 3)
##  X5  X6 
## 0.8 0.2
#Gráfico de los pesos normalizados (por default p=2)
barplot(pesos_ranking_suma$w_normalizados,
        main = "Ponderadores Ranking Exponencial",
        ylim = c(0,0.9),col = "orange")

#Comparación de valores de "p"

par(mfrow=c(1,3)) 
for(p in 2:4){  
  
  pesos<-ponderadores_subjetivos_rank_exponencial(vector_jerarquias = rj,p = p)
  barplot(pesos$w_normalizados,main = paste0("p=",p),ylim = c(0,1.5),col = "orange",cex.main=3,cex.axis = 3)
}