Customer Churn Prediction in Python: Complete End-to-End Classification Project

In subscription businesses — from SaaS platforms and telecommunications networks to consumer streaming services and fitness clubs — acquiring a new customer costs between 5 to 25 times more than retaining an existing one. When a customer quietly cancels their subscription, the business suffers not only the immediate loss of recurring monthly revenue (MRR), but also forfeits the multi-year Customer Lifetime Value (CLV) required to amortize initial acquisition expenditures.

Customer Churn Prediction is the application of supervised binary classification to identify at-risk subscribers before they depart. By detecting subtle behavioral signals — such as declining session frequency, spikes in customer support tickets, or approaching contract renewal milestones — machine learning models empower customer success teams to initiate high-touch retention workflows, offer targeted incentives, and preserve enterprise revenue.

1. Key Concepts & Mathematical Notation Glossary

Before establishing loss functions and code pipelines, review the core mathematical symbols and statistical definitions governing imbalanced binary classification:

2. Mathematical Foundations & Cost-Sensitive Optimization

Binary classification algorithms optimize parameterized probability distributions over categorical outputs. However, real-world churn datasets present asymmetric base rates — typically only 10% to 30% of customers churn in any given evaluation window. Standard unweighted optimization objective functions fail under these skewed distributions.

2.1 Binary Cross-Entropy & Balanced Loss Functions

Standard binary logistic regression minimizes unweighted Binary Cross-Entropy (Log-Loss):

When negative examples ($y = 0$, retained subscribers) outnumber positive examples ($y = 1$, churners) by 4-to-1, the non-churn sum dominates the objective gradient. The optimizer easily minimizes overall empirical loss by pushing all probabilities toward zero, leaving the minority class severely underfitted. To correct this distortion, we introduce cost-sensitive sample weighting:

Where Scikit-Learn's class_weight='balanced' computes inverse-frequency weights:

Here $n_1$ denotes the total number of churners and $n_0$ represents non-churners. If only 25% of customers churn ($n_1 = 0.25n$), then $w_1 = 2.0$ while $w_0 = 0.67$. Every misclassified churner penalizes the gradient three times more heavily than a misclassified loyal customer, forcing the decision surface to balance minority class sensitivity.

2.2 The Confusion Matrix & Metric Equations

Evaluating a binary classifier requires dissecting predictions across four fundamental quadrants:

• Precision (Positive Predictive Value): Out of all customers our model flagged as churners, what fraction actually canceled? High precision avoids wasting intervention dollars.
  $$\text{Precision} = \frac{\text{TP}}{\text{TP} + \text{FP}}$$
• Recall (Sensitivity / True Positive Rate): Out of all customers who truly canceled, what fraction did our model catch? High recall prevents customers from slipping through undetected.
  $$\text{Recall} = \frac{\text{TP}}{\text{TP} + \text{FN}}$$
• Specificity (True Negative Rate): Out of all customers who remained active, what fraction were correctly left uninterrupted?
  $$\text{Specificity} = \frac{\text{TN}}{\text{TN} + \text{FP}}$$
• $F_1$ Score: The harmonic mean of precision and recall. Because the harmonic mean punishes extreme imbalances, a model cannot achieve a strong $F_1$ by sacrificing one metric entirely.
  $$F_1 = 2 \cdot \frac{\text{Precision} \cdot \text{Recall}}{\text{Precision} + \text{Recall}} = \frac{2\text{TP}}{2\text{TP} + \text{FP} + \text{FN}}$$

2.3 The Accuracy Paradox

Why is raw classification accuracy dangerously misleading for churn modeling? Consider a telecommunications provider with a 15% annual churn rate. A naive, brainless "dummy classifier" that predicts churn = 0 for every single customer achieves 85% accuracy without running a single calculation:

Yet in production, this 85% accurate model is completely catastrophic: its Recall is exactly 0.0%. It fails to identify a single departing customer, forfeiting 100% of preventable churn revenue. In imbalanced domains, accuracy is a vanity metric; Precision, Recall, PR-AUC, and financial ROI are the ground truth.

2.4 The Business Expected Value Equation

Machine learning models in business do not operate in a vacuum of dimensionless metrics; they operate on a financial ledger. Let each prediction outcome correspond to an empirical dollar value:

• $V_{\text{retain}}$: Customer Lifetime Value saved if an at-risk customer is successfully retained (e.g., 300 USD annual recurring margin).
• $r_{\text{success}}$: Probability that a customer success outreach or discount offer successfully convinces an at-risk subscriber to stay (e.g., 50%).
• $C_{\text{offer}}$: Cost of the intervention campaign, including staff time and promotional discounts (e.g., 20 USD per contacted customer).
• $L_{\text{churn}}$: Net financial loss from an undetected churner, equivalent to the customer acquisition cost (CAC) required to replace them (e.g., 200 USD).

The total expected financial net benefit of operating the model across a customer cohort is given by:

Because missing a true churner ($\text{FN}$) inflicts a 200 USD replacement loss, while needlessly contacting a loyal subscriber ($\text{FP}$) only incurs a 20 USD coupon expense, a False Negative is ten times more expensive than a False Positive. This mathematical asymmetry dictates that the optimal classification threshold $\tau^*$ must be lower than the default 0.50 cutoff.

3. The 6-Stage Machine Learning Pipeline Architecture

A robust churn modeling system follows an end-to-end engineering architecture ensuring data hygiene, reproducible transformations, and automated business activation:

4. Hand-Worked Trace: 20-Customer Confusion Matrix & Metric Calculation

To understand how raw probability predictions translate into confusion matrix counts, work through this complete trace of 20 representative customer accounts evaluated at threshold $\tau = 0.50$:

Summing the outcomes across our 20 customer cohort:

• Actual Churners ($y = 1$): 5 customers (25.0% baseline churn rate).
• Actual Retained ($y = 0$): 15 customers (75.0% retention rate).
• True Positives (TP): 4 (C01, C02, C03, C04 correctly flagged).
• False Negatives (FN): 1 (C05 missed; canceled without intervention).
• False Positives (FP): 3 (C06, C07, C08 incorrectly targeted with retention promotions).
• True Negatives (TN): 12 (C09 through C20 correctly recognized as loyal).

Now, compute all standard performance metrics with step-by-step arithmetic:

5. Complete, Self-Contained Python Implementation

Below is the production-ready Python code. It generates a realistic synthetic Telco churn dataset ($n = 2,000$), handles missing data, engineers behavioral features, executes leak-free ColumnTransformer pipelines, benchmarks 3 model families with StratifiedKFold, optimizes decision thresholds, and scores an unseen at-risk account.

import numpy as np
import pandas as pd
from sklearn.model_selection import StratifiedKFold, cross_val_score, train_test_split
from sklearn.linear_model import LogisticRegression
from sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifier
from sklearn.preprocessing import StandardScaler, OneHotEncoder
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.metrics import classification_report, confusion_matrix, precision_score, recall_score, f1_score

# =============================================================================
# STAGE 1 & 2: SYNTHETIC TELCO DATASET GENERATION WITH MISSINGNESS
# =============================================================================
np.random.seed(42)
n = 2000

data = pd.DataFrame({
    "tenure_months": np.random.randint(0, 72, n),
    "monthly_charges": np.random.normal(65, 25, n).clip(15, 150),
    "contract_type": np.random.choice(["month-to-month", "one-year", "two-year"], n, p=[0.55, 0.25, 0.20]),
    "has_tech_support": np.random.choice([0, 1], n, p=[0.6, 0.4]),
    "num_support_calls": np.random.poisson(1.5, n),
    "payment_method": np.random.choice(["electronic_check", "credit_card", "bank_transfer"], n),
    "internet_service": np.random.choice(["dsl", "fiber", "none"], n, p=[0.35, 0.45, 0.20])
})

# Non-linear probability logit with realistic churn dynamics
churn_logit = (
    -1.5
    - 0.04 * data["tenure_months"]
    + 0.015 * data["monthly_charges"]
    + data["contract_type"].map({"month-to-month": 1.2, "one-year": -0.3, "two-year": -1.0})
    - data["has_tech_support"] * 0.5
    + data["num_support_calls"] * 0.35
    + data["payment_method"].map({"electronic_check": 0.4, "credit_card": -0.2, "bank_transfer": -0.2})
)
churn_probability = 1 / (1 + np.exp(-churn_logit))
data["churned"] = np.random.binomial(1, churn_probability)

# Simulate 15 real-world missing values in monthly_charges (e.g. billing sync latency)
data.loc[np.random.choice(n, 15, replace=False), "monthly_charges"] = np.nan

print(f"Total Customers: {len(data)}")
print(f"Churn Count:    {data['churned'].sum()} ({data['churned'].mean():.1%} baseline rate)")

# =============================================================================
# STAGE 3: DOMAIN FEATURE ENGINEERING
# =============================================================================
# Ratio of monthly bill relative to lifetime tenure
data["avg_charge_per_tenure_month"] = data["monthly_charges"] / (data["tenure_months"] + 1)
# Binary flag for onboarding vulnerability (< 6 months tenure)
data["is_new_customer"] = (data["tenure_months"] < 6).astype(int)
# Frequent support friction indicator (>= 3 support calls)
data["high_support_contact"] = (data["num_support_calls"] >= 3).astype(int)

X = data.drop(columns=["churned"])
y = data["churned"]

# =============================================================================
# STAGE 4: LEAK-FREE PREPROCESSING PIPELINE WITH COLUMNTRANSFORMER
# =============================================================================
num_cols = [
    "tenure_months", "monthly_charges", "has_tech_support",
    "num_support_calls", "avg_charge_per_tenure_month",
    "is_new_customer", "high_support_contact"
]
cat_cols = ["contract_type", "payment_method", "internet_service"]

num_pipe = Pipeline([
    ("imputer", SimpleImputer(strategy="median")),
    ("scaler", StandardScaler())
])
cat_pipe = Pipeline([
    ("imputer", SimpleImputer(strategy="most_frequent")),
    ("encoder", OneHotEncoder(drop="first", handle_unknown="ignore"))
])

preprocessor = ColumnTransformer([
    ("num", num_pipe, num_cols),
    ("cat", cat_pipe, cat_cols)
])

# =============================================================================
# STAGE 5: 5-FOLD STRATIFIED CROSS-VALIDATION BENCHMARK
# =============================================================================
skf = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)

models = {
    "Logistic Regression": Pipeline([
        ("prep", preprocessor),
        ("model", LogisticRegression(class_weight="balanced", max_iter=1000))
    ]),
    "Random Forest": Pipeline([
        ("prep", preprocessor),
        ("model", RandomForestClassifier(n_estimators=100, max_depth=8, class_weight="balanced", random_state=42))
    ]),
    "Gradient Boosting": Pipeline([
        ("prep", preprocessor),
        ("model", GradientBoostingClassifier(n_estimators=100, max_depth=3, random_state=42))
    ])
}

print("\n" + "="*65)
print("5-FOLD STRATIFIED CROSS-VALIDATION BENCHMARK RESULTS")
print("="*65)
for name, model in models.items():
    f1 = cross_val_score(model, X, y, cv=skf, scoring="f1")
    auc = cross_val_score(model, X, y, cv=skf, scoring="roc_auc")
    rec = cross_val_score(model, X, y, cv=skf, scoring="recall")
    print(f"{name:20s} | F1: {f1.mean():.3f} (+/- {f1.std():.3f}) | AUC: {auc.mean():.3f} (+/- {auc.std():.3f}) | Recall: {rec.mean():.3f}")

# =============================================================================
# STAGE 6: DECISION THRESHOLD TUNING & PRODUCTION INFERENCE
# =============================================================================
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, stratify=y, random_state=42)

# Train production model (Gradient Boosting)
prod_pipe = models["Gradient Boosting"].fit(X_train, y_train)
probs_test = prod_pipe.predict_proba(X_test)[:, 1]

print("\n" + "="*65)
print("DECISION THRESHOLD SWEEP ON TEST SET (n=400 customers)")
print("="*65)
for tau in [0.20, 0.35, 0.50, 0.65, 0.80]:
    preds = (probs_test >= tau).astype(int)
    cm = confusion_matrix(y_test, preds)
    tn, fp, fn, tp = cm.ravel()
    p = precision_score(y_test, preds, zero_division=0)
    r = recall_score(y_test, preds)
    f = f1_score(y_test, preds)
    print(f"tau={tau:.2f} -> Prec: {p:.3f} | Rec: {r:.3f} | F1: {f:.3f} | TP={tp:2d} FN={fn:2d} FP={fp:2d} TN={tn:3d}")

# Live scoring on an unseen at-risk subscriber account
new_account = pd.DataFrame([{
    "tenure_months": 2,
    "monthly_charges": 89.50,
    "contract_type": "month-to-month",
    "has_tech_support": 0,
    "num_support_calls": 4,
    "payment_method": "electronic_check",
    "internet_service": "fiber",
    "avg_charge_per_tenure_month": 29.83,
    "is_new_customer": 1,
    "high_support_contact": 1
}])

churn_risk = prod_pipe.predict_proba(new_account)[0, 1]
print("\n" + "="*65)
print(f"LIVE ACCOUNT CHURN RISK SCORE: {churn_risk:.1%}")
if churn_risk >= 0.35:
    print(">> ACTION TRIGGERED: High Churn Probability (>= 35% threshold).")
    print(">> DISPATCH: Trigger automated 20 USD monthly discount & priority CS phone outreach.")
else:
    print(">> ACTION: Low risk. Maintain standard automated newsletter cadences.")
print("="*65)

6. Empirical Results & Business Impact Analysis

Let us analyze the verified empirical results produced across 5-Fold Stratified Cross-Validation on the 2,000 customer dataset:

Why does balanced Logistic Regression achieve the highest recall (74.4%) and top ROC-AUC (0.810) over unweighted Gradient Boosting? Because Gradient Boosting defaults to standard unweighted loss; under a 28.1% minority class distribution, its trees focus predominantly on correctly partitioning the 71.9% non-churners. Conversely, class_weight='balanced' forces the model to treat minority misclassifications with triple priority.

6.1 Business ROI of Threshold Tuning (tau = 0.35 vs. tau = 0.50)

Examine the confusion matrix transitions across our holdout test set ($n = 400$ subscribers, containing 112 true churners and 288 loyal customers):

Now calculate the economic bottom line using our corporate cost parameters ($V_{\text{retain}} = 300\text{ USD}$, $r_{\text{success}} = 0.50$, $C_{\text{offer}} = 20\text{ USD}$, $L_{\text{churn}} = 200\text{ USD}$):

• At Default Threshold $\tau = 0.50$:
  $$\mathbb{E}[\text{Net}] = 50 \times (150 - 20) - 32 \times 20 - 62 \times 200 = 6{,}500 - 640 - 12{,}400 = -6{,}540\text{ USD}$$
• At Optimized Threshold $\tau = 0.35$:
  $$\mathbb{E}[\text{Net}] = 71 \times (150 - 20) - 59 \times 20 - 41 \times 200 = 9{,}230 - 1{,}180 - 8{,}200 = -150\text{ USD}$$
• Net Economic Gain: Shifting the decision threshold down to 0.35 captures 21 additional churners (+18.8% recall gain). While it causes 27 additional false alarms (costing 540 USD in promotions), it prevents 21 customers from walking away (saving 4,200 USD in acquisition replacement costs) — yielding a net bottom-line gain of +6,390 USD on just 400 customers!

7. Production Gotchas, Data Leakage & Engineering Pitfalls

In production environments, churn models fail not due to mathematical nuances, but due to insidious temporal and behavioral data leakage. Watch for these four common failure modes:

• Post-Decision Feature Leakage: Never include features generated after the subscriber made their decision to churn. Examples include date_of_cancellation_request, exit_survey_rating, or retention_call_duration. A production churn model must strictly use data available 30 to 60 days prior to the cancellation event.
• Survivorship Bias in Tenure: If you compute average monthly usage by dividing total usage by tenure without an offset ($x / t$), newborn subscribers ($t = 0$) create zero-division errors or extreme spikes. Always add Laplace smoothing ($t + 1$).
• Collinearity Between Tenure and Total Charges: In subscription billing, $\text{TotalCharges} \approx \text{Tenure} \times \text{MonthlyCharges}$. Including all three raw features creates severe multicollinearity in linear models. Engineer ratio features or use tree ensembles that naturally handle correlated inputs.
• Concept Drift & Seasonality: Customer churn patterns evolve. A model trained in November will miscalculate holiday gift subscriptions expiring in January. Establish automated quarterly retraining pipelines with data drift alerts using Population Stability Index (PSI).

8. Practice Exercises & Next Steps

Consolidate your applied machine learning mastery with these three production challenges:

1. The Real Kaggle Telco Churn Migration: Download the real-world Kaggle IBM Telco Churn dataset (7,043 rows, 21 features). Replace the synthetic generator with real data, inspect the 11 empty space strings in TotalCharges, and benchmark your pipeline.
2. Explainability with SHAP Values: Install shap and calculate TreeExplainer values for the Gradient Boosting model. Generate a SHAP summary plot to verify whether contract type or tenure is the primary driver of customer cancellations.
3. Survival Analysis Transition: Binary classification predicts if a customer will churn in a fixed 30-day window. Implement a Kaplan-Meier survival curve using the lifelines library to model when churn will occur across the entire customer lifecycle.