Unsupervised Learning

[1]:

library(tidyverse)

Registered S3 methods overwritten by 'ggplot2':
  method         from
  [.quosures     rlang
  c.quosures     rlang
  print.quosures rlang
── Attaching packages ─────────────────────────────────────── tidyverse 1.2.1 ──
✔ ggplot2 3.1.1     ✔ purrr   0.3.2
✔ tibble  2.1.2     ✔ dplyr   0.8.1
✔ tidyr   0.8.3     ✔ stringr 1.4.0
✔ readr   1.3.1     ✔ forcats 0.4.0
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()

[2]:

head(iris)
A data.frame: 6 × 5
Sepal.LengthSepal.WidthPetal.LengthPetal.WidthSpecies
<dbl><dbl><dbl><dbl><fct>
5.13.51.40.2setosa
4.93.01.40.2setosa
4.73.21.30.2setosa
4.63.11.50.2setosa
5.03.61.40.2setosa
5.43.91.70.4setosa

Review of tidyverse

1. How many rows and columns are there in the DataFrame?

[ ]:

2. Show just the Species column and any columns with information about Petal.

[ ]:

3. Show a DataFrame with just two columns showing the Length: Width ratio for Sepal and Petal.

[ ]:

4. Count the number of each species where Sepal.Length is less than 6.

[ ]:

5. Summarize the mean and standard deviation for each measurement grouped by Species.

[ ]:

6. Convert iris to “tall” form and assign the resulting DataFrame to iris_t with three columns Species, Measurement and Value.

[ ]:

Pairwise scatter plot

[3]:

pairs(iris)
../_images/cliburn_R04_Unsupervised_Learning_17_0.png 
[4]:

pairs(iris, col=iris$Species)
../_images/cliburn_R04_Unsupervised_Learning_18_0.png

Calculating distances

Mathematical conditions for a distance function - positivity, symmetry, triangle inequality.

[5]:

small <- iris %>% select(-Species) %>% sample_n(6)
small
A data.frame: 6 × 4
Sepal.LengthSepal.WidthPetal.LengthPetal.Width
<dbl><dbl><dbl><dbl>
5.93.24.81.8
5.03.41.60.4
4.83.01.40.3
5.52.54.01.3
7.13.05.92.1
6.73.14.71.5

Euclidean distance

[6]:

iris %>% select(-Species) -> df
x <- df[1,]
y <- df[2,]
sqrt(sum((y-x)^2))
0.53851648071345

Distance matrix

[7]:

dist(small)

          1         2         3         4         5
2 3.6124784
3 3.8807216 0.5000000
4 1.2409674 2.7622455 2.9154759
5 1.6673332 5.0941143 5.3646994 2.6570661
6 0.8660254 3.7148351 3.9937451 1.5264338 1.4035669

[8]:

dist(small, upper = T, diag=T)

          1         2         3         4         5         6
1 0.0000000 3.6124784 3.8807216 1.2409674 1.6673332 0.8660254
2 3.6124784 0.0000000 0.5000000 2.7622455 5.0941143 3.7148351
3 3.8807216 0.5000000 0.0000000 2.9154759 5.3646994 3.9937451
4 1.2409674 2.7622455 2.9154759 0.0000000 2.6570661 1.5264338
5 1.6673332 5.0941143 5.3646994 2.6570661 0.0000000 1.4035669
6 0.8660254 3.7148351 3.9937451 1.5264338 1.4035669 0.0000000

Scaling before distance

[9]:

dist(scale(small), upper=T, diag=T)

         1        2        3        4        5        6
1 0.000000 2.837472 3.074112 2.498835 1.633162 1.017358
2 2.837472 0.000000 1.357298 3.525992 4.217938 3.080620
3 3.074112 1.357298 0.000000 2.680814 4.270058 3.203481
4 2.498835 3.525992 2.680814 0.000000 2.834187 2.427062
5 1.633162 4.217938 4.270058 2.834187 0.000000 1.179145
6 1.017358 3.080620 3.203481 2.427062 1.179145 0.000000

What is the distance matrix is showing?

[10]:

dist(small, method = "maximum", diag=T, upper=T)

    1   2   3   4   5   6
1 0.0 3.2 3.4 0.8 1.2 0.8
2 3.2 0.0 0.4 2.4 4.3 3.1
3 3.4 0.4 0.0 2.6 4.5 3.3
4 0.8 2.4 2.6 0.0 1.9 1.2
5 1.2 4.3 4.5 1.9 0.0 1.2
6 0.8 3.1 3.3 1.2 1.2 0.0

Agglomerative hierarchical clustering

[11]:

iris %>% select(-Species) -> df

[12]:

c1 <- hclust(dist(df))

[13]:

plot(c1)
../_images/cliburn_R04_Unsupervised_Learning_33_0.png 
[14]:

z1 <- cutree(c1, 3)

[15]:

z1
  1. 1
  2. 1
  3. 1
  4. 1
  5. 1
  6. 1
  7. 1
  8. 1
  9. 1
  10. 1
  11. 1
  12. 1
  13. 1
  14. 1
  15. 1
  16. 1
  17. 1
  18. 1
  19. 1
  20. 1
  21. 1
  22. 1
  23. 1
  24. 1
  25. 1
  26. 1
  27. 1
  28. 1
  29. 1
  30. 1
  31. 1
  32. 1
  33. 1
  34. 1
  35. 1
  36. 1
  37. 1
  38. 1
  39. 1
  40. 1
  41. 1
  42. 1
  43. 1
  44. 1
  45. 1
  46. 1
  47. 1
  48. 1
  49. 1
  50. 1
  51. 2
  52. 2
  53. 2
  54. 3
  55. 2
  56. 3
  57. 2
  58. 3
  59. 2
  60. 3
  61. 3
  62. 3
  63. 3
  64. 2
  65. 3
  66. 2
  67. 3
  68. 3
  69. 2
  70. 3
  71. 2
  72. 3
  73. 2
  74. 2
  75. 2
  76. 2
  77. 2
  78. 2
  79. 2
  80. 3
  81. 3
  82. 3
  83. 3
  84. 2
  85. 3
  86. 2
  87. 2
  88. 2
  89. 3
  90. 3
  91. 3
  92. 2
  93. 3
  94. 3
  95. 3
  96. 3
  97. 3
  98. 2
  99. 3
  100. 3
  101. 2
  102. 2
  103. 2
  104. 2
  105. 2
  106. 2
  107. 3
  108. 2
  109. 2
  110. 2
  111. 2
  112. 2
  113. 2
  114. 2
  115. 2
  116. 2
  117. 2
  118. 2
  119. 2
  120. 2
  121. 2
  122. 2
  123. 2
  124. 2
  125. 2
  126. 2
  127. 2
  128. 2
  129. 2
  130. 2
  131. 2
  132. 2
  133. 2
  134. 2
  135. 2
  136. 2
  137. 2
  138. 2
  139. 2
  140. 2
  141. 2
  142. 2
  143. 2
  144. 2
  145. 2
  146. 2
  147. 2
  148. 2
  149. 2
  150. 2

Note that the label values are arbitrary - all we know is that the 1s belong to the same cluster, the 2s belong a another cluster, and the 3s belong to the final cluster. We have no idea what Species thee cluster labels represent, or even if the assignment is “correct” compared to the ground truth.

[16]:

par(mfrow=c(1,2))
plot(iris[, 1], iris[,2], col=iris$Species)
plot(iris[, 1], iris[,2], col=z1)
../_images/cliburn_R04_Unsupervised_Learning_37_0.png

Different linkage methods can give different cluster assignments

[17]:

c2 <- hclust(dist(df), method="average")

[18]:

plot(c2)
../_images/cliburn_R04_Unsupervised_Learning_40_0.png 
[19]:

z2 <- cutree(c2, k=3)

[20]:

par(mfrow=c(1,2))
plot(iris[, 1], iris[,2], col=iris$Species)
plot(iris[, 1], iris[,2], col=z2)
../_images/cliburn_R04_Unsupervised_Learning_42_0.png

But in reality, we would not know the true number of clusters!

[21]:

plot(c2)
abline(h=1.5, col='red')
../_images/cliburn_R04_Unsupervised_Learning_44_0.png

Slightly more informative plot

[22]:

plot(c2)
rect.hclust(c2, h=1.5, border="red")
../_images/cliburn_R04_Unsupervised_Learning_46_0.png 
[23]:

z3 <- cutree(c2, h=1.5)

[24]:

par(mfrow=c(1,2))
plot(iris[, 1], iris[,2], col=iris$Species)
plot(iris[, 1], iris[,2], col=z3)
../_images/cliburn_R04_Unsupervised_Learning_48_0.png

K-means clustering

[25]:

k1 <- kmeans(dist(df), centers=3)

[26]:

str(k1)

List of 9
 $ cluster     : Named int [1:150] 2 2 2 2 2 2 2 2 2 2 ...
  ..- attr(*, "names")= chr [1:150] "1" "2" "3" "4" ...
 $ centers     : num [1:3, 1:150] 5.095 0.5 3.464 5.146 0.622 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:3] "1" "2" "3"
  .. ..$ : chr [1:150] "1" "2" "3" "4" ...
 $ totss       : num 55479
 $ withinss    : num [1:3] 1745 495 2593
 $ tot.withinss: num 4833
 $ betweenss   : num 50646
 $ size        : int [1:3] 38 50 62
 $ iter        : int 2
 $ ifault      : int 0
 - attr(*, "class")= chr "kmeans"

[27]:

as.vector(k1$cluster)
  1. 2
  2. 2
  3. 2
  4. 2
  5. 2
  6. 2
  7. 2
  8. 2
  9. 2
  10. 2
  11. 2
  12. 2
  13. 2
  14. 2
  15. 2
  16. 2
  17. 2
  18. 2
  19. 2
  20. 2
  21. 2
  22. 2
  23. 2
  24. 2
  25. 2
  26. 2
  27. 2
  28. 2
  29. 2
  30. 2
  31. 2
  32. 2
  33. 2
  34. 2
  35. 2
  36. 2
  37. 2
  38. 2
  39. 2
  40. 2
  41. 2
  42. 2
  43. 2
  44. 2
  45. 2
  46. 2
  47. 2
  48. 2
  49. 2
  50. 2
  51. 3
  52. 3
  53. 3
  54. 3
  55. 3
  56. 3
  57. 3
  58. 3
  59. 3
  60. 3
  61. 3
  62. 3
  63. 3
  64. 3
  65. 3
  66. 3
  67. 3
  68. 3
  69. 3
  70. 3
  71. 3
  72. 3
  73. 3
  74. 3
  75. 3
  76. 3
  77. 3
  78. 1
  79. 3
  80. 3
  81. 3
  82. 3
  83. 3
  84. 3
  85. 3
  86. 3
  87. 3
  88. 3
  89. 3
  90. 3
  91. 3
  92. 3
  93. 3
  94. 3
  95. 3
  96. 3
  97. 3
  98. 3
  99. 3
  100. 3
  101. 1
  102. 3
  103. 1
  104. 1
  105. 1
  106. 1
  107. 3
  108. 1
  109. 1
  110. 1
  111. 1
  112. 1
  113. 1
  114. 3
  115. 1
  116. 1
  117. 1
  118. 1
  119. 1
  120. 3
  121. 1
  122. 3
  123. 1
  124. 3
  125. 1
  126. 1
  127. 3
  128. 3
  129. 1
  130. 1
  131. 1
  132. 1
  133. 1
  134. 3
  135. 1
  136. 1
  137. 1
  138. 1
  139. 3
  140. 1
  141. 1
  142. 1
  143. 3
  144. 1
  145. 1
  146. 1
  147. 3
  148. 1
  149. 1
  150. 3

Finding cluster means

Base R

These incantations are hard to remember. I suggest you stick to tidyverse methods.

[28]:

aggregate(df,by=list(k1$cluster),FUN=mean)
A data.frame: 3 × 5
Group.1Sepal.LengthSepal.WidthPetal.LengthPetal.Width
<int><dbl><dbl><dbl><dbl>
16.8210533.0657895.7473682.094737
25.0060003.4280001.4620000.246000
35.9193552.7532264.3903231.419355

dplyr

[29]:

df %>%
mutate(z = k1$cluster)  %>%
group_by(z) %>%
summarise_all(mean) -> centroids

[30]:

centroids
A tibble: 3 × 5
zSepal.LengthSepal.WidthPetal.LengthPetal.Width
<int><dbl><dbl><dbl><dbl>
16.8210533.0657895.7473682.094737
25.0060003.4280001.4620000.246000
35.9193552.7532264.3903231.419355
 
[31]:

par(mfrow=c(1,2))
plot(iris[,"Sepal.Length"], iris[,"Sepal.Width"])
points(centroids[["Sepal.Length"]], centroids[["Sepal.Width"]], col="red",  pch="x")
plot(iris[,"Sepal.Length"], iris[,"Sepal.Width"], col=k1$cluster)
../_images/cliburn_R04_Unsupervised_Learning_59_0.png

Dimension reduction and ordination

PCA

[32]:

pc <- prcomp(df, center=T, scale=T)

[33]:

str(pc)

List of 5
 $ sdev    : num [1:4] 1.708 0.956 0.383 0.144
 $ rotation: num [1:4, 1:4] 0.521 -0.269 0.58 0.565 -0.377 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
  .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
 $ center  : Named num [1:4] 5.84 3.06 3.76 1.2
  ..- attr(*, "names")= chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
 $ scale   : Named num [1:4] 0.828 0.436 1.765 0.762
  ..- attr(*, "names")= chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width"
 $ x       : num [1:150, 1:4] -2.26 -2.07 -2.36 -2.29 -2.38 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : NULL
  .. ..$ : chr [1:4] "PC1" "PC2" "PC3" "PC4"
 - attr(*, "class")= chr "prcomp"

[34]:

summary(pc)

Importance of components:
                          PC1    PC2     PC3     PC4
Standard deviation     1.7084 0.9560 0.38309 0.14393
Proportion of Variance 0.7296 0.2285 0.03669 0.00518
Cumulative Proportion  0.7296 0.9581 0.99482 1.00000

[35]:

plot(pc$x[,1], pc$x[,2])
../_images/cliburn_R04_Unsupervised_Learning_65_0.png

Classical MDS

[36]:

mds <- cmdscale(dist(df), k = 2)

[37]:

str(mds)

 num [1:150, 1:2] -2.68 -2.71 -2.89 -2.75 -2.73 ...
 - attr(*, "dimnames")=List of 2
  ..$ : NULL
  ..$ : NULL

[38]:

summary(mds)

      V1                V2
Min.   :-3.2238   Min.   :-1.26597
1st Qu.:-2.5303   1st Qu.:-0.32542
Median : 0.5546   Median :-0.02216
Mean   : 0.0000   Mean   : 0.00000
3rd Qu.: 1.5501   3rd Qu.: 0.32492
Max.   : 3.7956   Max.   : 1.37417

[39]:

plot(mds[,1], mds[,2])
../_images/cliburn_R04_Unsupervised_Learning_70_0.png

Heatmaps

[40]:

library(pheatmap)

[41]:

pheatmap(scale(df))
../_images/cliburn_R04_Unsupervised_Learning_73_0.png 
[42]:

pheatmap(df, kmeans_k = 3)
../_images/cliburn_R04_Unsupervised_Learning_74_0.png 
[43]:

pheatmap(dist(df))
../_images/cliburn_R04_Unsupervised_Learning_75_0.png 
[ ]: