Environnement et package¶

In [66]:
# activer l'environnement
using Pkg
Pkg.activate("env_julia_acp")
  Activating project at `c:\Users\ricco\Desktop\demo\env_julia_acp`
In [67]:
# liste des packages installés
Pkg.status()
Status `C:\Users\ricco\Desktop\demo\env_julia_acp\Project.toml`
  [a93c6f00] DataFrames v1.8.2
  [be5e9834] FactorRotations v0.5.2
  [7073ff75] IJulia v1.34.4
  [add582a8] MLJ v0.23.2
  [1b6a4a23] MLJMultivariateStatsInterface v0.5.4
  [23777cdb] MLJTransforms v0.1.5
  [6f286f6a] MultivariateStats v0.10.4
  [91a5bcdd] Plots v1.41.6
  [30f210dd] ScientificTypesBase v3.1.0
  [f3b207a7] StatsPlots v0.15.8
  [fdbf4ff8] XLSX v0.11.6

Importation et préparation des données¶

Importation, variables actives¶

In [68]:
# packages
import DataFrames as DFR
import XLSX

# lecture des données
df = DFR.DataFrame(XLSX.readtable("./autos_acp_livre.xlsx"))

# premières lignes
println(DFR.size(df))
(18, 10)
In [69]:
# premières lignes
DFR.first(df,5)
5×10 DataFrame
RowModeleFINITIONPRIXR_POID_PUISCYLPUISSLONGLARGPOIDSV_MAX
StringStringInt64Float64Int64Int64Int64Int64Int64Int64
1Alfasud TI2_B3057011.0127135079393161870165
2Audi 1003_TB3999013.05881588854681771110160
3Simca 13001_M2960015.44121294684241681050152
4Citroen GS Club1_M2825015.7627122259412161930151
5Fiat 1322_B3490011.27551585984391641105165
In [70]:
# description
DFR.describe(df)
10×7 DataFrame
Rowvariablemeanminmedianmaxnmissingeltype
SymbolUnion…AnyUnion…AnyInt64DataType
1ModeleAlfasud TIToyota Corolla0String
2FINITION1_M3_TB0String
3PRIX34158.62210033345.0477000Int64
4R_POID_PUIS13.18079.7247713.181818.36360Float64
5CYL1631.6711661577.526640Int64
6PUISS84.61115582.01280Int64
7LONG433.5393434.54690Int64
8LARG166.667157167.01770Int64
9POIDS1078.838151087.513700Int64
10V_MAX158.278140160.01800Int64
In [71]:
# récupération des variables actives
X = df[:,[:CYL,:PUISS,:LONG,:LARG,:POIDS,:V_MAX]]
DFR.describe(X)
6×7 DataFrame
Rowvariablemeanminmedianmaxnmissingeltype
SymbolFloat64Int64Float64Int64Int64DataType
1CYL1631.6711661577.526640Int64
2PUISS84.61115582.01280Int64
3LONG433.5393434.54690Int64
4LARG166.667157167.01770Int64
5POIDS1078.838151087.513700Int64
6V_MAX158.278140160.01800Int64

Types des variables¶

In [72]:
# types reconnus pas MLJ
import MLJ
MLJ.schema(X)
┌───────┬──────────┬───────┐
│ names │ scitypes │ types │
├───────┼──────────┼───────┤
│ CYL   │ Count    │ Int64 │
│ PUISS │ Count    │ Int64 │
│ LONG  │ Count    │ Int64 │
│ LARG  │ Count    │ Int64 │
│ POIDS │ Count    │ Int64 │
│ V_MAX │ Count    │ Int64 │
└───────┴──────────┴───────┘
In [73]:
# transformer les X count en variables continues
# sinon elles sont ignorées par les outils en aval
# utilisation du package ScientificTypesBase
import ScientificTypesBase
X = MLJ.coerce(X,ScientificTypesBase.Count => ScientificTypesBase.Continuous)

# schéma
MLJ.schema(X)
┌───────┬────────────┬─────────┐
│ names │ scitypes   │ types   │
├───────┼────────────┼─────────┤
│ CYL   │ Continuous │ Float64 │
│ PUISS │ Continuous │ Float64 │
│ LONG  │ Continuous │ Float64 │
│ LARG  │ Continuous │ Float64 │
│ POIDS │ Continuous │ Float64 │
│ V_MAX │ Continuous │ Float64 │
└───────┴────────────┴─────────┘

Standardisation pour ACP normée¶

In [74]:
# attention, l'ACP ici peut centrer automatiquement les données
# mais ne les standardise pas
# si on veut une ACP normée, on doit le faire en amont
Standardizer = @MLJ.load Standardizer pkg=MLJTransforms
┌ Info: For silent loading, specify `verbosity=0`. 
└ @ Main C:\Users\ricco\.julia\packages\MLJModels\9LbNu\src\loading.jl:159
import MLJTransforms ✔
MLJTransforms.Standardizer
In [75]:
# instancier -- ATTENTION, utilise écart-type << 1/(n-1) >> ici
# cela aura des répercussions sur les coordonnées des individus !!!
# qu'il faudra corriger plus bas
std = MLJ.machine(Standardizer(),X)

# entraîner i.e. calculer moyennes et écarts-type
MLJ.fit!(std)

# transformer
Z = MLJ.transform(std,X)

# description
DFR.describe(Z)
┌ Info: Training machine(Standardizer(features = Symbol[], …), …).
└ @ MLJBase C:\Users\ricco\.julia\packages\MLJBase\krfwA\src\machines.jl:499
6×7 DataFrame
Rowvariablemeanminmedianmaxnmissingeltype
SymbolFloat64Float64Float64Float64Int64DataType
1CYL-2.22045e-16-1.24533-0.1448582.760770Float64
2PUISS-1.4803e-16-1.45321-0.1281452.129380Float64
3LONG2.46716e-17-1.831970.04523381.60580Float64
4LARG1.77636e-15-1.81920.06273111.944660Float64
5POIDS5.67447e-16-1.926380.06327982.125960Float64
6V_MAX4.93432e-16-1.505540.1418591.789250Float64

ACP avec MultivariateStats sous Julia¶

Instanciation, entraînement¶

In [76]:
# chargement de la classe de calcul
PCA = @MLJ.load PCA pkg=MultivariateStats
import MLJMultivariateStatsInterface ✔
┌ Info: For silent loading, specify `verbosity=0`. 
└ @ Main C:\Users\ricco\.julia\packages\MLJModels\9LbNu\src\loading.jl:159
MLJMultivariateStatsInterface.PCA
In [77]:
# instancier
# mean = 0 pour dire que les variables sont déjà centrées
# par défaut, s'arrête dès que 99% de l'info est captée
acp = PCA(mean=0)
PCA(
  maxoutdim = 0, 
  method = :auto, 
  variance_ratio = 0.99, 
  mean = 0)
In [78]:
# associer l'instance aux données
mach = MLJ.machine(acp,Z)

# lancer l'entraînement
MLJ.fit!(mach)
┌ Info: Training machine(PCA(maxoutdim = 0, …), …).
└ @ MLJBase C:\Users\ricco\.julia\packages\MLJBase\krfwA\src\machines.jl:499
trained Machine; caches model-specific representations of data
  model: PCA(maxoutdim = 0, …)
  args: 
    1:	Source @966 ⏎ ScientificTypesBase.Table{AbstractVector{ScientificTypesBase.Continuous}}

Accès aux résultats¶

In [79]:
# accès aux résultats
res = MLJ.report(mach)
keys(res)
(:indim, :outdim, :tprincipalvar, :tresidualvar, :tvar, :mean, :principalvars, :loadings)
In [80]:
# dimensions en entrée et en sortie
# nombre de variables vs. nb de composantes
println("Entree = $(res.indim) ; sortie = $(res.outdim)")
Entree = 6 ; sortie = 5

Variance restituée - Scree Plot¶

In [81]:
# info captée sur les composantes
res.tprincipalvar
5.956709727324061
In [82]:
# soit en proportion sur les 5 premières composantes
# 99.28% => on est bien au-delà de 99%
100.0*(res.tprincipalvar/res.indim)
99.278495455401
In [83]:
# valeurs propres
res.principalvars
5-element Vector{Float64}:
 4.420858059608889
 0.8560622892459593
 0.37306607743997394
 0.21392208902330967
 0.09280121200592834
In [84]:
# sous forme de scree plot
import StatsPlots as STP
STP.plot(res.principalvars,
        xlabel = "Composante",
        ylabel = "Valeur propre",
        title = "Scree Plot",
        legend=false)
No description has been provided for this image
In [85]:
# en proportion de variances cumulées
explained_var = 100.0*(res.principalvars / res.indim)
cum_var = cumsum(explained_var)
println("Proportion variance cumulées : $(round.(cum_var,digits=2))")

# graphique
STP.plot(cum_var,
        xlabel = "Composante",
        ylabel = "% variance",
        title = "Variance cumulée",
        legend=false)
Proportion variance cumulées : [73.68, 87.95, 94.17, 97.73, 99.28]
No description has been provided for this image

Corrélations, contributions (CTR), cosinus^2 (COS2)¶

In [86]:
# les loadings
res.loadings
6×5 Matrix{Float64}:
 -0.893464   0.114906   0.215983    0.373615    0.0461763
 -0.886858   0.384689   0.112948   -0.165485   -0.0894812
 -0.886155  -0.381029  -0.0413102  -0.12939     0.222555
 -0.813536  -0.412736  -0.369448    0.0978545  -0.145672
 -0.905187  -0.224532   0.295865   -0.139547   -0.0927785
 -0.75471    0.573519  -0.296522   -0.0340294   0.0574706
In [87]:
# somme des carrés des loadings en colonne = valeurs propres
sum(res.loadings .^ 2,dims=1)
1×5 Matrix{Float64}:
 4.42086  0.856062  0.373066  0.213922  0.0928012
In [88]:
# les loadings sur les deux premières composantes
# qui sont aussi les corrélations en réalité
corr_fact = DFR.DataFrame(var = names(X))
corr_fact = DFR.hcat(corr_fact,DFR.DataFrame(res.loadings[:,1:2],["F1","F2"]))
corr_fact
6×3 DataFrame
RowvarF1F2
StringFloat64Float64
1CYL-0.8934640.114906
2PUISS-0.8868580.384689
3LONG-0.886155-0.381029
4LARG-0.813536-0.412736
5POIDS-0.905187-0.224532
6V_MAX-0.754710.573519

Cercle des corrélations¶

In [89]:
# coordonnées = corrélation
coord = corr_fact

# dessiner le cercle
theta = range(0, 2*pi, length=200)

STP.plot(cos.(theta), sin.(theta),
     aspect_ratio=:equal,
     label="",
     linestyle=:dash,
     xlabel="PC1",
     ylabel="PC2",
     title="Cercle des corrélations",
     xlims=(-1,+1),
     ylims=(-1,+1),
     framestyle=:box)

# lignes centrales
STP.hline!([0], linestyle = :dash, color = :gray,label="")
STP.vline!([0], linestyle = :dash, color = :gray,label="")

# placer les noms des variables
# les points sont masqués avec markersize=0
STP.scatter!(coord[:,:F1], coord[:,:F2], markersize=0,
            series_annotations=(coord[:,:var],STP.font(9,:green)),
            label="")
No description has been provided for this image

Contributions¶

In [90]:
# contributions en %
# diviser le carré des corrélations par les valeurs propres
CTR = 100.0.*((corr_fact[:,[:F1,:F2]] .^2) ./ res.principalvars[1:2]')

# arrondir à 1 décimale
CTR = round.(CTR,digits=1)

# ajouter les noms des variables
CTR = DFR.hcat(DFR.DataFrame(var=names(X)),CTR)
CTR
6×3 DataFrame
RowvarF1F2
StringFloat64Float64
1CYL18.11.5
2PUISS17.817.3
3LONG17.817.0
4LARG15.019.9
5POIDS18.55.9
6V_MAX12.938.4

Cosinus carré¶

In [91]:
# cosinus^2 en %
COS2 = 100.0.*((corr_fact[:,[:F1,:F2]] .^2))

# arrondir
COS2 = round.(COS2,digits=1)

# ajouter les noms des variables
COS2 = DFR.hcat(DFR.DataFrame(var=names(X)),COS2)
COS2
6×3 DataFrame
RowvarF1F2
StringFloat64Float64
1CYL79.81.3
2PUISS78.714.8
3LONG78.514.5
4LARG66.217.0
5POIDS81.95.0
6V_MAX57.032.9
In [92]:
# qui peuvent s'additionner d'une colonne à l'autre
# pour mesurer la qualité de représentation dans le repère
DFR.DataFrame(var=names(X),Cos2=vec(sum(Matrix(COS2[:,[:F1,:F2]]),dims=2)))
6×2 DataFrame
RowvarCos2
StringFloat64
1CYL81.1
2PUISS93.5
3LONG93.0
4LARG83.2
5POIDS86.9
6V_MAX89.9

Résultats additionnels - Représentation des individus¶

In [93]:
# paramètres supplémentaires
fp = MLJ.fitted_params(mach)
keys(fp)
(:projection,)
In [94]:
# matrice de projection
# à appliquer aux valeurs centrées et réduites
fp.projection
6×5 Matrix{Float64}:
 -0.424936   0.124191   0.353613    0.807786    0.15158
 -0.421794   0.415774   0.18492    -0.357792   -0.293735
 -0.42146   -0.411818  -0.0676339  -0.279752    0.730569
 -0.386922  -0.446087  -0.604868    0.211569   -0.47819
 -0.430512  -0.242676   0.484396   -0.301711   -0.304558
 -0.358944   0.619863  -0.485472   -0.0735743   0.188655
In [95]:
# pour obtenir les coordonnées des individus
(Matrix(Z) * fp.projection)[1:8,:]
8×5 Matrix{Float64}:
  2.07866    1.73537   -0.55575    0.196238   -0.292866
 -1.51746   -1.48402   -1.27821   -0.205397    0.144332
  1.08785   -0.655501  -0.443724  -0.162903   -0.364788
  2.50123    0.109704  -0.144384  -0.0168546   0.220462
 -0.415801   0.67597    0.18784   -0.610068    0.256081
  0.295666  -0.190622  -0.656763  -0.540278   -0.432546
 -0.664659  -0.906768   0.249587   0.197411    0.202773
  1.89359   -0.952824   0.602119   0.61256     0.284905
In [96]:
# mais que l'on peut aussi obtenir directement avec la méthode transform()
MLJ.transform(mach,Z)[1:8,:]
8×5 DataFrame
Rowx1x2x3x4x5
Float64Float64Float64Float64Float64
12.078661.73537-0.555750.196238-0.292866
2-1.51746-1.48402-1.27821-0.2053970.144332
31.08785-0.655501-0.443724-0.162903-0.364788
42.501230.109704-0.144384-0.01685460.220462
5-0.4158010.675970.18784-0.6100680.256081
60.295666-0.190622-0.656763-0.540278-0.432546
7-0.664659-0.9067680.2495870.1974110.202773
81.89359-0.9528240.6021190.612560.284905

ATTENTION !!! En réalité, ces coordonnées sont biaisées parce que l'écart-type << 1/(n-1) >> a été utilisé pour la standardisation. On doit appliquer un facteur de correction pour retrouver les valeurs produites par les librairies qui font référence !!! (factorminer, sklearn, fanalysis, scientisttools, etc.)

In [97]:
# effectifs
n = size(X,1)
print(n)
18
In [98]:
# appliquer la correction
dfFact = sqrt(n/(n-1)) .* MLJ.transform(mach,Z)
DFR.first(dfFact,8)
8×5 DataFrame
Rowx1x2x3x4x5
Float64Float64Float64Float64Float64
12.138921.78568-0.5718620.201927-0.301357
2-1.56146-1.52704-1.31527-0.2113520.148516
31.11939-0.674505-0.456588-0.167626-0.375364
42.573740.112884-0.14857-0.01734320.226853
5-0.4278550.6955670.193286-0.6277540.263505
60.304238-0.196149-0.675803-0.555941-0.445086
7-0.683928-0.9330570.2568230.2031340.208652
81.94849-0.9804480.6195750.6303190.293165
In [99]:
# position des individus dans le premier plan factoriel
STP.scatter(dfFact[:,:x1], dfFact[:,:x2],
     aspect_ratio=:equal,
     label="",
     xlabel="PC1",
     ylabel="PC2",
     title="Graphique des individus",
     xlims=(-5,+5),
     ylims=(-5,+5),
     framestyle=:box,
     markersize=0,
     series_annotations=(df.Modele,STP.font(6,:green)),
     size=(600,600)
     )

# lignes centrales
STP.hline!([0], linestyle = :dash, color = :gray,label="")
STP.vline!([0], linestyle = :dash, color = :gray,label="")
No description has been provided for this image

Contributions des individus aux facteurs¶

In [100]:
# contributions des individus aux composantes
ctr_ind = 100.0.*((dfFact[:,[:x1,:x2]] .^2) ./ (n .* res.principalvars[1:2]'))
DFR.hcat(DFR.DataFrame(Modele=df.Modele),round.(ctr_ind,digits=1))
18×3 DataFrame
RowModelex1x2
StringFloat64Float64
1Alfasud TI5.720.7
2Audi 1003.115.1
3Simca 13001.63.0
4Citroen GS Club8.30.1
5Fiat 1320.23.1
6Lancia Beta0.10.2
7Peugeot 5040.65.6
8Renault 16 TL4.86.2
9Renault 3024.47.3
10Toyota Corolla20.00.4
11Alfetta-1.660.223.7
12Princess-18001.34.6
13Datsun-200L10.92.0
14Taunus-20002.21.5
15Rancho0.65.2
16Mazda-92950.20.8
17Opel-Rekord6.60.1
18Lada-13009.20.1
In [101]:
# vérification de la somme en colonne
# égale à 100% effectivement
sum(Matrix(ctr_ind),dims=1)
1×2 Matrix{Float64}:
 100.0  100.0

COS2 - Qualité de représentation des individus¶

In [102]:
# distance à l'origine des individus dans l'espace initial
# sachant que les individus sont centrés et résuits
# ATTENTION encore au facteur de correction de la réduction
disto_ind = vec(sum(n/(n-1) .* Matrix(Z) .^ 2,dims=2))
18-element Vector{Float64}:
  8.22517550693188
  6.6737551583965145
  2.159326874204685
  6.78014507568542
  1.169124371607794
  1.1349503123144318
  1.5127929344721436
  5.636825753612391
 21.789657262877334
 16.29014276553667
  4.456770360239878
  1.952512623718658
 11.112624034411022
  2.452985724307104
  1.9633731069400353
  0.6845212812064279
  6.08311856016591
  7.922198293371691
In [103]:
# COS2 des individus
cos2_ind = 100.0.*((dfFact[:,[:x1,:x2]] .^2) ./ disto_ind)
DFR.hcat(DFR.DataFrame(Modele=df.Modele),round.(cos2_ind,digits=1))
18×3 DataFrame
RowModelex1x2
StringFloat64Float64
1Alfasud TI55.638.8
2Audi 10036.534.9
3Simca 130058.021.1
4Citroen GS Club97.70.2
5Fiat 13215.741.4
6Lancia Beta8.23.4
7Peugeot 50430.957.5
8Renault 16 TL67.417.1
9Renault 3089.25.2
10Toyota Corolla97.50.3
11Alfetta-1.664.382.1
12Princess-180053.136.3
13Datsun-200L77.82.8
14Taunus-200070.59.6
15Rancho24.341.0
16Mazda-929521.718.5
17Opel-Rekord86.20.2
18Lada-130092.60.3
In [104]:
# qualité de représentation des idnvidus dans le plan factoriel
DFR.DataFrame(Modele=df.Modele,Cos2=vec(sum(Matrix(cos2_ind),dims=2)))
18×2 DataFrame
RowModeleCos2
StringFloat64
1Alfasud TI94.3889
2Audi 10071.4741
3Simca 130079.0978
4Citroen GS Club97.8871
5Fiat 13257.0405
6Lancia Beta11.5454
7Peugeot 50488.469
8Renault 16 TL84.4075
9Renault 3094.4351
10Toyota Corolla97.8645
11Alfetta-1.6686.363
12Princess-180089.3802
13Datsun-200L80.6527
14Taunus-200080.1315
15Rancho65.3742
16Mazda-929540.2674
17Opel-Rekord86.369
18Lada-130092.8659

Variables supplémentaires¶

Variables quantitatives - Corrélation¶

In [105]:
# rappel des variables disponibles
DFR.describe(df)
10×7 DataFrame
Rowvariablemeanminmedianmaxnmissingeltype
SymbolUnion…AnyUnion…AnyInt64DataType
1ModeleAlfasud TIToyota Corolla0String
2FINITION1_M3_TB0String
3PRIX34158.62210033345.0477000Int64
4R_POID_PUIS13.18079.7247713.181818.36360Float64
5CYL1631.6711661577.526640Int64
6PUISS84.61115582.01280Int64
7LONG433.5393434.54690Int64
8LARG166.667157167.01770Int64
9POIDS1078.838151087.513700Int64
10V_MAX158.278140160.01800Int64
In [106]:
# rappel des coordonnées factorielles
DFR.first(dfFact,5)
5×5 DataFrame
Rowx1x2x3x4x5
Float64Float64Float64Float64Float64
12.138921.78568-0.5718620.201927-0.301357
2-1.56146-1.52704-1.31527-0.2113520.148516
31.11939-0.674505-0.456588-0.167626-0.375364
42.573740.112884-0.14857-0.01734320.226853
5-0.4278550.6955670.193286-0.6277540.263505
In [107]:
# corrélation de PRIX avec les 2 premières composantes
import Statistics
corr_prix = [Statistics.cor(df.PRIX,col) for col in DFR.eachcol(dfFact[:,[:x1,:x2]])]
print(corr_prix)
[-0.7724752409002815, 0.08670844446517484]
In [108]:
# placer PRIX dans le cercle des corrélations
# dessiner le cercle
theta = range(0, 2*pi, length=200)

STP.plot(cos.(theta), sin.(theta),
     aspect_ratio=:equal,
     label="",
     linestyle=:dash,
     xlabel="PC1",
     ylabel="PC2",
     title="Cercle des corrélations",
     xlims=(-1,+1),
     ylims=(-1,+1),
     framestyle=:box)

# lignes centrales
STP.hline!([0], linestyle = :dash, color = :gray,label="")
STP.vline!([0], linestyle = :dash, color = :gray,label="")

# placer les noms des variables actives
# les points sont masqués avec markersize=0
STP.scatter!(coord[:,:F1], coord[:,:F2], markersize=0,
            series_annotations=(coord[:,:var],STP.font(7,:green)),
            label="")

# placer la variable illustrative
STP.scatter!([corr_prix[1]],[corr_prix[2]], markersize=0,
            series_annotations=(["PRIX"],STP.font(9,:blue)),
            label="")
No description has been provided for this image

Variable qualitative - Coordonnées factorielles¶

In [109]:
# dataframe avec les coordonnées dans le plan
dfTemp = dfFact[:,1:2]
dfTemp.FINITION = df.FINITION

DFR.first(dfTemp,5)
5×3 DataFrame
Rowx1x2FINITION
Float64Float64String
12.138921.785682_B
2-1.56146-1.527043_TB
31.11939-0.6745051_M
42.573740.1128841_M
5-0.4278550.6955672_B
In [110]:
# calcul des moyennes conditionnelles
coord_finition = DFR.combine(DFR.groupby(dfTemp,:FINITION),[:x1,:x2] .=> Statistics.mean)
coord_finition
3×3 DataFrame
RowFINITIONx1_meanx2_mean
StringFloat64Float64
12_B-0.2353130.0452712
23_TB-1.39243-0.0340006
31_M2.00035-0.022579
In [111]:
# position des individus dans le premier plan factoriel
STP.scatter(dfFact[:,:x1], dfFact[:,:x2],
     aspect_ratio=:equal,
     label="",
     xlabel="PC1",
     ylabel="PC2",
     title="Graphique des individus",
     xlims=(-5,+5),
     ylims=(-5,+5),
     framestyle=:box,
     markersize=0,
     series_annotations=(df.Modele,STP.font(7,:silver)),
     size=(600,600)
     )

# lignes centrales
STP.hline!([0], linestyle = :dash, color = :gray,label="")
STP.vline!([0], linestyle = :dash, color = :gray,label="")

# placer les points additionnels -- moyennes conditionnelles
STP.scatter!(coord_finition.x1_mean, coord_finition.x2_mean,
            label="",
            markersize=0,
            series_annotations=(coord_finition.FINITION,STP.font(9,:blue)))
No description has been provided for this image

Traitements des individus supplémentaires¶

In [112]:
# chargement des données
dfSupp = DFR.DataFrame(XLSX.readtable("./autos_acp_livre.xlsx","ind_illust"))
dfSupp
2×7 DataFrame
RowModeleCYLPUISSLONGLARGPOIDSV_MAX
StringInt64Int64Int64Int64Int64Int64
1Peugeot 60426641364721771410180
2Peugeot 304 S128874414157915160
In [113]:
# variables actives
XSupp = DFR.select(dfSupp, DFR.Not(:Modele))
XSupp
2×6 DataFrame
RowCYLPUISSLONGLARGPOIDSV_MAX
Int64Int64Int64Int64Int64Int64
126641364721771410180
2128874414157915160
In [114]:
# appliquer la standardisation
ZSupp = MLJ.transform(std,XSupp)
ZSupp
2×6 DataFrame
RowCYLPUISSLONGLARGPOIDSV_MAX
Float64Float64Float64Float64Float64Float64
12.760772.5221.74151.944662.418021.78925
2-0.919067-0.520758-0.882059-1.8192-1.196230.141859
In [115]:
# obtenir les coordonnées factorielles
# ne pas oublier d'appliquer le même facteur de correction
# que sur les individus actifs
dfFactSupp = sqrt(n/(n-1)) .* MLJ.transform(mach,ZSupp)
dfFactSupp
2×5 DataFrame
Rowx1x2x3x4x5
Float64Float64Float64Float64Float64
1-5.563290.3386090.4642890.402146-0.389811
22.212241.257780.0930439-0.3537020.648528
In [116]:
# position des individus dans le premier plan factoriel
STP.scatter(dfFact[:,:x1], dfFact[:,:x2],
     aspect_ratio=:equal,
     label="",
     xlabel="PC1",
     ylabel="PC2",
     title="Graphique des individus",
     xlims=(-6,+6),
     ylims=(-6,+6),
     framestyle=:box,
     markersize=0,
     series_annotations=(df.Modele,STP.font(7,:gray)),
     size=(600,600)
     )

# lignes centrales
STP.hline!([0], linestyle = :dash, color = :gray,label="")
STP.vline!([0], linestyle = :dash, color = :gray,label="")

# placer les points additionnels -- individus illustratifs
STP.scatter!(dfFactSupp.x1, dfFactSupp.x2,
            label="",
            markersize=0,
            series_annotations=(dfSupp.Modele,STP.font(9,:red)))
No description has been provided for this image

Rotation des facteurs¶

In [117]:
# récupérations des loadings pour les 2 premières composantes
loadings_init = res.loadings[:,1:2]
loadings_init
6×2 Matrix{Float64}:
 -0.893464   0.114906
 -0.886858   0.384689
 -0.886155  -0.381029
 -0.813536  -0.412736
 -0.905187  -0.224532
 -0.75471    0.573519
In [118]:
# vérification de la variance restituée
sum(loadings_init .^ 2, dims = 1)
1×2 Matrix{Float64}:
 4.42086  0.856062
In [119]:
# et la somme des deux
sum(sum(loadings_init .^ 2, dims = 1))
5.276920348854848
In [120]:
# rotation des facteurs -- VARIMAX ici
# mais d'autres sont possibles, ex. Quartimax
import FactorRotations as FR
L_rot = FR.rotate(loadings_init, FR.Varimax())
L_rot
FactorRotations.FactorRotation{Float64} with loading matrix:
6×2 Matrix{Float64}:
 0.593779  0.677427
 0.410288  0.87531
 0.916503  0.300792
 0.883045  0.228963
 0.827204  0.430711
 0.186252  0.92942
In [121]:
# accès aux loadings
loadings_rot = FR.loadings(L_rot)
DFR.hcat(DFR.DataFrame(var=names(X)),DFR.DataFrame(loadings_rot,["F1","F2"]))
6×3 DataFrame
RowvarF1F2
StringFloat64Float64
1CYL0.5937790.677427
2PUISS0.4102880.87531
3LONG0.9165030.300792
4LARG0.8830450.228963
5POIDS0.8272040.430711
6V_MAX0.1862520.92942
In [122]:
# carré de la somme en colonne
# répartition différente
sum(loadings_rot .^ 2, dims = 1)
1×2 Matrix{Float64}:
 2.85961  2.41731
In [123]:
# la la somme totale reste la même
# i.e. l'information restituée sur les 2 composantes
sum(sum(loadings_rot .^ 2, dims = 1))
5.2769203488548495
In [124]:
# cercle des corrélations << après >> rotation varimax
# dessiner le cercle
theta = range(0, 2*pi, length=200)

STP.plot(cos.(theta), sin.(theta),
     aspect_ratio=:equal,
     label="",
     linestyle=:dash,
     xlabel="PC1",
     ylabel="PC2",
     title="Cercle des corrélations",
     xlims=(-1,+1),
     ylims=(-1,+1),
     framestyle=:box)

# lignes centrales
STP.hline!([0], linestyle = :dash, color = :gray,label="")
STP.vline!([0], linestyle = :dash, color = :gray,label="")

# placer les noms des variables actives
# les points sont masqués avec markersize=0
STP.scatter!(loadings_rot[:,1], loadings_rot[:,2], markersize=0,
            series_annotations=(names(X),STP.font(7,:green)),
            label="")
No description has been provided for this image
In [125]:
# matrice de rotation orthogonale
# permettant de calculer les coordoonnées
# des individus après rotation des facteurs
mat_rot = FR.rotation(L_rot)
mat_rot
2×2 Matrix{Float64}:
 -0.749692  -0.661787
 -0.661787   0.749692
In [126]:
# anciennes coordonnées
old_coord = Matrix(dfFact[:,1:2])
old_coord[1:5,:]
5×2 Matrix{Float64}:
  2.13892    1.78568
 -1.56146   -1.52704
  1.11939   -0.674505
  2.57374    0.112884
 -0.427855   0.695567
In [127]:
# appliquer la matrice pour les nouvelles coordonnées
new_coord = old_coord * mat_rot
new_coord[1:5,:]
5×2 Matrix{Float64}:
 -2.78527   -0.0768024
  2.18119   -0.111456
 -0.392815  -1.24647
 -2.00422   -1.61864
 -0.139558   0.80461
In [128]:
# vérifier si coordonnées correctement calculés
# somme des carrés en colonne / n = valeur propre => oui
sum(new_coord .^ 2, dims = 1)/n
1×2 Matrix{Float64}:
 2.85961  2.41731
In [129]:
# nouvelles positions des individus dans le premier plan factoriel
STP.scatter(new_coord[:,1], new_coord[:,2],
     aspect_ratio=:equal,
     label="",
     xlabel="PC1",
     ylabel="PC2",
     title="Graphique des individus (après rotation VARIMAX)",
     xlims=(-4,+4),
     ylims=(-4,+4),
     framestyle=:box,
     markersize=0,
     series_annotations=(df.Modele,STP.font(7,:green)),
     size=(600,600)
     )

# lignes centrales
STP.hline!([0], linestyle = :dash, color = :gray,label="")
STP.vline!([0], linestyle = :dash, color = :gray,label="")
No description has been provided for this image