Dans ce TP, toutes les comparaison entre Arbre et AFD/ACP seront à la fin (cas DATA TEST et cas DATA REDUITES par ACP notamment)

Fonction AFD FROM SCRATCH

AFD<-function(X,y){

n<-nrow(X)
p <- ncol(X)
K<-length(levels(factor(y)))
nk<-as.vector(table(y))
names(nk) <- levels(y)
pk<-nk/n

Xk <-list()
Xk <- split(X,y)

V<-var(X)*(n-1)/n
R<-cor(X)

g<-matrix(apply(X,2,mean),ncol=1)
rownames(g) <- colnames(X)

Vk<-list()
Rk<-list()
gk<-list()
for (k in 1:K){
    gk[[k]]<-matrix(apply(Xk[[k]],2,mean),ncol=1)
    rownames(gk[[k]]) <- colnames(X)
    Vk[[k]]<-var(Xk[[k]])*(nk[k]-1)/nk[k]
    Rk[[k]]<-cor(Xk[[k]])
}

names(Vk) <- levels(y)
names(Rk) <- levels(y)
names(gk) <- levels(y)

Xcent<-X-matrix(rep(1,n),ncol=1)%*%t(g)
W<-matrix(0,nrow=p,ncol=p)
B<-matrix(0,nrow=p,ncol=p)
for (k in 1:K){
    W<-W+pk[k]*Vk[[k]]
    B<-B+pk[k]*(gk[[k]]-g)%*%t(gk[[k]]-g)
}

# Calcul des FD :

res<-eigen(solve(V)%*%B)

eig <- Re(res$values[1:(K-1)])


UU<-Re(res$vectors[,1:(K-1)])
normes <- sqrt(diag(t(UU)%*%V%*%UU))

if (K>2) {
U <- sweep(UU,2,normes,"/")
} else {
U <- UU/normes
}

if (K>2) {
    rownames(U) <- colnames(X)
    colnames(U)<-paste("u",1:ncol(U),sep="")
} else
{
    names(U) <- colnames(X)

}

# Calcul des VD :
S<-as.matrix(Xcent)%*%U
if (K>2) {
    colnames(S) <-paste("s",1:ncol(S),sep="")
    rownames(S)<-rownames(X)
} else names(S)<-rownames(X)




list(eig=eig,U=U,S=S,nk=nk,K=K,y=y,X=X,Xk=Xk,W=W,B=B,g=data.frame(g=g),gk=data.frame(gk),V=V,R=R,Vk=Vk,Rk=Rk)
}

Fonction qui plot nos resulats AFD

plotAFD<-function(res)
  {
dim = c(1)
nbaxes <- length(res$eig)
if ((length(dim) > 1) && (nbaxes >1)) {
min.x<-min(res$S[,dim[1]])
max.x<-max(res$S[,dim[1]])
min.y<-min(res$S[,dim[2]])
max.y<-max(res$S[,dim[2]])
dim1 <- dim[1]
dim2 <- dim[2]
quali <- as.factor(res$y)
plot(res$S[,dim],xlab=paste("Dim",dim1),ylab=paste("Dim",dim2), xlim=c(min.x,max.x),ylim=c(min.y,max.y),pch=" ",main="Individus")
abline(h = 0, lty = 2,)
abline(v = 0, lty = 2)
legend("topleft",legend = levels(quali), text.col = c(1:length(levels(quali))),cex=0.8)
text(res$S[,dim],labels=rownames(res$S),col=as.numeric(quali))

mat.corr <- res$mat.corr
xlim <- c(-1, 1)*1.3
ylim <- c(-1, 1)*1.3
plot(0, 0, type="n",xlab=paste("Dim",dim1),ylab=paste("Dim",dim2), xlim = xlim, ylim = ylim, main="Variables")
x.cercle <- seq(-1, 1, by = 0.01)
y.cercle <- sqrt(1 - x.cercle^2)
lines(x.cercle, y = y.cercle)
lines(x.cercle, y = -y.cercle)
abline(v = 0, lty = 2)
abline(h = 0, lty = 2)
for (j in 1:nrow(mat.corr)) {
        arrows(0, 0, mat.corr[j,dim1], mat.corr[j,dim2], length = 0.1, angle = 15, code = 2)
        if (abs(mat.corr[j,dim1]) > abs(mat.corr[j,dim2])) {
                        if (mat.corr[j,dim1] >= 0) pos <- 4
                else pos <- 2
             } else {
                        if (mat.corr[j,dim2] >= 0) pos <- 3
                        else pos <- 1
                    }
              text(mat.corr[j,dim1],mat.corr[j,dim2],labels = rownames(mat.corr)[j], pos = pos)
}
} else {
    if (nbaxes==1) {
        quali <- as.factor(res$y)
        plot(res$S,rep(0,length(res$S)),xlab="Axe 1",ylab="",col=as.numeric(quali))
        legend("topleft",legend = levels(quali), text.col = c(1:length(levels(quali))),cex=0.8)
    } else {
        quali <- as.factor(res$y)
        plot(res$S[,dim],rep(0,length(res$S[,dim])),xlab=paste("Axe",dim),ylab="",col=as.numeric(quali))
        legend("topleft",legend = levels(quali), text.col = c(1:length(levels(quali))),cex=0.8)

    }
}   
}

Fonction acp_from_scratch du TD ACP

acp_from_scratch <- function(X, nombre_composants=NULL, reduire = FALSE) {
  
  n <- nrow(X)
  X <- as.matrix(X)
  data <- centrage(X)
  if (reduire) {
    data <- centrer_reduire(X)
  }

  #calculer la matrice de covariance matriciellement
  matrice_covariance <- (1/n) * t(data) %*% data
  
  # Calcul des valeurs propres et vecteurs propres
  eig_result <- eigen(matrice_covariance)
  valeurs_propres <- eig_result$values
  vecteurs_propres <- eig_result$vectors
  
  # Sélectionner le nombre spécifié de composants principaux
  if (!is.null(nombre_composants)) {
    vecteurs_propres <- vecteurs_propres[, 1:nombre_composants, drop = FALSE]
  }
  
  # Projection des données sur les composants principaux
  nouvelles_coordonnees <- data %*% vecteurs_propres
  
  return(list(matrice_covariance=matrice_covariance,
              nouvelles_coordonnees=nouvelles_coordonnees, 
              valeurs_propres=valeurs_propres, 
              vecteurs_propres=vecteurs_propres))
  
}

Fonction qui centre les données

centrage <- function(data) {
  n <- nrow(data)
  
  # Centrage des données
  # Définir la matrice des poids
  matrice_de_poids <- (1/n) * diag(n)
  
  # Calculer la moyenne pondérée matriciellement
  mat1_n <- matrix(rep(1, n), nrow=1)
  moyenne_ponderee <- mat1_n %*% data / n
  
  #calculer data_centree
  data_centree <- sweep(data, 2, moyenne_ponderee)
  
  return (data_centree)
}

Fonction qui normalise les données

centrer_reduire <- function(donnees) {
  n <- nrow(donnees)  # Nombre d'observations
  p <- ncol(donnees)  # Nombre de variables
  
  # Centrage des données
  donnees_centrees = centrage(donnees)
  
  # Réduction des données (calcul de l'écart-type de façon matricielle avec biais)
  # Calculer la matrice de covariance avec biais (diviser par n!)
  matrice_covariance <- (1 / n) * (t(donnees_centrees) %*% donnees_centrees)
  # Extraire l'écart-type à partir de la diagonale de la matrice de covariance
  ecarts_types <- sqrt(diag(matrice_covariance))
  
  # Réduire les données : diviser par l'écart-type
  donnees_reduites <- sweep(donnees_centrees, 2, ecarts_types, FUN = "/")
  
  return(donnees_reduites)
}

Importer les données

library(readr)
# Importer les données de l'énoncé
data <- read.csv("Income_Inequality.csv", sep=";")

# Il y a t il des NA ? NON
print(sum(is.na(data)))
## [1] 0

On importe les data du binôme (qui est un sur Python et l’autre R) pour avoir les mêmes données TRAIN (et même indices surtout) pour comparer nos 2 méthodes.

# Importer les données de notre binôme
# X correspond a l'indice initale des data sans train test split

data_train_X <- read.csv("data_train-2.csv")
data_train <- data_train_X
data_train$X <- NULL
data_test_X <- read.csv("data_test-2.csv")
data_test <- data_test_X
data_test$X <- NULL

#print(data_train)
#print(data_test)

Sélectionner uniquement les colonnes numériques pour l’ACP

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

—- QUESTION 4A —-

###Pour savoir si on norme ou non les données, on regarde les variances de nos variables, mais on reviendra sur cette question plus tard.

On represente les box plot

# Parcourir chaque colonne et créer un boxplot pour celle-ci
for (i in 1:ncol(data_train_numeric)) {
    # Créer un nouveau graphique pour chaque colonne
    plot.new()
    
    # Créer le boxplot pour la colonne courante
    boxplot(data_train_numeric[[i]], 
            main = paste("Boxplot pour", colnames(data_train_numeric)[i]), 
            ylab = "Valeurs",
            col = "lightblue", 
            border = "brown")
}

##################################################################################8

On remarque qu’il y a beaucoup de ce qu’on peut considérer de manière statistique comme des valeurs aberrantes (variables Eco 1, Eco3, Energy 2, Energy 3, Health 1, Health 2, Finance 2, Finance 3, Env, Other 1, Other 3)

Or pour ces variables (comme Health ou Energy) on peut s’attendre à une grande variance dans les données entre pays H et pays L.

On n’a pas plus d’information sur la signification des données donc on ne peut pas savoir s’il s’agit de valeurs aberrantes, d’anomalies ou si il n’y a effectivement pas d’erreurs dans les données.

##################################################################################8

Matrice des corrélations

library(corrplot)
## corrplot 0.92 loaded
library(dplyr)

# Définir la mise en page pour afficher deux graphiques côte à côte
par(mfrow = c(1, 2))

# Calculer la matrice de corrélation et créer un graphique de corrélation
cor(data_train_numeric) %>% 
  corrplot::corrplot()

# Calculer la matrice de corrélation et créer un second graphique de corrélation
# avec une méthode de coloration différente et une commande de regroupement hiérarchique
cor(data_train_numeric) %>% 
  corrplot::corrplot(method = "color", order = "hclust")

###########################################################################8

Une forte corrélation négative est observable entre la variable Poverty et les variables économiques, énergétiques, financières, gouvernementales et environnementales. En revanche, une corrélation positive prononcée existe avec Health 1 et 2. Cela suggère que Health 1 et 2 tendent à augmenter dans les pays plus pauvres. Cette observation est également corroborée par la corrélation négative de Health 1 et 2 avec les variables financières, économiques, gouvernementales, énergétiques et environnementales. Par ailleurs, les variables économiques affichent une forte corrélation positive entre elles, tout comme les variables Health. Finalement, à l’exception des variables Other 1, 2 et 3, qui montrent peu de corrélation avec les autres variables, la plupart des corrélations observées sont remarquablement fortes.

###########################################################################8

Calcul de la variabilité pour l’ensemble des données numeriques

data_train_numeric <- data.frame(data_train_numeric)
library(readxl)
library(dplyr)
library(ggplot2)
library(pander)

# Calcul de la variabilité pour l'ensemble des données numeriques
variability_overall <- data_train_numeric %>%
  summarise(across(where(is.numeric), list(moyenne = mean, ecart_type = sd, mediane = median)))

panderOptions('digits', 3)
pander(variability_overall)
Table continues below
Eco1_moyenne Eco1_ecart_type Eco1_mediane Eco2_moyenne Eco2_ecart_type
16046 19703 6291 15641 11885
Table continues below
Eco2_mediane Eco3_moyenne Eco3_ecart_type Eco3_mediane Energy1_moyenne
11855 21580 18411 14834 51.3
Table continues below
Energy1_ecart_type Energy1_mediane Energy2_moyenne Energy2_ecart_type
24 55 85 26.9
Table continues below
Energy2_mediane Energy3_moyenne Energy3_ecart_type Energy3_mediane
99.7 90.6 50.3 80.9
Table continues below
Health1_moyenne Health1_ecart_type Health1_mediane Health2_moyenne
27.3 33.4 13.8 2.65
Table continues below
Health2_ecart_type Health2_mediane Finan1_moyenne Finan1_ecart_type
1.42 2.07 0.404 0.256
Table continues below
Finan1_mediane Finan2_moyenne Finan2_ecart_type Finan2_mediane
0.353 4495758 10473592 1092634
Table continues below
Finan3_moyenne Finan3_ecart_type Finan3_mediane Finan4_moyenne
0.604 0.288 0.502 2.73
Table continues below
Finan4_ecart_type Finan4_mediane Finan5_moyenne Finan5_ecart_type
0.704 2.83 0.641 0.236
Table continues below
Finan5_mediane Governance_moyenne Governance_ecart_type
0.635 0.197 0.928
Table continues below
Governance_mediane Poverty_moyenne Poverty_ecart_type Poverty_mediane
0.0469 26.3 19.9 23.6
Table continues below
Env_moyenne Env_ecart_type Env_mediane Other1_moyenne
4.9 5.44 3.8 5.18
Table continues below
Other1_ecart_type Other1_mediane Other2_moyenne Other2_ecart_type
6.83 2.51 15.1 15.5
Other2_mediane Other3_moyenne Other3_ecart_type Other3_mediane
11.8 4.58 1.46 4.48

L’analyse précédente confirme nos observations initiales : pour certaines variables (comme Eco) il existe une différence notable entre les médianes et les moyennes, indiquant la présence de valeurs extrêmes qui contribuent à une variance élevée. Cette constatation souligne l’importance de normaliser nos données avant de procéder à l’ACP et à l’AFD. De plus, nous observons des écarts significatifs d’échelles entre les variables – certaines, comme Finance 2, se chiffrent en millions, tandis que d’autres, telles que Finance 1 et Health 2, présentent des valeurs proches de 1. Cette hétérogénéité des échelles renforce la nécessité d’une normalisation soignée pour assurer la cohérence et l’efficacité de nos analyses.

Calcul de la variabilité pour chaque type de pays H et L

data_train2 <- data_train
data_train2$Country <- NULL
data_train2$Year <- NULL
# Calcul de la variabilité pour chaque campagne
type_variability <- data_train2 %>%
  group_by(Income_Inequality) %>%
  summarise(across(where(is.numeric), list(moyenne = mean, ecart_type = sd, mediane = median)))

pander(type_variability)
Table continues below
Income_Inequality Eco1_moyenne Eco1_ecart_type Eco1_mediane
H 7306 9270 5025
L 29668 23577 28414
Table continues below
Eco2_moyenne Eco2_ecart_type Eco2_mediane Eco3_moyenne Eco3_ecart_type
10264 8773 8419 13359 11658
24024 11242 26820 34396 19664
Table continues below
Eco3_mediane Energy1_moyenne Energy1_ecart_type Energy1_mediane
11419 42.4 23.2 46
37311 65.2 17.9 67.5
Table continues below
Energy2_moyenne Energy2_ecart_type Energy2_mediane Energy3_moyenne
76.4 30.9 95.1 87.8
98.5 8.48 100 95
Table continues below
Energy3_ecart_type Energy3_mediane Health1_moyenne Health1_ecart_type
41.6 82 38.4 37.2
61.3 79 9.97 14.2
Table continues below
Health1_mediane Health2_moyenne Health2_ecart_type Health2_mediane
19.1 3.12 1.53 2.46
4.45 1.9 0.778 1.69
Table continues below
Finan1_moyenne Finan1_ecart_type Finan1_mediane Finan2_moyenne
0.306 0.202 0.246 4332084
0.557 0.257 0.614 4750898
Table continues below
Finan2_ecart_type Finan2_mediane Finan3_moyenne Finan3_ecart_type
12284089 480434 0.503 0.159
6762158 1902460 0.762 0.364
Table continues below
Finan3_mediane Finan4_moyenne Finan4_ecart_type Finan4_mediane
0.463 2.37 0.621 2.48
0.726 3.29 0.391 3.36
Table continues below
Finan5_moyenne Finan5_ecart_type Finan5_mediane Governance_moyenne
0.573 0.208 0.534 -0.202
0.747 0.239 0.777 0.82
Table continues below
Governance_ecart_type Governance_mediane Poverty_moyenne
0.665 -0.225 36.7
0.939 0.998 10.2
Table continues below
Poverty_ecart_type Poverty_mediane Env_moyenne Env_ecart_type
17.4 34.1 3.17 4.1
10.5 6.21 7.6 6.13
Table continues below
Env_mediane Other1_moyenne Other1_ecart_type Other1_mediane
1.72 6.1 6.69 3.71
6.28 3.74 6.81 0.708
Table continues below
Other2_moyenne Other2_ecart_type Other2_mediane Other3_moyenne
17.9 16.8 21 4.27
10.9 12.2 9.1 5.06
Other3_ecart_type Other3_mediane
1.42 4.22
1.4 4.93 
data_train2 <- data_train

##################################################################################8

Il apparaît que les pays classés dans la catégorie ‘L’ correspondent aux nations ‘riches’, en supposant que les indicateurs économiques augmentent avec la richesse d’un pays. Inversement, les pays classés comme ‘H’ semblent être relativement ‘pauvres’. Cette répartition s’aligne logiquement avec la signification des labels ‘H’ et ‘L’ en termes d’inégalités de revenus : on observe généralement des inégalités de revenus moins marquées dans les pays plus aisés.

##################################################################################8

Test de Mann-Whitney U pour chaque variable entre H et L:

##################################################################################8

On veut tester si les pays classés comme “L” et “H” sont statistiquement similaires ou non sur les 19 variables quantitatives, on va utiliser un test statistique. L’approche la plus courante serait d’effectuer des tests t de Student pour chaque variable quantitative. Mais la normalité des données est en question, il n’y a pas de raison que cela soit le cas, on va plutôt opter pour des tests non paramétriques comme le test de Mann-Whitney U.

##################################################################################8

library(dplyr)
library(stringr)

wilcox.test(data_train$Eco1 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Eco1 by data_train$Income_Inequality
## W = 18488, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Eco2 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Eco2 by data_train$Income_Inequality
## W = 15341, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Eco3 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Eco3 by data_train$Income_Inequality
## W = 16092, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Energy1 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Energy1 by data_train$Income_Inequality
## W = 19456, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Energy2 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Energy2 by data_train$Income_Inequality
## W = 15180, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Energy3 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Energy3 by data_train$Income_Inequality
## W = 43992, p-value = 0.9409
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Health1 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Health1 by data_train$Income_Inequality
## W = 74035, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Health2 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Health2 by data_train$Income_Inequality
## W = 69216, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Finan1 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Finan1 by data_train$Income_Inequality
## W = 20357, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Finan2 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Finan2 by data_train$Income_Inequality
## W = 30176, p-value = 4.248e-11
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Finan3 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Finan3 by data_train$Income_Inequality
## W = 27206, p-value = 1.275e-15
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Finan4 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Finan4 by data_train$Income_Inequality
## W = 8340, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Finan5 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Finan5 by data_train$Income_Inequality
## W = 25705, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Governance ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Governance by data_train$Income_Inequality
## W = 18101, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Poverty ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Poverty by data_train$Income_Inequality
## W = 80170, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Other1 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Other1 by data_train$Income_Inequality
## W = 61838, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Other2 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Other2 by data_train$Income_Inequality
## W = 52848, p-value = 2.047e-05
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(data_train$Other3 ~ data_train$Income_Inequality)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  data_train$Other3 by data_train$Income_Inequality
## W = 30108, p-value = 3.42e-11
## alternative hypothesis: true location shift is not equal to 0

###########################################################################8

CCL : Toutes les p-valeurs sont extrêmement faibles (sauf Energy 3 supérieur au seuil 5%).

On peut rejeter l’hypothèse nulle (qui serait que les deux échantillons proviennent de la même distribution ou que leur différence médiane est égale à zéro). Il y a donc des preuves statistiques suffisantes pour conclure que les valeurs des variables (sauf Energy 3) diffèrent de manière significative entre les classes h et L de Income_Inequality.

Dans le cas de Energy 3 on ne peut rien dire.

###########################################################################8

Calcul de la variabilité pour chaque pays

data_train2 <- data_train
data_train2$Year <- NULL

# Calcul de la variabilité pour chaque campagne
pays_variability <- data_train2 %>%
  group_by(Country) %>%
  summarise(across(where(is.numeric), list(moyenne = mean, ecart_type = sd, mediane = median)))

print(pays_variability)
## # A tibble: 87 × 58
##    Country Eco1_moyenne Eco1_ecart_type Eco1_mediane Eco2_moyenne
##      <int>        <dbl>           <dbl>        <dbl>        <dbl>
##  1       1        3064.           101.         3101.        4470.
##  2       2       13517.           414.        13568.       18059.
##  3       3        3485.           382.         3534.       11000.
##  4       4       55994.          1668.        56025.       36352.
##  5       5       45089.           973.        44590.       35713.
##  6       6         292.            13.7         294.         859.
##  7       7       41176.          1706.        41176.       32649.
##  8       8        1055.            61.4        1050.        2481.
##  9       9         628.            44.3         629.        1583.
## 10      10        6997.           586.         6745.       14717.
## # ℹ 77 more rows
## # ℹ 53 more variables: Eco2_ecart_type <dbl>, Eco2_mediane <dbl>,
## #   Eco3_moyenne <dbl>, Eco3_ecart_type <dbl>, Eco3_mediane <dbl>,
## #   Energy1_moyenne <dbl>, Energy1_ecart_type <dbl>, Energy1_mediane <dbl>,
## #   Energy2_moyenne <dbl>, Energy2_ecart_type <dbl>, Energy2_mediane <dbl>,
## #   Energy3_moyenne <dbl>, Energy3_ecart_type <dbl>, Energy3_mediane <dbl>,
## #   Health1_moyenne <dbl>, Health1_ecart_type <dbl>, Health1_mediane <dbl>, …
data_train2 <- data_train

###########################################################################8

L’analyse des données, une fois regroupées par pays, révèle un phénomène intéressant : pour l’ensemble des variables considérées, les moyennes et les médianes sont remarquablement proches les unes des autres. Cette constatation suggère que la variabilité et les valeurs “extrêmes” des données est significativement plus importante entre les différents pays (inter-pays) qu’à l’intérieur d’un même pays au fil des années intra-pays).

En d’autres termes, il semble que les fluctuations observées au début de l’analyse dans les variables semblent être plus fortement influencées par les différences nationales que par les variations annuelles au sein d’un même pays.

Cette tendance s’explique logiquement par la nature même des indicateurs économiques, politiques etc… En effet, des événements majeurs tels que l’effondrement d’une économie (variable Eco) ou un changement de gouvernement suite à un coup d’État par exemple (variable Governance) ne se produisent pas fréquemment. Par conséquent, il est naturel que les données relatives à ces variables ne varient pas considérablement d’une année à l’autre. Ainsi, cette relative stabilité annuelle au sein d’un même pays souligne la rareté des événements dramatiques impactant significativement l’économie ou la gouvernance d’une nation, d’autant plus qu’une periode de 10 ans (ce qui est énorme) est prise en compte dans l’analyse.

###########################################################################8

Calcul de la variabilité pour chaque année

data_train2 <- data_train
data_train2$Country <- NULL

# Calcul de la variabilité pour chaque campagne
year_variability <- data_train2 %>%
  group_by(Year) %>%
  summarise(across(where(is.numeric), list(moyenne = mean, ecart_type = sd, mediane = median)))

print(year_variability)
## # A tibble: 10 × 58
##     Year Eco1_moyenne Eco1_ecart_type Eco1_mediane Eco2_moyenne Eco2_ecart_type
##    <int>        <dbl>           <dbl>        <dbl>        <dbl>           <dbl>
##  1  2010       13544.          18668.        5208.       13040.          11040.
##  2  2011       16933.          21465.        5511.       15385.          12675.
##  3  2012       14725.          17983.        6643.       15029.          11448.
##  4  2013       16483.          18791.        9216.       16681.          11593.
##  5  2014       15641.          19730.        6363.       15078.          11258.
##  6  2015       15967.          20537.        5709.       15346.          12071.
##  7  2016       15888.          19571.        6186.       15553.          12127.
##  8  2017       15616.          19903.        6340.       15672.          12062.
##  9  2018       18489.          21136.        8739.       17266.          12868.
## 10  2019       17403.          20213.        8592.       17543.          12011.
## # ℹ 52 more variables: Eco2_mediane <dbl>, Eco3_moyenne <dbl>,
## #   Eco3_ecart_type <dbl>, Eco3_mediane <dbl>, Energy1_moyenne <dbl>,
## #   Energy1_ecart_type <dbl>, Energy1_mediane <dbl>, Energy2_moyenne <dbl>,
## #   Energy2_ecart_type <dbl>, Energy2_mediane <dbl>, Energy3_moyenne <dbl>,
## #   Energy3_ecart_type <dbl>, Energy3_mediane <dbl>, Health1_moyenne <dbl>,
## #   Health1_ecart_type <dbl>, Health1_mediane <dbl>, Health2_moyenne <dbl>,
## #   Health2_ecart_type <dbl>, Health2_mediane <dbl>, Finan1_moyenne <dbl>, …
data_train2 <- data_train

###########################################################################8

On retrouve que les moyennes et médiannes sont différentes à une année fixée et que les moyennes d’une variable fixée restent trés proches, donc ceci confirme que les fluctuations observées au début de l’analyse dans les variables semblent être plus fortement influencées par les différences nationales que par les variations annuelles au sein d’un même pays, qui se produirait lors de catastrophe majeur.

###########################################################################8

Correlations pour les pays H et Correlation pour les pays L

library(dplyr)
library(ggplot2)
library(reshape2)



data_H <- data_train %>%
  filter(Income_Inequality == 'H') %>%
  select_if(is.numeric)

data_L <- data_train %>%
  filter(Income_Inequality == 'L') %>%
  select_if(is.numeric)

correlation_matrix_H <- cor(data_H)
correlation_matrix_L <- cor(data_L)




correlation_matrix_melted_H <- melt(correlation_matrix_H)
correlation_matrix_melted_L <- melt(correlation_matrix_L)


heatmap_H <- ggplot(correlation_matrix_melted_H, aes(Var1, Var2, fill = value)) +
  geom_tile() +
  scale_fill_gradient2(low = "blue", high = "red", mid = "white", 
                       midpoint = 0, limit = c(-1,1), space = "Lab", 
                       name="Correlation") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1),
        axis.title = element_blank()) +
  labs(title = "Matrice de Corrélation pour H")

heatmap_L <- ggplot(correlation_matrix_melted_L, aes(Var1, Var2, fill = value)) +
  geom_tile() +
  scale_fill_gradient2(low = "blue", high = "red", mid = "white", 
                       midpoint = 0, limit = c(-1,1), space = "Lab", 
                       name="Correlation") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1),
        axis.title = element_blank()) +
  labs(title = "Matrice de Corrélation pour L")


grid.arrange(heatmap_H, heatmap_L, ncol = 2)

###########################################################################8

À première vue, les deux matrices de corrélation apparaissent trés similaires. Pour confirmer cette impression de manière quantitative, nous pouvons calculer la matrice des différences absolues entre les éléments correspondants des deux matrices case à case. Des valeurs proches de zéro dans cette matrice de différences indiqueront une grande similarité, tandis que des valeurs plus élevées signaleront des divergences notables entre les matrices comparées H et L.

###########################################################################8

Calcul de la différence absolue entre les matrices de corrélation pour reperer les différence entre pays H et L

correlation_diff <- abs(correlation_matrix_H - correlation_matrix_L)
correlation_diff_melted <- melt(correlation_diff)

heatmap_diff <- ggplot(correlation_diff_melted, aes(Var1, Var2, fill = value)) +
  geom_tile() +
  scale_fill_gradient2(low = "blue", high = "red", mid = "white", 
                       midpoint = 1, limit = c(0,2), space = "Lab", 
                       name="Différence Absolue") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1),
        axis.title = element_blank()) +
  labs(title = "Différence Absolue des Matrices de Corrélation H et L")
print(heatmap_diff)

###########################################################################8

Interprétation : • Stabilité des Relations entre Variables : La similitude des matrices de corrélation suggère que les relations entre les différentes variables étudiées sont stables et cohérentes, quel que soit le classement des pays en classes “H” ou “L”. Cela peut indiquer que les dynamiques sous-jacentes ou les interactions entre ces variables sont similaires dans les deux groupes et indépendantes des inégalitées de revenus.

• Indépendance par rapport au classement Géo-Politique : La correspondance étroite entre les matrices signifie que les relations entre les variables ne sont pas fortement influencées par les facteurs économiques, sociaux ou politiques qui distinguent les pays de classe “H” des pays de classe “L”. En d’autres termes, les variables interagissent de manière similaire indépendamment du niveau de développement ou de richesse des pays.

Exploitation : • Modélisation et Prédictions : La constance des relations entre les variables peut permettre de développer des modèles prédictifs ou des analyses qui seront valables pour les deux classes de pays. Cela pourrait être utile dans les domaines de la planification économique, de la politique publique…

• Facteurs Communs : Cette similarité peut également motiver une recherche de facteurs communs ou universels qui influencent ces variables, suggérant que certaines forces ou tendances agissent de manière uniforme à travers différentes économies, ceci ne sera pas fait ici.

###########################################################################8

—- QUESTION 4B —-

ACP avec notre fonction FROM SCRATCH sur les 70%

# Exécution de l'ACP
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=2, reduire = TRUE)

# Préparation des données pour le graphique
individus <- as.data.frame(res_pca$nouvelles_coordonnees)

# Graphique des individus
plot(individus[,1], individus[,2], xlab="CP 1", ylab="CP 2", main="ACP - Graphique des Individus sur les 70%", col = "blue",ylim=c(-7, 4))

Plot en 3D sur les 70%

# Chargement des packages nécessaires
library(plotly)
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:MASS':
## 
##     select
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
# Exécution de l'ACP pour obtenir au moins trois composantes principales
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=3, reduire = TRUE)

# Préparation des données pour le graphique
individus_3d <- as.data.frame(res_pca$nouvelles_coordonnees)

# Création du graphique 3D avec plotly
fig <- plot_ly(individus_3d, x = ~V1, y = ~V2, z = ~V3, type = 'scatter3d', mode = 'markers',
               marker = list(size = 2, color = 'blue'))

# Ajout des labels pour les axes
fig <- fig %>% layout(scene = list(xaxis = list(title = 'CP 1'),
                                   yaxis = list(title = 'CP 2'),
                                   zaxis = list(title = 'CP 3')))

# Affichage du graphique
fig

###########################################################################8

On voit un groupe en bas à gauche de 7 points et un deuxieme de 8 points aussi en bas à gauche. On se demande pourquoi et on y répondra dans la suite

Il semble y avoir les 2 classes à droite et à gauche et un mélange au centre, on le vérifira quand on mettra en couleur les classses H et L.

###########################################################################8

Barplot de l’inertie expliquée et cumulée :

# Exécution de l'ACP
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=NULL, reduire = TRUE)

# Préparation pour l'inertie expliquée et cumulée
inertie_expliquee <- res_pca$valeurs_propres / sum(res_pca$valeurs_propres)
cumul_inertie <- cumsum(inertie_expliquee)
ylim_max <- max(cumul_inertie) + .1

# Barplot de l'inertie expliquée
midpoints <- barplot(inertie_expliquee, main="ACP - Inertie expliquée et cumulée sur les 70%", xlab="Composantes", ylab="Inertie", ylim=c(0, ylim_max), col = "blue")
lines(midpoints, cumul_inertie, type="b", col="red")
abline(h = 1, col="grey", lty=2)  # Ligne pointillée pour l'asymptote

Règle de Karlis - Saporta - Spinaki

sueil = mean(VP) + 2sd(VP)

# Seuil calculé à partir de la moyenne et de l'écart-type des valeurs propres
karlis_threshold_value <- mean(res_pca$valeurs_propres) + 2 * sd(res_pca$valeurs_propres)
cat("Valeur de karlis_sueil :", karlis_threshold_value)
## Valeur de karlis_sueil : 5.673652

Règle de Kaiser - Guttman NORMEE

sueil = 1

# Barplot de l'inertie expliquée pour la normalisée
midpoints <- barplot(inertie_expliquee, main="ACP NORMEE + règles classiques sur les 70%", xlab="Composantes", ylab="Inertie expliquée", ylim=c(0, ylim_max), col = "blue")
lines(midpoints, cumul_inertie, type="b", pch=19, col="red")

# Ajouter les valeurs des valeurs propres en abscisse, écriture verticale
text(x=midpoints, y=-0.05, labels=round(res_pca$valeurs_propres, 2), srt=90, adj=1, xpd=TRUE, cex=0.8)

# Règle de Kaiser - Guttman
# Identifier la dernière composante avec une valeur propre >= 1 et la position de la barre correspondante
kaiser_comps <- which(res_pca$valeurs_propres >= 1)
if(length(kaiser_comps) > 0) {
  kaiser_limit <- max(kaiser_comps)
  # Utilisez les midpoints pour aligner la ligne verticale avec la barre correspondante
  abline(v=midpoints[kaiser_limit]-.1, col="green", lty=2)
  abline(v=midpoints[kaiser_limit]+.1, col="violet", lty=2)
}
abline(v=.75, col="orange", lty=2)
abline(h = .8, col="grey", lty=2)  # Ligne pointillée pour l'asymptote
# Ajouter une légende pour les seuils
legend("topright",inset=c(0, 0.2), legend=c("Karlis-Saporta-Spinaki","Kaiser-Guttman seuil 1", "règle des 80% inertie cumulée"), col=c("orange", "green", "violet"), lty=2, cex=0.8)

On revient sur la question : faut il NORMEE les données ?

# Exécution de l'ACP
result_normee <- acp_from_scratch(data_train_numeric, nombre_composants=NULL, reduire = TRUE)
result_centree <-acp_from_scratch(data_train_numeric, nombre_composants=NULL, reduire = FALSE)

# Définissez d'abord vos données pour la centrée
inertie_expliquee_centree <- result_centree$valeurs_propres / sum(result_centree$valeurs_propres)
cumul_inertie_centree <- cumsum(inertie_expliquee_centree)
ylim_max_centree <- max(cumul_inertie_centree) + .1

# Puis pour la normalisée
inertie_expliquee_normee <- result_normee$valeurs_propres / sum(result_normee$valeurs_propres)
cumul_inertie_normee <- cumsum(inertie_expliquee_normee)
ylim_max_normee <- max(cumul_inertie_normee) + .1

# Définir les paramètres de la fenêtre graphique pour les subplots
par(mfrow=c(1, 2))  # 1 ligne, 2 colonnes

# Premier subplot pour la centrée
midpoints_centree <- barplot(inertie_expliquee_centree, main="CENTREE", xlab="Composantes", ylab="Inertie expliquée et cumulée", ylim=c(0, ylim_max_centree), col = "blue")
lines(midpoints_centree, cumul_inertie_centree, type="b", col="red")
abline(h = 1, col="grey", lty=2)  # Ligne pointillée pour l'asymptote

# Second subplot pour la normalisée
midpoints_normee <- barplot(inertie_expliquee_normee, main="NORMEE", xlab="Composantes", ylab="Inertie expliquée et cumulée", ylim=c(0, ylim_max_normee), col = "blue")
lines(midpoints_normee, cumul_inertie_normee, type="b", col="red")
abline(h = 1, col="grey", lty=2)  # Ligne pointillée pour l'asymptote

# Réinitialiser les paramètres de la fenêtre graphique
par(mfrow=c(1, 1))

###########################################################################8

On voit bien que oui car on a un cas classique d’ACP CENTREE : la premiere composante écrase toutes les autres en captant largement la variance.

###########################################################################8

CCL : Il faut NORMEE nos données

PLOT 3D des classes en couleur.

library(plotly)
library(FactoMineR)
library(dplyr)
library(stringr)

# Exécution de l'ACP
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=NULL, reduire = TRUE)

# Coordonnées des individus obtenues à partir de votre fonction ACP
ind_coords <- res_pca$nouvelles_coordonnees

# On va visualiser la classe "Income_Inequality"
data_train_numeric$Income_Inequality <- data_train$Income_Inequality

# Filtration des données pour chaque catégorie
data_L <- data_train_numeric[data_train_numeric$Income_Inequality == "L", ]
data_H <- data_train_numeric[data_train_numeric$Income_Inequality == "H", ]

# Coordonnées pour chaque catégorie
coords_L <- ind_coords[data_train_numeric$Income_Inequality == "L", ]
coords_H <- ind_coords[data_train_numeric$Income_Inequality == "H", ]

# Création du graphique Plotly
plot_ly() %>%
  add_trace(data = data_L, x = ~coords_L[, 1], y = ~coords_L[, 2], z = ~coords_L[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'red', opacity = 0.7), name = "L") %>%
  add_trace(data = data_H, x = ~coords_H[, 1], y = ~coords_H[, 2], z = ~coords_H[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'green', opacity = 0.7), name = "H") %>%
  layout(title = "ACP: Pays L (rouge) vs H (vert) sur les 70%",
         scene = list(xaxis = list(title = 'PC1'),
                      yaxis = list(title = 'PC2'),
                      zaxis = list(title = 'PC3')),
         legend = list(x = 0.1, y = 0.9))
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

###########################################################################8

On voit apparaite nos 2 classes : en vert pays H d’un côté. En rouge les pays H de l’autre. Au centre, il semble que les pays H soient un peu entouré de pays L. Le mélange n’est pas homogène mais plutôt hétérogène. On revient sur nos 2 groupes de 7 et 8 points. L’un est totalement pays L (groupe de 8) et l’autre pays H (groupe de 7). On y reviendra plus tard.

###########################################################################8

Representation dans les différents plan PC1-PC2, PC1-PC3 PC1-PC4, PC1-PC5, PC2-PC3

library(plotly)
library(FactoMineR)
library(dplyr)
library(stringr)

# Exécution de l'ACP
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=NULL, reduire = TRUE)

# Coordonnées des individus obtenues à partir de votre fonction ACP
ind_coords <- res_pca$nouvelles_coordonnees

# On va visualiser la classe "Income_Inequality"
data_train_numeric$Income_Inequality <- data_train$Income_Inequality

# Créer une fonction pour générer un graphique Plotly
generate_plot <- function(x, y, title) {
  data_L <- ind_coords[data_train_numeric$Income_Inequality == "L", ]
  data_H <- ind_coords[data_train_numeric$Income_Inequality == "H", ]
  
  plot_ly() %>%
    add_trace(x = ~data_L[, x], y = ~data_L[, y], 
              type = 'scatter', mode = 'markers',
              marker = list(size = 6, color = 'red'), name = "L") %>%
    add_trace(x = ~data_H[, x], y = ~data_H[, y], 
              type = 'scatter', mode = 'markers',
              marker = list(size = 6, color = 'green'), name = "H") %>%
    layout(title = title,
           xaxis = list(title = paste("PC", x)),
           yaxis = list(title = paste("PC", y)))
}

# Créer et afficher les graphiques pour chaque paire de composantes
generate_plot(1, 2, "ACP: PC1 vs PC2")
generate_plot(1, 3, "ACP: PC1 vs PC3")
generate_plot(1, 4, "ACP: PC1 vs PC4")
generate_plot(1, 5, "ACP: PC1 vs PC5")
generate_plot(2, 3, "ACP: PC2 vs PC3")
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

###########################################################################8

On voit que en augmentant les composantes, la qualité des 2 classes distinctes diminue à l’oeuil.

###########################################################################8

Qualité des projections des individus Qik=5 selon la formule du cours NORMEE

library(pander)

k <- 5  # Nombre de composantes à considérer

# Initialisation d'un data frame pour stocker les résultats
results_df <- data.frame(
moyenne = numeric(),
ecart_type = numeric(),
medianne= numeric()
)

# ACP
res_pca <- acp_from_scratch(as.matrix(data_train_numeric), nombre_composants = NULL, reduire = TRUE)
nouvelles_coordonnees <- as.matrix(res_pca$nouvelles_coordonnees)

# Initialisation du vecteur pour stocker la qualité des projections
qualite_projections_k5 <- numeric(nrow(res_pca$nouvelles_coordonnees))

# Boucle sur chaque individu
for (i in 1:nrow(res_pca$nouvelles_coordonnees)) {
    Ci <- res_pca$nouvelles_coordonnees[i, ]
    
    # Calcul de la qualité de la projection pour l'individu i
    # Vérification que le dénominateur n'est pas nul
    if(sum(Ci^2) != 0) {
        qualite_projections_k5[i] <- sum(Ci[1:k]^2) / sum(Ci^2)
        #print(sum(Ci[1:k]^2))
        #print(sum(Ci^2))
    } else {
        qualite_projections_k5[i] <- NA  # Assigner NA si le dénominateur nul
    }
}


# Calcul des moyennes et des écarts types 
moyenne <- mean(qualite_projections_k5)
ecart_type <- sd(qualite_projections_k5)
medianne <- median(qualite_projections_k5)

# Ajout des résultats dans le data frame 
results_df_1 <- rbind(results_df, data.frame(moyenne, ecart_type, medianne))


# Renommage des colonnes pour la clarté
names(results_df_1) <- c("Moyenne", "Ecart_type", "Médianne")

# Affichage du data frame avec pander pour un joli formatage
#print(qualite_projections)
pander(results_df_1)
Moyenne Ecart_type Médianne
0.748 0.195 0.782

Qualité des projections des individus Qik selon la formule du cours NORMEE SELON LE NBR K CHOISIT

library(pander)
library(ggplot2)



# ACP sur toutes les composantes
res_pca <- acp_from_scratch(as.matrix(data_train_numeric), nombre_composants = NULL, reduire = TRUE)
nouvelles_coordonnees <- as.matrix(res_pca$nouvelles_coordonnees)

# Initialisation d'un data frame pour stocker les résultats
results_df <- data.frame(
  K = integer(),
  Moyenne = numeric(),
  Ecart_type = numeric(),
  Médianne = numeric()
)

# Calculer les statistiques pour chaque valeur de k
for (k in 1:19) {
  qualite_projections <- apply(nouvelles_coordonnees, 1, function(Ci) {
    if(sum(Ci^2) != 0) {
      return(sum(Ci[1:k]^2) / sum(Ci^2))
    } else {
      return(NA)  # Assigner NA si le dénominateur est nul
    }
  })
  
  # Calcul des moyennes, écarts types et médianes
  moyenne <- mean(qualite_projections, na.rm = TRUE)
  ecart_type <- sd(qualite_projections, na.rm = TRUE)
  medianne <- median(qualite_projections, na.rm = TRUE)

  # Ajout des résultats dans le data frame 
  results_df <- rbind(results_df, data.frame(K = k, Moyenne = moyenne, Ecart_type = ecart_type, Médianne = medianne))
}


# Transformer le dataframe pour la heatmap
results_long <- melt(results_df, id.vars = 'K')

# Création de la heatmap avec valeurs dans les cases
ggplot(results_long, aes(x = variable, y = factor(K), fill = value)) +
  geom_tile() + 
  geom_text(aes(label = round(value, 2)), size = 3) +  # Ajouter des textes
  scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = median(results_df$Moyenne, na.rm = TRUE)) +
  labs(x = "Qualité des projections des individus selon K", y = "Nombre de Composantes (K)", fill = "Valeur") +
  theme_minimal()

###########################################################################8

On voit bien que la qualitée de projection tend vers 1 pour K se rapprochant de 19, c’est évident. Pour notre choix de réduction de k=5, on a une qualitée de projection de 75%

###########################################################################8

Tri des qualités de projection pour trouver les plus faibles (mal projetés) et récupération des indices (numéros des individus)

  # Tri des qualités de projection pour trouver les plus faibles (mal projetés)
  # et récupération des indices (numéros des individus)
sorted_indices <- order(qualite_projections_k5, decreasing = FALSE)

# Sélection des 10 individus les plus mal projetés
worst_projected_indices <- head(sorted_indices, 10)
worst_projected_qualities <- qualite_projections_k5[worst_projected_indices]

# Création d'un data frame pour les afficher
worst_projected_df <- data.frame(
  Numero_Individu = worst_projected_indices,
  Qualite_proj = worst_projected_qualities
)

# Affichage du data frame avec pander
pander(worst_projected_df)
Numero_Individu Qualite_proj
413 0.0645
445 0.101
92 0.101
142 0.106
246 0.111
416 0.111
429 0.112
247 0.13
577 0.142
507 0.164

Voici les individus les moins bien projetées

# Calculer les moyennes pour chaque variable
mean_data <- colMeans(data_train_numeric)

# Calculer l'écart type pour chaque variable
sd_data <- apply(data_train_numeric, 2, sd)

# Extraire les valeurs des individus mal projetés
worst_projected_data <- data_train_numeric[worst_projected_indices, ]

# Créer un tableau pour la comparaison
comparison_df <- data.frame(
  Variable = colnames(data_train_numeric),
  Mean_All_Data = mean_data,
  SD_All_Data = sd_data,  # Ajout de la colonne pour l'écart type
  Worst_Projected = colMeans(worst_projected_data)
)

# Utilisation de pander pour afficher le tableau
library(pander)
pander(comparison_df)
  Variable Mean_All_Data SD_All_Data Worst_Projected
Eco1 Eco1 16046 19703 7064
Eco2 Eco2 15641 11885 12650
Eco3 Eco3 21580 18411 16731
Energy1 Energy1 51.3 24 51.9
Energy2 Energy2 85 26.9 95.1
Energy3 Energy3 90.6 50.3 67.2
Health1 Health1 27.3 33.4 19.9
Health2 Health2 2.65 1.42 2.02
Finan1 Finan1 0.404 0.256 0.592
Finan2 Finan2 4495758 10473592 3795663
Finan3 Finan3 0.604 0.288 0.473
Finan4 Finan4 2.73 0.704 2.61
Finan5 Finan5 0.641 0.236 0.667
Governance Governance 0.197 0.928 0.333
Poverty Poverty 26.3 19.9 32.2
Env Env 4.9 5.44 5.4
Other1 Other1 5.18 6.83 3.15
Other2 Other2 15.1 15.5 25.3
Other3 Other3 4.58 1.46 4.6

#########################################8

Il nous semble que notre code ne marche pas bien pour gerer les indices des individus et nous n’arrivons pas à le fixer, ont le voit car les 10 individus les plus mal projetées ne sont pas statistiquement trés différents des autres individus, pas de signature que l’on peut constater.

D’autant plus que nous n’avons pas les mêmes individus entre ACP et la méthode des arbres. Nous aurions pu espérer quelques individues que l’on retrouve entre les 2 méthodes.

#########################################8

library(plotly)
library(FactoMineR)
library(dplyr)
library(stringr)

# Exécution de l'ACP
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=3, reduire = TRUE)

# Coordonnées des individus obtenues à partir de votre fonction ACP
ind_coords <- res_pca$nouvelles_coordonnees
#print(dim(ind_coords))

# On va visualiser la classe "Income_Inequality"
data_train_numeric$Income_Inequality <- data_train$Income_Inequality
data_train_numeric$X <- data_train_X$X
#print(data_train_numeric)
#print(data_train)

# Filtration des données pour chaque catégorie
data_L <- data_train_numeric[data_train_numeric$Income_Inequality == "L", ]
data_H <- data_train_numeric[data_train_numeric$Income_Inequality == "H", ]
#print(data_L)
#print(data_H)


# Coordonnées pour chaque catégorie
coords_L <- ind_coords[data_train_numeric$Income_Inequality == "L", ]
coords_H <- ind_coords[data_train_numeric$Income_Inequality == "H", ]





# Création du graphique Plotly
p <- plot_ly() %>%
  add_trace(data = data_L, x = ~coords_L[, 1], y = ~coords_L[, 2], z = ~coords_L[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'red', opacity = 0.7), name = "L") %>%
  add_trace(data = data_H, x = ~coords_H[, 1], y = ~coords_H[, 2], z = ~coords_H[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'green', opacity = 0.7), name = "H") %>%
  # Ajouter les 10 individus les plus mal projetés en noir
  add_trace(x = ind_coords[worst_projected_indices, 1], y = ind_coords[worst_projected_indices, 2], z = ind_coords[worst_projected_indices, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 4, color = 'black', opacity = 1), name = "Mal Projetés") %>%
  layout(title = "ACP: Pays L (rouge) vs H (vert) SUR LES 70%",
         scene = list(xaxis = list(title = 'PC1'),
                      yaxis = list(title = 'PC2'),
                      zaxis = list(title = 'PC3')),
         legend = list(x = 0.1, y = 0.9))

# Afficher le graphique
p
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL
data_train_numeric$X <- NULL

####################################################8

Ce critère n’a pas de signification pour les individus proches de l’origine.Or ils sont proches de l’origine. Quand on détecte un individu pour lequel le cosinus carré est faible, on doit tenir compte de sa distance à l’origine avant d’indiquer qu’il est mal représenté, donc ici ils sont proches et on ne peut rien dire pour les 10 points.

####################################################8

library(plotly)
library(FactoMineR)
library(dplyr)
library(stringr)

# Exécution de l'ACP
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=3, reduire = TRUE)

# Coordonnées des individus obtenues à partir de votre fonction ACP
ind_coords <- res_pca$nouvelles_coordonnees

# On va visualiser la classe "Income_Inequality"
data_train_numeric$Income_Inequality <- data_train$Income_Inequality
data_train_numeric$X <- data_train_X$X
#print(data_train_numeric)
#print(data_train)

# Filtration des données pour chaque catégorie
data_L <- data_train_numeric[data_train_numeric$Income_Inequality == "L", ]
data_H <- data_train_numeric[data_train_numeric$Income_Inequality == "H", ]

# Coordonnées pour chaque catégorie
coords_L <- ind_coords[data_train_numeric$Income_Inequality == "L", ]
coords_H <- ind_coords[data_train_numeric$Income_Inequality == "H", ]

# Calculer la distance de chaque point à l'origine
distances <- apply(ind_coords, 1, function(x) sqrt(sum(x^2)))
#print(ind_coords)
#print(distances)

distance = 4
# Filtrer les individus éloignés de l'origine
indices_far_enough <- which(distances > distance)

# Tri des qualités de projection pour trouver les plus faibles parmi ceux éloignés
sorted_indices_far_enough <- order(qualite_projections_k5[indices_far_enough], decreasing = FALSE)

# Sélection des 10 individus les plus mal projetés parmi ceux éloignés
worst_projected_indices_far_enough <- indices_far_enough[head(sorted_indices_far_enough, 10)]

# Ajout d'une sphère
theta <- seq(0, 2*pi, length.out = 100)
phi <- seq(0, pi, length.out = 100)
sphere_x <- distance * cos(theta) %o% sin(phi)
sphere_y <- distance * sin(theta) %o% sin(phi)
sphere_z <- distance * rep(1, length(theta)) %o% cos(phi)


# Création du graphique Plotly
p <- plot_ly() %>%
  # Ajouter la sphère en premier
  add_trace(x = c(sphere_x), y = c(sphere_y), z = c(sphere_z), 
            type = 'mesh3d', 
            opacity = 0.01,
            color = 'blue', 
            name = 'Sphère de Distance') %>%
  # Ajouter les autres éléments du graphique
  add_trace(data = data_L, x = ~coords_L[, 1], y = ~coords_L[, 2], z = ~coords_L[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'red', opacity = 0.7), name = "L") %>%
  add_trace(data = data_H, x = ~coords_H[, 1], y = ~coords_H[, 2], z = ~coords_H[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'green', opacity = 0.7), name = "H") %>%
  add_trace(x = ind_coords[worst_projected_indices_far_enough, 1], y = ind_coords[worst_projected_indices_far_enough, 2], z = ind_coords[worst_projected_indices_far_enough, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 4, color = 'black', opacity = 1), name = "Mal Projetés") %>%
  layout(title = "ACP: Pays L (rouge) vs H (vert) SUR LES 70%",
         scene = list(xaxis = list(title = 'PC1'),
                      yaxis = list(title = 'PC2'),
                      zaxis = list(title = 'PC3')),
         legend = list(x = 0.1, y = 0.9))

# Afficher le graphique
p
## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels

## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

Comparaison avec PCA :

library(caret)
## Loading required package: lattice
## 
## Attaching package: 'caret'
## The following object is masked from 'package:purrr':
## 
##     lift
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

# Effectuer l'ACP
pca_result <- PCA(data_train_numeric, graph = TRUE, scale.unit = TRUE)

# Exécution de l'ACP
res_pca <- acp_from_scratch(data_train_numeric, nombre_composants=2, reduire = TRUE)

ind_coords <- res_pca$nouvelles_coordonnees

data_train_numeric$Income_Inequality <- data_train$Income_Inequality

# Créer une fonction pour générer un graphique Plotly
generate_plot <- function(x, y, title) {
  data_L <- ind_coords[data_train_numeric$Income_Inequality == "L", ]
  data_H <- ind_coords[data_train_numeric$Income_Inequality == "H", ]
  
  plot_ly() %>%
    add_trace(x = ~-data_L[, x], y = ~-data_L[, y], 
              type = 'scatter', mode = 'markers',
              marker = list(size = 6, color = 'red'), name = "L") %>%
    add_trace(x = ~-data_H[, x], y = ~-data_H[, y], 
              type = 'scatter', mode = 'markers',
              marker = list(size = 6, color = 'green'), name = "H") %>%
    layout(title = title,
           xaxis = list(title = paste("PC", x)),
           yaxis = list(title = paste("PC", y)))
}

# Créer et afficher les graphiques pour chaque paire de composantes
generate_plot(1, 2, "ACP: PC1 vs PC2")
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

####################################################8

Ce sont les mêmes resultats, ouf !

####################################################8

PLOT 3D des résultats de PCA

library(plotly)
library(FactoMineR)
library(dplyr)
library(stringr)

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric <- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

# Exécution de l'ACP
pca_result <- PCA(data_train_numeric, graph = FALSE, scale.unit = TRUE)

# Coordonnées des individus
ind_coords <- pca_result$ind$coord

# Création de vecteurs de couleurs pour les individus en fonction de 'Income Inequality'
data_train_numeric$Income_Inequality <- data_train$Income_Inequality

# Séparation des individus en deux groupes
ind_L <- ind_coords[data_train_numeric$Income_Inequality == "L", ]
ind_H <- ind_coords[data_train_numeric$Income_Inequality == "H", ]

# Création du graphique Plotly
p <- plot_ly() %>%
  add_trace(x = -ind_L[, 1], y = -ind_L[, 2], z = -ind_L[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 3, color = 'red', opacity = 0.7), name = "L") %>%
  add_trace(x = -ind_H[, 1], y = -ind_H[, 2], z = -ind_H[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 3, color = 'green', opacity = 0.7), name = "H") %>%
  layout(title = "ACP avec fonction PCA pays H vs L",
         scene = list(xaxis = list(title = 'PC1'),
                      yaxis = list(title = 'PC2'),
                      zaxis = list(title = 'PC3')),
         legend = list(x = 0.1, y = 0.9))

# Afficher le graphique
p
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric <- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL
library(plotly)
library(FactoMineR)
library(dplyr)
library(stringr)

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric <- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

# Exécution de l'ACP
res_pca <- acp_from_scratch(as.matrix(data_train_numeric), nombre_composants=NULL, reduire = TRUE)

# Coordonnées des individus obtenues à partir de votre fonction ACP
ind_coords <- res_pca$nouvelles_coordonnees

# On va visualiser la classe "Income_Inequality"
data_train_numeric$Income_Inequality <- data_train$Income_Inequality

# Filtration des données pour chaque catégorie
data_L <- data_train_numeric[data_train_numeric$Income_Inequality == "L", ]
data_H <- data_train_numeric[data_train_numeric$Income_Inequality == "H", ]

# Coordonnées pour chaque catégorie
coords_L <- ind_coords[data_train_numeric$Income_Inequality == "L", ]
coords_H <- ind_coords[data_train_numeric$Income_Inequality == "H", ]

# Création du graphique Plotly
plot_ly() %>%
  add_trace(data = data_L, x = ~coords_L[, 1], y = ~coords_L[, 2], z = ~coords_L[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'red', opacity = 0.7), name = "L") %>%
  add_trace(data = data_H, x = ~coords_H[, 1], y = ~coords_H[, 2], z = ~coords_H[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'green', opacity = 0.7), name = "H") %>%
  layout(title = "ACP FROM SCRATCH: Pays L (rouge) vs H (vert) sur les 70%",
         scene = list(xaxis = list(title = 'PC1'),
                      yaxis = list(title = 'PC2'),
                      zaxis = list(title = 'PC3')),
         legend = list(x = 0.1, y = 0.9))
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL

####################################################8

On obtient les mêmes résultats avec FROM SCRATCH et avec PCA, encore une fois

####################################################8

—- QUESTION 4C —-

Projection des 30% par ACP avec uniquement les VP des 70%

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL


#Appliquer l'ACP sur les données d'entraînement
k <- 19
res_pca_train <- acp_from_scratch(data_train_numeric, nombre_composants=k, reduire=TRUE)

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL

# Normaliser les données de test
data_test_normalise <- centrer_reduire(as.matrix(data_test_numeric))

# Projection des données de test sur les vecteurs propres de l'ACP
projections_test <- data_test_normalise %*% res_pca_train$vecteurs_propres[, 1:k]

#print(dim(data_test_normalise))
#print(dim(res_pca_train$vecteurs_propres[, 1:k]))

# Coordonnées des individus obtenues à partir de votre fonction ACP
ind_coords <- projections_test

# On va visualiser la classe "Income_Inequality"
data_test_numeric$Income_Inequality <- data_test$Income_Inequality

# Créer une fonction pour générer un graphique Plotly
generate_plot <- function(x, y, title) {
  data_L <- ind_coords[data_test_numeric$Income_Inequality == "L", ]
  data_H <- ind_coords[data_test_numeric$Income_Inequality == "H", ]
  
  plot_ly() %>%
    add_trace(x = ~data_L[, x], y = ~data_L[, y], 
              type = 'scatter', mode = 'markers',
              marker = list(size = 6, color = 'red'), name = "L") %>%
    add_trace(x = ~data_H[, x], y = ~data_H[, y], 
              type = 'scatter', mode = 'markers',
              marker = list(size = 6, color = 'green'), name = "H") %>%
    layout(title = title,
           xaxis = list(title = paste("PC", x)),
           yaxis = list(title = paste("PC", y)))
}

# Créer et afficher les graphiques pour chaque paire de composantes
generate_plot(1, 2, "ACP 30% : PC1 vs PC2")
generate_plot(1, 3, "ACP 30% : PC1 vs PC3")
generate_plot(1, 4, "ACP 30% : PC1 vs PC4")
generate_plot(1, 5, "ACP 30% : PC1 vs PC5")
generate_plot(2, 3, "ACP 30% : PC2 vs PC3")
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL

####################################################8

On retrouve des résultats trés proche de ceux des 70%: on à nos 2 classes H et F. Et en augmentant les composante, la distinction des 2 classes est aussi moindre.

On a aussi, il semble, nos 2 petits groups de pays H et L en bas à droite, qui reviennent.

####################################################8

Qualité des projections des individus Qik=5 selon la formule du cours NORMEE CAS 30% TEST DATA

k <- 5  # Nombre de composantes à considérer

# Initialisation d'un data frame pour stocker les résultats
results_df <- data.frame(
moyenne = numeric(),
ecart_type = numeric(),
medianne= numeric()
)

# PAS ACP ICI
nouvelles_coordonnees <- projections_test

# Initialisation du vecteur pour stocker la qualité des projections
qualite_projections_k5 <- numeric(nrow(nouvelles_coordonnees))

# Boucle sur chaque individu
for (i in 1:nrow(nouvelles_coordonnees)) {
    Ci <- nouvelles_coordonnees[i, ]
    
    # Calcul de la qualité de la projection pour l'individu i
    # Vérification que le dénominateur n'est pas nul
    if(sum(Ci^2) != 0) {
        qualite_projections_k5[i] <- sum(Ci[1:k]^2) / sum(Ci^2)
        #print(sum(Ci[1:k]^2))
        #print(sum(Ci^2))
    } else {
        qualite_projections_k5[i] <- NA  # Assigner NA si le dénominateur nul
    }
}


# Calcul des moyennes et des écarts types 
moyenne <- mean(qualite_projections_k5)
ecart_type <- sd(qualite_projections_k5)
medianne <- median(qualite_projections_k5)

# Ajout des résultats dans le data frame 
results_df <- rbind(results_df, data.frame(moyenne, ecart_type, medianne))


# Renommage des colonnes pour la clarté
names(results_df) <- c("Moyenne", "Ecart_type", "Médianne")

# Affichage du data frame avec pander pour un joli formatage
cat("qualite projections 70%")
## qualite projections 70%
pander(results_df)
Moyenne Ecart_type Médianne
0.752 0.176 0.795 
cat("qualite projections des 30% sur les 70%")
## qualite projections des 30% sur les 70%
pander(results_df_1)
Moyenne Ecart_type Médianne
0.748 0.195 0.782

Qualité des projections des individus Qik selon la formule du cours NORMEE SELON LE NBR K CHOISIT des 30% projetés

# PAS de ACP ici
nouvelles_coordonnees <- projections_test

# Initialisation d'un data frame pour stocker les résultats
results_df <- data.frame(
  K = integer(),
  Moyenne = numeric(),
  Ecart_type = numeric(),
  Médianne = numeric()
)

# Initialisation du vecteur pour stocker la qualité des projections
qualite_projections <- numeric(nrow(nouvelles_coordonnees))

# Calculer les statistiques pour chaque valeur de k
for (k in 1:19) {
  qualite_projections <- apply(nouvelles_coordonnees, 1, function(Ci) {
    if(sum(Ci^2) != 0) {
      return(sum(Ci[1:k]^2) / sum(Ci^2))
    } else {
      return(NA)  # Assigner NA si le dénominateur est nul
    }
  })
  
  # Calcul des moyennes, écarts types et médianes
  moyenne <- mean(qualite_projections, na.rm = TRUE)
  ecart_type <- sd(qualite_projections, na.rm = TRUE)
  medianne <- median(qualite_projections, na.rm = TRUE)

  # Ajout des résultats dans le data frame 
  results_df <- rbind(results_df, data.frame(K = k, Moyenne = moyenne, Ecart_type = ecart_type, Médianne = medianne))
}


# Transformer le dataframe pour la heatmap
results_long <- melt(results_df, id.vars = 'K')

# Création de la heatmap avec valeurs dans les cases
ggplot(results_long, aes(x = variable, y = factor(K), fill = value)) +
  geom_tile() + 
  geom_text(aes(label = round(value, 2)), size = 3) +  # Ajouter des textes
  scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = median(results_df$Moyenne, na.rm = TRUE)) +
  labs(x = "Qualité de projection selon K", y = "Nombre de Composantes (K)", fill = "Valeur") +
  theme_minimal()

# Tri des qualités de projection pour trouver les plus faibles (mal projetés)
# et récupération des indices (numéros des individus)
sorted_indices <- order(qualite_projections_k5, decreasing = FALSE)

# Sélection des 10 individus les plus mal projetés
worst_projected_indices <- head(sorted_indices, 10)
worst_projected_qualities <- qualite_projections_k5[worst_projected_indices]

# Création d'un data frame pour les afficher
worst_projected_df <- data.frame(
  Numero_Individu = worst_projected_indices,
  Qualite_proj = worst_projected_qualities
)

# Affichage du data frame avec pander
pander(worst_projected_df)
Numero_Individu Qualite_proj
68 0.101
53 0.112
11 0.129
203 0.225
214 0.244
187 0.274
218 0.343
85 0.347
60 0.369
108 0.376 
ind_coords <- projections_test[,1:3]
#print(dim(ind_coords))

# On va visualiser la classe "Income_Inequality"
data_test_numeric$Income_Inequality <- data_test$Income_Inequality

# Filtration des données pour chaque catégorie
data_L <- data_test_numeric[data_test_numeric$Income_Inequality == "L", ]
data_H <- data_test_numeric[data_test_numeric$Income_Inequality == "H", ]

# Coordonnées pour chaque catégorie
coords_L <- ind_coords[data_test_numeric$Income_Inequality == "L", ]
coords_H <- ind_coords[data_test_numeric$Income_Inequality == "H", ]

# Identifiez les 10 individus les plus mal projetés (voir code précédent)

# Création du graphique Plotly
p <- plot_ly() %>%
  add_trace(data = data_L, x = ~coords_L[, 1], y = ~coords_L[, 2], z = ~coords_L[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'red', opacity = 0.7), name = "L") %>%
  add_trace(data = data_H, x = ~coords_H[, 1], y = ~coords_H[, 2], z = ~coords_H[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'green', opacity = 0.7), name = "H") %>%
  # Ajouter les 10 individus les plus mal projetés en noir
  add_trace(x = ind_coords[worst_projected_indices, 1], y = ind_coords[worst_projected_indices, 2], z = ind_coords[worst_projected_indices, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 4, color = 'black', opacity = 1), name = "Mal Projetés") %>%
  layout(title = "ACP: Pays L (rouge) vs H (vert) SUR LES 30%",
         scene = list(xaxis = list(title = 'PC1'),
                      yaxis = list(title = 'PC2'),
                      zaxis = list(title = 'PC3')),
         legend = list(x = 0.1, y = 0.9))

# Afficher le graphique
p
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL

####################################################8

Même commentaire que le cas 70% :

Ce critère n’a pas de signification pour les individus proches de l’origine. Quand on détecte un individu pour lequel le cosinus carré est faible, on doit tenir compte de sa distance à l’origine avant d’indiquer qu’il est mal représenté, donc ici ils sont proches et on ne peut rien dire pour les 10 points.

####################################################8

ind_coords <- projections_test[,1:3]
#print(dim(ind_coords))

# On va visualiser la classe "Income_Inequality"
data_test_numeric$Income_Inequality <- data_test$Income_Inequality

# Filtration des données pour chaque catégorie
data_L <- data_test_numeric[data_test_numeric$Income_Inequality == "L", ]
data_H <- data_test_numeric[data_test_numeric$Income_Inequality == "H", ]

# Coordonnées pour chaque catégorie
coords_L <- ind_coords[data_test_numeric$Income_Inequality == "L", ]
coords_H <- ind_coords[data_test_numeric$Income_Inequality == "H", ]

# Calculer la distance de chaque point à l'origine
distances <- apply(ind_coords, 1, function(x) sqrt(sum(x^2)))

# Filtrer les individus éloignés de l'origine
indices_far_enough <- which(distances > distance)

# Tri des qualités de projection pour trouver les plus faibles parmi ceux éloignés
sorted_indices_far_enough <- order(qualite_projections_k5[indices_far_enough], decreasing = FALSE)

# Sélection des 10 individus les plus mal projetés parmi ceux éloignés
worst_projected_indices_far_enough <- indices_far_enough[head(sorted_indices_far_enough, 10)]

# Ajout d'une sphère
theta <- seq(0, 2*pi, length.out = 100)
phi <- seq(0, pi, length.out = 100)
sphere_x <- distance * cos(theta) %o% sin(phi)
sphere_y <- distance * sin(theta) %o% sin(phi)
sphere_z <- distance * rep(1, length(theta)) %o% cos(phi)

# Création du graphique Plotly avec les individus mal projetés en noir
p <- plot_ly() %>%
  # Ajouter la sphère en premier
  add_trace(x = c(sphere_x), y = c(sphere_y), z = c(sphere_z), 
            type = 'mesh3d', 
            opacity = 0.01,
            color = 'blue', 
            name = 'Sphère de Distance') %>%
  add_trace(data = data_L, x = ~coords_L[, 1], y = ~coords_L[, 2], z = ~coords_L[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'red', opacity = 0.7), name = "L") %>%
  add_trace(data = data_H, x = ~coords_H[, 1], y = ~coords_H[, 2], z = ~coords_H[, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 2, color = 'green', opacity = 0.7), name = "H") %>%
  add_trace(x = ind_coords[worst_projected_indices_far_enough, 1], y = ind_coords[worst_projected_indices_far_enough, 2], z = ind_coords[worst_projected_indices_far_enough, 3],
            type = 'scatter3d', mode = 'markers',
            marker = list(size = 4, color = 'black', opacity = 1), name = "Mal Projetés") %>%
  layout(title = "ACP: Pays L (rouge) vs H (vert) SUR LES 30%",
         scene = list(xaxis = list(title = 'PC1'),
                      yaxis = list(title = 'PC2'),
                      zaxis = list(title = 'PC3')),
         legend = list(x = 0.1, y = 0.9))

# Afficher le graphique
p
## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels

## Warning in RColorBrewer::brewer.pal(N, "Set2"): minimal value for n is 3, returning requested palette with 3 different levels
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL

——- QUESTION 6 (et aussi 5B sur la projection plus loin quandon arrive au cas DATA TEST) ——-

AFD

data_train2 <- as.data.frame(data_train)
#head(data_train2)
data_train2$Country <- NULL
data_train2$Year <- NULL
data_train2 <- as.data.frame(data_train2)
#head(data_train2)
data_train_numeric2 <- data_train2[,2:20]
data_train_numeric2 <- as.data.frame(centrer_reduire(as.matrix(data_train_numeric2)))
y <- data_train2[,1]


res <- AFD(data_train_numeric2, y)

##############################################################8

2 groupes donc au plus 1 axes discriminant

##############################################################8

res$nk #nombre de pays L et H dans chaque groupes
## [1] 371 238

##############################################################8

Il y a 371 pays H contre 238 pays L pour les data TRAIN

##############################################################8

Afficher les matrices V, B et W obtenues par calcul dans le code du début:

head(res$V)
##               Eco1       Eco2       Eco3    Energy1    Energy2    Energy3
## Eco1     1.0000000  0.9080034  0.9425410  0.5650859  0.4133341 -0.2087527
## Eco2     0.9080034  1.0000000  0.9508460  0.7311746  0.6011981 -0.1682211
## Eco3     0.9425410  0.9508460  1.0000000  0.6626844  0.5441478 -0.1794532
## Energy1  0.5650859  0.7311746  0.6626844  1.0000000  0.7579802 -0.1311375
## Energy2  0.4133341  0.6011981  0.5441478  0.7579802  1.0000000 -0.1461514
## Energy3 -0.2087527 -0.1682211 -0.1794532 -0.1311375 -0.1461514  1.0000000
##            Health1    Health2     Finan1      Finan2     Finan3      Finan4
## Eco1    -0.4849284 -0.4716964  0.7920762 0.200062501  0.8691177  0.60902181
## Eco2    -0.6631237 -0.6673156  0.8646535 0.260548519  0.8145357  0.75823816
## Eco3    -0.6062391 -0.6045611  0.8051368 0.217012383  0.7661956  0.67795475
## Energy1 -0.7805115 -0.8073815  0.6917169 0.267382265  0.4561979  0.77038120
## Energy2 -0.8868724 -0.9093691  0.5727215 0.207286300  0.3334483  0.80503547
## Energy3  0.1524014  0.1096786 -0.1924987 0.009182499 -0.1888878 -0.04649566
##             Finan5 Governance    Poverty         Env     Other1     Other2
## Eco1     0.6616285  0.8462453 -0.6980418  0.61847191 -0.1989375 -0.3221600
## Eco2     0.7483751  0.9009978 -0.8415264  0.67745684 -0.3003948 -0.3674652
## Eco3     0.7092772  0.8512549 -0.8018838  0.76947906 -0.1775236 -0.3572246
## Energy1  0.5888335  0.7355785 -0.7775939  0.41999488 -0.5040555 -0.5037842
## Energy2  0.5846327  0.5732584 -0.8046352  0.44240140 -0.3921128 -0.2637850
## Energy3 -0.1941946 -0.2352324  0.1242897 -0.01740092  0.1596215  0.2009650
##              Other3
## Eco1     0.28906391
## Eco2     0.26475296
## Eco3     0.20587837
## Energy1  0.26117060
## Energy2  0.22378564
## Energy3 -0.03769449
head(res$W)
##               Eco1       Eco2       Eco3    Energy1    Energy2    Energy3
## Eco1     0.6928377  0.5946579  0.6332922  0.3083264  0.1917161 -0.2476118
## Eco2     0.5946579  0.6803468  0.6353720  0.4692466  0.3751189 -0.2078625
## Eco3     0.6332922  0.6353720  0.6886506  0.4041809  0.3210243 -0.2185762
## Energy1  0.3083264  0.4692466  0.4041809  0.7853728  0.5727279 -0.1636202
## Energy2  0.1917161  0.3751189  0.3210243  0.5727279  0.8401023 -0.1741883
## Energy3 -0.2476118 -0.2078625 -0.2185762 -0.1636202 -0.1741883  0.9950839
##            Health1    Health2     Finan1      Finan2     Finan3      Finan4
## Eco1    -0.2542488 -0.2397177  0.5270885 0.189240067  0.6256729  0.25577748
## Eco2    -0.4278005 -0.4306672  0.5943316 0.249508228  0.5661903  0.39788298
## Eco3    -0.3739925 -0.3710067  0.5383491 0.206116435  0.5210971  0.32231091
## Energy1 -0.5876846 -0.6134686  0.4702115 0.258335706  0.2527005  0.47510145
## Energy2 -0.7204367 -0.7419961  0.3815325 0.199477900  0.1578025  0.55016921
## Energy3  0.1815847  0.1390263 -0.2260223 0.007813352 -0.2196860 -0.09118463
##             Finan5 Governance    Poverty         Env      Other1     Other2
## Eco1     0.4619193  0.5481098 -0.3379302  0.39792353 -0.10533499 -0.2003847
## Eco2     0.5446458  0.5968608 -0.4741657  0.45246880 -0.20490807 -0.2432386
## Eco3     0.5082115  0.5510943 -0.4393260  0.54743255 -0.08328524 -0.2346221
## Energy1  0.4218951  0.4863647 -0.4765738  0.23563673 -0.42581243 -0.4019913
## Energy2  0.4405423  0.3581533 -0.5448142  0.28327545 -0.32457848 -0.1759239
## Energy3 -0.2194598 -0.2729495  0.1698475 -0.04530252  0.17146318  0.2163709
##              Other3
## Eco1     0.14148430
## Eco2     0.11420256
## Eco3     0.05729629
## Energy1  0.13780763
## Energy2  0.11730676
## Energy3 -0.05636480
head(res$B)
##               Eco1       Eco2      Eco3    Energy1    Energy2     Energy3
## Eco1    0.30716230 0.31334551 0.3092488 0.25675941 0.22161801 0.038859134
## Eco2    0.31334551 0.31965319 0.3154740 0.26192801 0.22607921 0.039641373
## Eco3    0.30924878 0.31547399 0.3113494 0.25850352 0.22312341 0.039123095
## Energy1 0.25675941 0.26192801 0.2585035 0.21462724 0.18525226 0.032482660
## Energy2 0.22161801 0.22607921 0.2231234 0.18525226 0.15989769 0.028036918
## Energy3 0.03885913 0.03964137 0.0391231 0.03248266 0.02803692 0.004916073
##            Health1     Health2     Finan1      Finan2     Finan3     Finan4
## Eco1    -0.2306796 -0.23197867 0.26498769 0.010822434 0.24344480 0.35324433
## Eco2    -0.2353232 -0.23664842 0.27032192 0.011040291 0.24834537 0.36035518
## Eco3    -0.2322465 -0.23355444 0.26678769 0.010895949 0.24509846 0.35564384
## Energy1 -0.1928269 -0.19391282 0.22150532 0.009046559 0.20349744 0.29527975
## Energy2 -0.1664356 -0.16737292 0.19118897 0.007808401 0.17564575 0.25486626
## Energy3 -0.0291833 -0.02934765 0.03352362 0.001369147 0.03079823 0.04468898
##             Finan5 Governance     Poverty       Env      Other1      Other2
## Eco1    0.19970913 0.29813546 -0.36011162 0.2205484 -0.09360254 -0.12177530
## Eco2    0.20372930 0.30413696 -0.36736071 0.2249880 -0.09548677 -0.12422665
## Eco3    0.20106571 0.30016062 -0.36255778 0.2220465 -0.09423836 -0.12260249
## Energy1 0.16693845 0.24921381 -0.30102018 0.1843581 -0.07824311 -0.10179295
## Energy2 0.14409041 0.21510514 -0.25982102 0.1591260 -0.06753436 -0.08786104
## Energy3 0.02526522 0.03771715 -0.04555776 0.0279016 -0.01184167 -0.01540581
##             Other3
## Eco1    0.14757960
## Eco2    0.15055040
## Eco3    0.14858208
## Energy1 0.12336297
## Energy2 0.10647888
## Energy3 0.01867031

ACP sur matrice des correlations

library(FactoMineR)

pca<-PCA(data_train2,quali.sup=1,graph=FALSE) 
plot(pca,habillage=1)

plot(pca,choix="var")

plot(pca,invisible="ind")

##############################################################8

On voit bien que nos 2 classes sont sur quasi la même dimension et de signe contraire, on a 1 seul axe discriminant en AFD d’ailleurs. On peut en déduire que pour notre AFD on aura les classes L d’un côté et celles H de l’autre sur une droite 1D car 1 seul axe car 2 classes.

##############################################################8

La valeur propres qu’on va utiliser avec son vecteur propre associée en AFD:

res$eig # valeur propre
## [1] 0.6141451
res$U # vecteur propre
##         Eco1         Eco2         Eco3      Energy1      Energy2      Energy3 
## -0.473554214  0.802178114 -0.159233798  0.048485851  0.366497876 -0.164142115 
##      Health1      Health2       Finan1       Finan2       Finan3       Finan4 
## -0.227342897 -0.220339393 -0.203279001  0.085541954  0.246794613 -0.760776222 
##       Finan5   Governance      Poverty          Env       Other1       Other2 
##  0.001066804 -0.204325150  1.034689307  0.112008863 -0.022390423 -0.012909606 
##       Other3 
## -0.103691643

On utilise notre fonction plotAFD pour representer nos 2 classes TRAIN

plotAFD(res)

1er facteurs discriminants

u1<-res$U #1er facteurs discriminants
print(u1)
##         Eco1         Eco2         Eco3      Energy1      Energy2      Energy3 
## -0.473554214  0.802178114 -0.159233798  0.048485851  0.366497876 -0.164142115 
##      Health1      Health2       Finan1       Finan2       Finan3       Finan4 
## -0.227342897 -0.220339393 -0.203279001  0.085541954  0.246794613 -0.760776222 
##       Finan5   Governance      Poverty          Env       Other1       Other2 
##  0.001066804 -0.204325150  1.034689307  0.112008863 -0.022390423 -0.012909606 
##       Other3 
## -0.103691643

Coordonnéees (1D) des centres des groupes H et L en AFD /// DATA TRAIN

gL <- res$gk[,1]
gH <- res$gk[,2]
t(gL)%*%u1
##           [,1]
## [1,] 0.6276781
t(gH)%*%u1
##            [,1]
## [1,] -0.9784395

Vérification des coordonnées des 2 centres: On obtient bien les mêmes resultats !

S1 <- res$S[,1] # VD
#print(S1)
Sk <- split(S1,y)
lapply(Sk,mean)
## $H
## [1] 0.6276781
## 
## $L
## [1] -0.9784395

On calcul le seuil entre les 2 classes :

on est en 1D donc le seuil entre L et H c’est simplement (gH + gL)/2

S1 <- res$S[,1]
#print(dim(S1))
#print(S1)
Sk <- split(S1,y)
#print(dim(Sk))
#print(Sk)
smoy <- lapply(Sk,mean)
seuil1 <- (as.numeric(smoy$L)+as.numeric(smoy$H))/2 #seuil entre L et H
print(seuil1)
## [1] -0.1753807

Prédition de nos DATA TRAIN avec ce seuil :

predict <- cut(S1,breaks=c(-4,seuil1,4),labels=c("L","H"))
y <- factor(y, levels = c("L", "H"))
table(y,predict)
##    predict
## y     L   H
##   L 212  26
##   H  35 336

On verra l’accuracy plus tard avec les autres accuracy

—- QUESTION 6 et 5B —-

Préparation des données pur le cas DATA TEST

# Préparation des données TRAIN
data_train2 <- as.data.frame(data_train)
#head(data_train2)
data_train2$Country <- NULL
data_train2$Year <- NULL
data_train2 <- as.data.frame(data_train2)
#head(data_train2)
data_train_numeric2 <- data_train2[,2:20]
data_train_numeric2 <- as.data.frame(centrer_reduire(as.matrix(data_train_numeric2)))
y_train2 <- data_train2[,1]


# Préparation des données TEST
data_test2 <- as.data.frame(data_test)
#head(data_test2)
data_test2$Country <- NULL
data_test2$Year <- NULL
data_test2 <- as.data.frame(data_test2)
#head(data_test2)
data_test_numeric2 <- data_test2[,2:20]
data_test_numeric2 <- as.data.frame(centrer_reduire(as.matrix(data_test_numeric2)))
y_test2 <- data_test2[,1]


# DATA SPLIT
Xtrain <- data_train_numeric2
ytrain <- y_train2
Xtest <- data_test_numeric2
ytest <- y_test2

Construction de la règle seuil sur l’ensemble d’apprentissage (C’est le même qui précédemment en fait!)

res <- AFD(Xtrain,ytrain)
S <- res$S
#print(S)
Sk <- split(S,ytrain)
smoy <- lapply(Sk,mean)
#print(smoy)
seuil1 <- (as.numeric(smoy$L)+as.numeric(smoy$H))/2
print(seuil1)
## [1] -0.1753807

Application du seuil sur DATA TEST

Calcul de la matrice de confusion sur DATA TEST

Stest <- as.matrix(Xtest)%*%u1 

#print(Stest)

predict <- cut(Stest,breaks=c(-4,seuil1,4),labels=c("L","H"))

#print(predict)

#levels(ytest)
#levels(predict)

ytest <- factor(ytest, levels = c("L", "H"))
conf_matrix <- table(ytest, predict)
pander(conf_matrix)
  L H
L 81 9
H 32 139

—– QUESTION 5B ——-

On cherche les 20 individus les plus mal projetées : 10 pour la classe L et 10 pour la classe H:

Stest <- as.matrix(Xtest)%*%u1 
mean_H <- mean(Stest[ytest == "H"])
mean_L <- mean(Stest[ytest == "L"])
deviation_H <- abs(Stest[ytest == "H"] - mean_H)
deviation_L <- abs(Stest[ytest == "L"] - mean_L)

# Identifier les 10 plus grands écarts pour chaque groupe
top10_deviation_indices_H <- order(deviation_H, decreasing = TRUE)[1:10]
top10_deviation_indices_L <- order(deviation_L, decreasing = TRUE)[1:10]


global_indices_H <- which(ytest == "H")[top10_deviation_indices_H]
global_indices_L <- which(ytest == "L")[top10_deviation_indices_L]

# Créer un data.frame pour afficher les résultats
results_H <- data.frame(
  Index = global_indices_H,
  Group = "H",
  Value = Stest[global_indices_H],
  Deviation = deviation_H[top10_deviation_indices_H]
)

results_L <- data.frame(
  Index = global_indices_L,
  Group = "L",
  Value = Stest[global_indices_L],
  Deviation = deviation_L[top10_deviation_indices_L]
)

# Combiner les résultats pour H et L
combined_results <- rbind(results_H, results_L)

# Affichage du tableau
pander(combined_results)
Index Group Value Deviation
197 H -1.27 1.83
241 H -1.26 1.82
186 H -1.25 1.81
120 H -1.07 1.63
10 H -1.02 1.58
101 H 2.04 1.48
18 H 2 1.44
245 H -0.783 1.34
157 H 1.87 1.31
83 H 1.85 1.29
108 L 0.64 1.7
121 L 0.634 1.7
203 L 0.555 1.62
67 L 0.524 1.59
74 L 0.337 1.4
46 L -2.41 1.35
212 L -2.37 1.31
256 L 0.25 1.31
145 L -2.34 1.28
92 L 0.148 1.21

############################################8

On retrouve l’individus 203 seulement en comparant avec la partie précédente en ACP.

On ne retrouve pas les mêmes individus qu’en ACP globalement, c’est surement encore du à notre code qui a du mal avec les indices.

############################################8

Calcul des métriques de performance DATA TEST

accuracy <- sum(diag(conf_matrix)) / sum(conf_matrix)
sensitivity <- diag(conf_matrix) / rowSums(conf_matrix)
specificity <- diag(conf_matrix) / colSums(conf_matrix)
precision <- diag(conf_matrix) / colSums(conf_matrix)
recall <- sensitivity
F1 <- 2 * (precision * recall) / (precision + recall)

# Data frame pour stocker les métriques avec descriptions
performance_metrics <- data.frame(
  Metric = c("Accuracy", "Sensitivity", "Specificity", "Precision", "Recall", "F1 Score"),
  Value = c(accuracy, mean(sensitivity), mean(specificity), mean(precision), mean(recall), mean(F1)),
  Description = c(
    "Proportion de toutes les prédictions correctes",
    "Proportion de positifs réels correctement identifiés",
    "Proportion de négatifs réels correctement identifiés",
    "Proportion de prédictions positives correctes",
    "Identique à la sensibilité",
    "Moyenne harmonique de la précision et du rappel"
  )
)

pander(performance_metrics)
Metric Value Description
Accuracy 0.843 Proportion de toutes les prédictions correctes
Sensitivity 0.856 Proportion de positifs réels correctement identifiés
Specificity 0.828 Proportion de négatifs réels correctement identifiés
Precision 0.828 Proportion de prédictions positives correctes
Recall 0.856 Identique à la sensibilité
F1 Score 0.835 Moyenne harmonique de la précision et du rappel

Comparaison entre notre fonction AFD et la fonction LDA de R:

Facteur discriminant

library(MASS)
library(readr)
Xtrain2 <- Xtrain
Xtrain2$Income_Inequality <- data_train$Income_Inequality


res3 <- lda(Income_Inequality~., Xtrain2)
u1 <- -res$U

# Extraire les données
scaling <- res3$scaling
u1_data <- u1

relative_difference <- ifelse(scaling != 0, (u1_data - scaling) / scaling, NA)


comparison_df <- data.frame(Scaling = scaling, U1 = u1_data, RelativeDifference = relative_difference)
colnames(comparison_df) <- c("Facteur Discriminant", "U1", "Différence Relative")


print(comparison_df)
##            Facteur Discriminant           U1 Différence Relative
## Eco1                0.761102907  0.473554214          -0.3778053
## Eco2               -1.289271801 -0.802178114          -0.3778053
## Eco3                0.255922770  0.159233798          -0.3778053
## Energy1            -0.077927132 -0.048485851          -0.3778053
## Energy2            -0.589040475 -0.366497876          -0.3778053
## Energy3             0.263811486  0.164142115          -0.3778053
## Health1             0.365388660  0.227342897          -0.3778053
## Health2             0.354132531  0.220339393          -0.3778053
## Finan1              0.326712832  0.203279001          -0.3778053
## Finan2             -0.137484217 -0.085541954          -0.3778053
## Finan3             -0.396651728 -0.246794613          -0.3778053
## Finan4              1.222730106  0.760776222          -0.3778053
## Finan5             -0.001714582 -0.001066804          -0.3778053
## Governance          0.328394218  0.204325150          -0.3778053
## Poverty            -1.662967019 -1.034689307          -0.3778053
## Env                -0.180022199 -0.112008863          -0.3778053
## Other1              0.035986198  0.022390423          -0.3778053
## Other2              0.020748498  0.012909606          -0.3778053
## Other3              0.166654649  0.103691643          -0.3778053

#print((0.761102907 - 0.473554214)/0.761102907) #print((1.289271801 - 0.802178114)/1.289271801)

##############################################################8

En comparant le facteur disciminant, on remarque que les signes sont tous bon. Et que les coefficient sont les même a leur différence relative prés. Surement la manière dont code et interprete la fonction LDA ses coefficients.

##############################################################8

Histogramme : histogramme “empilé” pour les valeurs de la fonction discriminante DATA TRAIN

data_train2 <- as.data.frame(data_train)
#head(data_train2)
data_train2$Country <- NULL
data_train2$Year <- NULL
data_train2 <- as.data.frame(data_train2)
#head(data_train2)
data_train_numeric2 <- data_train2[,2:20]
data_train_numeric2 <- as.data.frame(centrer_reduire(as.matrix(data_train_numeric2)))
y_train2 <- data_train2[,1]

data_test2 <- as.data.frame(data_test)
#head(data_test2)
data_test2$Country <- NULL
data_test2$Year <- NULL
data_test2 <- as.data.frame(data_test2)
#head(data_test2)
data_test_numeric2 <- data_test2[,2:20]
data_test_numeric2 <- as.data.frame(centrer_reduire(as.matrix(data_test_numeric2)))
y_test2 <- data_test2[,1]


lda_model <- lda(y_train2 ~ ., data = data_train_numeric2)
p <- predict(lda_model, data_train_numeric2)
ldahist(data = p$x[,1], g = y_train2)

Histogramme : histogramme “empilé” pour les valeurs de la fonction discriminante DATA TEST:

p <- predict(lda_model, data_test_numeric2)
ldahist(data = p$x[,1], g = y_test2)

##############################################################8

Ces histogrammes sont basés sur LD1. On observe un certain chevauchement entre les groupes H et L de pays, ce qui n’est pas vraiment bon pour un LD1, on aurait pu esperer que ce ne soit pas le cas pour bien différencier les 2 classes. C’était un résultat prévisible au vu de notre analyse ACP.

##############################################################8

Confusion matrix et accuracy – DATA TEST

p2 <- predict(lda_model, data_test_numeric2)$class
tab1 <- table(Predicted = p2, Actual = y_test2)
tab1
##          Actual
## Predicted   H   L
##         H 144  11
##         L  27  79

DATA TEST : Accuracy du modele est de 0.8544061

sum(diag(tab1))/sum(tab1)
## [1] 0.8544061

##############################################################8

Pour conclure sur l’accuracy de notre fonction AFD et de la fonction LDA : La notre donne 0.839 contre 0.854 pour celle de LDA qui est donc meileur. Mais la notre s’en rapproche fortement. Il doit nous manquer des otpimisation caar nosu avons coder seulement l’AFD “basique”.

##############################################################8

6.b

Réduction de X_train par ACP

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL


#Appliquer l'ACP sur les données d'entraînement
k <- 5
res_pca_train <- acp_from_scratch(data_train_numeric, nombre_composants=k, reduire=TRUE)

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL


data_test_normalise <- centrer_reduire(as.matrix(data_test_numeric))
# Projection des données de test sur les vecteurs propres de l'ACP
projections_test <- data_test_normalise %*% res_pca_train$vecteurs_propres[, 1:k]

ind_coords_train <- as.data.frame(res_pca_train$nouvelles_coordonnees)
ind_coords_train$target <- data_train$Income_Inequality
#print(dim(ind_coords_train))
#write.csv(ind_coords_train, "train.csv", row.names = TRUE)


ind_coords_test <- as.data.frame(projections_test)
ind_coords_test$target <- data_test$Income_Inequality
#print(dim(ind_coords_test))
#write.csv(ind_coords_test, "test.csv", row.names = TRUE)
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL


#Appliquer l'ACP sur les données d'entraînement
k <- 3
res_pca_train <- acp_from_scratch(data_train_numeric, nombre_composants=k, reduire=TRUE)

#print(res_pca_train$nouvelles_coordonnees)

# Sélectionner uniquement les colonnes numériques pour l'ACP
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL

# Normaliser les données de test
data_test_normalise <- centrer_reduire(as.matrix(data_test_numeric))

# Projection des données de test sur les vecteurs propres de l'ACP
projections_test <- data_test_normalise %*% res_pca_train$vecteurs_propres[, 1:k]

ind_coords_train <- as.data.frame(res_pca_train$nouvelles_coordonnees)
#print(dim(ind_coords_train))
#print(res_pca_train$ind$coord)
ind_coords_test <- as.data.frame(projections_test)
colnames(ind_coords_train) <- c('V1', 'V2', 'V3')
colnames(ind_coords_test) <- c('V1', 'V2', 'V3')




fig <- plot_ly() %>%
    add_trace(data = ind_coords_train, x = ~V1, y = ~V2, z = ~V3, 
              type = 'scatter3d', mode = 'markers',
              marker = list(size = 2, color = 'blue'), name = 'Train Data') %>%
    add_trace(data = ind_coords_test, x = ~V1, y = ~V2, z = ~V3, 
              type = 'scatter3d', mode = 'markers',
              marker = list(size = 2, color = 'red'), name = 'Test Data') %>%
    layout(title = '3D  Plot PCA Coordinates reduites ',
           scene = list(xaxis = list(title = 'PC1'),
                        yaxis = list(title = 'PC2'),
                        zaxis = list(title = 'PC3')))

fig
# Sélectionner uniquement les colonnes numériques pour l'ACP
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL

On refait notre AFD sur donnes réduites de 19 à 5 dimensions.

Préparation des données

data_train_numeric<- select_if(data_train, is.numeric)
data_train_numeric$Country <- NULL
data_train_numeric$Year <- NULL
k <- 5
res_pca_train <- acp_from_scratch(data_train_numeric, nombre_composants=k, reduire=TRUE)
data_test_numeric<- select_if(data_test, is.numeric)
data_test_numeric$Country <- NULL
data_test_numeric$Year <- NULL
data_test_normalise <- centrer_reduire(as.matrix(data_test_numeric))
projections_test <- data_test_normalise %*% res_pca_train$vecteurs_propres[, 1:k]
ind_coords_train <- as.data.frame(res_pca_train$nouvelles_coordonnees)
ind_coords_train$target <- data_train$Income_Inequality
ind_coords_test <- as.data.frame(projections_test)
ind_coords_test$target <- data_test$Income_Inequality





Xtrain <- ind_coords_train[,1:5]
ytrain <- ind_coords_train[,6]

Xtest <- ind_coords_test[,1:5]
ytest <- ind_coords_test[,6]

On re calcul le modele sur les données TRAIN

res1 <- AFD(Xtrain, ytrain)

On utilise notre fonction plotAFD pour representer nos 2 classes TRAIN

plotAFD(res1)

1er facteurs discriminants

u11<-res1$U #1er facteurs discriminants pour 5 dimension
pander(u11)
V1 V2 V3 V4 V5
0.287 0.0694 0.121 0.318 0.0129

Comparaison des coordonnées des centres

S1 <- res$S[,1] # VD
#print(S1)
Sk <- split(S1,y)
lapply(Sk,mean)
## $L
## [1] -0.9784395
## 
## $H
## [1] 0.6276781
S11 <- res1$S[,1] # VD
#print(S1)
Sk1 <- split(S11,ytrain)
lapply(Sk1,mean)
## $H
## [1] 0.502422
## 
## $L
## [1] -0.7831872

Les coordonnées des centres sont assez différents, on peut s’attendre à une accuracy moins bonne

Calcul du seuil pour data REDUITES en 5D

smoy <- lapply(Sk,mean)
seuil1 <- (as.numeric(smoy$L)+as.numeric(smoy$H))/2 #seuil entre L et H
print(seuil1)
## [1] -0.1753807
smoy <- lapply(Sk1,mean)
seuil11 <- (as.numeric(smoy$L)+as.numeric(smoy$H))/2 #seuil entre L et H
print(seuil11)
## [1] -0.1403826

Les 2 seuils sont tout de même assez différent pour nos données, on s’attend à une accuracy assez inférieur donc encore une fois

Calcul de la matrice de confusion sur DATA TEST

Stest1 <- as.matrix(Xtest)%*%u11

predict1 <- cut(Stest1,breaks=c(-4,seuil11,4),labels=c("L","H"))

#levels(ytest)
#levels(predict)

ytest <- factor(ytest, levels = c("L", "H"))
conf_matrix1 <- table(ytest, predict1)
pander(conf_matrix1)
  L H
L 77 13
H 44 127 
pander(conf_matrix)
  L H
L 81 9
H 32 139

Calcul des métriques de performance DATA TEST

accuracy1 <- sum(diag(conf_matrix1)) / sum(conf_matrix1)
sensitivity1 <- diag(conf_matrix1) / rowSums(conf_matrix1)
specificity1 <- diag(conf_matrix1) / colSums(conf_matrix1)
precision1 <- diag(conf_matrix1) / colSums(conf_matrix1)
recall1 <- sensitivity1
F11 <- 2 * (precision1 * recall1) / (precision1 + recall1)

# Data frame pour stocker les métriques avec descriptions
performance_metrics1 <- data.frame(
  Metric = c("Accuracy", "Sensitivity", "Specificity", "Precision", "Recall", "F1 Score"),
  Value = c(accuracy1, mean(sensitivity1), mean(specificity1), mean(precision1), mean(recall1), mean(F11)),
  Description = c(
    "Proportion de toutes les prédictions correctes",
    "Proportion de positifs réels correctement identifiés",
    "Proportion de négatifs réels correctement identifiés",
    "Proportion de prédictions positives correctes",
    "Identique à la sensibilité",
    "Moyenne harmonique de la précision et du rappel"
  )
)

pander(performance_metrics1)
Metric Value Description
Accuracy 0.782 Proportion de toutes les prédictions correctes
Sensitivity 0.799 Proportion de positifs réels correctement identifiés
Specificity 0.772 Proportion de négatifs réels correctement identifiés
Precision 0.772 Proportion de prédictions positives correctes
Recall 0.799 Identique à la sensibilité
F1 Score 0.773 Moyenne harmonique de la précision et du rappel 
pander(performance_metrics)
Metric Value Description
Accuracy 0.843 Proportion de toutes les prédictions correctes
Sensitivity 0.856 Proportion de positifs réels correctement identifiés
Specificity 0.828 Proportion de négatifs réels correctement identifiés
Precision 0.828 Proportion de prédictions positives correctes
Recall 0.856 Identique à la sensibilité
F1 Score 0.835 Moyenne harmonique de la précision et du rappel

CCL : on obtient une accuracy de 78% pour nos données réduites en 5D contre 84% pour nos données initiales en 19D.

On voit donc un avantage a la réduction des données : on perd certe de l’accuracy dans notre modèle, mais on gagne des ressources en temps et en energies monstrueuses en appliquant nos predictions à des données de dimension bien inférieurs.

Il faut donc trouver un bon compromie entre ressource et accuracy, ie entre ACP et modèle de prédiction.

Comparaison entre Foret et AFD :

Cas où l’on prend les DATA TRAIN en 19D :

AFD : 0.843 \ Matrice de confusion : \[\begin{array}{c|cc} & H & L \\ \hline H & 139 & 32 \\ L & 9 & 81 \\ \end{array}\] Arbre de décision : 0.95 \ Matrice de confusion :} \[\begin{array}{c|cc} & H & L \\ \hline H & 166 & 5 \\ L & 8 & 82 \\ \end{array}\]

Cas où l’on prend les DATA TRAIN REDUITES en 5D :

AFD : 0.782 \ Matrice de confusion : \[\begin{array}{c|cc} & H & L \\ \hline H & 127 & 44 \\ L & 13 & 77 \\ \end{array}\] Arbre de décision : 0.896 \ Matrice de confusion : \[\begin{array}{c|cc} & H & L \\ \hline H & 160 & 11 \\ L & 16 & 74 \\ \end{array}\]

En termes de prédiction, il apparaît que les modèles basés sur les arbres de décision surpassent l’Analyse Factorielle Discriminante (AFD) en termes de précision. Cette efficacité accrue peut être attribuée à la capacité de ces modèles à s’adapter finement grâce à l’ajustement des hyperparamètres. Contrairement à l’AFD, où les options d’optimisation sont relativement limitées, les modèles d’arbres de décision offrent une plus grande latitude pour affiner et ajuster les paramètres, permettant ainsi une personnalisation et une optimisation plus poussées. Cette flexibilité se traduit par une meilleure adaptation du modèle aux spécificités des données, conduisant potentiellement à une amélioration des performances de prédiction, en particulier dans des contextes où les relations entre les variables sont complexes, commes les notres.

Concernant les matrices de confusion entre les deux méthodes, une différence notable est observée au niveau des faux négatifs (situés en haut à droite dans notre cas). Avec l’AFD, il y a manifestement plus de faux négatifs que de faux positifs. En revanche, avec les arbres de décision, la répartition des faux positifs et des faux négatifs est plus équilibrée, avec peut-être une légère prédominance des faux positifs.

Ainsi, selon que l’on souhaite minimiser les faux positifs ou les faux négatifs, l’une ou l’autre approche peut être privilégiée. Cette considération est particulièrement pertinente dans le secteur médical, où la question des faux négatifs et positifs est cruciale. Cette observation conclut notre travail pratique.