Kapitel 21 Klyngeanalyse Premium


Vi har et datasæt bestående af 32 bilmodeller med 11 variable:

  1. mpg Miles/(US) gallon
  2. cyl Number of cylinders
  3. disp Displacement (cu.in.)
  4. hp Gross horsepower
  5. drat Rear axle ratio
  6. wt Weight (1000 lbs)
  7. qsec 1/4 mile time
  8. vs Engine (0 = V-shaped, 1 = straight)
  9. am Transmission (0 = automatic, 1 = manual)
  10. gear Number of forward gears
  11. carb Number of carburetors

Vi vil gerne gruppere de forskellige bilmodeller i forskellige grupper eller klynger udfra deres specifikationer. For at undersøge, om vi på baggrund af tekniske karakteristika, kan gruppere bilmodellerne, benytter vi klyngeanalyse. Bemærk i faktoranalysen grupperer vi variablene, i klyngeanalysen grupperer vi respondenterne eller observationerne, her altså bilerne.

Der findes overordnet 2 typer af klyngeudvælgelse:

Ikke-hierarkisk, k-means metoden benyttes, hvis vi har store datasæt hvor der kræves mange observationer, man vælger på forhånd hvor mange klynger man vil have.

Hierarkisk klyngedannelse, agglomerative metode hvor man starter med at hver respondent har sin egen klynge og man derefter sammenhober disse trin for trin kaldes den sammenhobede eller agglomerative metode. Vi benytter til bildatasættet den agglomerative metode, da vi ikke har en stort datasæt, er denne klart at foretrække.

pacman::p_load("datasets")#Vi henter pakken datasets der indeholder en del datasæt
head(mtcars) #Vi kan se starten af datasættet med Head kommandoen
##                    mpg cyl disp  hp drat    wt  qsec vs am gear carb
## Mazda RX4         21.0   6  160 110 3.90 2.620 16.46  0  1    4    4
## Mazda RX4 Wag     21.0   6  160 110 3.90 2.875 17.02  0  1    4    4
## Datsun 710        22.8   4  108  93 3.85 2.320 18.61  1  1    4    1
## Hornet 4 Drive    21.4   6  258 110 3.08 3.215 19.44  1  0    3    1
## Hornet Sportabout 18.7   8  360 175 3.15 3.440 17.02  0  0    3    2
## Valiant           18.1   6  225 105 2.76 3.460 20.22  1  0    3    1

Vi kan se at når vi sammenligner forskellige variable ser det ud til at der er forskellige grupper. Nedenfor ser vi fx de 32 biler plottet i et diagram efter hestekræfter og miles per gallon. Bilerne er farvekodet med antal cylindre.

#Her angiver vi modelnavne for hver bil
#med label rownames vi bruger geom_text i stedet for point her.
pacman::p_load("ggplot2")
ggplot(mtcars, aes(hp, mpg, color = cyl)) +
  geom_point() #Plot med kun punkter

ggplot(mtcars, aes(hp, mpg, color = cyl,label=rownames(mtcars)))+
  geom_text(size=3,check_overlap = TRUE)

Vi kan benytte t til at transponere data matricen, så kan vi tegne et corrplot, det er meget mørkt da alle bilerne er positivt korrelerede, men man kan ane nogle sammenhænge. t(mtcars) betyder vi transponerer (vender) matricen, så ser vi i stedet på grupper af respondenter, som vi netop analyserer i klyngeanalysen. Hvis ikke vi vender matricen ser vi på korrelationsmatricen mellem de 11 variable i stedet, ligesom vi tidligere har gjort med faktoranalysen.

#korrelationsmatricen for transponeret mtcars data, hclust betyder vi ordner efter variable der passer sammen
pacman::p_load(corrplot)
corrplot(cor(t(mtcars)), order = "hclust", tl.col='black', tl.cex=.5)

#Her er korrelationsmatricen ikke transponeret, hvilket svarer til en matrice baseret på variable som ved faktoranalysen vi tidligere så på.
corrplot(cor(mtcars), order = "hclust", tl.col='black', tl.cex=.5)

21.1 Hierakisk klyngeanalyse hclust kommandoen Premium


Vi får nu R til at danne klynger vha. af hclust (hierarcical cluster) kommandoen, denne benytter default en metode der hedder complete til at finde ens klynger, der findes mange andre metoder. For at benytte hclust skal R først beregne afstandene mellem bilerne dette gøres med dist kommandoen. Algoritmen beregner afstandene mellem de forskellige bilmodeller vha. af den euklidiske metrik. Biler med kort afstand kommer i klynger sammen, biler med lang afstand kommer i forskellige klynger.

Nedenfor ses et udsnit af afstandende mellem hver af de 32 biler, det er en meget stor 32 \(\times\) 32 matrice, derfor har vi benyttet head for kun at vise noget af matricen. Fx er afstanden mellem to forskellige biler som en Mazda RX4 og en Lincoln Continental 318.05 hvilket er en stor afstand i forhold til fx. Mazda RX4 og Mazda RX4 Wag på kun 0.62. Bemærk hvordan supersportsvognen Maserati Bora har store afstande til de fleste af de øvrige biler, Maseratien var en komfortabel, rummeligere og kraftigere og tungere sportsvogn end fx. Ferrari Dino.

head(as.matrix(dist(mtcars)))
##                     Mazda RX4 Mazda RX4 Wag Datsun 710 Hornet 4 Drive
## Mazda RX4           0.0000000     0.6153251   54.90861       98.11252
## Mazda RX4 Wag       0.6153251     0.0000000   54.89152       98.09589
## Datsun 710         54.9086059    54.8915169    0.00000      150.99352
## Hornet 4 Drive     98.1125212    98.0958939  150.99352        0.00000
## Hornet Sportabout 210.3374396   210.3358546  265.08316      121.02976
## Valiant            65.4717710    65.4392224  117.75470       33.55087
##                   Hornet Sportabout   Valiant Duster 360 Merc 240D  Merc 230
## Mazda RX4                  210.3374  65.47177  241.40765  50.15327  25.46831
## Mazda RX4 Wag              210.3359  65.43922  241.40887  50.11461  25.32845
## Datsun 710                 265.0832 117.75470  294.47902  49.65848  33.18038
## Hornet 4 Drive             121.0298  33.55087  169.42996 121.27397 118.24331
## Hornet Sportabout            0.0000 152.12414   70.17673 241.50697 233.49240
## Valiant                    152.1241   0.00000  194.60945  89.59111  85.00796
##                    Merc 280 Merc 280C Merc 450SE Merc 450SL Merc 450SLC
## Mazda RX4          15.36419  15.67247  135.43070  135.40144   135.47947
## Mazda RX4 Wag      15.29569  15.58377  135.42548  135.39604   135.47232
## Datsun 710         66.93635  67.02614  189.19549  189.16317   189.23454
## Hornet 4 Drive     91.42240  91.46129   72.49643   72.43135    72.57185
## Hornet Sportabout 199.33450 199.34066   84.38885   84.36840    84.43324
## Valiant            60.29098  60.26557   90.69703   90.67697    90.70930
##                   Cadillac Fleetwood Lincoln Continental Chrysler Imperial
## Mazda RX4                   326.3396            318.0470         304.72034
## Mazda RX4 Wag               326.3355            318.0429         304.71692
## Datsun 710                  381.0926            372.8012         359.30149
## Hornet 4 Drive              234.4404            227.9726         218.15483
## Hornet Sportabout           116.2804            108.0624          97.20491
## Valiant                     266.6281            259.6304         248.77133
##                    Fiat 128 Honda Civic Toyota Corolla Toyota Corona
## Mazda RX4          93.26800   102.83076       100.6040      42.30752
## Mazda RX4 Wag      93.25310   102.82387       100.5888      42.26592
## Datsun 710         40.99338    52.77046        47.6535      12.96547
## Hornet 4 Drive    184.96897   191.55187       192.6714     138.53047
## Hornet Sportabout 302.03772   310.03246       309.5582     252.33320
## Valiant           152.11533   158.96158       159.8303     105.28764
##                   Dodge Challenger AMC Javelin Camaro Z28 Pontiac Firebird
## Mazda RX4                163.11508   149.60472  233.22288        248.67803
## Mazda RX4 Wag            163.11342   149.60145  233.22487        248.67620
## Datsun 710               217.77958   204.31889  286.00492        303.35839
## Hornet 4 Drive            72.44039    61.36019  163.66326        156.22403
## Hornet Sportabout         48.98389    61.42742   70.96653         40.00525
## Valiant                  103.43107    91.04443  187.84638        188.52721
##                   Fiat X1-9 Porsche 914-2 Lotus Europa Ford Pantera L
## Mazda RX4          92.50484      44.40337     65.73284       245.4247
## Mazda RX4 Wag      92.49400      44.40736     65.73626       245.4294
## Datsun 710         39.88151      13.13571     25.09486       297.2940
## Hornet 4 Drive    184.44712     139.15795    163.23674       180.1140
## Hornet Sportabout 301.56695     254.14526    272.35824        89.5934
## Valiant           151.43794     106.05858    130.82482       203.0178
##                   Ferrari Dino Maserati Bora Volvo 142E
## Mazda RX4             66.76610      265.6454   39.18940
## Mazda RX4 Wag         66.77642      265.6491   39.16260
## Datsun 710            90.24155      309.7718   20.69394
## Hornet 4 Drive       130.55230      229.3419  137.03633
## Hornet Sportabout    215.06739      170.7094  248.00634
## Valiant              106.56948      242.4393  104.18637

Vi gemmer hclust data i clusters variablen, vi så kan benytte til at tegne en oversigt over klyngerne. Vi kan nu plotte en grafisk oversigt over bilerne. I nederste linje er den fineste inddeling, hvor samtlige biler er i deres egen klynge. Den blå linje med 3 skæringer i dendogrammet indikerer der er 3 klynger, den røde 4 klynger.

clusters <- hclust(dist(mtcars))
plot(clusters,cex=0.5,main = "Dendogram af mtcars",xlab = "Klyngetræ",sub="Bilmodeller")
abline(h = 190, col="red") #Tegn rød vandret linje h betyder horisontal
abline(h = 230, col="blue")

Hvis vi ønsker at undersøge en indeling med et bestemt antal klynger, kan vi bruge cutree i R, til at undersøge klyngerne i en skæring med fx. 4 klynger nærmere. Her ser vi som nævnt, Maserati Bora skiller sig ud ved at have sin egen klynge. Nummeret ved hver af de 32 biler angiver hvilken klynge bilen tilhører.

clusterCut <- cutree(clusters, 4) #Opdeling i 4 klynger.
clusterCut
##           Mazda RX4       Mazda RX4 Wag          Datsun 710      Hornet 4 Drive 
##                   1                   1                   1                   2 
##   Hornet Sportabout             Valiant          Duster 360           Merc 240D 
##                   3                   2                   3                   1 
##            Merc 230            Merc 280           Merc 280C          Merc 450SE 
##                   1                   1                   1                   2 
##          Merc 450SL         Merc 450SLC  Cadillac Fleetwood Lincoln Continental 
##                   2                   2                   3                   3 
##   Chrysler Imperial            Fiat 128         Honda Civic      Toyota Corolla 
##                   3                   1                   1                   1 
##       Toyota Corona    Dodge Challenger         AMC Javelin          Camaro Z28 
##                   1                   2                   2                   3 
##    Pontiac Firebird           Fiat X1-9       Porsche 914-2        Lotus Europa 
##                   3                   1                   1                   1 
##      Ford Pantera L        Ferrari Dino       Maserati Bora          Volvo 142E 
##                   3                   1                   4                   1

Vi kan benytte subset kommandoen til at se på hvilke variable der er i hver klynge, herunder ser vi på klynge 3.

subset(clusterCut,clusterCut==3)
##   Hornet Sportabout          Duster 360  Cadillac Fleetwood Lincoln Continental 
##                   3                   3                   3                   3 
##   Chrysler Imperial          Camaro Z28    Pontiac Firebird      Ford Pantera L 
##                   3                   3                   3                   3

Vi kan ligeledes sammenligne klyngeinddelingen med de enkelte variable og se om disse passer sammen. Fx. passer hp meget fint med inddelingen i klynger.

table(clusterCut, mtcars$cyl)
##           
## clusterCut  4  6  8
##          1 11  5  0
##          2  0  2  5
##          3  0  0  8
##          4  0  0  1
table(clusterCut, mtcars$mpg)
##           
## clusterCut 10.4 13.3 14.3 14.7 15 15.2 15.5 15.8 16.4 17.3 17.8 18.1 18.7 19.2
##          1    0    0    0    0  0    0    0    0    0    0    1    0    0    1
##          2    0    0    0    0  0    2    1    0    1    1    0    1    0    0
##          3    2    1    1    1  0    0    0    1    0    0    0    0    1    1
##          4    0    0    0    0  1    0    0    0    0    0    0    0    0    0
##           
## clusterCut 19.7 21 21.4 21.5 22.8 24.4 26 27.3 30.4 32.4 33.9
##          1    1  2    1    1    2    1  1    1    2    1    1
##          2    0  0    1    0    0    0  0    0    0    0    0
##          3    0  0    0    0    0    0  0    0    0    0    0
##          4    0  0    0    0    0    0  0    0    0    0    0
table(clusterCut, mtcars$hp)
##           
## clusterCut 52 62 65 66 91 93 95 97 105 109 110 113 123 150 175 180 205 215 230
##          1  1  1  1  2  1  1  1  1   0   1   2   1   2   0   1   0   0   0   0
##          2  0  0  0  0  0  0  0  0   1   0   1   0   0   2   0   3   0   0   0
##          3  0  0  0  0  0  0  0  0   0   0   0   0   0   0   2   0   1   1   1
##          4  0  0  0  0  0  0  0  0   0   0   0   0   0   0   0   0   0   0   0
##           
## clusterCut 245 264 335
##          1   0   0   0
##          2   0   0   0
##          3   2   1   0
##          4   0   0   1
table(clusterCut, mtcars$carb)
##           
## clusterCut 1 2 3 4 6 8
##          1 5 6 0 4 1 0
##          2 2 2 3 0 0 0
##          3 0 2 0 6 0 0
##          4 0 0 0 0 0 1
table(clusterCut, mtcars$wt)
##           
## clusterCut 1.513 1.615 1.835 1.935 2.14 2.2 2.32 2.465 2.62 2.77 2.78 2.875
##          1     1     1     1     1    1   1    1     1    1    1    1     1
##          2     0     0     0     0    0   0    0     0    0    0    0     0
##          3     0     0     0     0    0   0    0     0    0    0    0     0
##          4     0     0     0     0    0   0    0     0    0    0    0     0
##           
## clusterCut 3.15 3.17 3.19 3.215 3.435 3.44 3.46 3.52 3.57 3.73 3.78 3.84 3.845
##          1    1    0    1     0     0    2    0    0    0    0    0    0     0
##          2    0    0    0     1     1    0    1    1    0    1    1    0     0
##          3    0    1    0     0     0    1    0    0    1    0    0    1     1
##          4    0    0    0     0     0    0    0    0    1    0    0    0     0
##           
## clusterCut 4.07 5.25 5.345 5.424
##          1    0    0     0     0
##          2    1    0     0     0
##          3    0    1     1     1
##          4    0    0     0     0

Hvis dendogrammet virker lidt uoverskueligt, kan man vælge ape pakken for at lave mere fancy plots, her er rigtig mange muligheder.

#install.packages("ape")
library("ape")
colors = c("red", "blue", "green", "pink")
clus4 = cutree(clusters, 4)
plot(as.phylo(clusters), type = "fan", tip.color = colors[clus4],
     label.offset = 0, cex = 0.5)

Herunder er et plot, hvor farvekoden er baseret på klyngerne.

ggplot(mtcars, aes(hp, mpg)) +
  geom_point(alpha = 0.4, size = 3.5) + geom_point(col = clusterCut)

21.2 Ikke hierakisk klyngeanalyse kmeans kommandoen Premium


Vi kunne også have brugt kmeans metoden, her skal vi så angive hvor mange klynger, vi ønsker i analysen. Her benytter vi K-means og får 4 klynger med 7, 6, 9, 10 biler. I output fra R kan vi under Cluster means se gennemsnit, for de 4 klynger for alle 11 variable.

mtcarsCluster4 <- kmeans(mtcars, 4)
mtcarsCluster4
## K-means clustering with 4 clusters of sizes 9, 7, 6, 10
## 
## Cluster means:
##        mpg      cyl     disp       hp     drat       wt     qsec        vs
## 1 14.64444 8.000000 388.2222 232.1111 3.343333 4.161556 16.40444 0.0000000
## 2 19.94286 5.714286 166.5714 120.1429 3.705714 3.107857 18.47143 0.5714286
## 3 16.83333 7.666667 284.5667 158.3333 3.033333 3.625000 17.76833 0.1666667
## 4 27.05000 4.000000 101.5700  81.4000 4.086000 2.199300 18.76100 0.9000000
##          am     gear     carb
## 1 0.2222222 3.444444 4.000000
## 2 0.4285714 4.000000 3.571429
## 3 0.0000000 3.000000 2.333333
## 4 0.8000000 4.100000 1.500000
## 
## Clustering vector:
##           Mazda RX4       Mazda RX4 Wag          Datsun 710      Hornet 4 Drive 
##                   2                   2                   4                   3 
##   Hornet Sportabout             Valiant          Duster 360           Merc 240D 
##                   1                   2                   1                   4 
##            Merc 230            Merc 280           Merc 280C          Merc 450SE 
##                   2                   2                   2                   3 
##          Merc 450SL         Merc 450SLC  Cadillac Fleetwood Lincoln Continental 
##                   3                   3                   1                   1 
##   Chrysler Imperial            Fiat 128         Honda Civic      Toyota Corolla 
##                   1                   4                   4                   4 
##       Toyota Corona    Dodge Challenger         AMC Javelin          Camaro Z28 
##                   4                   3                   3                   1 
##    Pontiac Firebird           Fiat X1-9       Porsche 914-2        Lotus Europa 
##                   1                   4                   4                   4 
##      Ford Pantera L        Ferrari Dino       Maserati Bora          Volvo 142E 
##                   1                   2                   1                   4 
## 
## Within cluster sum of squares by cluster:
## [1] 46659.317  8808.032  6355.581 10247.471
##  (between_SS / total_SS =  88.4 %)
## 
## Available components:
## 
## [1] "cluster"      "centers"      "totss"        "withinss"     "tot.withinss"
## [6] "betweenss"    "size"         "iter"         "ifault"

Nedenfor ser vi på en tabel inddeling med 4 klynger med de 32 biler, sorteret efter antallet af cylindre.

table(mtcarsCluster4$cluster, mtcars$cyl)
##    
##      4  6  8
##   1  0  0  9
##   2  1  6  0
##   3  0  1  5
##   4 10  0  0

Nedenfor ser vi på en tabel inddeling med 3 klynger med de 32 biler, sorteret efter antallet af cylindre. Kører vi kmeans analysen igen, falder bilerne ikke nødvendigvis i samme klynger som tidligere, det skyldes algoritmen kan give forskellig optimale inddelinger.

mtcarsCluster3 <- kmeans(mtcars, 3)
table(mtcarsCluster3$cluster, mtcars$cyl)
##    
##      4  6  8
##   1 11  0  0
##   2  0  7  0
##   3  0  0 14

Vi kan lave et klyngeplot der viser forskellene på de 32 biler, herunder ses plottet med 3 klynger. Bemærk vi skal hente pakken factoextra, der indeholder plotfunktionen fviz_cluster(). Vi har benyttet funktionen scale(), det er en rigtig god ide at benytte hvis, der er stor forskel på måleenhederne i en data.frame. Funktionen scale() bringer variablene i samme skale. Der er fx. stor forskel på enhederne i carb og disp, prøv at sammenligne nedenstående plot med et tilsvarende plot uden scale.

km3.res <- kmeans(scale(mtcars), 3, nstart = 25)
pacman::p_load(factoextra)
fviz_cluster(km3.res, data = mtcars, main = "Klyngeplot biler opdelt i 3 klynger",repel = TRUE)

Vi kan lave et klyngeplot der viser forskellene på de 32 biler, herunder ses plottet med 4 klynger.

km4.res <- kmeans(scale(mtcars), 4, nstart = 25)
fviz_cluster(km4.res, data = mtcars, main = "Klyngeplot biler opdelt i 4 klynger",repel = TRUE)

Vi kan lave et klyngeplot der viser forskellene på de 32 biler, herunder ses plottet med 5 klynger.

km5.res <- kmeans(scale(mtcars), 5, nstart = 25)
fviz_cluster(km5.res, data = mtcars, main = "Klyngeplot biler opdelt i 5 klynger",repel = TRUE)

21.3 Validering med ANOVA Premium

Skal man undersøge om grupperne/klyngerne er forskellige med hensyn til de forskellige variable, kan man benytte Anova, hvor klyngerne er den uafhængige variabel. Man skal da gerne nå frem til at klyngegennemsnittede er signifikant forskellige mht. flere af variablene der indgår i analysen.

Spørgsmål US arrestationer samt urbaniseringsgrad

Lav en klyngeanalyse for datasæt med 50 observationer for amerikanske stater på 4 variable, data stammer fra World Almanac and Book of facts 1975. (Crime rates).

  1. Mord, antal arrestationer (pr 100,000)
  2. Overfald, antal arrestationer (pr 100,000)
  3. Urbaniseringsgrad andel af bybefolkning.
  4. Voldtægt, antal arrestationer (pr 100,000)
pacman::p_load(datasets)
arrest <- USArrests
head(arrest)
##            Murder Assault UrbanPop Rape
## Alabama      13.2     236       58 21.2
## Alaska       10.0     263       48 44.5
## Arizona       8.1     294       80 31.0
## Arkansas      8.8     190       50 19.5
## California    9.0     276       91 40.6
## Colorado      7.9     204       78 38.7

Undersøg om der kan dannes klynger og hvorledes disse kan karakteriseres. Illustrer grafisk og kommenter på karakteristika for klyngerne.


Svar US arrestationer samt urbaniseringsgrad kort version

Spørgsmål GDP

Download filen om GDP og GDP per capita her. Importer denne i R. Sørg for at få navnene for de enkelte lande som label, dette kan du fx. gøre ved nedenstående kommandoer:

Hvad betyder: GDP[1:50,] og GDP[,1]

GDP <- GDP[1:50,]
row.names(GDP) <- as.matrix(GDP[,1])

plot(clusters,cex=0.5,main = "GDP",xlab = "Klyngetræ",sub="GDP")
Undersøg om der kan dannes klynger og hvorledes disse kan karakteriseres, der er rigtig mange observationer dvs. lande, se i stedet på deldatasæt der giver mening. Illustrer grafisk og kommenter på karakteristika for klyngerne.

Spørgsmål Forbes 100 US

Download filen om de 100 rigeste i USA her. Filen er tilrettet, dvs. binære variable er kodet om til .

Undersøg om der kan dannes klynger og hvorledes disse kan karakteriseres. Illustrer grafisk og kommenter på karakteristika for klyngerne.

Bemærk det er godt at sætte navnene på de velhavende som rækkenavne i din data.frame, for at det er nemmere at få et overblik i klyngetræ diagrammet, dette kan du fx. gøre som nedenfor:

row.names(Forbes100) <- as.matrix(Forbes100[,1])
clusters <- hclust(dist(Forbes100))
plot(clusters,cex=0.5,main = "US 100 Rigeste",xlab = "Klyngetræ",sub="US top 100")

Spørgsmål Valgfri datasæt Find et valgfrit datasæt fx. på nettet, gennemfør en klyngeanalyse på dette datasæt.