# 2.3 Least Squares and Nearest Neighbors

### 2.3.3 From Least Squares to Nearest Neighbors

1. Generates 10 means $m_k$ from a bivariate Gaussian distrubition for each color:
   - $N((1, 0)^T, \textbf{I})$ for <span style="color: blue">BLUE</span>
   - $N((0, 1)^T, \textbf{I})$ for <span style="color: orange">ORANGE</span>
2. For each color generates 100 observations as following:
   - For each observation it picks $m_k$ at random with probability 1/10.
   - Then generates a $N(m_k,\textbf{I}/5)$



```python
%matplotlib inline

import random
import numpy as np
import matplotlib.pyplot as plt

sample_size = 100

def generate_data(size, mean):
    identity = np.identity(2)
    m = np.random.multivariate_normal(mean, identity, 10)
    return np.array([
        np.random.multivariate_normal(random.choice(m), identity / 5)
        for _ in range(size)
    ])

def plot_data(orange_data, blue_data): 
    axes.plot(orange_data[:, 0], orange_data[:, 1], 'o', color='orange')
    axes.plot(blue_data[:, 0], blue_data[:, 1], 'o', color='blue')
    
blue_data = generate_data(sample_size, [1, 0])
orange_data = generate_data(sample_size, [0, 1])

data_x = np.r_[blue_data, orange_data]
data_y = np.r_[np.zeros(sample_size), np.ones(sample_size)]

# plotting
fig = plt.figure(figsize = (8, 8))
axes = fig.add_subplot(1, 1, 1)
plot_data(orange_data, blue_data)

plt.show()
```

### 2.3.1 Linear Models and Least Squares

$$\hat{Y} = \hat{\beta_0} + \sum_{j=1}^{p} X_j\hat{\beta_j}$$

where $\hat{\beta_0}$ is the intercept, also know as the *bias*. It is convenient to include the constant variable 1 in X and $\hat{\beta_0}$ in the vector of coefficients $\hat{\beta}$, and then write as: 

$$\hat{Y} = X^T\hat{\beta} $$

#### Residual sum of squares

How to fit the linear model to a set of training data? Pick the coefficients $\beta$ to minimize the *residual sum of squares*:

$$RSS(\beta) = \sum_{i=1}^{N} (y_i - x_i^T\beta) ^ 2 = (\textbf{y} - \textbf{X}\beta)^T (\textbf{y} - \textbf{X}\beta)$$

where $\textbf{X}$ is an $N \times p$ matrix with each row an input vector, and $\textbf{y}$ is an N-vector of the outputs in the training set. Differentiating w.r.t. β we get the normal equations:

$$\mathbf{X}^T(\mathbf{y} - \mathbf{X}\beta) = 0$$

If $\mathbf{X}^T\mathbf{X}$ is nonsingular, then the unique solution is given by: 

$$\hat{\beta} = (\mathbf{X}^T\mathbf{X})^{-1}\mathbf{X}^T\mathbf{y}$$




```python
class LinearRegression:
    def fit(self, X, y):
        X = np.c_[np.ones((X.shape[0], 1)), X]
        self.beta = np.linalg.inv(X.T @ X) @ X.T @ y

        return self
    
    def predict(self, x):
        return np.dot(self.beta, np.r_[1, x])

model = LinearRegression().fit(data_x, data_y)
print("beta = ", model.beta)
```

    beta =  [ 0.52677771 -0.15145005  0.15818643]


#### Example of the linear model in a classification context

The fitted values $\hat{Y}$ are converted to a fitted class variable $\hat{G}$ according to the rule:

$$
\begin{equation}
\hat{G} = \begin{cases}
\text{ORANGE} & \text{ if } \hat{Y} \gt 0.5 \\ 
\text{BLUE    } & \text{ if } \hat{Y} \leq 0.5 
\end{cases}
\end{equation}
$$


```python
from itertools import filterfalse, product

def plot_grid(orange_grid, blue_grid):
    axes.plot(orange_grid[:, 0], orange_grid[:, 1], '.', zorder = 0.001,
              color='orange', alpha = 0.3, scalex = False, scaley = False)

    axes.plot(blue_grid[:, 0], blue_grid[:, 1], '.', zorder = 0.001,
          color='blue', alpha = 0.3, scalex = False, scaley = False)

plot_xlim = axes.get_xlim()
plot_ylim = axes.get_ylim()

grid = np.array([*product(np.linspace(*plot_xlim, 50), np.linspace(*plot_ylim, 50))])

is_orange = lambda x: model.predict(x) > 0.5

orange_grid = np.array([*filter(is_orange, grid)])
blue_grid = np.array([*filterfalse(is_orange, grid)])

axes.clear()
axes.set_title("Linear Regression of 0/1 Response")
plot_data(orange_data, blue_data)
plot_grid(orange_grid, blue_grid)

find_y = lambda x: (0.5 - model.beta[0] - x * model.beta[1]) / model.beta[2]
axes.plot(plot_xlim, [*map(find_y, plot_xlim)], color = 'black', 
          scalex = False, scaley = False)


fig
```

### 2.3.2 Nearest-Neighbor Methods

$$\hat{Y}(x) = \frac{1}{k} \sum_{x_i \in N_k(x)} y_i$$

where $N_k(x)$ is the neighborhood of $x$ defined by the $k$ closest points $x_i$ in the training sample.


```python
class KNeighborsRegressor:
    def __init__(self, k):
        self._k = k

    def fit(self, X, y):
        self._X = X
        self._y = y
        return self
    
    def predict(self, x):
        X, y, k = self._X, self._y, self._k
        distances = ((X - x) ** 2).sum(axis=1)
      
        return np.mean(y[distances.argpartition(k)[:k]])
```


```python
def plot_k_nearest_neighbors(k):
    model = KNeighborsRegressor(k).fit(data_x, data_y)
    is_orange = lambda x: model.predict(x) > 0.5
    orange_grid = np.array([*filter(is_orange, grid)])
    blue_grid = np.array([*filterfalse(is_orange, grid)])

    axes.clear()
    axes.set_title(str(k) + "-Nearest Neighbor Classifier")

    plot_data(orange_data, blue_data)
    plot_grid(orange_grid, blue_grid)

plot_k_nearest_neighbors(1)
fig

```

It appears that k-nearest-neighbor have a single parameter (*k*), however the effective number of parameters is N/k and is generally bigger than the p parameters in least-squares fits. **Note:** if the neighborhoods
were nonoverlapping, there would be N/k neighborhoods and we would fit one parameter (a mean) in each neighborhood.


```python
plot_k_nearest_neighbors(15)

fig
```

e9b7a65be312d9e9c032aa28d4fa5726c2dedef6
ipynb
Jupyter Notebook
chapter-02/2.3-least-squares-and-nearest-neighbors.ipynb
debajitd/the-elements-of-statistical-learning
67fae7a4f938b5375dd0d5b52dc159b01add7803
MIT
chapter-02/2.3-least-squares-and-nearest-neighbors.ipynb
debajitd/the-elements-of-statistical-learning
67fae7a4f938b5375dd0d5b52dc159b01add7803
MIT
chapter-02/2.3-least-squares-and-nearest-neighbors.ipynb
debajitd/the-elements-of-statistical-learning
67fae7a4f938b5375dd0d5b52dc159b01add7803
MIT
Qwen/Qwen-72B
1. YES
2. YES