In [1]:

```
# Set up packages for lecture. Don't worry about understanding this code,
# but sure to run it if you're following along.
import numpy as np
import babypandas as bpd
import pandas as pd
from matplotlib_inline.backend_inline import set_matplotlib_formats
import matplotlib.pyplot as plt
from scipy import stats
set_matplotlib_formats("svg")
plt.style.use('ggplot')
np.set_printoptions(threshold=20, precision=2, suppress=True)
pd.set_option("display.max_rows", 7)
pd.set_option("display.max_columns", 8)
pd.set_option("display.precision", 2)
import warnings
warnings.filterwarnings('ignore')
# New minimize function (wrapper around scipy.optimize.minimize)
from inspect import signature
from scipy import optimize
def minimize(function):
n_args = len(signature(function).parameters)
initial = np.zeros(n_args)
return optimize.minimize(lambda x: function(*x), initial).x
# All of the following code is for visualization.
def plot_regression_line(df, x, y, margin=.02):
'''Computes the slope and intercept of the regression line between columns x and y in df (in original units) and plots it.'''
m = slope(df, x, y)
b = intercept(df, x, y)
df.plot(kind='scatter', x=x, y=y, s=100, figsize=(10, 5), label='original data')
left = df.get(x).min()*(1 - margin)
right = df.get(x).max()*(1 + margin)
domain = np.linspace(left, right, 10)
plt.plot(domain, m*domain + b, color='orange', label='regression line', lw=4)
plt.suptitle(format_equation(m, b), fontsize=18)
plt.legend();
def format_equation(m, b):
if b > 0:
return r'$y = %.2fx + %.2f$' % (m, b)
elif b == 0:
return r'$y = %.2fx' % m
else:
return r'$y = %.2fx %.2f$' % (m, b)
def plot_errors(df, m, b, ax=None):
x = df.get('x')
y = m * x + b
df.plot(kind='scatter', x='x', y='y', s=100, label='original data', ax=ax, figsize=(10, 5) if ax is None else None)
if ax:
plotter = ax
else:
plotter = plt
plotter.plot(x, y, color='orange', lw=4)
for k in np.arange(df.shape[0]):
xk = df.get('x').iloc[k]
yk = np.asarray(y)[k]
if k == df.shape[0] - 1:
plotter.plot([xk, xk], [yk, df.get('y').iloc[k]], linestyle=(0, (1, 1)), c='r', lw=4, label='errors')
else:
plotter.plot([xk, xk], [yk, df.get('y').iloc[k]], linestyle=(0, (1, 1)), c='r', lw=4)
plt.title(format_equation(m, b), fontsize=18)
plt.xlim(50, 90)
plt.ylim(40, 100)
plt.legend();
```

- Quiz 4 scores are released, and the Grade Report has been updated – it reflects your scores on everything other than Homework 6, Lab 7, the Final Project, and the Final Exam.
- Lab 7 is due
**Thursday 12/7 at 11:59 PM**.- It is true that your lowest lab is dropped.
- However, it's a bad idea to simply ignore this lab, because it's the only assignment on regression, which will be tested on the Final Exam.

- Homework 6 is due
**today at 11:59PM**. - The Final Project is due
**Tuesday 12/5 at 11:59PM**. - The Final Exam is on
**Saturday 12/9 from 7-10PM**.- We'll send out more details early next week, including your room and seat assignment.
- Start preparing by working through old Final Exams at practice.dsc10.com. Take them as if they are real exams – time yourself, and don't use resources other than the reference sheet.
- We'll work through the solutions to the Spring 2023 Final Exam in class on Wednesday.

- The regression line in standard units.
- The regression line in original units.
- Outliers.
- Errors in prediction.

Recall, in the last lecture, we aimed to use a mother's height to predict her adult son's height.

In [2]:

```
galton = bpd.read_csv('data/galton.csv')
male_children = galton[galton.get('gender') == 'male']
mom_son = bpd.DataFrame().assign(mom = male_children.get('mother'),
son = male_children.get('childHeight'))
mom_son
```

Out[2]:

mom | son | |
---|---|---|

0 | 67.0 | 73.2 |

4 | 66.5 | 73.5 |

5 | 66.5 | 72.5 |

... | ... | ... |

925 | 60.0 | 66.0 |

929 | 66.0 | 64.0 |

932 | 63.0 | 66.5 |

481 rows × 2 columns

In [3]:

```
mom_son.plot(kind='scatter', x='mom', y='son', figsize=(10, 5));
```

Recall, the correlation coefficient $r$ of two variables $x$ and $y$ is defined as the

**average**value of the**product**of $x$ and $y$- when both are measured in
**standard units**.

In [4]:

```
def standard_units(col):
return (col - col.mean()) / np.std(col)
def calculate_r(df, x, y):
'''Returns the average value of the product of x and y,
when both are measured in standard units.'''
x_su = standard_units(df.get(x))
y_su = standard_units(df.get(y))
return (x_su * y_su).mean()
```

In [5]:

```
r_mom_son = calculate_r(mom_son, 'mom', 'son')
r_mom_son
```

Out[5]:

0.32300498368490554

- The regression line is the line through $(0,0)$ with slope $r$, when both variables are measured in
**standard units**.

- We use the regression line to make predictions!

**Example**: If a mother's height is 0.5 SDs above the average mother's height, and $r = 0.32$, our prediction is that her son's height will be $0.5 \cdot 0.32 = 0.16$ SDs above the average son's height.

**Issue**: To use this form of the regression line, we need to know mothers' heights in standard units, but it would be more convenient to think in terms of inches.

A course has a midterm (mean 80, standard deviation 15) and a really hard final (mean 50, standard deviation 12).

If the scatter plot comparing midterm & final scores for students looks linearly associated with correlation 0.75, then what is the predicted final exam score for a student who received a 90 on the midterm?

- A. 54
- B. 56
- C. 58
- D. 60
- E. 62

Each time we want to predict the height of an adult son given the height of his mother, we have to:

- Convert the mother's height from inches to standard units.

- Multiply by the correlation coefficient to predict the son's height in standard units.

- Convert the son's predicted height from standard units back to inches.

When $x$ and $y$ are in standard units, the regression line is given by

What is the regression line when $x$ and $y$ are in their original units (e.g. inches)?

$$\frac{\text{predicted } y - \text{mean of }y}{\text{SD of }y} = r \cdot \frac{x - \text{mean of } x}{\text{SD of }x}$$

- Re-arranging the above equation into the form $\text{predicted } y = mx + b$ yields the formulas:

- $m$ is the slope of the regression line and $b$ is the intercept.

Let's implement these formulas in code and try them out.

In [6]:

```
def slope(df, x, y):
"Returns the slope of the regression line between columns x and y in df (in original units)."
r = calculate_r(df, x, y)
return r * np.std(df.get(y)) / np.std(df.get(x))
def intercept(df, x, y):
"Returns the intercept of the regression line between columns x and y in df (in original units)."
return df.get(y).mean() - slope(df, x, y) * df.get(x).mean()
```

In [7]:

```
m_heights = slope(mom_son, 'mom', 'son')
m_heights
```

Out[7]:

0.36506116024257595

In [8]:

```
b_heights = intercept(mom_son, 'mom', 'son')
b_heights
```

Out[8]:

45.85803797199311

So, the regression line is

$$\text{predicted son's height in inches} = 0.365 \cdot \text{mother's height in inches} + 45.858$$In [9]:

```
def predict_son(mom):
return m_heights * mom + b_heights
```

What's the predicted height of a son whose mother is 62 inches tall?

In [10]:

```
predict_son(62)
```

Out[10]:

68.49182990703282

What if the mother is 55 inches tall? 73 inches tall?

In [11]:

```
predict_son(55)
```

Out[11]:

65.9364017853348

In [12]:

```
predict_son(73)
```

Out[12]:

72.50750266970115

In [13]:

```
xs = np.arange(57, 72)
ys = predict_son(xs)
mom_son.plot(kind='scatter', x='mom', y='son', figsize=(10, 5), title='Regression line predictions, in original units', label='original data');
plt.plot(xs, ys, color='orange', lw=4, label='regression line')
plt.legend();
```

Consider the dataset below. What is the correlation between $x$ and $y$?

In [14]:

```
outlier = bpd.read_csv('data/outlier.csv')
outlier.plot(kind='scatter', x='x', y='y', s=100, figsize=(10, 5));
```

In [15]:

```
calculate_r(outlier, 'x', 'y')
```

Out[15]:

-0.02793982443854457

In [16]:

```
plot_regression_line(outlier, 'x', 'y')
```

In [17]:

```
without_outlier = outlier[outlier.get('y') > 40]
```

In [18]:

```
calculate_r(without_outlier, 'x', 'y')
```

Out[18]:

0.9851437295364016

In [19]:

```
plot_regression_line(without_outlier, 'x', 'y')
```

**Takeaway**: Even a single outlier can have a massive impact on the correlation, and hence the regression line. Look for these before performing regression. **Always visualize first!**

- In examples we've seen so far, the regression line seems to fit our data pretty well.
- But how well?
- What makes the regression line good?
- Would another line be better?

In [20]:

```
outlier.plot(kind='scatter', x='x', y='y', s=100, figsize=(10, 5));
```

In [21]:

```
m_no_outlier = slope(without_outlier, 'x', 'y')
b_no_outlier = intercept(without_outlier, 'x', 'y')
m_no_outlier, b_no_outlier
```

Out[21]:

(0.975927715724588, 3.042337135297444)

In [22]:

```
plot_errors(without_outlier, m_no_outlier, b_no_outlier)
```

- A good prediction line is one where the errors tend to be small.

- To measure the rough size of the errors, for a particular set of predictions:
- Square the errors so that they don't cancel each other out.
- Take the mean of the squared errors.
- Take the square root to fix the units.

- This is called
**root mean square error**(RMSE).- Notice the similarities to computing the SD!

In [23]:

```
without_outlier
```

Out[23]:

x | y | |
---|---|---|

0 | 50 | 53.53 |

1 | 55 | 54.21 |

2 | 60 | 65.65 |

... | ... | ... |

6 | 80 | 79.61 |

7 | 85 | 88.17 |

8 | 90 | 91.05 |

9 rows × 2 columns

First, let's compute the regression line's predictions for the entire dataset.

In [24]:

```
predicted_y = m_no_outlier * without_outlier.get('x') + b_no_outlier
predicted_y
```

Out[24]:

array([51.84, 56.72, 61.6 , 66.48, 71.36, 76.24, 81.12, 86. , 90.88])

To find the RMSE, we need to start by finding the errors and squaring them.

In [25]:

```
# Errors.
without_outlier.get('y') - predicted_y
```

Out[25]:

0 1.69 1 -2.51 2 4.06 ... 6 -1.51 7 2.18 8 0.18 Name: y, Length: 9, dtype: float64

In [26]:

```
# Squared errors.
(without_outlier.get('y') - predicted_y) ** 2
```

Out[26]:

0 2.86 1 6.31 2 16.45 ... 6 2.27 7 4.74 8 0.03 Name: y, Length: 9, dtype: float64

In [27]:

```
# Mean squared error.
((without_outlier.get('y') - predicted_y) ** 2).mean()
```

Out[27]:

4.823770221019621

In [28]:

```
# Root mean squared error.
np.sqrt(((without_outlier.get('y') - predicted_y) ** 2).mean())
```

Out[28]:

2.19630831647554

- We've been using the regression line to make predictions. But we could use a different line!
- To make a prediction for
`x`

using an arbitrary line defined by`slope`

and`intercept`

, compute`x * slope + intercept`

.

- To make a prediction for

- For this dataset, if we choose a
**different line**, we will end up with different predictions, and hence a**different RMSE**.

In [29]:

```
def rmse(slope, intercept):
'''Calculates the RMSE of the line with the given slope and intercept,
using the 'x' and 'y' columns of without_outlier.'''
# The true values of y.
true = without_outlier.get('y')
# The predicted values of y, from plugging the x values from the
# given DataFrame into the line with the given slope and intercept.
predicted = slope * without_outlier.get('x') + intercept
return np.sqrt(((true - predicted) ** 2).mean())
```

In [30]:

```
# Check that our function works on the regression line.
rmse(m_no_outlier, b_no_outlier)
```

Out[30]:

2.19630831647554

Let's compute the RMSEs of several different lines on the same dataset.

In [31]:

```
# Experiment by changing one of these!
lines = [(1.2, -15), (0.75, 11.5), (-0.4, 100)]
fig, ax = plt.subplots(1, 3, figsize=(14, 4))
for i, line in enumerate(lines):
plt.subplot(1, 3, i + 1)
m, b = line
plot_errors(without_outlier, m, b, ax=ax[i])
ax[i].set_title(format_equation(m, b) + f'\nRMSE={np.round(rmse(m, b), 2)}')
```

- RMSE describes how well a line fits the data.
**The lower the RMSE of a line is, the better it fits the data**.

- If you take DSC 40A, you'll learn how to do this using calculus. For now, we'll use a Python function that can do it automatically –
`minimize`

.

`minimize`

¶- The function
`minimize`

takes in a function as an argument, and returns the inputs to that function that produce the smallest output.

- For instance, we know that the minimizing input to the function $f(x) = (x - 5)^2 + 4$ is $x = 5$.
`minimize`

can find this, too:

In [32]:

```
def f(x):
return (x - 5) ** 2 + 4
# Plot of f(x).
x = np.linspace(0, 10)
y = f(x)
plt.plot(x, y)
plt.title(r'$f(x) = (x - 5)^2 + 4$');
```