---
title: "KNN & trees: Idea day (Notes)"
subtitle: "Stat 253"
author: "Your Name"
format:
  html:
    toc: true
    toc-depth: 2
    embed-resources: true
---


```{r include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE, 
  warning = FALSE,
  message = FALSE,
  fig.height = 2.75, 
  fig.width = 4.25,
  fig.env='figure',
  fig.pos = 'h',
  fig.align = 'center')
```

# Learning Goals {.unnumbered}

-   Clearly describe the recursive binary splitting algorithm for tree building for both regression and classification
-   Compute the weighted average Gini index to measure the quality of a classification tree split
-   Compute the sum of squared residuals to measure the quality of a regression tree split
-   Explain how recursive binary splitting is a greedy algorithm
-   Explain how different tree parameters relate to the bias-variance tradeoff





\
\

# Notes: Nonparametric Classification {.unnumbered}

## Where are we? {.unnumbered .smaller}

![](https://mac-stat.github.io/STAT253/images/MLdiagram3.jpg){width="90%"}

**CONTEXT**

-   **world = supervised learning**\
    We want to model some output variable $y$ using a set of potential predictors ($x_1, x_2, ..., x_p$).

-   **task = CLASSIFICATION**\
    $y$ is categorical

-   **algorithm = NONparametric**



\

**GOAL**

Just as least squares and LASSO in the regression setting, the **parametric** logistic regression model makes very specific assumptions about the relationship of a binary categorical outcome $y$ with predictors $x$.

Specifically, it assumes this relationship can be written as a specific formula with *parameters* $\beta$:

$$\text{log(odds that y is 1)} = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + ... + \beta_k x_k$$

**NONparametric** algorithms will be necessary when this model is too rigid to capture more complicated relationships.







\

## Motivating Example {.unnumbered .smaller}

![Image Source: Google Maps](https://mac-stat.github.io/STAT253/images/Macalester_AerialPhoto.jpg){width=100%}


Aerial photography studies of land cover are important to land conservation, land management, and understanding environmental impact of land use.

> IMPORTANT: Other aerial photography studies focused on people and movement can be used for surveillance, raising major ethical questions.



Let's load data on a sample of aerial images from the [UCI Machine Learning Repository](https://archive.ics.uci.edu/ml/datasets/Urban+Land+Cover):

```{r}
#| eval: true
#| code-fold: true
# Load packages
library(tidymodels)
library(tidyverse)

# Load & process data
# There are 9 types of land use. For now, we'll only consider 3 types.
# There are 147 possible predictors of land use. For now, we'll only consider 4 predictors.
land_3 <- read.csv("https://mac-stat.github.io/data/land_cover.csv") %>% 
  rename(type = class) %>% 
  filter(type %in% c("asphalt ","grass ","tree ")) %>% 
  mutate(type = as.factor(type)) %>% 
  select(type, NDVI, Mean_G, Bright_100, SD_NIR)
```

```{r}
#| eval: true
# Check it out
head(land_3)
```

```{r}
#| eval: true
# Table of land types
land_3 %>% 
  count(type)
```



Thus we have the following variables:

-   `type` = observed type of land cover, hand-labelled by a *human* (asphalt, grass, or tree)

-   factors *computed* from the image\
    Though the data includes measurements of size and shape, we'll focus on texture and "spectral" measurements, i.e. how the land interacts with sun radiation.

    -   `NDVI` = vegetation index
    -   `Mean_G` = the green-ness of the image
    -   `Bright_100` = the brightness of the image
    -   `SD_NIR` = a texture measurement (calculated by the standard deviation of "near infrared")







\

## Limits of Logistic Regression {.unnumbered}

Why *can't* we use logistic regression to model land `type` (y) by the possible predictors `NDVI`, `Mean_G`, `Bright_100`, `SD_NIR` (x)?









\

## Parametric v. Nonparametric Classification {.unnumbered .smaller}

There *are* parametric classifications algorithms that can model y outcomes with more than 2 categories. But they're complicated and not very common.

We'll consider two **nonparametric** algorithms today: K Nearest Neighbors & classification trees.







\

## Nonparametric Classification {.unnumbered .smaller}

-   Pro: **flexibility**
    -   Assume there *is* a relationship between y and x, but don't make any assumptions about the "shape" of this relationship.
    -   This is good when our relationships are complicated!
    
-   Cons:
    -   **lack of insights:** Can be useful for classification, but provide fewer insights into the *relationships* we're modeling (eg: no coefficients, p-values, etc).
    -   **ignoring information about relationships**: When the assumptions of a parametric model are appropriate, i.e. when the shape of a relationship is "known", nonparametric algorithms will typically provide worse classifications by not utilizing that shape. (When the shape of a relationship is "known", we should use that shape.)
    -   more computationally intense








\
\

# Exercises: Part 1 (KNN) {.unnumbered}

**GOAL**

Build *intuition* for the K Nearest Neighbors classification algorithm.





\

## Exercise 1: Check out the data {.unnumbered .smaller}

We'll start by classifying land `type` using vegetation index (`NDVI`) and green-ness (`Mean_G`).

First, plot and describe this relationship.

```{r}
#| eval: true
#| code-fold: true
# Store the plot because we'll use it later
veg_green_plot <- land_3 %>%
  ggplot(aes(x = Mean_G, y = NDVI, color = type)) + 
  geom_point() +
  scale_color_manual(values = c("asphalt " = "black",
                            "grass " = "#E69F00",
                            "tree " = "#56B4E9")) +
  theme_minimal()


# Check it out
veg_green_plot
```







\

## Exercise 2: Intuition {.unnumbered .smaller}

The red dot below represents a *new* image with `NDVI` = 0.335 and `Mean_G` = 110:

```{r}
#| eval: true
#| code-fold: true
veg_green_plot + 
  geom_point(aes(y = 0.335, x = 110), color = "red")
```

How would you classify this image (asphalt, grass, or tree) using...

a.  the **1** nearest neighbor

b.  the **3** nearest neighbors

c.  all **277** neighbors in the sample






\

## Exercise 3: Details {.unnumbered .smaller}

Just as with KNN in the regression setting, it will be important to *standardize* quantitative x predictors before using them in the algorithm.
*Why*?








\

## Exercise 4: Check your work {.unnumbered .smaller}

In exercise two, your answers should be grass (a), tree (b), and grass (c).
<!--If any of your answers are different, revisit!--> Let's confirm these results by doing this "by hand" ...

Try writing R code to find the closest neighbors (without `tidymodels`).







\

## Exercise 5: Tuning the KNN algorithm {.unnumbered .smaller}

The KNN algorithm depends upon the tuning parameter $K$, the number of neighbors we consider in our classifications.
Play around with the shiny app (code below) to build some intuition here.
After, answer the following questions.

a.  In general, how would you describe the **classification regions** defined by the KNN algorithm?

b.  Let's explore the goldilocks problem here.
    When K is too small:

    -   the classification regions are very \_\_\_\_
    -   the model is (overfit or underfit?)
    -   the model will have \_\_\_\_\_ bias and \_\_\_\_ variance

c.  When K is too big:

    -   the classification regions are very \_\_\_\_
    -   the model is (overfit or underfit?)
    -   the model will have \_\_\_\_\_ bias and \_\_\_\_ variance

d.  What K value would you pick, based solely on the plots?
    (We'll eventually tune the KNN using classification accuracy as a guide.)

```{r}
#| code-fold: true
# Define KNN plot
library(gridExtra)
library(FNN)
knnplot <- function(x1, x2, y, k, lab_1, lab_2){
 x1 <- (x1 - mean(x1)) / sd(x1)
 x2 <- (x2 - mean(x2)) / sd(x2)
 x1s <- seq(min(x1), max(x1), len = 100)
 x2s <- seq(min(x2), max(x2), len = 100) 
 testdata <- expand.grid(x1s,x2s)
 knnmod <- knn(train = data.frame(x1, x2), test = testdata, cl = y, k = k, prob = TRUE)
 testdata <- testdata %>% mutate(class = knnmod)
 g1 <- ggplot(testdata, aes(x = Var1, y = Var2, color = class)) + 
     geom_point() + 
     labs(x = paste(lab_1), y = paste(lab_2), title = "KNN classification regions") + 
     theme_minimal() +
     scale_color_manual(values = c("asphalt " = "black",
                            "grass " = "#E69F00",
                            "tree " = "#56B4E9")) 
 g2 <- ggplot(NULL, aes(x = x1, y = x2, color = y)) + 
     geom_point() + 
     labs(x = paste(lab_1), y = paste(lab_2), title = "Raw data") + 
     theme_minimal() +
     scale_color_manual(values = c("asphalt " = "black",
                            "grass " = "#E69F00",
                            "tree " = "#56B4E9")) 
 grid.arrange(g1, g2)
}
```

```{r}
#| code-fold: true
#| eval: false
library(shiny)
# Build the shiny server
server_KNN <- function(input, output) {
  output$model_plot <- renderPlot({
    knnplot(x1 = land_3$Mean_G, x2 = land_3$NDVI, y = land_3$type, k = input$k_pick, lab_1 = "Mean_G (standardized)", lab_2 = "NDVI (standardized)")
  })
}

# Build the shiny user interface
ui_KNN <- fluidPage(
  sidebarLayout(
    sidebarPanel(
      h4("Pick K:"), 
      sliderInput("k_pick", "K", min = 1, max = 277, value = 1)
    ),
    mainPanel(
      plotOutput("model_plot")
    )
  )
)


# Run the shiny app!
shinyApp(ui = ui_KNN, server = server_KNN)
```








\

## PAUSE TO REFLECT: KNN pros & cons {.unnumbered .smaller}

Pros:

-   flexible

-   intuitive

-   can model y variables that have more than 2 categories

Cons:

-   "memory-based" aka **lazy learner** (technical term!)\
    In logistic regression, we run the algorithm *one time* and get a formula that we can use to classify all future data points.
    In KNN, we have to re-calculate distances and identify neighbors, hence start the algorithm over, *each time* we want to classify a new data point.

-   computationally expensive (given that we have to start over each time, and calculating the distance between every pair of points gets very time-consuming as the sample size increases)

-   provides classifications, but no real sense of the relationship of y with predictors x








\

# Exercises: Part 2 (Trees) {.unnumbered}

For the rest of class, work together to develop intuition on a new method: Classification Trees

Work together on Exercises 6--11.



\

## Exercise 6: Intuition {.unnumbered .smaller}

Classification trees are another nonparametric algorithm.

Let's build intuition for trees here.

Ask your instructor for the [paper handout](https://docs.google.com/document/d/10T6EEO73TUUZP3GbRxCkEXfoeVpG-hEc5gBmkNncwUQ/edit?usp=sharing).

Complete it using only your intuition!








\

## Exercise 7: Use the tree {.unnumbered .smaller}

Build the classification tree in R.

We'll use this tree for prediction in this exercise, and dig into the tree details in the next exercise.

```{r}
#| eval: true
#| code-fold: true
# Build the tree in R
# This is only demo code! It will change in the future.
library(rpart)
library(rpart.plot)

demo_model <- rpart(type ~ NDVI + SD_NIR, land_3, maxdepth = 2)
rpart.plot(demo_model) 
```

Predict the land type of images with the following properties:

a.  `NDVI` = 0.05 and `SD_NIR` = 7

b.  `NDVI` = -0.02 and `SD_NIR` = 7








\

## Exercise 8: Understand the tree {.unnumbered .smaller}

The classification tree is like an upside down real world tree.

The **root node** at the top starts with all 277 sample images, i.e. 100% of the data.

Among the images in this node, 21% are asphalt, 40% are grass, and 38% are tree.
Thus if we stopped our tree here, we'd classify any new image as "grass" since it's the most common category.

Let's explore the *terminal* or **leaf nodes** at the bottom of the tree.

a.  What percent of images would be classified as "tree"?

b.  Among images classified as tree:

    -   What percent are actually trees?
    -   What percent are actually grass (thus misclassified as trees)?
    -   What percent are actually asphalt (thus misclassified as trees)?







\

## Exercise 9: tuning trees: `min_n` {.unnumbered .smaller}

In the above tree, we only made 2 splits.
But we can keep growing our tree!
Run the shiny app below.

Keep `cost_complexity` at 0 and change `min_n`.
This tuning parameter controls the minimum size, or number of images, that can fall into any (terminal) leaf node.

a.  Let's explore the goldilocks problem here.
    When `min_n` is too *small* (i.e. the tree is too *big*):

    -   the classification regions are very \_\_\_\_
    -   the model is (overfit or underfit?)
    -   the model will have \_\_\_\_\_ bias and \_\_\_\_ variance

b.  When `min_n` is too *big* (i.e. the tree is too *small*):

    -   the classification regions are very \_\_\_\_
    -   the model is (overfit or underfit?)
    -   the model will have \_\_\_\_\_ bias and \_\_\_\_ variance

c.  What `min_n` value would you pick, based solely on the plots?

d.  Tuning classification trees is also referred to as "pruning".
    Why does this make sense?

e.  Check out the classification regions.
    In what way do these differ from the KNN classification regions?
    What feature of these regions reflects *binary* splitting?

```{r}
#| code-fold: true
# Define tree plot functions
library(gridExtra)
tree_plot <- function(x1, x2, y, lab_1, lab_2, cp = 0, minbucket = 1){
  model <- rpart(y ~ x1 + x2, cp = cp, minbucket = minbucket)
  x1s <- seq(min(x1), max(x1), len = 100)
  x2s <- seq(min(x2), max(x2), len = 100) 
  testdata <- expand.grid(x1s,x2s) %>% 
    mutate(type = predict(model, newdata = data.frame(x1 = Var1, x2 = Var2), type = "class"))
  g1 <- ggplot(testdata, aes(x = Var1, y = Var2, color = type)) + 
    geom_point() + 
    labs(x = paste(lab_1), y = paste(lab_2), title = "tree classification regions")  + 
     theme_minimal() +
     scale_color_manual(values = c("asphalt " = "black",
                            "grass " = "#E69F00",
                            "tree " = "#56B4E9")) + 
    theme(legend.position = "bottom")
  g2 <- ggplot(NULL, aes(x = x1, y = x2, color = y)) + 
    geom_point() + 
    labs(x = paste(lab_1), y = paste(lab_2), title = "Raw data") + 
     theme_minimal() +
     scale_color_manual(values = c("asphalt " = "black",
                            "grass " = "#E69F00",
                            "tree " = "#56B4E9")) + 
    theme(legend.position = "bottom")
  grid.arrange(g1, g2, ncol = 2)
  
}

tree_plot_2 <- function(cp, minbucket){
  model <- rpart(type ~ NDVI + SD_NIR, land_3, cp = cp, minbucket = minbucket)
  rpart.plot(model,
    box.palette = 0,
    extra = 0
  )  
}
```

```{r}
#| eval: false
#| code-fold: true
library(shiny)
# Build the shiny server
server_tree <- function(input, output) {
  output$model_plot <- renderPlot({
    tree_plot(x1 = land_3$SD_NIR, x2 = land_3$NDVI, y = land_3$type, cp = input$cp_pick, lab_1 = "SD_NIR", lab_2 = "NDVI", minbucket = input$bucket)
  })
  output$trees <- renderPlot({
    tree_plot_2(cp = input$cp_pick, minbucket = input$bucket)
  })
}

# Build the shiny user interface
ui_tree <- fluidPage(
  sidebarLayout(
    sidebarPanel(
      sliderInput("bucket", "min_n", min = 1, max = 100, value = 1),
      sliderInput("cp_pick", "cost_complexity", min = 0, max = 0.36, value = -1)
    ),
    mainPanel(
      plotOutput("model_plot"),
      plotOutput("trees")
    )
  )
)


# Run the shiny app!
shinyApp(ui = ui_tree, server = server_tree)
```








\

## Exercise 10: tuning trees: cost-complexity parameter {.unnumbered .smaller}

There's another tuning parameter to consider: `cost_complexity`!
Like the LASSO $\lambda$ *penalty* parameter penalizes the inclusion of more predictors, `cost_complexity` penalizes the introduction of new splits.
When `cost_complexity` is 0, there's no penalty -- we can make a split even if it doesn't improve our classification accuracy.
But the bigger (more positive) the `cost_complexity`, the greater the penalty -- we can only make a split if the complexity it adds to the tree is offset by its improvement to the classification accuracy.

a.  Let's explore the goldilocks problem here. When `cost_complexity` is too *small* (i.e. the tree is too *big*):
    -   the classification regions are very \_\_\_\_
    -   the model is (overfit or underfit?)
    -   the model will have \_\_\_\_\_ bias and \_\_\_\_ variance
b.  When `cost_complexity` is too *big* (i.e. the tree is too *small*):
    -   the classification regions are very \_\_\_\_
    -   the model is (overfit or underfit?)
    -   the model will have \_\_\_\_\_ bias and \_\_\_\_ variance








\

## OPTIONAL: cost complexity details {.unnumbered .smaller}

Define some notation:

-   Greek letter $\alpha$ ("alpha") = cost complexity parameter
-   T = number of terminal or leaf nodes in the tree
-   R(T) = total *mis*classification rate corresponding to the tree with T leaf nodes

Then our final tree is that which *minimizes* the combined misclassification rate and penalized number of nodes:

  $$R(T) + \alpha T $$


Hence, $\alpha$ plays a role similar to $\lambda$ in the case of the LASSO algorithm. 





\

## Exercise 11: tree properties {.unnumbered .smaller}

NOTE: If you don't get to this during class, no big deal.
These ideas will also be in our next video!

a.  We called the KNN algorithm **lazy** since we have to re-build the algorithm every time we want to make a new prediction / classification.
    Are classification trees lazy?

b.  We called the backward stepwise algorithm **greedy** -- it makes the best (local) decision in each step, but these might not end up being *globally* optimal.
    Mainly, once we kick out a predictor, we can't bring it back in even if it would be useful later.
    Explain why classification trees are greedy.

c.  In the KNN algorithm, it's important to **standardize** our quantitative predictors to the same scale.
    Is this necessary for classification trees?








\
\




# Done!

- Render your notes.
- Check the solutions in the course website (Solution drop downs).
- If you finish all that during class, start your homework!


