Application de DBSCAN

Données

Nous allons continuer à utiliser le jeu de données “Spotify” (csv à télécharger si nécessaire). Pour cela, nous allons les charger, et utiliser les packages readr (chargement de données), tidyr (rotation de données) et dplyr (manipulation de données). Ces trois packages proviennent du tidyverse (en plus de ggplot2 - réalisation de graphiques), ensemble de packages très puissants pour l’analyse de données.

library(readr)
library(tidyr)
library(dplyr)
library(ggplot2)

songs = read_csv("spotify-2017songs.csv")
songs_num = songs %>%
  select(acousticness, danceability, energy, instrumentalness, liveness, speechiness, valence) %>%
  # mutate_all(scale) %>%
  tibble()

Recherche des paramètres

Nous sommes en dimension 7, la proposition de base est donc de choisir \(MinPts = 14\). Nous allons nous basé sur cette option.

Pour \(\varepsilon\), nous devons créer le graphique des \(k\text{-}dist\) (avec \(k=14\)) de chaque point, triées dans l’ordre décroissant. Ici, nous avons 2017 individus, ce qui reste raisonnable. Dans le cas de données de tailles plus importantes, nous aurions dû procédé à un échantillonage (voire plusieurs) pour chercher cette valeur.

Nous calculons d’abord la matrice de distance :

d = dist(songs_num)

##         1      2      3      4      5      6      7      8      9
## 2  0.5150                                                        
## 3  0.1858 0.5056                                                 
## 4  0.9420 0.7832 0.8860                                          
## 5  0.9471 0.7018 0.9959 0.9164                                   
## 6  0.2797 0.4458 0.2075 0.8922 0.8985                            
## 7  0.2977 0.3688 0.2281 0.8495 0.8472 0.1292                     
## 8  0.9528 0.8560 0.9402 0.6693 0.7770 0.9022 0.8388              
## 9  0.6049 0.7310 0.6526 1.0686 0.8800 0.5154 0.5525 1.0196       
## 10 0.4915 0.6008 0.5053 1.0525 0.7965 0.4314 0.4146 1.0198 0.4872

Ensuite, nous la transformons en matrice (pour la manipuler ensuite)

m = as.matrix(d)

##         1      2      3      4      5      6      7      8      9     10
## 1  0.0000 0.5150 0.1858 0.9420 0.9471 0.2797 0.2977 0.9528 0.6049 0.4915
## 2  0.5150 0.0000 0.5056 0.7832 0.7018 0.4458 0.3688 0.8560 0.7310 0.6008
## 3  0.1858 0.5056 0.0000 0.8860 0.9959 0.2075 0.2281 0.9402 0.6526 0.5053
## 4  0.9420 0.7832 0.8860 0.0000 0.9164 0.8922 0.8495 0.6693 1.0686 1.0525
## 5  0.9471 0.7018 0.9959 0.9164 0.0000 0.8985 0.8472 0.7770 0.8800 0.7965
## 6  0.2797 0.4458 0.2075 0.8922 0.8985 0.0000 0.1292 0.9022 0.5154 0.4314
## 7  0.2977 0.3688 0.2281 0.8495 0.8472 0.1292 0.0000 0.8388 0.5525 0.4146
## 8  0.9528 0.8560 0.9402 0.6693 0.7770 0.9022 0.8388 0.0000 1.0196 1.0198
## 9  0.6049 0.7310 0.6526 1.0686 0.8800 0.5154 0.5525 1.0196 0.0000 0.4872
## 10 0.4915 0.6008 0.5053 1.0525 0.7965 0.4314 0.4146 1.0198 0.4872 0.0000

Maintenant, pour chaque ligne (ou colonne, résultat identique), on trie les valeurs par ordre croissants.

m_tri = apply(m, 2, sort)

##            1      2      3      4      5      6      7      8      9     10
##  [1,] 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
##  [2,] 0.0000 0.1163 0.1749 0.2487 0.2604 0.1278 0.1172 0.2245 0.1514 0.1071
##  [3,] 0.1569 0.1228 0.1815 0.2517 0.3353 0.1292 0.1268 0.3206 0.1745 0.2522
##  [4,] 0.1754 0.1761 0.1858 0.2795 0.3929 0.1396 0.1292 0.3332 0.2140 0.2532
##  [5,] 0.1763 0.1796 0.1858 0.2826 0.3970 0.1434 0.1363 0.3934 0.2165 0.2579
##  [6,] 0.1858 0.1937 0.1995 0.2837 0.3980 0.1448 0.1520 0.4198 0.2231 0.2606
##  [7,] 0.2306 0.2020 0.2075 0.2926 0.4058 0.1451 0.1522 0.4452 0.2366 0.2654
##  [8,] 0.2321 0.2032 0.2086 0.3401 0.4092 0.1460 0.1541 0.4591 0.2425 0.2694
##  [9,] 0.2360 0.2045 0.2206 0.3405 0.4246 0.1473 0.1603 0.4608 0.2491 0.2789
## [10,] 0.2423 0.2056 0.2281 0.3615 0.4312 0.1504 0.1620 0.4662 0.2546 0.2907
## [11,] 0.2543 0.2066 0.2390 0.3717 0.4313 0.1605 0.1669 0.4764 0.2557 0.3009
## [12,] 0.2605 0.2129 0.2403 0.3801 0.4350 0.1628 0.1683 0.4786 0.2575 0.3019
## [13,] 0.2632 0.2136 0.2451 0.3826 0.4383 0.1652 0.1696 0.4836 0.2607 0.3077
## [14,] 0.2647 0.2166 0.2482 0.3878 0.4418 0.1759 0.1783 0.4865 0.2616 0.3088
## [15,] 0.2797 0.2173 0.2527 0.3904 0.4453 0.1769 0.1795 0.4883 0.2626 0.3128
## [16,] 0.2881 0.2191 0.2549 0.3973 0.4525 0.1781 0.1827 0.4974 0.2712 0.3129
## [17,] 0.2912 0.2212 0.2559 0.3996 0.4527 0.1818 0.1836 0.4979 0.2726 0.3145
## [18,] 0.2975 0.2220 0.2607 0.4005 0.4620 0.1820 0.1857 0.4990 0.2729 0.3168
## [19,] 0.2977 0.2236 0.2618 0.4098 0.4641 0.1854 0.1917 0.5098 0.2755 0.3226
## [20,] 0.3035 0.2329 0.2636 0.4147 0.4645 0.1870 0.2004 0.5111 0.2761 0.3254

Ensuite, on récupère la \(k+1\)ième valeur la plus petite (\(+1\) car on aura le 0), donc la \(k+1\)ième ligne.

v = m_tri[15,]

##      1      2      3      4      5      6      7      8      9     10 
## 0.2797 0.2173 0.2527 0.3904 0.4453 0.1769 0.1795 0.4883 0.2626 0.3128

On créé un data.frame (plus exactement un tibble pour rester dans l’environnement tidyverse) contenant ces valeurs triées dans l’ordre décroissant.

df = tibble(n = 1:2017, kdist = sort(v, decreasing = TRUE)) %>%
  mutate(diff = c(NA, -diff(kdist)))

On va pouvoir afficher le graphique des \(k\text{-}dist\) triées par ordre décroissant, et chercher un coude.

coeff = 50
ggplot(df) +
  geom_vline(xintercept = c(170, 247), alpha = 0.25, lty = 2) +
  geom_hline(yintercept = df$kdist[c(170, 247)], alpha = 0.25, lty = 2) +
  geom_line(aes(n, diff * coeff), color = "red", alpha = 0.5) +
  geom_line(aes(n, kdist)) +
  scale_y_continuous(name = "k-dist (k = 14)", sec.axis = sec_axis( trans=~. / coeff, name="Différence (en rouge)")) +
  theme_minimal()

On remarque 2 pics en terme de différence, vers \(\varepsilon=0.343\) et \(\varepsilon=0.314\) (on testera probablement d’autres valeurs).

Nous pouvons chercher avec une valeur différente de \(MinPts\) (\(k\)). Testez différentes valeurs (ici \(k=5\)).

k = 5
ggplot(tibble(n = 1:2017, kdist = sort(m_tri[k+1,], decreasing = TRUE)) %>%
         mutate(diff = c(NA, -diff(kdist)))) +
  geom_line(aes(n, kdist)) +
  labs(y = paste0("k-dist avec k=",k)) +
  theme_minimal()

Recherche de la partition

On applique donc la mthode DBSCAN avec \(\varepsilon=0.343\) et \(MinPts=14\) dans un premier temps. On utilise ici la fonction dbscan du package dbscan.

library(dbscan)
res = dbscan(songs_num, eps = .343, minPts = 14)
res

## DBSCAN clustering for 2017 objects.
## Parameters: eps = 0.343, minPts = 14
## The clustering contains 1 cluster(s) and 35 noise points.
## 
##    0    1 
##   35 1982 
## 
## Available fields: cluster, eps, minPts

On voit tout de suite que notre choix de \(\varepsilon=2\) était mauvais. Nous obtenons une seule classe, avec 35 individus considérés comme outliers. Nous allons tester l’autre valeur \(\varepsilon=0.314\).

res = dbscan(songs_num, eps = .314, minPts = 14)
res

## DBSCAN clustering for 2017 objects.
## Parameters: eps = 0.314, minPts = 14
## The clustering contains 1 cluster(s) and 62 noise points.
## 
##    0    1 
##   62 1955 
## 
## Available fields: cluster, eps, minPts

On trouve encore une seule classe. On peut modifier aussi \(MinPts\) pour avoir plus de classes. On va le doubler.

res = dbscan(songs_num, eps = .314, minPts = 28)
res

## DBSCAN clustering for 2017 objects.
## Parameters: eps = 0.314, minPts = 28
## The clustering contains 2 cluster(s) and 108 noise points.
## 
##    0    1    2 
##  108 1844   65 
## 
## Available fields: cluster, eps, minPts

Après quelques tests, on peut affiner un peu en choisissant \(\varepsilon=0.32\) et \(MinPts=42\).

res = dbscan(songs_num, eps = 0.32, minPts = 42)
res

## DBSCAN clustering for 2017 objects.
## Parameters: eps = 0.32, minPts = 42
## The clustering contains 3 cluster(s) and 202 noise points.
## 
##    0    1    2    3 
##  202 1679   79   57 
## 
## Available fields: cluster, eps, minPts

Nous obtenons ici 3 classes (une très grande et 2 petites), ainsi que 202 points outliers. Voici la caractérisation des classes :

centres = songs_num %>%
  mutate(classe = res$cluster) %>%
  group_by(classe) %>%
  summarise_all(mean)

Représentation sur l’ACP

En clustering, il est très courant d’utiliser l’ACP pour représenter les individus et la partition obtenue.

library(FactoMineR)
acp = PCA(songs_num %>%
  mutate(classe = res$cluster), scale.unit = FALSE, quali.sup = 8)

Les deux premiers facteurs représentent presque 50% de l’information totale du nuage de points, ce qui est beaucoup dans un cadre réel. Les variables ne sont pour autant pas très corrélées aux axes. Toutefois, on peut dire les élémetns suivants :

l’axe 1 oppose les chansons énergiques et joyeuses (à droite du graphique) aux chansons acoustiques et plus instrumentales (à gauche donc) ;
l’axe 2 oppose les chansons acoustiques (en bas sur le graphique) aux chansons instrumentales (en haut).

acp_ind = tibble(data.frame(acp$ind$coord)) %>% 
         mutate(classe = res$cluster)
ggplot(acp_ind, aes(Dim.1, Dim.2, color = factor(classe))) +
  geom_point() +
  theme_minimal()

La représentation graphique nous montre deux phénomènes intéressants dans l’étude de nos données :

Elles sont cantonnées dans un carré, ce qui est logique car comprise entre 0 et 1 (et nous avons fait le choix de ne pas les standardiser ni pour DBSCAN, ni pour l’ACP) ;
Il y une masse de chansons très similaires (en vert à droite donc), et d’autres éparpillées dans l’espace sans séparation claire (il y a un continuum).

C’est ce dernier phénomène qui tend à regrouper tous les individus dans une même classe avec DBSCAN.

Application de DBSCAN

FX Jollois

11 mars 2022

Données

Recherche des paramètres

Recherche de la partition

Représentation sur l’ACP