Learning Algorithms

Inspired by Sam Lau, who co-authored the Learning Data Science book.

In this lesson, we'll introduce machine learning from the ground up. By the end of this lesson, students will be able to:

• Describe the difference between traditional algorithms and machine learning algorithms.
• Identify the components of a machine learning model and dataset features and labels.
• Apply sklearn to train a decision tree for classification and regression tasks.

A while back, we discussed data visualization using the Puget Sound Clean Air Agency's EPA-grade air quality sensors (AQS). However, these sensors are typically expensive, costing anywhere between \$15,000 and \$40,000 each, making it hard to deploy a large number of these sensors. Furthermore, EPA-grade AQS measurements also undergo calibration and accuracy checks that lead to delays of one or two hours, leading to data that is very accurate but not necessarily timely.

In contrast, "PurpleAir makes sensors that empower Community Scientists who collect hyper-local air quality data and share it with the public." In this lesson, we'll learn how we can use more accurate but less timely AQS measurements to calibrate the less accurate but more timely PurpleAir sensor (PAS) measurements so that we can provide the best information to the general public. The concepts in this lesson are actually used in the real world when you visit the EPA AirNow Fire and Smoke Map: the PAS data in this map are calibrated using the approach we will learn today.

import pandas as pd
import seaborn as sns

sns.set_theme()

Our dataset includes over 12,000 matched observations where we've paired each AQS measurement with a nearby PAS measurement, along with 3 other variables that experts have identified as potentially impacting PAS measurement quality. The dataset includes 5 columns:

• The very accurate EPA-grade air quality sensor (AQS) measurement of the PM2.5 value.
• The temperature in degrees celsius.
• The relative humidity as a percentage between 0% and 100%.
• The dew point, where a higher dew point means more moisture in the air.
• The less accurate but more timely PurpleAir sensor (PAS) measurement of the PM2.5 value.

How can we use the PAS measurement together with the temperature, humidity, and dew point to predict the AQS measurement?

sensor_data = pd.read_csv("sensor_data.csv")
sensor_data

	AQS	temp	humidity	dew	PAS
0	6.7	18.027263	38.564815	3.629662	8.616954
1	3.8	16.115280	49.404315	5.442318	3.493916
2	4.0	19.897634	29.972222	1.734051	3.799601
3	4.7	21.378334	32.474513	4.165624	4.369691
4	3.2	18.443822	43.898226	5.867611	3.191071
...	...	...	...	...	...
12092	5.5	-12.101337	54.188889	-19.555834	2.386120
12093	16.8	4.159967	56.256030	-3.870659	32.444987
12094	15.6	1.707895	65.779221	-4.083768	25.297018
12095	14.0	-14.380144	48.206481	-23.015378	8.213208
12096	5.8	5.081813	52.200000	-4.016401	9.436011

12097 rows × 5 columns

Let's use data visualization to understand this dataset by scatter-plotting PAS to AQS. (Why do we place PAS on the x-axis and AQS on the y-axis?)

sns.relplot(sensor_data, x="PAS", y="AQS")

<seaborn.axisgrid.FacetGrid at 0x7ac52847c110>

Guessing game

Let's start with a simpler task: How can we use the PAS measurement alone to predict the AQS measurement? To do this, let's choose a line that best describes the trend, or the regression line that is produced by calling lmplot.

def plot_line(slope, intercept=0):
    grid = sns.relplot(sensor_data, x="PAS", y="AQS")
    grid.facet_axis(0, 0).plot([0, 260], [intercept, slope * 260 + intercept], c="orange")
    grid.set(title=f"Slope = {slope:.2f}, intercept = {intercept:.2f}")
    return grid


plot_line(1)

<seaborn.axisgrid.FacetGrid at 0x7ac462905f10>

plot_line(.8)

<seaborn.axisgrid.FacetGrid at 0x7d6db28fbc50>

plot_line(.5)

<seaborn.axisgrid.FacetGrid at 0x7d6db2776cd0>

sensor_data["AQS"]

0         6.7
1         3.8
2         4.0
3         4.7
4         3.2
         ... 
12092     5.5
12093    16.8
12094    15.6
12095    14.0
12096     5.8
Name: AQS, Length: 12097, dtype: float64

sensor_data["PAS"]

0         8.616954
1         3.493916
2         3.799601
3         4.369691
4         3.191071
           ...    
12092     2.386120
12093    32.444987
12094    25.297018
12095     8.213208
12096     9.436011
Name: PAS, Length: 12097, dtype: float64

(sensor_data["AQS"] / sensor_data["PAS"]).mean()

0.9044468847554963

plot_line(0.9)

<seaborn.axisgrid.FacetGrid at 0x7ac4609162d0>

sensor_data["AQS"].max() - sensor_data["AQS"].min()

107.8

sensor_data["PAS"].max() - sensor_data["PAS"].min()

256.12481068004786

(sensor_data["AQS"] / sensor_data["PAS"]).mode()

0    0.493501
dtype: float64

107.8 / 256.12

0.42089645478681864

plot_line(0.42)

<seaborn.axisgrid.FacetGrid at 0x7ac45fd980d0>

err42 = sensor_data['PAS'] * 0.42 - sensor_data['AQS']
err42

0        -3.080880
1        -2.332555
2        -2.404168
3        -2.864730
4        -1.859750
           ...    
12092    -4.497829
12093    -3.173106
12094    -4.975252
12095   -10.550453
12096    -1.836875
Length: 12097, dtype: float64

err42.mean()

-2.9652289987410287

err42.abs().mean()

3.1648363175956207

(err42 ** 2).mean()

15.170646931726855

What differentiates machine learning from just human guessing is the use of an algorithm to find the best line, which requires a metric for the quality of a trend.

We can visualize all of our guesses so far by plotting them against their mean squared errors on what's called a loss surface.

def plot_loss(slopes):
    from sklearn.metrics import mean_squared_error as mse
    losses = [mse(sensor_data['PAS'] * s, sensor_data['AQS']) for s in slopes]
    grid = sns.relplot(x=slopes, y=losses)
    grid.set(title="Loss surface", xlabel="Slope", ylabel="MSE", xlim=[0, 1], ylim=[0, None])
    return grid


plot_loss([0.45, 0.9, 0.8, 0.2, 0.42, 0.5, 0.6, 0.575])

<seaborn.axisgrid.FacetGrid at 0x7ac45ea3e710>

Gradient descent

So how do we write a machine learning algorithm that can optimize this metric and find the minimum mean squared error in the loss surface so that it selects the best possible line? Machine learning scientists can apply concepts from linear algebra to solve this system by selecting a random initial theta (slope) value and then rolling down the hill toward the minimum value at the bottom of the bowl. We can express this using numpy, a numeric computation module for Python that is a building block for pandas and sklearn (as we'll see later).

$$ \nabla_{\!\theta}\; \text{MSE}(\boldsymbol{\theta}, \mathbf{X}, \mathbf{y}) = -\frac{2}{n}(\mathbf{X}^\top \mathbf{y} - \mathbf{X}^\top \mathbf{X} \boldsymbol{\theta}) $$

import numpy as np


def grad_mse(theta, X, y):
    return np.array(- 2 / len(X) * (X.T @ y - X.T @ X * theta))


thetas = [np.random.random()]
print("Random initial theta value:", thetas[-1])

We can then take a small step in the opposite direction of the gradient to roll down the hill until we converge on a good guess. To make this a machine learning algorithm, we simply put the update step in a loop until the theta values no longer noticeably change.

plot_line(thetas[-1])
plot_loss(thetas)

# Take a small step in the opposite direction of the gradient to roll downhill
thetas.append(thetas[-1] - 0.002 * grad_mse(thetas[-1], sensor_data["PAS"], sensor_data["AQS"]))

Linear regression models

What we've just described is the gradient descent algorithm for fitting a linear regression model. A linear regression model is a machine learning model that is used to predict a numeric value (like AQS measurements) using a linear combination of coefficients and features (columns from the training dataset). scikit-learn provides an easy way to do define a linear regression model, fit our training dataset X to the target values y, and examine the coefficients to look inside the model.

from sklearn.linear_model import LinearRegression
X = sensor_data[["PAS"]]
y = sensor_data["AQS"]
reg = LinearRegression().fit(X, y)

print("Model:", " + ".join([f"{reg.intercept_:.2f}"] + [f"{coef:.2f}({X.columns[i]})" for i, coef in enumerate(reg.coef_)]))
print("Error:", mean_squared_error(y, reg.predict(X)))
plot_line(reg.coef_[0], reg.intercept_)

This procedure is more or less how lmplot works!

sns.lmplot(sensor_data, x="PAS", y="AQS")

But the advantage of designing our own model is that we can combine other variables to reduce the mean squared error loss. The final model that the EPA uses only takes into account the sensor measurement and the relative humidity, but not any other variables. Later, we'll learn why they made this decision.

Classification versus regression

Everything we've seen so far fall under the category of regression, where we aim to predict a numeric target value (one column) from one or more features (one or more other columns). The other main category of problems is classification, which is just like regression except we aim to predict a categorical target value. For example, we might want to answer the question: How can we predict whether a house belongs in NY or SF from its beds, baths, price, year of construction, square footage, price per square foot, and elevation?

homes = pd.read_csv("homes.csv")
homes

Let's learn about decision trees, a machine learning algorithm that can be used for classification (and also, as it turns out, regression too). Decision trees learn a nested if-then-else logical hierarchy to fit a training dataset. In the following visualization, the color and opacity of each box represents whether that subset of homes is more likely to be in NY or SF. The notes inside each node indicate information about the values:

1. The first line shows the condition. If the condition is true, go left; if not, go right.
2. The second line shows the percentage of samples represented by that node.
3. The third line shows the proportion of homes within that sample that belong in ["NY", "SF"].
4. The fourth line shows the majority class in that category, corresponding to the bigger number on line 3.

from sklearn.tree import DecisionTreeClassifier, plot_tree
X = homes.drop("city", axis=1)
y = homes["city"]
clf = DecisionTreeClassifier(max_depth=2).fit(X, y)

import matplotlib.pyplot as plt
plt.figure(dpi=300)
plot_tree(
    clf,
    feature_names=X.columns,
    class_names=["NY", "SF"],
    label="root",
    filled=True,
    impurity=False,
    proportion=True,
    rounded=False
);

We can also use this dataset for regression too. Write a one-line expression to train a DecisionTreeRegressor model to predict the price of a home in this dataset from all other variables.

from sklearn.tree import DecisionTreeRegressor

reg = ...
reg

Consider each of the following tasks and answer whether they would be best modeled as classification or regression.

Predict whether an email is spam or not spam.

Classification, since the target value is the category "spam" or "not spam".

Predict the number of views a video will receive based on subscriber count.

Regression, since the target value is a number.

Predict the next word to appear in a sentence.

Classification, since the target value is to choose from the dictionary of all possible next words.

If you're curious about how ChatGPT works, Jay Mody has a good introduction to GPT in 60 Lines of NumPy.