Environnement et packages¶

In [63]:
# activer l'environnement
using Pkg
Pkg.activate("env_julia_kmeans")
  Activating project at `c:\Users\ricco\Desktop\demo\env_julia_kmeans`
In [64]:
# liste des packages installés
Pkg.status()
Status `C:\Users\ricco\Desktop\demo\env_julia_kmeans\Project.toml`
  [aaaa29a8] Clustering v0.15.8
  [a93c6f00] DataFrames v1.8.2
  [da1fdf0e] FreqTables v1.0.0
  [7073ff75] IJulia v1.34.4
  [6f286f6a] MultivariateStats v0.10.4
  [2913bbd2] StatsBase v0.34.10
  [f3b207a7] StatsPlots v0.15.8
  [fdbf4ff8] XLSX v0.11.7

Importation et préparation des données¶

Importation et inspection des données¶

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

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

# premières lignes
println(DFR.size(df))
(26, 10)
In [66]:
# premières lignes
DFR.first(df,5)
5×10 DataFrame
RowModelePaysco2puissancecylindreevitesselongueurlargeurhauteurpoids
StringStringInt64Int64Int64Int64Int64Int64Int64Int64
1PANDA Autre135541108150354159154860
2TWINGO FR143601149151344163143840
3YARIS Autre13465998155364166150880
4CITRONC2 FR141611124158367166147932
5CORSA Autre1277012481653841651441035
In [67]:
# description
DFR.describe(df)
10×7 DataFrame
Rowvariablemeanminmedianmaxnmissingeltype
SymbolUnion…AnyUnion…AnyInt64DataType
1ModeleALFA 156 YARIS 0String
2PaysAutreFR0String
3co2174.846113161.02870Int64
4puissance135.96254139.02500Int64
5cylindree1885.699981939.032220Int64
6vitesse198.077150200.02500Int64
7longueur423.962344427.54860Int64
8largeur174.423159176.01940Int64
9hauteur147.923134146.51690Int64
10poids1288.588401350.017350Int64
In [68]:
# récupération des variables actives
X = df[:,4:end]
DFR.describe(X)
7×7 DataFrame
Rowvariablemeanminmedianmaxnmissingeltype
SymbolFloat64Int64Float64Int64Int64DataType
1puissance135.96254139.02500Int64
2cylindree1885.699981939.032220Int64
3vitesse198.077150200.02500Int64
4longueur423.962344427.54860Int64
5largeur174.423159176.01940Int64
6hauteur147.923134146.51690Int64
7poids1288.588401350.017350Int64
In [69]:
# dimensions
n, p = size(X)
println("Obsservation = $n, variables actives = $p")
Obsservation = 26, variables actives = 7

Standardisation¶

In [70]:
# moyennes par colonne
import Statistics
mean_actifs = Statistics.mean(Matrix(X),dims=1)
println(mean_actifs)
[135.96153846153845 1885.6923076923076 198.07692307692307 423.96153846153845 174.42307692307693 147.92307692307693 1288.576923076923]
In [71]:
# écarts-type par colonne -- avec 1/n (corrected = false)
etype_actifs = Statistics.std(Matrix(X),dims=1,corrected=false)
println(etype_actifs)
[59.37849958543708 599.542054034755 30.679388065458635 43.775195618972276 7.869967632242485 7.231996622505364 248.92185568303637]
In [72]:
# centrage-réduction
Z = (X .- mean_actifs) ./ etype_actifs
DFR.first(Z,5)
5×7 DataFrame
Rowpuissancecylindreevitesselongueurlargeurhauteurpoids
Float64Float64Float64Float64Float64Float64Float64
1-1.38032-1.29714-1.56708-1.5982-1.959740.840283-1.72173
2-1.27928-1.22876-1.53448-1.82664-1.45148-0.680736-1.80208
3-1.19507-1.48062-1.4041-1.36976-1.070280.287185-1.64139
4-1.26244-1.27046-1.30631-1.30123-1.07028-0.127638-1.43249
5-1.11087-1.06363-1.07815-0.912881-1.19735-0.542461-1.0187
In [73]:
# vérifions -- pour le fun : autre manière avec combine()
DFR.combine(Z, names(Z) .=> Statistics.mean)
1×7 DataFrame
Rowpuissance_meancylindree_meanvitesse_meanlongueur_meanlargeur_meanhauteur_meanpoids_mean
Float64Float64Float64Float64Float64Float64Float64
11.96424e-168.54018e-183.75768e-161.87884e-16-1.31519e-15-1.5383e-15-1.15292e-16
In [74]:
# et pour l'écart-type
DFR.combine(Z, names(Z) .=> x -> Statistics.std(x,corrected=false))
1×7 DataFrame
Rowpuissance_functioncylindree_functionvitesse_functionlongueur_functionlargeur_functionhauteur_functionpoids_function
Float64Float64Float64Float64Float64Float64Float64
11.01.01.01.01.01.01.0

K-Means avec le package Clustering¶

Transposition (parce que le package le veut ainsi)¶

In [75]:
# transformation de Z en Matrix pour la suite
# et transposition parce que Clustering le manipule ainsi
ZT = Matrix(Z)'
size(ZT)
(7, 26)

Clustering avec la méthode des K-Means¶

In [76]:
# importation de Clustering
import Clustering

# définir le générateur de nombre aléatoires
# pour disposer de résultats reproductibles
import Random
random_gen = Random.Xoshiro(42)

# lancement des K-Means
# 4 clusters pour commencer
# initialisation un peu plus "intelligente" que aléatoire
# cf. la doc. pour les références
k = 4
result = Clustering.kmeans(ZT,k;init=:kmpp,rng=random_gen,display=:final)
K-means converged with 7 iterations (objv = 45.819306801648224)
Clustering.KmeansResult{Matrix{Float64}, Float64, Int64}([0.2570068941995722 -1.139736178049572 -0.48774453150111985 1.26934488742589; -0.07232015513704199 -1.122963321309747 -0.03159462720726769 1.43966719237178; … ; -0.61928022265599 -0.2264053376722154 1.5662235020512365 0.1489108932340721; 0.3034453838009675 -1.3515649921976434 -0.1529673759358768 1.2236359991530503], [2, 2, 2, 2, 2, 2, 2, 3, 3, 3  …  1, 1, 4, 1, 1, 4, 4, 4, 4, 4], [2.0973187724565108, 0.9235303797633918, 0.5749094015636622, 0.10692731595061389, 0.3447622122972245, 1.231360088875098, 1.5091460813552207, 0.6503400181715344, 2.0365409210767034, 3.2085130412803062  …  0.23751437073437032, 0.5571407635109935, 3.63761170030517, 3.9109076524646786, 4.931668905002222, 5.541282472994869, 0.2802785936448995, 1.6501647519066793, 3.8238842523296945, 1.5997427752913183], [9, 7, 4, 6], [9, 7, 4, 6], 45.819306801648224, 7, true)

Inspection des résultats¶

In [77]:
# informations portées par les résultats
dump(result)
Clustering.KmeansResult{Matrix{Float64}, Float64, Int64}
  centers: Array{Float64}((7, 4)) [0.2570068941995722 -1.139736178049572 -0.48774453150111985 1.26934488742589; -0.07232015513704199 -1.122963321309747 -0.03159462720726769 1.43966719237178; … ; -0.61928022265599 -0.2264053376722154 1.5662235020512365 0.1489108932340721; 0.3034453838009675 -1.3515649921976434 -0.1529673759358768 1.2236359991530503]
  assignments: Array{Int64}((26,)) [2, 2, 2, 2, 2, 2, 2, 3, 3, 3  …  1, 1, 4, 1, 1, 4, 4, 4, 4, 4]
  costs: Array{Float64}((26,)) [2.0973187724565108, 0.9235303797633918, 0.5749094015636622, 0.10692731595061389, 0.3447622122972245, 1.231360088875098, 1.5091460813552207, 0.6503400181715344, 2.0365409210767034, 3.2085130412803062  …  0.23751437073437032, 0.5571407635109935, 3.63761170030517, 3.9109076524646786, 4.931668905002222, 5.541282472994869, 0.2802785936448995, 1.6501647519066793, 3.8238842523296945, 1.5997427752913183]
  counts: Array{Int64}((4,)) [9, 7, 4, 6]
  wcounts: Array{Int64}((4,)) [9, 7, 4, 6]
  totalcost: Float64 45.819306801648224
  iterations: Int64 7
  converged: Bool true

Assignation aux classes¶

In [78]:
# classes assignées
# par numéro de cluster
result.assignments
26-element Vector{Int64}:
 2
 2
 2
 2
 2
 2
 2
 3
 3
 3
 ⋮
 1
 4
 1
 1
 4
 4
 4
 4
 4
In [79]:
# nombre d'observations par cluster
# par ordre de numéro de cluster
result.counts
4-element Vector{Int64}:
 9
 7
 4
 6
In [80]:
# vérifions
import StatsBase
StatsBase.countmap(result.assignments)
Dict{Int64, Int64} with 4 entries:
  4 => 6
  2 => 7
  3 => 4
  1 => 9

Distances et inertie intra-classe¶

In [81]:
# carré des distances de
# chaque individu au centre de classe
# qui lui a été assignée
result.costs
26-element Vector{Float64}:
 2.0973187724565108
 0.9235303797633918
 0.5749094015636622
 0.10692731595061389
 0.3447622122972245
 1.231360088875098
 1.5091460813552207
 0.6503400181715344
 2.0365409210767034
 3.2085130412803062
 ⋮
 0.5571407635109935
 3.63761170030517
 3.9109076524646786
 4.931668905002222
 5.541282472994869
 0.2802785936448995
 1.6501647519066793
 3.8238842523296945
 1.5997427752913183
In [82]:
# somme
Statistics.sum(result.costs)
45.819306801648224
In [83]:
# que l'on peut avoir directement
# et qui correspond à l'inertie intra-classes
result.totalcost
45.819306801648224

Centres de classes - Moyennes conditionnelles¶

In [84]:
# coordonnées des centres de classes
# p = 7 variables, k = 4 clusters
result.centers
7×4 Matrix{Float64}:
  0.257007   -1.13974   -0.487745   1.26934
 -0.0723202  -1.12296   -0.0315946  1.43967
  0.490042   -1.20387   -0.499584   1.00251
  0.599899   -1.24249   -0.638753   0.975555
  0.609805   -1.26996   -0.149057   0.666278
 -0.61928    -0.226405   1.56622    0.148911
  0.303445   -1.35156   -0.152967   1.22364
In [85]:
# avec les noms des variables
DFR.hcat(DFR.DataFrame(var=names(Z)),DFR.DataFrame(result.centers,string.("cluster_",1:k)))
7×5 DataFrame
Rowvarcluster_1cluster_2cluster_3cluster_4
StringFloat64Float64Float64Float64
1puissance0.257007-1.13974-0.4877451.26934
2cylindree-0.0723202-1.12296-0.03159461.43967
3vitesse0.490042-1.20387-0.4995841.00251
4longueur0.599899-1.24249-0.6387530.975555
5largeur0.609805-1.26996-0.1490570.666278
6hauteur-0.61928-0.2264051.566220.148911
7poids0.303445-1.35156-0.1529671.22364
In [86]:
# sous forme de heatmap
import StatsPlots as SPS
SPS.heatmap(result.centers,
            xticks = (1:k,string.("cluster_",1:k)),
            yticks = (1:p,names(Z)),
            c = :RdBu,
            clims = (-1.5, 1.5),
            yflip=true) #pour ne pas inverser les lignes
No description has been provided for this image

Rapport de corrélation - Importance des variables¶

In [87]:
# calcul des carrés de rapport de corrélation
eta2 = map(x -> sum(result.counts .* x.^2),eachrow(result.centers)) ./ n
eta2
7-element Vector{Float64}:
 0.7810177058970296
 0.819778493123995
 0.7436511603725736
 0.8225998033501063
 0.6687959459341628
 0.5290637901651591
 0.8728120703243497
In [88]:
# associés aux noms de variables et triés de manière décroissante
dfEta2 = DFR.DataFrame(var = names(Z), eta2 = round.(eta2,digits=3))
DFR.sort(dfEta2,:eta2,rev=true)
7×2 DataFrame
Rowvareta2
StringFloat64
1poids0.873
2longueur0.823
3cylindree0.82
4puissance0.781
5vitesse0.744
6largeur0.669
7hauteur0.529

Visualisation dans un plan factoriel (ACP)¶

ACP avec le package MultivatiateStats.jl¶

In [89]:
# ACP de MultivariateStats (utilisation directe sans passer par MLJ)
# nous utilisons les données centrées et réduites pour une ACP normée
# ici aussi les données utilisées doivent être transposées
import MultivariateStats as MVS

# instanciation - entraînement
# on ne conserve que le premier plan factoriel
acp = MVS.fit(MVS.PCA,ZT; maxoutdim=2)
PCA(indim = 7, outdim = 2, principalratio = 0.8614816741559015)

Pattern matrix (unstandardized loadings):
────────────────────────
         PC1         PC2
────────────────────────
1   0.918891  -0.103308
2   0.889152   0.219829
3   0.957785  -0.168914
4   0.947379  -0.0319105
5   0.865393   0.0334297
6  -0.166778   0.996773
7   0.950739   0.240646
────────────────────────

Importance of components:
─────────────────────────────────────────────
                                PC1       PC2
─────────────────────────────────────────────
SS Loadings (Eigenvalues)  5.13045   1.14113
Variance explained         0.704733  0.156749
Cumulative variance        0.704733  0.861482
Proportion explained       0.818047  0.181953
Cumulative proportion      0.818047  1.0
─────────────────────────────────────────────
In [90]:
# coordonnées factorielles des individus
acp_coord = MVS.predict(acp,ZT)
acp_coord
2×26 Matrix{Float64}:
 3.93358    3.67505    3.3517     …  -3.70226   -2.82673   -3.53498
 0.496966  -0.918521  -0.0614611      0.266714  -0.803486   0.0209832
In [91]:
# sous la forme de dataframe
# en transposant les coordonnées
dfFact = DFR.hcat(DFR.DataFrame(modele=df.Modele),DFR.DataFrame(acp_coord',["F1","F2"]))
dfFact
26×3 DataFrame
RowmodeleF1F2
StringFloat64Float64
1PANDA 3.933580.496966
2TWINGO 3.67505-0.918521
3YARIS 3.3517-0.0614611
4CITRONC2 3.10829-0.369219
5CORSA 2.55102-0.686828
6FIESTA 2.08384-0.526902
7CLIO 2.07097-1.07609
8MODUS 1.431631.32487
9MUSA 1.183432.87182
10GOLF 0.9294730.390791
11MERC_A 0.1561361.65677
12AUDIA3 0.596604-0.675629
13CITRONC4 -0.578346-0.162753
14AVENSIS -0.5226080.177911
15VECTRA -1.24241-0.236413
16PASSAT -0.950595-0.261212
17LAGUNA -1.14852-0.727192
18MEGANECC -1.21535-0.92033
19PTCRUISER -1.313691.07915
20MONDEO -1.99861-0.573619
21MAZDARX8 -1.39977-2.25075
22VELSATIS -2.125321.77886
23CITRONC5 -2.512550.185565
24MERC_E -3.702260.266714
25ALFA 156 -2.82673-0.803486
26BMW530 -3.534980.0209832

Points dans le premier plan factoriel¶

In [92]:
# représentation graphique
# position des individus dans le premier plan factoriel
SPS.scatter(dfFact[:,:F1], dfFact[:,:F2],
     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,SPS.font(:black,7)),
     size=(600,600)
     )

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

Points selon leur cluster d'appartenance¶

In [93]:
# graphique qui sert de maquette
SPS.scatter(dfFact[:,:F1], dfFact[:,:F2],
     aspect_ratio=:equal,
     label="", markersize=0,
     xlabel="PC1", ylabel="PC2",
     title="Graphique des individus",
     xlims=(-5,+5), ylims=(-5,+5),
     framestyle=:box, size=(600,600)
     )

# couleurs selon le cluster d'appartenance
couleurs = [:red,:green,:blue,:orange]

# pour chaque cluster
for i in 1:k
     # filtrer : indices des individus du cluster n)i
     idx = findall(result.assignments .== i)
     # afficher les Modèles avec la couleur du cluster
     SPS.scatter!(dfFact[idx,:F1], dfFact[idx,:F2],    
          aspect_ratio=:equal, label="", markersize=0,
          series_annotations=(df.Modele[idx],SPS.font(couleurs[i],7))
     )
     # faire apparaître la légende avec la couleur
     SPS.scatter!([NaN],[NaN],
          markersize=6,label="Cluster_$i",markercolor=couleurs[i])     
end

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

Traitement des variables supplémentaires¶

In [94]:
# rappel des variables disponibles
DFR.describe(df)
10×7 DataFrame
Rowvariablemeanminmedianmaxnmissingeltype
SymbolUnion…AnyUnion…AnyInt64DataType
1ModeleALFA 156 YARIS 0String
2PaysAutreFR0String
3co2174.846113161.02870Int64
4puissance135.96254139.02500Int64
5cylindree1885.699981939.032220Int64
6vitesse198.077150200.02500Int64
7longueur423.962344427.54860Int64
8largeur174.423159176.01940Int64
9hauteur147.923134146.51690Int64
10poids1288.588401350.017350Int64

Variable supplémentaire quantitative vs. clusters¶

In [95]:
# dataframe avec clusters et co2
dfTemp = DFR.DataFrame(co2=df.co2,cluster=result.assignments)
dfTemp
26×2 DataFrame
Rowco2cluster
Int64Int64
11352
21432
31342
41412
51272
61172
71132
81633
91463
101433
111413
121681
131421
141551
151591
161971
171961
181911
192354
201891
212841
221884
232384
241834
252874
262314
In [96]:
# moyennes conditionnelles
round.(DFR.combine(DFR.groupby(dfTemp,:cluster),:co2 => StatsBase.mean),digits=1)
4×2 DataFrame
Rowclusterco2_mean
Float64Float64
11.0186.8
22.0130.0
33.0148.2
44.0227.0

Variable supplémentaire qualitative vs. clusters¶

In [97]:
# associer cluster et Pays
dfTemp = DFR.DataFrame(Pays=df.Pays,cluster=result.assignments)
dfTemp
26×2 DataFrame
RowPayscluster
StringInt64
1Autre2
2FR2
3Autre2
4FR2
5Autre2
6Autre2
7FR2
8FR3
9Autre3
10Autre3
11Autre3
12Autre1
13FR1
14Autre1
15Autre1
16Autre1
17FR1
18FR1
19Autre4
20Autre1
21Autre1
22FR4
23FR4
24Autre4
25Autre4
26Autre4
In [98]:
# tableau croisé
import FreqTables as FT
tbl = FT.freqtable(dfTemp,:cluster,:Pays)
tbl
4×2 Named Matrix{Int64}
cluster ╲ Pays │ Autre     FR
───────────────┼─────────────
1              │     6      3
2              │     4      3
3              │     3      1
4              │     4      2
In [99]:
# profil ligne
FT.prop(tbl,margins=1)
4×2 Named Matrix{Float64}
cluster ╲ Pays │    Autre        FR
───────────────┼───────────────────
1              │ 0.666667  0.333333
2              │ 0.571429  0.428571
3              │     0.75      0.25
4              │ 0.666667  0.333333
In [100]:
# profil colonne
FT.prop(tbl,margins=2)
4×2 Named Matrix{Float64}
cluster ╲ Pays │    Autre        FR
───────────────┼───────────────────
1              │ 0.352941  0.333333
2              │ 0.235294  0.333333
3              │ 0.176471  0.111111
4              │ 0.235294  0.222222

Individus supplémentaires¶

Importation, préparation¶

In [101]:
# chargement des données
dfSupp = DFR.DataFrame(XLSX.readtable("./autos_kmeans.xlsx","illustratifs"))
dfSupp
4×10 DataFrame
RowModelePaysco2puissancecylindreevitesselongueurlargeurhauteurpoids
StringStringInt64Int64Int64Int64Int64Int64Int64Int64
1P1007 FR1537513601653741691611181
2P407 FR19413619972124681821451415
3P307CC FR21018019972254351761431490
4P607 FR22320427212304911841451723
In [102]:
# centrage réduction avec les moyennes et écarts-type
# des individus actifs
# on ne traite que les variables actives bien sûr
ZSupp = (dfSupp[:,4:end] .- mean_actifs) ./ etype_actifs
ZSupp
4×7 DataFrame
Rowpuissancecylindreevitesselongueurlargeurhauteurpoids
Float64Float64Float64Float64Float64Float64Float64
1-1.02666-0.876823-1.07815-1.14132-0.6890851.8082-0.432171
20.0006477350.1856550.4538251.006010.962764-0.4041870.507883
30.7416570.1856550.8775620.2521620.200372-0.6807360.809182
41.145841.393241.040541.531431.21689-0.4041871.74522
In [103]:
# que l'on va transposer
ZTSupp = Matrix(ZSupp)'
ZTSupp
7×4 adjoint(::Matrix{Float64}) with eltype Float64:
 -1.02666    0.000647735   0.741657   1.14584
 -0.876823   0.185655      0.185655   1.39324
 -1.07815    0.453825      0.877562   1.04054
 -1.14132    1.00601       0.252162   1.53143
 -0.689085   0.962764      0.200372   1.21689
  1.8082    -0.404187     -0.680736  -0.404187
 -0.432171   0.507883      0.809182   1.74522

Coordonnées factorielles¶

In [104]:
# coordonnées factorielles
ZFactSupp = MVS.predict(acp,ZTSupp)
2×4 Matrix{Float64}:
 2.27178  -1.2966    -1.31663   -3.31953
 1.69174  -0.296274  -0.626452   0.0197047
In [105]:
# graphique qui sert de maquette
SPS.scatter(dfFact[:,:F1], dfFact[:,:F2],
     aspect_ratio=:equal,
     label="", markersize=0,
     xlabel="PC1", ylabel="PC2",
     title="Graphique des individus",
     xlims=(-5,+5), ylims=(-5,+5),
     framestyle=:box, size=(600,600)
     )

# couleurs selon le cluster d'appartenance
couleurs = [:red,:green,:blue,:orange]

# pour chaque cluster
for i in 1:k
     idx = findall(result.assignments .== i)
     SPS.scatter!(dfFact[idx,:F1], dfFact[idx,:F2],
          label="", markersize=0,
          series_annotations=(df.Modele[idx],SPS.font(couleurs[i],4))
     )     
end

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

# individus supplémentaires
SPS.scatter!(ZFactSupp[1,:], ZFactSupp[2,:],
     label="", markersize=0,
     series_annotations=(dfSupp.Modele,SPS.font(:black,9))
     )
No description has been provided for this image

Rattachement aux clusters via distances aux centres de classes¶

In [106]:
# pour rappel, voici les centroïdes
result.centers
7×4 Matrix{Float64}:
  0.257007   -1.13974   -0.487745   1.26934
 -0.0723202  -1.12296   -0.0315946  1.43967
  0.490042   -1.20387   -0.499584   1.00251
  0.599899   -1.24249   -0.638753   0.975555
  0.609805   -1.26996   -0.149057   0.666278
 -0.61928    -0.226405   1.56622    0.148911
  0.303445   -1.35156   -0.152967   1.22364
In [107]:
# calculer les distances aux centroïdes
distances = []
for y in eachcol(ZTSupp)
    push!(distances,map(x -> sum((x .- y) .^2), eachcol(result.centers)))
end

# contrôle
# chaque ligne corresp. à un individu supplémentaire
# pour chaque individu => distance aux 4 centroïdes des clusters
distances
4-element Vector{Any}:
 [15.907016795352956, 5.421740837986061, 2.0202913517138494, 26.780032217023596]
 [0.5111522835842484, 19.290860656426496, 9.455357675944878, 4.390247387163716]
 [0.999710313991838, 18.855572941073603, 10.345543636003631, 3.4670658623765127]
 [6.602227699238375, 40.09323619128301, 21.131875334523173, 1.2089896922400951]
In [108]:
# transformer en matrice
distances = stack(distances;dims=1)
distances
4×4 Matrix{Float64}:
 15.907      5.42174   2.02029  26.78
  0.511152  19.2909    9.45536   4.39025
  0.99971   18.8556   10.3455    3.46707
  6.60223   40.0932   21.1319    1.20899
In [109]:
# assignation aux clusters
# min. par ligne
clus_supp = vec(mapslices(argmin,distances;dims=2))
clus_supp
4-element Vector{Int64}:
 3
 1
 1
 4
In [110]:
# graphique factoriel avec les groupes d'affectation

# graphique qui sert de maquette
SPS.scatter(dfFact[:,:F1], dfFact[:,:F2],
     aspect_ratio=:equal,
     label="", markersize=0,
     xlabel="PC1", ylabel="PC2",
     title="Graphique des individus",
     xlims=(-5,+5), ylims=(-5,+5),
     framestyle=:box, size=(600,600)
     )

# couleurs selon le cluster d'appartenance
couleurs = [:red,:green,:blue,:orange]

# pour chaque cluster
for i in 1:k
     idx = findall(result.assignments .== i)
     SPS.scatter!(dfFact[idx,:F1], dfFact[idx,:F2],    
          label="", markersize=0,          
          series_annotations=(df.Modele[idx],SPS.font(couleurs[i],4))
     )     
end

# pour chaque individu supplémentaire
for id in 1:size(dfSupp,1)
     SPS.scatter!([ZFactSupp[1,id]], [ZFactSupp[2,id]],
          label="", markersize=0,
          series_annotations=([dfSupp.Modele[id]],
                              SPS.font(couleurs[clus_supp[id]],9))
     )    
end

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

Identification du nombre de clusters¶

En passant par une boucle qui passe en revue les différentes solutions, courbe WSS puis repérage du "coude".

In [111]:
# valeurs de k à tester
ks = 1:7

# stockage des inerties intra-classes
wss = Float64[]

# définir de nouveau le générateur de nombres aléatoires
random_gen = Random.Xoshiro(42)

# boucler
for i in ks
    R = Clustering.kmeans(ZT,i;init=:kmpp,rng=random_gen)
    push!(wss,R.totalcost)
end

# affichage des valeurs
println(wss)
[182.0, 83.54649247605984, 66.45537108521084, 50.1163449111245, 46.27814717278716, 34.02141281398916, 30.91000480329152]
In [112]:
# courbe de décroissance du WSS
SPS.plot(ks, wss, marker=:o,
     xlabel="k", ylabel="WSS",
     title="Méthode du coude",
     label="")
No description has been provided for this image

Il y a d'autres pistes, comme l'utilisation de critères à optimiser (ex. silhouette statistic), ou encore basées sur des heuristiques (ex. gap statistic), etc...