Data Leakage: the Easiest Way to Fool Yourself

Chapters 1 and 2 built the master dataframe and showed six exploratory angles. Now I need to pick, for each of the four scenarios defined back in Chapter 1, which features actually make sense to use. Full executed notebook, public on Colab. This notebook also starts from scratch (reloads the 9 CSVs and rebuilds the master dataframe with the same code from Chapter 1) to stay self-contained.

A new feature first: customer-seller distance

Chapter 2 showed that customer state correlates with delivery time. Before diving into the scenarios, I computed a real distance, not just "same state or not": average latitude/longitude per zip code prefix (olist_geolocation_dataset, over a million rows, hence averaging first), and the Haversine formula between customer and seller.

geo_media = geolocation.groupby('geolocation_zip_code_prefix')[['geolocation_lat', 'geolocation_lng']].mean().reset_index()

def haversine(lat1, lon1, lat2, lon2):
    """Distância em km entre dois pontos lat/lng, aproximando a Terra como esfera de raio 6371km."""
    R = 6371
    phi1, phi2 = np.radians(lat1), np.radians(lat2)
    dphi = np.radians(lat2 - lat1)
    dlambda = np.radians(lon2 - lon1)
    a = np.sin(dphi / 2) ** 2 + np.cos(phi1) * np.cos(phi2) * np.sin(dlambda / 2) ** 2
    return 2 * R * np.arcsin(np.sqrt(a))

# força int64 dos dois lados da chave, alguns ambientes inferem tipos diferentes pro CEP e isso quebra o merge silenciosamente
m = mestre.copy()
m['customer_zip_code_prefix'] = m['customer_zip_code_prefix'].astype('int64')
m['seller_zip_code_prefix'] = m['seller_zip_code_prefix'].astype('int64')
geo_media['geolocation_zip_code_prefix'] = geo_media['geolocation_zip_code_prefix'].astype('int64')

m = m.merge(
    geo_media.rename(columns={'geolocation_zip_code_prefix': 'customer_zip_code_prefix', 'geolocation_lat': 'cust_lat', 'geolocation_lng': 'cust_lng'}),
    on='customer_zip_code_prefix', how='left',
)
m = m.merge(
    geo_media.rename(columns={'geolocation_zip_code_prefix': 'seller_zip_code_prefix', 'geolocation_lat': 'sell_lat', 'geolocation_lng': 'sell_lng'}),
    on='seller_zip_code_prefix', how='left',
)
m['distance_km'] = haversine(m['cust_lat'], m['cust_lng'], m['sell_lat'], m['sell_lng'])

groupby(...).mean() collapses the geolocation table (over a million rows, several coordinates per zip prefix) into a single average lat/lng per prefix, otherwise the next merge would multiply every order by every coordinate registered under that prefix. haversine is the classic distance formula between two points on a sphere: it converts degrees to radians, and the sine/cosine combination inside arcsin gives the angular distance, multiplied by Earth's radius (6,371 km) that becomes a real km distance, not a straight line on a flat map, which would distort quite a bit across a country the size of Brazil. The merge happens twice, once to bring in the customer's coordinate, once for the seller's, each renaming the generic columns (geolocation_lat/lng) to a specific name (cust_lat/sell_lat) before joining, otherwise the second merge call would overwrite the columns the first one had just brought in.

Found a real bug at this step before publishing: without the three .astype('int64') lines right at the start, the merge between customer_zip_code_prefix/seller_zip_code_prefix and the geolocation table lost 70% of the data (only 30% coverage), because the zip-code column's data type wasn't guaranteed to match on both sides across every environment (one side as int64, the other as object or float64, and pandas simply doesn't match keys when the types diverge, without raising any error). Forced explicit int64 on all three before merging, and coverage went back to what it should be: only 581 orders out of 118,310 ended up without a distance (99.5% coverage). Average customer-seller distance across all of Brazil is 597 km, with a max of 8,678 km, a number that proves just how huge the country is.

Proof in image form too: below, a sample of 500 orders, the customer placed on the plane (longitude and latitude) and the distance to the seller as height. Spin it and notice the distance peaks don't follow a simple "further north means further away" pattern, they show up scattered across very different latitudes, a sign that distance depends more on which specific seller handled that order than on the customer's region alone.

Loading real data...

Scenario 1: delivery delay, and the data leakage demonstration

Candidates: distance_km, product_weight_g, price, freight_value, approval time (approval_hours), purchase month. The atrasado target is defined by comparing order_delivered_customer_date against order_estimated_delivery_date. And this is where the classic trap of every ML course lives: if I include atraso_dias (the actual difference between those two dates) as a feature, the model isn't learning to predict delay, it's reading the answer straight off the exam.

m['approval_hours'] = (m['order_approved_at'] - m['order_purchase_timestamp']).dt.total_seconds() / 3600
m['month'] = m['order_purchase_timestamp'].dt.month
m['atraso_dias'] = (m['order_delivered_customer_date'] - m['order_estimated_delivery_date']).dt.days

delay_df = (
    m.drop_duplicates(subset='order_id')
    [['order_id', 'distance_km', 'product_weight_g', 'price', 'freight_value', 'approval_hours', 'month', 'atraso_dias']]
    .dropna()
)
delay_df['atrasado'] = (delay_df['atraso_dias'] > 0).astype(int)

approval_hours measures how long Olist took to approve payment after purchase, in hours (not days, because most orders approve within a few hours, and dt.days would round that down to zero almost every time). month pulls out just the purchase month, the seasonality feature. atraso_dias is the gap between actual and promised delivery (negative is early, positive is a real delay), and atrasado becomes the binary target. drop_duplicates(subset='order_id') brings the master dataframe back to one row per order (it had one row per item, remember Chapter 1), otherwise orders with many items would weigh more in training just for having more rows, not for being more delayed.

Output: 95,968 orders used, real delay rate of 6.76%.

With the dataframe ready, I trained two RandomForests, same seed, same train/test split, changing only the feature list:

features_honestas = ['distance_km', 'product_weight_g', 'price', 'freight_value', 'approval_hours', 'month']
features_vazadas = features_honestas + ['atraso_dias']

X = delay_df[features_honestas]
y = delay_df['atrasado']
Xtr, Xte, ytr, yte = train_test_split(X, y, test_size=0.2, random_state=42, stratify=y)

modelo_honesto = RandomForestClassifier(n_estimators=100, max_depth=8, random_state=42)
modelo_honesto.fit(Xtr, ytr)
acc_honesto = accuracy_score(yte, modelo_honesto.predict(Xte))
auc_honesto = roc_auc_score(yte, modelo_honesto.predict_proba(Xte)[:, 1])

Xl = delay_df[features_vazadas]
Xltr, Xlte, yltr, ylte = train_test_split(Xl, y, test_size=0.2, random_state=42, stratify=y)
modelo_vazado = RandomForestClassifier(n_estimators=100, max_depth=8, random_state=42)
modelo_vazado.fit(Xltr, yltr)
acc_vazado = accuracy_score(ylte, modelo_vazado.predict(Xlte))
auc_vazado = roc_auc_score(ylte, modelo_vazado.predict_proba(Xlte)[:, 1])

stratify=y guarantees the delay proportion (6.76%) stays the same in train and test, important with imbalance this size. Both blocks are structurally identical, RandomForestClassifier(n_estimators=100, max_depth=8, random_state=42) in both, same random_state=42 on both splits. The only real difference between the two entire training runs is a single extra column in the second feature list: atraso_dias.

Loading real data...

The honest model lands at AUC 0.7388. The leaky model lands at AUC 1.0, perfect, because atraso_dias > 0 is practically the target's own definition. One detail worth calling out: the honest model's raw accuracy (93.24%) looks great at first glance, but since only 6.76% of orders are late, a dumb model that always guesses "on time" already scores 93.24% without learning anything. Accuracy alone is misleading with classes this imbalanced, AUC is the metric that tells the real story here.

Also worth noting, the honest model's feature importance:

importancias = sorted(zip(features_honestas, modelo_honesto.feature_importances_), key=lambda x: -x[1])

zip pairs each feature name with its matching importance (same order features_honestas was passed into training), and sorted(..., key=lambda x: -x[1]) sorts from largest weight to smallest (the minus sign flips the default ascending order).

Loading real data...

month leads by a wide margin, at 0.36 importance, practically double the runner-up (distance_km, 0.18). This ties directly back to Chapter 2's finding: November (Black Friday) strains logistics in a way that physical distance alone can't explain. Seasonality matters more for predicting delay than how far the seller is from the customer.

Scenario 2: review score

Candidates: atraso_dias, delivery time, price, freight, number of installments. Since review_score is a bimodal ordinal scale (Chapter 2 already showed this), I used mutual information instead of plain Pearson correlation:

m['delivery_days'] = (m['order_delivered_customer_date'] - m['order_purchase_timestamp']).dt.days

s2 = (
    m.drop_duplicates(subset='order_id')
    [['order_id', 'atraso_dias', 'delivery_days', 'price', 'freight_value', 'payment_installments', 'review_score']]
    .dropna()
)

X2 = s2[['atraso_dias', 'delivery_days', 'price', 'freight_value', 'payment_installments']]
y2 = s2['review_score']
mi = mutual_info_classif(X2, y2, random_state=42)

delivery_days is the real delivery time (purchase to arrival), different from atraso_dias (which compares against the promised date, not the purchase date). Again drop_duplicates(subset='order_id') brings it back to one row per order. X2 gathers the five candidates, y2 is the score (review_score, the ordinal categorical target). mutual_info_classif measures how much information each continuous feature carries about a categorical target even when the relationship isn't linear, unlike Pearson correlation which I already used in Chapter 2, which assumes a linear relationship, and a bimodal scale like a review score doesn't behave that way.

atraso_dias and delivery_days dominate, almost 10 times more informative than price or freight (95,829 orders used). Confirms, with a different technique, exactly what Chapter 2's correlation already pointed at.

Scenario 3: freight value

Candidates: weight and the product's three dimensions. Plain Pearson correlation is enough here:

s3 = (
    m.drop_duplicates(subset='order_id')
    [['order_id', 'freight_value', 'product_weight_g', 'product_length_cm', 'product_height_cm', 'product_width_cm']]
    .dropna()
)

correlacao_frete = s3[['freight_value', 'product_weight_g', 'product_length_cm', 'product_height_cm', 'product_width_cm']].corr()['freight_value'].sort_values(ascending=False)

.corr() computes the full correlation matrix across the five columns, and ['freight_value'] grabs just the column for freight: each variable's correlation with itself (always 1.0, which is why freight_value doesn't show up in the table below) and with the other four.

Output: 98,650 orders used.

Weight moderately dominates over any single dimension, which makes physical sense: carriers charge by volumetric weight, but actual weight pulls the bill harder than any one side of the box alone.

Scenario 4: customer segmentation (RFM)

This isn't feature selection in the supervised sense, it's engineering the three variables that will feed Chapter 5's clustering: Recency, Frequency, and Monetary value per customer.

data_max = orders['order_purchase_timestamp'].max()

orders_validos = orders[orders['order_status'] != 'canceled']
oc = orders_validos.merge(customers, on='customer_id', how='left')
pagamento_por_pedido = payments.groupby('order_id')['payment_value'].sum().reset_index()
oc = oc.merge(pagamento_por_pedido, on='order_id', how='left')

rfm = oc.groupby('customer_unique_id').agg(
    recencia_dias=('order_purchase_timestamp', lambda x: (data_max - x.max()).days),
    frequencia=('order_id', 'nunique'),
    valor_monetario=('payment_value', 'sum'),
).reset_index()

data_max is the most recent purchase date in the whole dataset, the "today" that Recency uses as its reference (can't use today's real date, the dataset is from 2018, that would inflate everyone's recency equally and destroy the relative comparison between customers). orders_validos drops canceled orders before any aggregation, a customer who canceled shouldn't count as an "active customer" for RFM. pagamento_por_pedido sums the amount paid, which can come in several installments recorded as separate rows (the same row-multiplication I already saw in Chapter 1's merge), before joining it back onto the orders. With that ready, groupby('customer_unique_id') computes all three metrics at once: recencia_dias (days since the last purchase up to data_max), frequencia (unique orders), and valor_monetario (amount paid, summed).

95,560 unique customers. And the number that stood out the most: only 2,924 customers (3.06% of the total) placed more than one order. Average frequency is 1.034, meaning the entire marketplace is dominated by one-time purchases. That completely changes how to think about segmentation in Chapter 5: there won't be a large "loyal customer" category, RFM here will likely separate more by amount spent and recency than by purchase frequency.

Wrapping up the chapter

Customer-seller distance computed (and a real merge bug fixed along the way), data leakage demonstrated in practice (the gap between honest AUC 0.74 and leaky AUC 1.0 is the most concrete lesson I can give on the subject), and each scenario's features picked with a technique suited to its problem type. Seasonality turned out to matter more than distance for delay, delay dominates the review score, weight dominates freight, and the customer base is overwhelmingly one-time buyers. Next chapter I train the real models for each scenario, with experiment tracking and metric comparison.