---
title: "Vins italiens et K-means"
author: Mathieu Ribatet
output: html_notebook
---

## Introduction

In this lab we will try to classify italian wines from several wine makers. The dataset has 178 wines for which 13 quantitative variables have been recorded. You will perform a K--means classification for those wines. However you'll have to wait a bit and re-work on the Fisher's iris data set to see if you understand and are able to reproduce the results we discussed during the lecture.

## Re-analyzing the Fisher's iris dataset

The data set is already embbeded into the R software. To retrieve them, just invoke

```{r}
data(iris)
head(iris)
```

Recall that K--means classification is unsupervised so you should not use the last variable---which is even not quantitative so, unless you do not understand at all the K--means, you would not have used it, wouldn't you?

1.  Start with a (short) descriptive analysis

```{r}
plot(iris[,-5], col = iris[,5])
```

1.  Read the documentation of the *kmeans* function. Please pay particular attention to the arguments *x*, *centers*, *iter.max* et *nstart*.

2.  Perform a K--means classification

    ```{r}
    iris2 <- scale(iris[,-5])
    firstClassif <- kmeans(iris2, 5, nstart = 100, iter.max = 500)
    firstClassif
    ```

3.  Interpret results

4.  Vary the number of classes and get the "optimal" number of classes (if you're lost elbow rule)

    ```{r}
    metric <- rep(NA, 10)
    metric[1] <- 1
    for (nclust in 2:10){
      tmp <- kmeans(iris2, nclust, nstart = 100, iter.max = 50)
      metric[nclust] <- tmp$tot.withinss / tmp$totss
    }

    plot(1:10, metric, type = "h", xlab = "#clusters",
         ylab = expression(W(x,pi)/I(x)), ylim = c(0, 1))

    ## On garde 5 clusters
    fit <- kmeans(iris2, 5, nstart = 100, iter.max = 50)
    ```

5.  Create some fictious flowers and predict its species. (Eventually create a function dedicated to that purpose)

    ```{r}
    new_flowers <- matrix(rnorm(20 * 4, 2.5, 3), 20, 4)

    ## Je centre et reduis ces données avec les valeurs
    ## issues de iris !!!!
    center <- attributes(iris2)$`scaled:center`
    scale <- attributes(iris2)$`scaled:scale`
    new_flowers2 <- scale(new_flowers,
                          center = center, scale = scale)

    ## Je me fais une fonction qui trouve le centroie le plus proche d'une fleur
    whichCentroid <- function(newx, centroids){
      which.min(colSums((t(centroids) - newx)^2))
    }

    pred <- apply(new_flowers2, 1, whichCentroid, centroids = fit$centers)

    pairs(new_flowers, col = pred)
    ```

6.  Since we do have access to the *Species* for each flower (5th column), compute the confusion matrix. The *table* function might help.

    ```{r}
    table(iris$Species, fit$cluster)
    ```

    Il faudrait réorganiser les "étiquettes" du classifieur pour qu'on est une diagonale meilleure.

    ```{r}
    idx <- which(fit$cluster == 2)
    sub_data <- iris2[idx,]
    dim(sub_data)

    fit2 <- kmeans(sub_data, 2, nstart = 100, iter.max = 50)
    table(iris[idx,"Species"], fit2$cluster)
    ```

7.  Optional. I wrote an incomplete piece of code which is a (re)implementation of the LLoyd alogorithm. Fill in gaps.

```{r}
mykmeans <- function(data, n.class, iter.max = 50){
  ## data : matrix of the data
  ## n.class : number of classes
  ## iter.max : maximal number of iterations
  
  n.obs <- nrow(data)
  data <- as.matrix(data)## convert a data frame to matrix (if needed)
  
  ## Centroid initialization (taken randomly without replacement into individuals)
  centroids <- data[sample(n.obs, n.class),]
  
  for (i in 1:iter.max){
    ## Compute the distance for each individuals to each centroid
    dist <- matrix(NA, n.obs, n.class)
    for (class in 1:n.class)
      
      dist[,class] <- colSums((t(data) - centroids[class,])^2)## nasty piece of code that I can explain
      
    ## Assign labels to each individuals
    cluster <- 1#### To be filled
      
    ## Centroids update
    centroids <- 1## To be filled
    
    if (a_stopping_rule_to_define)
      break
  }
  
  ## Compute within and between inertia
  within <- 1## To be filled
  between <- 1## To be filled
  
  if (iter == iter.max)
    warning("Maximal number of iterations has been reached :-(")
  
  return(list(cluster = cluster, within = within, between = between, iter = iter, centroids = centroids))
}
```

## Application: Italian wines

The data can be retrieved from [here](https://mribatet.perso.math.cnrs.fr/CentraleNantes/BioSTIC/ItalianWine.RData). Start with a descriptive analysis of the data and, next, perform a K--means classification. Comment results. Perform predictions on the following new wines

```{r}
new.wines <- rbind(wine1 = c(12.37, 1.17, 1.92, 19.6, 78, 2.11, 2.00, 0.27, 1.04, 4.68, 1.12, 3.48, 510),
                  wine2 = c(13.45, 3.70, 2.60, 23.0, 111, 1.70, 0.92, 0.43, 1.46, 10.68, 0.85, 1.56, 695))
```

Good luck !
