Missing data

In lesson 4.1 we looked at the data before building a model. Today: what the data do not contain. Missing data is an empty cell in the table: someone did not measure, did not write it down or did not want to say. Most models will not accept an empty cell, so you have to do something with it. There are two simple ways out: drop the row, or write some value into the gap, which is called imputing missing values.

Try it

• Start with the penguins. How many penguins lack sex, and how many lack measurements? Are they the same rows?
• The model is to predict body mass from the three other measurements. Does it need sex? How many penguins stay if you drop only the rows that miss the model's columns?
• Switch to diabetes. Here BMI is missing for the heaviest patients. Fill the gap with the mean: does the link between BMI and disease progression weaken or strengthen?
• Add the “missing” column and compare R² in cross-validation. Then switch to random patients. What changes?

Once you have checked it in the lab, let's name what you see.

The idea in plain words

Dropping a row with a gap is simple, but it costs. In the penguins, the rows with any gap are 11 of 344, so 333 remain. A model that predicts body mass from the bill and flipper measurements never reads sex. If you drop only the rows that miss the model's columns, 342 remain. Drop rows because of the columns you actually use.

The more important question is who the data are missing for. In the diabetes dataset every one of the 442 patients has a BMI, so the lab hides it on purpose for 110 of them. Once at random, once for the 110 heaviest, that is, everyone with a BMI of at least 29.4. Such a gap is not random. The patients without BMI have an average disease progression of 209, the others 133. The mere fact that BMI is missing says something about the patient.

When you write the mean of the other patients (24.3) into the gap, the 110 people with the most progression land at an average BMI. The correlation of BMI with disease progression falls from 0.586 to 0.352: the link weakens. The model's line does not get flatter, because the patients with a BMI set it; the fit gets worse, because the 110 filled-in points lie far from the line. The R² of a linear regression in cross-validation falls from 0.328 to 0.106. An extra 0/1 column “BMI was missing” helps. The model then gives the gap itself a weight of 76 points of progression, and R² comes back to 0.284. With random gaps this column adds nothing: 0.264 against 0.267 with the mean alone.

Jak to policzono

• Penguins: Palmer Penguins (Gorman, Williams and Fraser, 2014; palmerpenguins package: Horst, Hill and Gorman, 2020), licence CC0. 344 rows, real gaps: sex for 11 (Adelie 6, Gentoo 5), the four measurements for the same 2 penguins.
• Diabetes: Diabetes (Efron, Hastie, Johnstone and Tibshirani, 2004), shipped with scikit-learn under the BSD-3-Clause licence. 442 patients, disease progression after one year. BMI is in the data for everyone; the gaps for 110 patients are ours: at random or for the 110 with the highest BMI.
• Model: linear regression of disease progression on BMI (with the missing column, on BMI and that column), cross-validation on 5 parts repeated 10 times. The mean to fill in is computed from the patients with a BMI; in cross-validation, in every split, from the training part only, like SimpleImputer in a pipeline.

Why it works

The filled-in mean sits in the middle of the BMI axis, but it belongs to patients from its right end. Their disease progression is high, so filling in creates a vertical stack of high points above an average BMI. Such a stack adds nothing to the joint variation of BMI and progression. Compared with the patients who have a BMI (correlation 0.459), it only adds variation to progression, so the correlation falls to 0.352. Compared with the full data (0.586), these 110 patients lose their share of the joint variation: they had a high BMI and high progression, and now they stand in the middle of the axis. The standard deviation of BMI shrinks from 4.4 to 2.3, because 110 values sit in one point.

Dropping those patients does not save the day either. What is left is a cloud without its right end, and the correlation on the remaining 332 is 0.459 instead of 0.586. On top of that, the model gives no answer at all for 110 patients. The missing column works because it gives the model a separate weight for the group whose BMI we do not know. That group differs in disease progression, so the weight carries the difference into the predictions. When the gap is random, the group without BMI does not differ from the rest (156 against 151), and the weight has nothing to carry.

A common mistake

True or false: “fill the gaps with the mean and you are done”. False. The mean shrinks the spread of the column and pretends the gap was never there. When the gap depends on the value itself, as here for the heaviest, the mean bends the relation: the correlation falls from 0.586 to 0.352 and R² from 0.328 to 0.106. The “missing” column brings R² back to 0.284. In the code below it beats the mean alone in 48 of 50 validation folds.

The mean and dropping rows are fine when gaps are few and random. But first check whether the rows with a gap differ from the rest.

import numpy as np
import pandas as pd
from sklearn.datasets import load_diabetes
from sklearn.model_selection import cross_val_score, RepeatedKFold
from sklearn.linear_model import LinearRegression
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline

df = pd.read_csv("data/raw/penguins.csv")
gaps = df.isna().sum()
print(f"{len(df)} penguins; missing: " + ", ".join(f"{c} {n}" for c, n in gaps[gaps > 0].items()))
measured = df.dropna(subset=["bill_length_mm", "bill_depth_mm", "flipper_length_mm", "body_mass_g"])  # drop only rows without the model's columns
print(f"dropna() leaves {len(df.dropna())}; dropping only rows without measurements leaves {len(measured)}")
X, y = load_diabetes(return_X_y=True, as_frame=True, scaled=False)
bmi, y = X.bmi.to_numpy(), y.to_numpy()
heavy = bmi >= np.sort(bmi)[-110]  # the 110 heaviest: BMI of at least 29.4
rnd = np.zeros(len(bmi), bool)
rnd[np.random.default_rng(0).choice(len(bmi), 110, replace=False)] = True  # 110 at random
cv = RepeatedKFold(n_splits=5, n_repeats=10, random_state=0)
r2 = lambda model, Z: cross_val_score(model, Z, y, cv=cv, scoring="r2")
print(f"{len(bmi)} patients, BMI vs progression: r {np.corrcoef(bmi, y)[0, 1]:.3f}, CV R2 {r2(LinearRegression(), bmi[:, None]).mean():.3f}")
for name, miss in [("random", rnd), ("heaviest", heavy)]:
    hidden = np.where(miss, np.nan, bmi)[:, None]
    filled = np.where(miss, bmi[~miss].mean(), bmi)  # the mean only of the patients with a BMI
    mean_r2 = r2(make_pipeline(SimpleImputer(), LinearRegression()), hidden)  # the mean is learnt in every split, from its training part only
    flag_r2 = r2(make_pipeline(SimpleImputer(add_indicator=True), LinearRegression()), hidden)  # second column: 1 = BMI was missing
    print(f"BMI hidden for 110 {name}: progression {y[miss].mean():.0f} vs {y[~miss].mean():.0f}; mean fill: r {np.corrcoef(filled, y)[0, 1]:.3f}, CV R2 {mean_r2.mean():.3f}; mean + flag: CV R2 {flag_r2.mean():.3f}")
print(f"heaviest: mean + flag beats the mean alone in {(flag_r2 > mean_r2).sum()} of 50 folds")

344 penguins; missing: bill_length_mm 2, bill_depth_mm 2, flipper_length_mm 2, body_mass_g 2, sex 11
dropna() leaves 333; dropping only rows without measurements leaves 342
442 patients, BMI vs progression: r 0.586, CV R2 0.332
BMI hidden for 110 random: progression 159 vs 150; mean fill: r 0.503, CV R2 0.241; mean + flag: CV R2 0.239
BMI hidden for 110 heaviest: progression 209 vs 133; mean fill: r 0.352, CV R2 0.112; mean + flag: CV R2 0.289
heaviest: mean + flag beats the mean alone in 48 of 50 folds