Workshop 1

Follow along with the slides, then work through the lab below.

Lab 1: Predict house prices

Getting started (2 minutes)

• Open the lab notebook in Colab and save your own copy.
• Click the play button beside a code cell to run it. Work from top to bottom.
• Double-click text to edit an answer. Click outside when done.

The first cell loads four libraries:

It then loads 20,640 neighbourhoods from the 1990 California census into a table called housing.

1
import matplotlib.pyplot as plt  # Draws graphs.
2
import numpy as np  # Does maths on many numbers at once.
3
import pandas as pd  # Stores data in tables.
4
from sklearn.datasets import fetch_california_housing
5
from sklearn.linear_model import LinearRegression
6
from sklearn.metrics import mean_absolute_error, mean_squared_error
7
from sklearn.model_selection import train_test_split
8
from sklearn.tree import DecisionTreeRegressor
9

10
housing = fetch_california_housing(as_frame=True).frame  # Download the data as a table.
11
housing.head()  # Show the first five rows.

1. Look before you fit (3 minutes)

Predict: do richer neighbourhoods have more expensive houses?

Replace each ... to plot income against house value.

1
sample = housing.sample(2000, random_state=0)
2

3
incomes = sample["MedInc"]
4
values = ...  # Select the MedHouseVal column.
5

6
plt.scatter(incomes, values, s=5)
7
plt.xlabel("Median income ($10,000s)")
8
plt.ylabel(...)  # Label the y-axis.
9
plt.show()

1
sample = housing.sample(2000, random_state=0)
2

3
incomes = sample["MedInc"]
4
values = sample["MedHouseVal"]
5

6
plt.scatter(incomes, values, s=5)
7
plt.xlabel("Median income ($10,000s)")
8
plt.ylabel("Median house value ($100,000s)")
9
plt.show()

Was your prediction right? Why is there a flat line of points at the top of the graph?

2. Fit a model and predict (5 minutes)

Train the model as in the flat-price demo. The LinearRegression documentation lists its methods: fit trains the model, and predict uses it on new inputs.

Replace each .... Expect predictions of about 1.70 and 3.79 ($170,000 and $379,000).

1
inputs = housing[["MedInc"]]
2
outputs = housing["MedHouseVal"]
3

4
model = LinearRegression()
5
# Fit the model here.
6

7
new_neighbourhoods = pd.DataFrame({"MedInc": [3, 8]})
8
model.predict(...)  # Predict values for the new neighbourhoods.

1
inputs = housing[["MedInc"]]
2
outputs = housing["MedHouseVal"]
3

4
model = LinearRegression()
5
model.fit(inputs, outputs)
6

7
new_neighbourhoods = pd.DataFrame({"MedInc": [3, 8]})
8
model.predict(new_neighbourhoods)

Run the cell below to draw the model’s line over your sample. The model learnt two parameters: a slope and an intercept. What do they mean in words?

1
line_ends = pd.DataFrame({"MedInc": [0, 15]})  # Predict at both ends of the x-axis.
2
plt.scatter(incomes, values, s=5)
3
plt.plot(line_ends["MedInc"], model.predict(line_ends), color="black")  # Join the two predictions.
4
plt.xlabel("Median income ($10,000s)")
5
plt.ylabel("Median house value ($100,000s)")
6
plt.show()
7
print("slope:", model.coef_[0], "intercept:", model.intercept_)

3. Measure the loss (4 minutes)

NumPy and pandas do arithmetic on whole columns at once.

1. Complete the MSE calculation. Expect about 0.70. The next line checks it with scikit-learn.
2. A baseline model always predicts the average house value. Calculate its MSE. Expect about 1.33.

1
actual = outputs
2
predicted = model.predict(inputs)
3

4
errors = actual - predicted
5
mse = ...  # Square the errors, then take the mean.
6
print(mse)
7
print(mean_squared_error(actual, predicted))

1
baseline = np.full(len(actual), actual.mean())  # The average, once per neighbourhood.
2
mean_squared_error(..., ...)  # Compare the actual values with the baseline.

1
actual = outputs
2
predicted = model.predict(inputs)
3

4
errors = actual - predicted
5
mse = np.mean(errors ** 2)
6
print(mse)
7
print(mean_squared_error(actual, predicted))

1
baseline = np.full(len(actual), actual.mean())  # The average, once per neighbourhood.
2
mean_squared_error(actual, baseline)

4. Give the model more inputs (3 minutes)

housing.drop(columns="MedHouseVal") makes a copy of the table without the MedHouseVal column. That leaves all 8 input columns.

Fit a new model on every input and calculate its MSE. Expect about 0.52.

1
all_inputs = housing.drop(columns="MedHouseVal")
2

3
big_model = LinearRegression()
4
...  # Fit big_model on all_inputs and actual.
5
mean_squared_error(actual, ...)  # Predict with big_model.

1
all_inputs = housing.drop(columns="MedHouseVal")
2

3
big_model = LinearRegression()
4
big_model.fit(all_inputs, actual)
5
mean_squared_error(actual, big_model.predict(all_inputs))

The loss fell. Does that prove the new model will predict better for neighbourhoods it has not seen?

Lab 1 check (1 minute)

Explain in one sentence each:

• What a loss function measures.
• Why we compare a model against a simple baseline.

Extension A: Errors in dollars

MSE is in squared units, which are hard to picture. Calculate the mean absolute error with mean_absolute_error, which works like mean_squared_error. Expect about 0.63. What is that in dollars?

1
mean_absolute_error(..., ...)

1
mean_absolute_error(actual, predicted)

Extension B: Train the model yourself

LinearRegression solves for the best slope and intercept directly. Neural networks cannot do this, so they use gradient descent: start with a guess, then take many small downhill steps.

.to_numpy() turns a column into a NumPy array: a plain list of numbers that supports fast maths.

Fill in the two update lines. After 1,000 steps, expect a slope of about 0.42, an intercept of about 0.44 and an MSE of about 0.70, close to LinearRegression.

Then try a learning rate of 0.1. What happens, and why?

1
x = housing["MedInc"].to_numpy()
2
y = housing["MedHouseVal"].to_numpy()
3
slope, intercept = 0.0, 0.0
4
learning_rate = 0.01
5

6
for step in range(1000):
7
    error = slope * x + intercept - y  # Forward pass, then compare with the answers.
8
    slope_gradient = 2 * np.mean(error * x)
9
    intercept_gradient = 2 * np.mean(error)
10
    slope = ...  # Step slope downhill.
11
    intercept = ...  # Step intercept downhill.
12

13
print(slope, intercept, np.mean((slope * x + intercept - y) ** 2))

1
x = housing["MedInc"].to_numpy()
2
y = housing["MedHouseVal"].to_numpy()
3
slope, intercept = 0.0, 0.0
4
learning_rate = 0.01
5

6
for step in range(1000):
7
    error = slope * x + intercept - y  # Forward pass, then compare with the answers.
8
    slope_gradient = 2 * np.mean(error * x)
9
    intercept_gradient = 2 * np.mean(error)
10
    slope = slope - learning_rate * slope_gradient
11
    intercept = intercept - learning_rate * intercept_gradient
12

13
print(slope, intercept, np.mean((slope * x + intercept - y) ** 2))

Take-home: Find the sweet spot (about 15 minutes)

You will choose how flexible a model should be, using validation data.

A decision tree predicts by asking a series of yes-or-no questions, such as “Is MedInc above 5?”. max_depth limits how many questions it can ask in a row. A deeper tree is more flexible.

T1. Split the data three ways (3 minutes)

train_test_split(table, test_size=0.3, random_state=0) shuffles the rows, then splits them into two tables. test_size=0.3 puts 30% of the rows in the second table.

Split rest in half to make the validation and test sets. Expect 14,448, 3,096 and 3,096 rows.

1
train, rest = train_test_split(housing, test_size=0.3, random_state=0)
2
val, test = ...  # Split rest in half.
3

4
input_columns = housing.columns.drop("MedHouseVal")  # The names of the 8 input columns.
5
print(len(train), len(val), len(test))  # len counts the rows.

1
train, rest = train_test_split(housing, test_size=0.3, random_state=0)
2
val, test = train_test_split(rest, test_size=0.5, random_state=0)
3

4
input_columns = housing.columns.drop("MedHouseVal")  # The names of the 8 input columns.
5
print(len(train), len(val), len(test))  # len counts the rows.

T2. Train trees of increasing depth (5 minutes)

train[input_columns] selects the input columns from the training table.

The loop trains one tree for each depth from 1 to 20 and records its training MSE. Complete the line that records the validation MSE, then run the cell to plot both.

Predict before you run: what will each line look like?

1
depths = range(1, 21)  # 1, 2, ..., 20
2
train_losses = []
3
val_losses = []
4
for depth in depths:
5
    tree = DecisionTreeRegressor(max_depth=depth, random_state=0)
6
    tree.fit(train[input_columns], train["MedHouseVal"])
7
    train_losses.append(mean_squared_error(train["MedHouseVal"], tree.predict(train[input_columns])))
8
    val_losses.append(...)  # The same, for the validation set.
9

10
plt.plot(depths, train_losses, marker="o", label="Training")
11
plt.plot(depths, val_losses, marker="o", label="Validation")
12
plt.xlabel("Tree depth")
13
plt.ylabel("MSE")
14
plt.legend()
15
plt.show()

1
depths = range(1, 21)  # 1, 2, ..., 20
2
train_losses = []
3
val_losses = []
4
for depth in depths:
5
    tree = DecisionTreeRegressor(max_depth=depth, random_state=0)
6
    tree.fit(train[input_columns], train["MedHouseVal"])
7
    train_losses.append(mean_squared_error(train["MedHouseVal"], tree.predict(train[input_columns])))
8
    val_losses.append(mean_squared_error(val["MedHouseVal"], tree.predict(val[input_columns])))
9

10
plt.plot(depths, train_losses, marker="o", label="Training")
11
plt.plot(depths, val_losses, marker="o", label="Validation")
12
plt.xlabel("Tree depth")
13
plt.ylabel("MSE")
14
plt.legend()
15
plt.show()

Which depths underfit? Which overfit? How can you tell from the two lines?

T3. Choose a depth, then test it once (3 minutes)

np.argmin(val_losses) gives the position of the smallest value in the list.

The cell picks the depth with the lowest validation MSE. Expect depth 9, with a validation MSE of about 0.43. Complete the last line to calculate its MSE on the test set. Do this only once, at the end.

1
best_position = np.argmin(val_losses)
2
best_depth = depths[best_position]
3
print("best depth:", best_depth, "validation MSE:", val_losses[best_position])
4

5
best_tree = DecisionTreeRegressor(max_depth=best_depth, random_state=0)
6
best_tree.fit(train[input_columns], train["MedHouseVal"])
7
...  # Calculate the MSE on the test set.

1
best_position = np.argmin(val_losses)
2
best_depth = depths[best_position]
3
print("best depth:", best_depth, "validation MSE:", val_losses[best_position])
4

5
best_tree = DecisionTreeRegressor(max_depth=best_depth, random_state=0)
6
best_tree.fit(train[input_columns], train["MedHouseVal"])
7
mean_squared_error(test["MedHouseVal"], best_tree.predict(test[input_columns]))

Why should you not use the test set to choose the depth?

Open-ended: Beat 0.40

Can you get the validation MSE below 0.40? Keep the test set for one final check.

Start from the cell below, and change the model or its settings. Ideas to try:

• Other tree settings, such as min_samples_leaf, the smallest number of examples allowed in each final group.
• Other models, such as RandomForestRegressor from sklearn.ensemble. Read its documentation first.

Record what you tried:

1
from sklearn.ensemble import RandomForestRegressor
2

3
candidate = DecisionTreeRegressor(max_depth=9, random_state=0)  # Change the model or its settings.
4
candidate.fit(train[input_columns], train["MedHouseVal"])
5
print("training MSE:", mean_squared_error(train["MedHouseVal"], candidate.predict(train[input_columns])))
6
print("validation MSE:", mean_squared_error(val["MedHouseVal"], candidate.predict(val[input_columns])))

1
from sklearn.ensemble import HistGradientBoostingRegressor, RandomForestRegressor
2

3
for candidate in [RandomForestRegressor(random_state=0, n_jobs=-1), HistGradientBoostingRegressor(random_state=0)]:
4
    candidate.fit(train[input_columns], train["MedHouseVal"])
5
    print(type(candidate).__name__)
6
    print("training MSE:", mean_squared_error(train["MedHouseVal"], candidate.predict(train[input_columns])))
7
    print("validation MSE:", mean_squared_error(val["MedHouseVal"], candidate.predict(val[input_columns])))