---
title: "TMA4268 Statistical Learning V2018"
subtitle: "Chapter 10: Solution sketches"
author: "Thiago G. Martins, Department of Mathematical Sciences, NTNU"
date: "Spring 2018"
output: #3rd letter intentation hierarchy
#  beamer_presentation:
###    incremental: true # or >* for one at a time
#  slidy_presentation:
#    font_adjustment: +1  
   prettydoc::html_pretty:
    theme: architect
    highlight: github
#  pdf_document:
  #  toc: true
  #  toc_depth: 2

---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

# Recommended exercise 1

## Get the data

The data can be downloaded from: https://www.math.ntnu.no/emner/TMA4268/2018v/data/pca-examples.rdata

```{r, message=FALSE, warning=FALSE}
# Modify the location based on your filesystem
load("../10_unsupervised_learning/datasets/pca-examples.Rdata")

# We will work with nyt.frame
nyt_data = nyt.frame
```

## Compute the PCA

```{r, message=FALSE, warning=FALSE}
nyt_pca = prcomp(nyt_data[,-1])
```

## Default biplot

Too much information on the graph, we should select only a few loading vectors to display.

```{r, message=FALSE, warning=FALSE}
biplot(nyt_pca, scale=0)
```

Lets pick some words with high PC1 weight and some words with high PC2 weight. Only looking at the graph we can see that PC1 is associated with music and PC2 with art.

```{r, message=FALSE, warning=FALSE}
nyt_loading = nyt_pca$rotation[, 1:2]
informative_loadings = rbind(
  head(nyt_loading[order(nyt_loading[,1], decreasing = TRUE),]),
  head(nyt_loading[order(nyt_loading[,2], decreasing = TRUE),])
)

biplot(x = nyt_pca$x[,1:2], y= informative_loadings, scale=0)
```

## PVE 

The numbers below show that although the graphs based on PC1 and PC2 give some insight about document types, the first two PCs explain only small portion of the variability contained in the data.

```{r, message=FALSE, warning=FALSE}
pr.var=nyt_pca$sdev^2
pr.var
```

We can then compute the proportion of variance explained by each principal component

```{r, message=FALSE, warning=FALSE}
pve=pr.var/sum(pr.var)
pve
```

We can plot the PVE explained by each component, as well
as the cumulative PVE, as follows:

```{r, message=FALSE, warning=FALSE}
plot(pve, 
     xlab="Principal Component", 
     ylab="Proportion of Variance Explained", 
     ylim=c(0,1),
     type='b')
plot(cumsum(pve), 
     xlab="Principal Component", 
     ylab="Cumulative Proportion of Variance Explained", 
     ylim=c(0,1),
     type='b')
```

# Recommended exercise 2

The answer is on page 388 of the book, with the explanation around Equation (10.12)

# Recommended exercise 3

## Perform k-means

```{r, message=FALSE, warning=FALSE}
km.out=kmeans(nyt_data[,-1],2,nstart=20)
```

## Cluster assignments

```{r, message=FALSE, warning=FALSE}
km.out$cluster
```

## Plot the data

To plot the data we need to use the PCA projections. Below we will use the plot based on PCA with true labels (A for art and M for music) and compare with the plot that color the points according to the k-means clustering.

```{r, message=FALSE, warning=FALSE}
# PCA with true labels
plot(nyt_pca$x[,1:2],type="n")
points(nyt_pca$x[nyt_data[,"class.labels"]=="art",1:2],pch="A")
points(nyt_pca$x[nyt_data[,"class.labels"]=="music",1:2],pch="M")

# PCA with true labels but colored by k-means clustering
plot(nyt_pca$x[,1:2],type="n")
points(nyt_pca$x[nyt_data[,"class.labels"]=="art",1:2],pch="A",col=(km.out$cluster+1)[nyt_data[,"class.labels"]=="art"])
points(nyt_pca$x[nyt_data[,"class.labels"]=="music",1:2],pch="M",col=(km.out$cluster+1)[nyt_data[,"class.labels"]=="music"])
```

# Recommended exercise 4

## Perform hierarchical clustering

I will use euclidean distance and complete linkage

```{r, message=FALSE, warning=FALSE}
hc.complete=hclust(dist(nyt_data[,-1]), method="complete")
str(hc.complete)
```

## Plot dendograms

```{r, message=FALSE, warning=FALSE}
plot(hc.complete,main="Complete Linkage", labels=as.character(nyt_data[,1]), xlab="", sub="", cex=.9)
```

## Cluster assignment

```{r, message=FALSE, warning=FALSE}
hc.clusters=cutree(hc.complete, 2)
```

## Plot the data

Based on the cluster assignment above and the plots below we see that hierarchical clustering performs worse than k-means for the same number of clusters (K=2)

```{r, message=FALSE, warning=FALSE}
# PCA with true labels
plot(nyt_pca$x[,1:2],type="n")
points(nyt_pca$x[nyt_data[,"class.labels"]=="art",1:2],pch="A")
points(nyt_pca$x[nyt_data[,"class.labels"]=="music",1:2],pch="M")

# PCA with true labels but colored by k-means clustering
plot(nyt_pca$x[,1:2],type="n")
points(nyt_pca$x[nyt_data[,"class.labels"]=="art",1:2],pch="A",col=(hc.clusters+1)[nyt_data[,"class.labels"]=="art"])
points(nyt_pca$x[nyt_data[,"class.labels"]=="music",1:2],pch="M",col=(hc.clusters+1)[nyt_data[,"class.labels"]=="music"])
```
