Customer Churn Analysis: Feature Identification, Visualization & Retention Recommendations

Overview

I'll build a complete, well-documented pipeline that handles data loading, feature analysis, visualization, and business recommendations.

---

Complete Python Script

"""
Customer Churn Analysis Pipeline
================================
Identifies top predictive features, visualizes correlations,
and generates actionable retention recommendations.
"""

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import seaborn as sns
from sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifier
from sklearn.preprocessing import LabelEncoder, StandardScaler
from sklearn.model_selection import train_test_split, cross_val_score
from sklearn.metrics import classification_report, roc_auc_score
from sklearn.inspection import permutation_importance
import warnings
warnings.filterwarnings('ignore')

# ─────────────────────────────────────────────
# 1.  DATA GENERATION  (replace with real data)
# ─────────────────────────────────────────────

def generate_churn_dataset(n: int = 5_000, seed: int = 42) -> pd.DataFrame:
    """
    Generates a realistic synthetic churn dataset.
    Replace this function with your actual data loader, e.g.:
        pd.read_csv('churn_data.csv')
    """
    rng = np.random.default_rng(seed)

    tenure          = rng.integers(1, 72, n)                          # months
    monthly_charges = rng.uniform(20, 120, n)                         # dollars
    total_charges   = tenure * monthly_charges * rng.uniform(0.9, 1.1, n)
    num_products    = rng.integers(1, 6, n)
    support_calls   = rng.integers(0, 15, n)
    contract_type   = rng.choice(['Month-to-Month', 'One Year', 'Two Year'], n,
                                  p=[0.55, 0.25, 0.20])
    payment_method  = rng.choice(['Electronic Check', 'Mailed Check',
                                  'Bank Transfer', 'Credit Card'], n)
    internet_service= rng.choice(['DSL', 'Fiber Optic', 'None'], n,
                                  p=[0.35, 0.45, 0.20])
    age             = rng.integers(18, 80, n)
    satisfaction    = rng.integers(1, 6, n)                           # 1-5 scale
    login_freq      = rng.integers(0, 30, n)                          # logins/month

    # Realistic churn probability driven by features
    churn_prob = (
        0.40 * (tenure < 12).astype(float)
      + 0.25 * (monthly_charges > 80).astype(float)
      + 0.20 * (support_calls > 8).astype(float)
      + 0.30 * (contract_type == 'Month-to-Month').astype(float)
      - 0.20 * (num_products > 3).astype(float)
      - 0.15 * (satisfaction > 3).astype(float)
      + 0.10 * rng.uniform(0, 1, n)             # noise
    )
    churn_prob = np.clip(churn_prob / churn_prob.max(), 0.05, 0.95)
    churn      = (rng.uniform(0, 1, n) < churn_prob).astype(int)

    df = pd.DataFrame({
        'CustomerID'      : [f'CUST-{i:05d}' for i in range(n)],
        'Tenure'          : tenure,
        'MonthlyCharges'  : monthly_charges.round(2),
        'TotalCharges'    : total_charges.round(2),
        'NumProducts'     : num_products,
        'SupportCalls'    : support_calls,
        'ContractType'    : contract_type,
        'PaymentMethod'   : payment_method,
        'InternetService' : internet_service,
        'Age'             : age,
        'Satisfaction'    : satisfaction,
        'LoginFrequency'  : login_freq,
        'Churn'           : churn
    })
    return df


# ─────────────────────────────────────────────
# 2.  PREPROCESSING
# ─────────────────────────────────────────────

def preprocess(df: pd.DataFrame) -> tuple[pd.DataFrame, pd.Series,
                                          list[str], dict]:
    """
    Encodes categoricals, scales numerics, returns feature matrix.
    Returns: X, y, feature_names, encoders
    """
    df = df.copy()

    # ── Drop ID column ────────────────────────
    df.drop(columns=['CustomerID'], errors='ignore', inplace=True)

    # ── Separate target ───────────────────────
    y = df.pop('Churn')

    # ── Identify column types ─────────────────
    cat_cols = df.select_dtypes(include=['object', 'category']).columns.tolist()
    num_cols = df.select_dtypes(include=[np.number]).columns.tolist()

    encoders = {}
    for col in cat_cols:
        le = LabelEncoder()
        df[col] = le.fit_transform(df[col].astype(str))
        encoders[col] = le

    scaler = StandardScaler()
    df[num_cols] = scaler.fit_transform(df[num_cols])
    encoders['scaler'] = scaler

    return df, y, df.columns.tolist(), encoders


# ─────────────────────────────────────────────
# 3.  FEATURE IMPORTANCE ANALYSIS
# ─────────────────────────────────────────────

def identify_top_features(X: pd.DataFrame, y: pd.Series,
                           top_n: int = 3) -> pd.DataFrame:
    """
    Combines three importance signals for robust feature ranking:
      • Random Forest impurity importance
      • Gradient Boosting importance
      • Permutation importance (model-agnostic)
    """
    X_tr, X_te, y_tr, y_te = train_test_split(
        X, y, test_size=0.2, random_state=42, stratify=y)

    # ── Model 1: Random Forest ─────────────────
    rf = RandomForestClassifier(n_estimators=300, random_state=42, n_jobs=-1)
    rf.fit(X_tr, y_tr)
    rf_imp = pd.Series(rf.feature_importances_, index=X.columns)

    # ── Model 2: Gradient Boosting ─────────────
    gb = GradientBoostingClassifier(n_estimators=200, random_state=42)
    gb.fit(X_tr, y_tr)
    gb_imp = pd.Series(gb.feature_importances_, index=X.columns)

    # ── Model 3: Permutation Importance ────────
    perm = permutation_importance(rf, X_te, y_te, n_repeats=15,
                                  random_state=42, n_jobs=-1)
    perm_imp = pd.Series(perm.importances_mean, index=X.columns)

    # ── Normalise & combine ────────────────────
    def norm(s): return (s - s.min()) / (s.max() - s.min() + 1e-9)

    importance_df = pd.DataFrame({
        'RF_Importance'   : norm(rf_imp),
        'GB_Importance'   : norm(gb_imp),
        'Perm_Importance' : norm(perm_imp),
    })
    importance_df['Composite_Score'] = importance_df.mean(axis=1)
    importance_df.sort_values('Composite_Score', ascending=False, inplace=True)

    # ── CV AUC scores ──────────────────────────
    auc_rf = cross_val_score(rf, X, y, cv=5, scoring='roc_auc').mean()
    auc_gb = cross_val_score(gb, X, y, cv=5, scoring='roc_auc').mean()
    print(f"\n{'='*55}")
    print(f"  Model Performance (5-Fold CV AUC)")
    print(f"  Random Forest     : {auc_rf:.4f}")
    print(f"  Gradient Boosting : {auc_gb:.4f}")
    print(f"{'='*55}")

    print(f"\n  Top {top_n} Predictive Features:")
    print(f"  {'Feature':<22} {'Score':>8}")
    print(f"  {'-'*32}")
    for feat, row in importance_df.head(top_n).iterrows():
        print(f"  {feat:<22} {row['Composite_Score']:>8.4f}")
    print(f"{'='*55}\n")

    return importance_df


# ─────────────────────────────────────────────
# 4.  VISUALISATION
# ─────────────────────────────────────────────

PALETTE   = {'primary': '#2E86AB', 'danger': '#E84855',
             'success': '#3BB273', 'warning': '#F4A261', 'neutral': '#6C757D'}
PLT_STYLE = {'figure.facecolor': '#0F1117', 'axes.facecolor': '#1A1D2E',
             'axes.edgecolor': '#3A3D52', 'text.color': '#E8E8F0',
             'axes.labelcolor': '#E8E8F0', 'xtick.color': '#A0A0C0',
             'ytick.color': '#A0A0C0', 'grid.color': '#2A2D3E',
             'grid.alpha': 0.6}
plt.rcParams.update(PLT_STYLE)


def plot_all(df_raw: pd.DataFrame, importance_df: pd.DataFrame) -> None:
    """Master visualisation: 6 coordinated panels."""

    # ── Prepare numeric-only frame ─────────────
    num_df = df_raw.select_dtypes(include=[np.number])

    fig = plt.figure(figsize=(22, 18))
    fig.suptitle('Customer Churn Analysis Dashboard',
                 fontsize=22, fontweight='bold',
                 color='#E8E8F0', y=0.98)

    gs = gridspec.GridSpec(3, 3, figure=fig,
                           hspace=0.45, wspace=0.40)

    # ╔══════════════════════════════════════════╗
    # ║  PANEL 1 – Full Correlation Matrix       ║
    # ╚══════════════════════════════════════════╝
    ax1 = fig.add_subplot(gs[0, :2])
    corr = num_df.corr()
    mask = np.triu(np.ones_like(corr, dtype=bool))

    cmap = sns.diverging_palette(220, 10, as_cmap=True)
    sns.heatmap(corr, mask=mask, cmap=cmap, vmax=1, vmin=-1, center=0,
                annot=True, fmt='.2f', linewidths=0.5,
                linecolor='#2A2D3E', ax=ax1,
                annot_kws={'size': 8, 'color': '#E8E8F0'},
                cbar_kws={'shrink': 0.8})

    ax1.set_title('Feature Correlation Matrix (Lower Triangle)',
                  fontsize=13, pad=10, color='#E8E8F0')
    ax1.tick_params(axis='x', rotation=45, labelsize=9)
    ax1.tick_params(axis='y', rotation=0, labelsize=9)

    # ╔══════════════════════════════════════════╗
    # ║  PANEL 2 – Churn Rate Donut              ║
    # ╚══════════════════════════════════════════╝
    ax2 = fig.add_subplot(gs[0, 2])
    churn_counts = df_raw['Churn'].value_counts()
    labels = [f"Retained\n{churn_counts[0]:,}", f"Churned\n{churn_counts[1]:,}"]
    colors = [PALETTE['success'], PALETTE['danger']]
    wedges, texts, autotexts = ax2.pie(
        churn_counts, labels=labels, autopct='%1.1f%%',
        colors=colors, startangle=90,
        wedgeprops=dict(width=0.55, edgecolor='#1A1D2E', linewidth=2),
        textprops={'color': '#E8E8F0', 'fontsize': 10})
    for at in autotexts:
        at.set_fontsize(11)
        at.set_fontweight('bold')
    ax2.set_title('Churn Distribution', fontsize=13,
                  pad=10, color='#E8E8F0')

    # ╔══════════════════════════════════════════╗
    # ║  PANEL 3 – Composite Feature Importance  ║
    # ╚══════════════════════════════════════════╝
    ax3 = fig.add_subplot(gs[1, :2])
    top10 = importance_df.head(10)
    colors_bar = [PALETTE['danger'] if i < 3 else PALETTE['primary']
                  for i in range(len(top10))]

    bars = ax3.barh(top10.index[::-1], top10['Composite_Score'][::-1],
                    color=colors_bar[::-1], edgecolor='none', height=0.65)

    for bar, score in zip(bars, top10['Composite_Score'][::-1]):
        ax3.text(bar.get_width() + 0.005, bar.get_y() + bar.get_height()/2,
                 f'{score:.3f}', va='center', ha='left',
                 color='#E8E8F0', fontsize=9)

    ax3.set_xlabel('Composite Importance Score', fontsize=10)
    ax3.set_title('Feature Importance (RF + GB + Permutation)',
                  fontsize=13, pad=10, color='#E8E8F0')
    ax3.axvline(x=0, color='#3A3D52', linewidth=0.8)
    ax3.set_xlim(0, top10['Composite_Score'].max() * 1.18)

    # Legend for top-3 highlight
    from matplotlib.patches import Patch
    legend_elems = [Patch(facecolor=PALETTE['danger'],  label='Top 3 Features'),
                    Patch(facecolor=PALETTE['primary'], label='Other Features')]
    ax3.legend(handles=legend_elems, loc='lower right',
               facecolor='#1A1D2E', edgecolor='#3A3D52',
               labelcolor='#E8E8F0', fontsize=9)

    # ╔══════════════════════════════════════════╗
    # ║  PANEL 4 – Churn by Contract Type        ║
    # ╚══════════════════════════════════════════╝
    ax4 = fig.add_subplot(gs[1, 2])
    contract_churn = (df_raw.groupby('ContractType')['Churn']
                            .agg(['mean', 'count'])
                            .rename(columns={'mean': 'ChurnRate', 'count': 'Count'}))
    contract_churn['ChurnRate'] *= 100

    bar_colors = [PALETTE['danger'] if r > 30 else PALETTE['warning']
                  if r > 15 else PALETTE['success']
                  for r in contract_churn['ChurnRate']]

    ax4.bar(contract_churn.index, contract_churn['ChurnRate'],
            color=bar_colors, edgecolor='none', width=0.55)
    ax4.axhline(y=df_raw['Churn'].mean()*100, color='#A0A0C0',
                linestyle='--', linewidth=1.2, label='Overall avg')

    for i, (idx, row) in enumerate(contract_churn.iterrows()):
        ax4.text(i, row['ChurnRate'] + 0.8, f"{row['ChurnRate']:.1f}%",
                 ha='center', va='bottom', color='#E8E8F0',
                 fontsize=10, fontweight='bold')
        ax4.text(i, -3.5, f"n={row['Count']:,}", ha='center',
                 color='#A0A0C0', fontsize=8)

    ax4.set_title('Churn Rate by Contract Type',
                  fontsize=13, pad=10, color='#E8E8F0')
    ax4.set_ylabel('Churn Rate (%)', fontsize=10)
    ax4.set_ylim(-6, contract_churn['ChurnRate'].max() * 1.2)
    ax4.tick_params(axis='x', rotation=10, labelsize=9)
    ax4.legend(facecolor='#1A1D2E', edgecolor='#3A3D52',
               labelcolor='#E8E8F0', fontsize=9)

    # ╔══════════════════════════════════════════╗
    # ║  PANEL 5 – Tenure vs Monthly Charges     ║
    # ╚══════════════════════════════════════════╝
    ax5 = fig.add_subplot(gs[2, :2])
    churned     = df_raw[df_raw['Churn'] == 1]
    not_churned = df_raw[df_raw['Churn'] == 0]

    ax5.scatter(not_churned['Tenure'], not_churned['MonthlyCharges'],
                c=PALETTE['success'], alpha=0.25, s=18, label='Retained',
                edgecolors='none')
    ax5.scatter(churned['Tenure'], churned['MonthlyCharges'],
                c=PALETTE['danger'],  alpha=0.45, s=22, label='Churned',
                edgecolors='none')

    # Risk zone annotation
    ax5.axvspan(0, 12, alpha=0.07, color=PALETTE['danger'],
                label='High-Risk Tenure Zone (<12 mo)')
    ax5.axhspan(80, df_raw['MonthlyCharges'].max(),
                alpha=0.06, color=PALETTE['warning'],
                label='High-Charge Zone (>$80)')
    ax5.axvline(x=12, color=PALETTE['danger'],  linestyle=':', linewidth=1.2)
    ax5.axhline(y=80, color=PALETTE['warning'], linestyle=':', linewidth=1.2)

    ax5.set_xlabel('Tenure (months)', fontsize=10)
    ax5.set_ylabel('Monthly Charges ($)', fontsize=10)
    ax5.set_title('Tenure vs Monthly Charges  (Churn Overlay)',
                  fontsize=13, pad=10, color='#E8E8F0')
    ax5.legend(facecolor='#1A1D2E', edgecolor='#3A3D52',
               labelcolor='#E8E8F0', fontsize=8, ncol=2)

    # ╔══════════════════════════════════════════╗
    # ║  PANEL 6 – Support Calls Distribution   ║
    # ╚══════════════════════════════════════════╝
    ax6 = fig.add_subplot(gs[2, 2])
    bins = range(0, df_raw['SupportCalls'].max() + 2)

    ax6.hist(not_churned['SupportCalls'], bins=bins, alpha=0.65,
             color=PALETTE['success'], label='Retained',
             edgecolor='#1A1D2E', linewidth=0.4, density=True)
    ax6.hist(churned['SupportCalls'],     bins=bins, alpha=0.65,
             color=PALETTE['danger'],  label='Churned',
             edgecolor='#1A1D2E', linewidth=0.4, density=True)

    ax6.axvline(x=not_churned['SupportCalls'].mean(),
                color=PALETTE['success'], linestyle='--', linewidth=1.5,
                label=f"Retained avg: {not_churned['SupportCalls'].mean():.1f}")
    ax6.axvline(x=churned['SupportCalls'].mean(),
                color=PALETTE['danger'],  linestyle='--', linewidth=1.5,
                label=f"Churned avg:  {churned['SupportCalls'].mean():.1f}")

    ax6.set_xlabel('Support Calls per Month', fontsize=10)
    ax6.set_ylabel('Density', fontsize=10)
    ax6.set_title('Support Calls Distribution', fontsize=13,
                  pad=10, color='#E8E8F0')
    ax6.legend(facecolor='#1A1D2E', edgecolor='#3A3D52',
               labelcolor='#E8E8F0', fontsize=8)

    plt.savefig('churn_analysis_dashboard.png', dpi=150,
                bbox_inches='tight', facecolor='#0F1117')
    plt.show()
    print("  Dashboard saved → churn_analysis_dashboard.png")


# ─────────────────────────────────────────────
# 5.  RETENTION RECOMMENDATIONS ENGINE
# ─────────────────────────────────────────────

def generate_recommendations(df: pd.DataFrame,
                              importance_df: pd.DataFrame) -> None:
    """
    Derives data-driven retention recommendations from the top features
    and prints a structured action plan with measurable KPIs.
    """

    top3 = importance_df.index[:3].tolist()

    # ── Statistical summaries ──────────────────
    churned     = df[df['Churn'] == 1]
    not_churned = df[df['Churn'] == 0]
    overall_churn_rate = df['Churn'].mean() * 100

    # Per-segment churn rates
    early_churn = df[df['Tenure'] <= 12]['Churn'].mean() * 100
    hi_charge_churn = df[df['MonthlyCharges'] > 80]['Churn'].mean() * 100
    hi_calls_churn  = df[df['SupportCalls'] > 8]['Churn'].mean() * 100
    mtm_churn  = df[df['ContractType'] == 'Month-to-Month']['Churn'].mean() * 100
    low_sat_churn = df[df['Satisfaction'] <= 2]['Churn'].mean() * 100

    banner = "═" * 60

    print(f"\n{banner}")
    print(f"  ACTIONABLE RETENTION RECOMMENDATIONS")
    print(f"  Based on top features: {', '.join(top3)}")
    print(f"{banner}")
    print(f"  Overall Churn Rate: {overall_churn_rate:.1f}%")
    print(f"{banner}\n")

    recs = [
        # ── REC 1 ────────────────────────────────
        {
            "title"   : "🚨 PRIORITY 1 — Early-Tenure Onboarding Program",
            "insight" : (
                f"Customers with tenure ≤ 12 months churn at {early_churn:.1f}% "
                f"vs {overall_churn_rate:.1f}% overall — "
                f"{early_churn/overall_churn_rate:.1f}× the baseline risk."
            ),
            "actions" : [
                "Assign a dedicated Success Manager for first 90 days",
                "Trigger automated check-in emails at Day 7, 30, 60, 90",
                "Offer a 'Loyalty Bonus' discount at the 6-month milestone",
                "Build an interactive onboarding wizard reducing time-to-value",
                "Create a peer community / user group for new customers",
            ],
            "kpis"    : [
                "Target: reduce <12-month churn by 25% within 2 quarters",
                "Track: 30/60/90-day activation rates",
                "Metric: Net Promoter Score at month 3",
            ],
            "effort"  : "Medium  |  Impact: HIGH",
        },
        # ── REC 2 ────────────────────────────────
        {
            "title"   : "💰 PRIORITY 2 — Value-Realignment for High-Charge Customers",
            "insight" : (
                f"Customers paying >$80/month churn at {hi_charge_churn:.1f}%. "
                "Price sensitivity compounds with low satisfaction — "
                "a toxic combination."
            ),
            "actions" : [
                "Proactively audit plans — migrate customers to better-fit tiers",
                "Introduce mid-tier 'Best Value' bundle to reduce sticker shock",
                "Send personalised 'Value Reports' showing ROI vs spend each quarter",
                "Offer loyalty discounts (5-15%) at 12-month anniversaries",
                "Create a transparent pricing page with competitor comparison",
            ],
            "kpis"    : [
                "Target: reduce high-charge segment churn by 20%",
                "Track: average revenue per churned account (ARPA)",
                "Metric: plan upgrade/downgrade ratio monthly",
            ],
            "effort"  : "Medium  |  Impact: HIGH",
        },
        # ── REC 3 ────────────────────────────────
        {
            "title"   : "📞 PRIORITY 3 — Support Experience Transformation",
            "insight" : (
                f"High-support customers (>8 calls/month) churn at "
                f"{hi_calls_churn:.1f}%. Each extra call signals "
                "recurring product friction and rising frustration."
            ),
            "actions" : [
                "Flag accounts with >5 calls/month for proactive outreach",
                "Implement root-cause tagging on every support ticket",
                "Build a self-service knowledge base targeting top 10 issue types",
                "Deploy AI chatbot for tier-1 deflection (target 40% deflection)",
                "Introduce 'Support Health Score' — alert CSMs when score drops",
                "Offer dedicated support queue for high-value at-risk accounts",
            ],
            "kpis"    : [
                "Target: reduce repeat-contact rate by 30% in 6 months",
                "Track: First-Contact Resolution (FCR) rate",
                "Metric: average calls per customer per month",
            ],
            "effort"  : "Low-Med  |  Impact: HIGH",
        },
        # ── REC 4 ────────────────────────────────
        {
            "title"   : "📋 PRIORITY 4 — Contract Conversion Campaign",
            "insight" : (
                f"Month-to-Month customers churn at {mtm_churn:.1f}%. "
                "Long-term contracts are the single strongest churn barrier."
            ),
            "actions" : [
                "Offer 10-20% discount for annual contract commitment",
                "Bundle premium features exclusively into longer contracts",
                "Run a 'Lock-In Loyalty' campaign targeting MTM customers >6 months",
                "Use in-app prompts with social proof ('87% of users on annual plan')",
            ],
            "kpis"    : [
                "Target: convert 15% of MTM base to annual contracts per quarter",
                "Track: Monthly Recurring Revenue (MRR) lock-in rate",
                "Metric: contract-length distribution shift",
            ],
            "effort"  : "Low  |  Impact: VERY HIGH",
        },
        # ── REC 5 ────────────────────────────────
        {
            "title"   : "⭐ PRIORITY 5 — Satisfaction-Driven Rescue Program",
            "insight" : (
                f"Customers with satisfaction score ≤ 2 churn at "
                f"{low_sat_churn:.1f}%. A Rescue Program can recover "
                "30-40% of at-risk accounts when triggered early."
            ),
            "actions" : [
                "Send CSAT / NPS surveys at key lifecycle events",
                "Auto-trigger 'Rescue Workflow' when CSAT drops below 3",
                "Empower frontline reps with retention offers (service credit, upgrade)",
                "Executive 'Save Call' for strategic accounts scoring ≤ 2",
                "Post-resolution follow-up at 30 days to confirm issue resolved",
            ],
            "kpis"    : [
                "Target: rescue 35% of flagged at-risk accounts",
                "Track: NPS trend month-over-month",
                "Metric: post-rescue 6-month retention rate",
            ],
            "effort"  : "Medium  |  Impact: HIGH",
        },
    ]

    # ── Print each recommendation ──────────────
    for rec in recs:
        print(f"  {rec['title']}")
        print(f"  {'─'*56}")
        print(f"  📊 Insight: {rec['insight']}")
        print(f"\n  ✅ Actions:")
        for a in rec['actions']:
            print(f"      • {a}")
        print(f"\n  📈 KPIs:")
        for k in rec['kpis']:
            print(f"      • {k}")
        print(f"\n  ⚙️  Effort / Impact: {rec['effort']}")
        print()

    # ── Churn cost model ───────────────────────
    avg_revenue = df['MonthlyCharges'].mean()
    avg_ltv     = avg_revenue * df['Tenure'].mean()
    n_churned   = int(df['Churn'].sum())

    print(f"{banner}")
    print(f"  BUSINESS IMPACT SNAPSHOT")
    print(f"{banner}")
    print(f"  Customers churned       : {n_churned:,}")
    print(f"  Avg monthly charge      : ${avg_revenue:.2f}")
    print(f"  Avg customer LTV        : ${avg_ltv:,.0f}")
    print(f"  Est. revenue at risk    : ${n_churned * avg_revenue:,.0f} / month")
    print(f"  5% retention improvement: ${n_churned * avg_revenue * 0.05:,.0f} saved/month")
    print(f"{banner}\n")


# ─────────────────────────────────────────────
# 6.  MAIN PIPELINE
# ─────────────────────────────────────────────

def main():
    print("\n🔄 Loading data...")
    df_raw = generate_churn_dataset()          # ← swap in your real DataFrame
    print(f"   Dataset shape : {df_raw.shape}")
    print(f"   Churn rate    : {df_raw['Churn'].mean()*100:.1f}%\n")

    print("🔄 Preprocessing...")
    X, y, features, encoders = preprocess(df_raw)

    print("🔄 Identifying top predictive features...")
    importance_df = identify_top_features(X, y, top_n=3)

    print("🔄 Generating visualisations...")

I'll analyze customer churn data and provide a complete solution. First, let me create a comprehensive Python script that includes data analysis, visualization, and actionable recommendations.

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.preprocessing import LabelEncoder
import warnings
warnings.filterwarnings('ignore')

# Set style for better visualizations
plt.style.use('seaborn-v0_8-darkgrid')
sns.set_palette("husl")

# ==================== 1. LOAD AND PREPARE SAMPLE DATA ====================
def create_sample_churn_data():
    """Create a realistic sample churn dataset if no data is provided"""
    np.random.seed(42)
    n_samples = 1000
    
    data = {
        'customer_id': [f'CUST{str(i).zfill(4)}' for i in range(1, n_samples + 1)],
        'tenure': np.random.randint(1, 72, n_samples),
        'monthly_charges': np.random.uniform(20, 120, n_samples).round(2),
        'total_charges': np.random.uniform(50, 5000, n_samples).round(2),
        'contract': np.random.choice(['Month-to-month', 'One year', 'Two year'], n_samples, p=[0.5, 0.3, 0.2]),
        'payment_method': np.random.choice(['Electronic check', 'Mailed check', 'Bank transfer', 'Credit card'], n_samples),
        'paperless_billing': np.random.choice(['Yes', 'No'], n_samples, p=[0.6, 0.4]),
        'tech_support': np.random.choice(['Yes', 'No', 'No internet service'], n_samples, p=[0.3, 0.6, 0.1]),
        'online_security': np.random.choice(['Yes', 'No', 'No internet service'], n_samples, p=[0.3, 0.6, 0.1]),
        'internet_service': np.random.choice(['Fiber optic', 'DSL', 'No'], n_samples, p=[0.5, 0.4, 0.1]),
        'senior_citizen': np.random.choice([0, 1], n_samples, p=[0.7, 0.3]),
        'dependents': np.random.choice(['Yes', 'No'], n_samples, p=[0.3, 0.7]),
        'partner': np.random.choice(['Yes', 'No'], n_samples, p=[0.5, 0.5]),
        'churn': np.random.choice(['Yes', 'No'], n_samples, p=[0.3, 0.7])
    }
    
    # Make correlations realistic
    df = pd.DataFrame(data)
    
    # Realistic correlations
    df.loc[df['contract'] == 'Month-to-month', 'churn'] = np.random.choice(['Yes', 'No'], 
                                                                          len(df[df['contract'] == 'Month-to-month']), 
                                                                          p=[0.6, 0.4])
    df.loc[df['contract'] == 'Two year', 'churn'] = np.random.choice(['Yes', 'No'], 
                                                                    len(df[df['contract'] == 'Two year']), 
                                                                    p=[0.1, 0.9])
    
    df.loc[df['tech_support'] == 'No', 'churn'] = np.random.choice(['Yes', 'No'], 
                                                                  len(df[df['tech_support'] == 'No']), 
                                                                  p=[0.5, 0.5])
    df.loc[df['tech_support'] == 'Yes', 'churn'] = np.random.choice(['Yes', 'No'], 
                                                                   len(df[df['tech_support'] == 'Yes']), 
                                                                   p=[0.1, 0.9])
    
    # Higher monthly charges correlate with higher churn
    high_charge_mask = df['monthly_charges'] > 80
    df.loc[high_charge_mask, 'churn'] = np.random.choice(['Yes', 'No'], 
                                                        len(df[high_charge_mask]), 
                                                        p=[0.6, 0.4])
    
    return df

# ==================== 2. IDENTIFY TOP PREDICTIVE FEATURES ====================
def analyze_churn_features(df):
    """Analyze dataset and identify top 3 predictive features"""
    
    # Create a copy for encoding
    df_encoded = df.copy()
    
    # Encode categorical variables
    le = LabelEncoder()
    categorical_cols = df.select_dtypes(include=['object']).columns
    categorical_cols = [col for col in categorical_cols if col != 'customer_id']
    
    for col in categorical_cols:
        df_encoded[col] = le.fit_transform(df_encoded[col])
    
    # Calculate correlation with churn
    churn_corr = df_encoded.corr()['churn'].abs().sort_values(ascending=False)
    
    # Drop churn itself from the list
    churn_corr = churn_corr.drop('churn')
    
    print("=" * 70)
    print("FEATURE IMPORTANCE ANALYSIS")
    print("=" * 70)
    
    print("\nTop 10 Features Correlated with Churn:")
    for i, (feature, corr) in enumerate(churn_corr.head(10).items(), 1):
        direction = df_encoded.corr()['churn'].loc[feature]
        direction_symbol = "(+)" if direction > 0 else "(-)"
        print(f"{i:2}. {feature:20} | Correlation: {corr:.4f} {direction_symbol}")
    
    top_features = churn_corr.head(3).index.tolist()
    
    print(f"\n🔍 TOP 3 PREDICTIVE FEATURES FOR CHURN:")
    for i, feature in enumerate(top_features, 1):
        corr_value = churn_corr[feature]
        direction = df_encoded.corr()['churn'].loc[feature]
        print(f"{i}. {feature} (Correlation: {corr_value:.3f})")
        
        # Provide interpretation
        if feature in ['tenure', 'total_charges']:
            if direction < 0:
                print(f"   → Lower {feature} values correlate with higher churn")
            else:
                print(f"   → Higher {feature} values correlate with higher churn")
        elif feature in df.columns:
            unique_vals = df[feature].unique()
            print(f"   → Feature values: {', '.join(map(str, unique_vals[:5]))}")
    
    return top_features, churn_corr

# ==================== 3. VISUALIZE CORRELATION MATRIX ====================
def visualize_correlations(df, top_features):
    """Create comprehensive correlation visualizations"""
    
    # Prepare data for correlation matrix
    df_vis = df.copy()
    
    # Encode categorical variables
    categorical_cols = df.select_dtypes(include=['object']).columns
    categorical_cols = [col for col in categorical_cols if col not in ['customer_id', 'churn']]
    
    le = LabelEncoder()
    for col in categorical_cols:
        df_vis[col] = le.fit_transform(df_vis[col])
    
    # Encode churn for correlation
    df_vis['churn_encoded'] = le.fit_transform(df['churn'])
    
    # Select top 10 features plus churn for visualization
    relevant_features = top_features + ['churn_encoded'] + ['tenure', 'monthly_charges']
    relevant_features = list(set(relevant_features))[:10]  # Get unique, limit to 10
    
    # Create correlation matrix
    corr_matrix = df_vis[relevant_features].corr()
    
    # Create visualization figure
    fig = plt.figure(figsize=(16, 12))
    
    # 1. Main Correlation Heatmap
    plt.subplot(2, 2, 1)
    mask = np.triu(np.ones_like(corr_matrix, dtype=bool))
    sns.heatmap(corr_matrix, mask=mask, annot=True, fmt='.2f', cmap='RdBu_r', 
                center=0, square=True, linewidths=0.5, cbar_kws={"shrink": 0.8})
    plt.title('Correlation Matrix (Top Features)', fontsize=14, fontweight='bold', pad=20)
    
    # 2. Churn Correlation Bar Plot
    plt.subplot(2, 2, 2)
    churn_corr = corr_matrix['churn_encoded'].drop('churn_encoded').sort_values()
    colors = ['red' if x > 0 else 'green' for x in churn_corr]
    churn_corr.plot(kind='barh', color=colors, edgecolor='black')
    plt.xlabel('Correlation with Churn', fontsize=12)
    plt.title('Feature Correlation with Churn', fontsize=14, fontweight='bold', pad=20)
    plt.axvline(x=0, color='black', linestyle='-', linewidth=0.5)
    
    # 3. Top Feature Distribution by Churn
    plt.subplot(2, 2, 3)
    top_feature = top_features[0]
    if df[top_feature].dtype == 'object':
        churn_rates = df.groupby(top_feature)['churn'].value_counts(normalize=True).unstack()
        churn_rates[['No', 'Yes']].plot(kind='bar', stacked=True, ax=plt.gca())
        plt.title(f'Churn Rate by {top_feature}', fontsize=14, fontweight='bold', pad=20)
        plt.ylabel('Proportion')
    else:
        for churn_status in ['Yes', 'No']:
            subset = df[df['churn'] == churn_status]
            sns.kdeplot(data=subset[top_feature], label=churn_status, fill=True)
        plt.title(f'Distribution of {top_feature} by Churn Status', fontsize=14, fontweight='bold', pad=20)
        plt.xlabel(top_feature)
        plt.ylabel('Density')
        plt.legend(title='Churned')
    
    # 4. Pair Plot of Top 3 Features
    plt.subplot(2, 2, 4)
    if len(top_features) >= 2:
        scatter_feature = top_features[1] if len(top_features) > 1 else 'tenure'
        if df[scatter_feature].dtype != 'object' and df[top_features[0]].dtype != 'object':
            scatter = plt.scatter(df[top_features[0]], df[scatter_feature], 
                                 c=df['churn'].map({'Yes': 1, 'No': 0}), 
                                 cmap='coolwarm', alpha=0.6, edgecolors='black')
            plt.xlabel(top_features[0])
            plt.ylabel(scatter_feature)
            plt.title(f'{top_features[0]} vs {scatter_feature}', fontsize=14, fontweight='bold', pad=20)
            plt.colorbar(scatter, label='Churn (1=Yes, 0=No)')
        else:
            # Alternative: Box plot if categorical
            if df[top_features[0]].dtype == 'object':
                sns.boxplot(x=top_features[0], y='tenure', hue='churn', data=df)
                plt.title(f'Tenure by {top_features[0]} and Churn', fontsize=14, fontweight='bold', pad=20)
    
    plt.tight_layout()
    plt.savefig('churn_analysis.png', dpi=300, bbox_inches='tight')
    plt.show()
    
    return corr_matrix

# ==================== 4. GENERATE ACTIONABLE RECOMMENDATIONS ====================
def generate_recommendations(df, top_features):
    """Generate actionable retention recommendations based on analysis"""
    
    print("\n" + "=" * 70)
    print("ACTIONABLE RETENTION RECOMMENDATIONS")
    print("=" * 70)
    
    recommendations = []
    
    # Analyze each top feature
    for feature in top_features[:3]:
        print(f"\n🔧 **Recommendations based on '{feature}':**")
        
        if feature == 'contract':
            churn_by_contract = df.groupby('contract')['churn'].value_counts(normalize=True).unstack()
            month_to_month_rate = churn_by_contract.loc['Month-to-month', 'Yes'] if 'Month-to-month' in churn_by_contract.index else 0
            
            recommendations.append({
                'feature': 'contract',
                'insight': f"Month-to-month contracts have {month_to_month_rate:.1%} churn rate",
                'actions': [
                    "Create incentives for annual contracts (10-15% discount)",
                    "Implement loyalty program for long-term customers",
                    "Add early termination fee for month-to-month contracts"
                ]
            })
            print(f"  • Month-to-month churn rate: {month_to_month_rate:.1%}")
            print("  → Offer annual contract discounts")
            print("  → Create loyalty program tiers")
            
        elif feature == 'tech_support':
            churn_with_support = df[df['tech_support'] == 'Yes']['churn'].value_counts(normalize=True).get('Yes', 0)
            churn_without_support = df[df['tech_support'] == 'No']['churn'].value_counts(normalize=True).get('Yes', 0)
            
            recommendations.append({
                'feature': 'tech_support',
                'insight': f"Customers without tech support churn at {churn_without_support:.1%} vs {churn_with_support:.1%} with support",
                'actions': [
                    "Bundle tech support with premium plans",
                    "Offer free tech support trial for 3 months",
                    "Implement proactive tech support outreach"
                ]
            })
            print(f"  • Churn without support: {churn_without_support:.1%}")
            print(f"  • Churn with support: {churn_with_support:.1%}")
            print("  → Bundle tech support automatically")
            print("  → Offer free trial of premium support")
            
        elif feature == 'monthly_charges':
            high_churn_threshold = df['monthly_charges'].quantile(0.75)
            high_charge_churn = df[df['monthly_charges'] > high_churn_threshold]['churn'].value_counts(normalize=True).get('Yes', 0)
            low_charge_churn = df[df['monthly_charges'] <= high_churn_threshold]['churn'].value_counts(normalize=True).get('Yes', 0)
            
            recommendations.append({
                'feature': 'monthly_charges',
                'insight': f"High-paying customers (>${high_churn_threshold:.0f}/month) churn at {high_charge_churn:.1%} rate",
                'actions': [
                    "Create value-added services for high-paying customers",
                    "Implement price-lock guarantees",
                    "Offer personalized discounts at renewal"
                ]
            })
            print(f"  • Churn for high charges (>${high_churn_threshold:.0f}): {high_charge_churn:.1%}")
            print("  → Add value services to justify cost")
            print("  → Offer price-lock guarantees")
            
        elif feature == 'tenure':
            low_tenure_threshold = 12  # Months
            low_tenure_churn = df[df['tenure'] < low_tenure_threshold]['churn'].value_counts(normalize=True).get('Yes', 0)
            
            recommendations.append({
                'feature': 'tenure',
                'insight': f"New customers (<{low_tenure_threshold} months) churn at {low_tenure_churn:.1%} rate",
                'actions': [
                    "Implement onboarding program for first 90 days",
                    "Assign dedicated account managers to new customers",
                    "Create 'welcome' discount that converts to loyalty program"
                ]
            })
            print(f"  • Churn in first {low_tenure_threshold} months: {low_tenure_churn:.1%}")
            print("  → Enhance onboarding process")
            print("  → Assign dedicated support for new customers")
    
    # Overall recommendations
    print("\n" + "=" * 70)
    print("🎯 STRATEGIC RETENTION INITIATIVES")
    print("=" * 70)
    
    strategic_initiatives = [
        {
            'name': "Proactive Retention Campaign",
            'description': "Identify at-risk customers using predictive model",
            'kpi': "Reduce churn by 25% in 6 months",
            'actions': [
                "Flag customers with month-to-month contracts",
                "Target customers without tech support",
                "Focus on high monthly charges with low tenure"
            ]
        },
        {
            'name': "Customer Success Program",
            'description': "Improve customer experience and perceived value",
            'kpi': "Increase NPS by 15 points",
            'actions': [
                "Regular check-ins at 30, 90, 180 days",
                "Personalized usage tips and optimizations",
                "Exclusive content and early access features"
            ]
        },
        {
            'name': "Win-back Campaign",
            'description': "Recapture customers who recently churned",
            'kpi': "Re-acquire 15% of churned customers",
            'actions': [
                "Special offer for return within 60 days",
                "Survey to understand churn reasons",
                "Address pain points before re-engagement"
            ]
        }
    ]
    
    for initiative in strategic_initiatives:
        print(f"\n📋 {initiative['name']}")
        print(f"   {initiative['description']}")
        print(f"   Target KPI: {initiative['kpi']}")
        print("   Key Actions:")
        for action in initiative['actions']:
            print(f"   • {action}")
    
    # Create implementation timeline
    print("\n" + "=" * 70)
    print("⏱️ RECOMMENDED IMPLEMENTATION TIMELINE")
    print("=" * 70)
    
    timeline = [
        ("Week 1-2", "Data analysis and segmentation", "Analytics Team"),
        ("Week 3-4", "Design retention offers and messaging", "Marketing + Product"),
        ("Month 2", "Pilot program with 500 at-risk customers", "Customer Success"),
        ("Month 3", "Full rollout and monitoring", "All Teams"),
        ("Month 4-6", "Optimize based on results", "Data Science + Business")
    ]
    
    for timeframe, task, team in timeline:
        print(f"{timeframe:10} | {task:40} | {team}")
    
    return recommendations

# ==================== 5. MAIN EXECUTION ====================
def main():
    """Main execution function"""
    
    print("CUSTOMER CHURN ANALYSIS & RETENTION STRATEGY")
    print("=" * 70)
    
    # Load or create data
    # In practice, you would load your actual data:
    # df = pd.read_csv('customer_churn.csv')
    
    # For demonstration, create sample data
    df = create_sample_churn_data()
    
    print(f"\nDataset Overview:")
    print(f"• Total customers: {len(df):,}")
    print(f"• Churn rate: {df['churn'].value_counts(normalize=True)['Yes']:.1%}")
    print(f"• Features available: {len(df.columns)}")
    print(f"• Data shape: {df.shape}")
    
    # Identify top features
    top_features, feature_correlations = analyze_churn_features(df)
    
    # Visualize correlations
    print("\n📊 Generating visualizations...")
    corr_matrix = visualize_correlations(df, top_features)
    
    # Generate recommendations
    recommendations = generate_recommendations(df, top_features)
    
    # Save analysis summary
    summary = {
        'total_customers': len(df),
        'churn_rate': df['churn'].value_counts(normalize=True)['Yes'],
        'top_features': top_features,
        'feature_correlations': feature_correlations.to_dict(),
        'recommendations': recommendations
    }
    
    print("\n" + "=" * 70)
    print("✅ ANALYSIS COMPLETE")
    print("=" * 70)
    print(f"\nResults saved to:")
    print("• churn_analysis.png (visualizations)")
    print("• Feature correlations and recommendations available in memory")
    print("\nNext steps:")
    print("1. Validate findings with business stakeholders")
    print("2. Implement pilot retention program")
    print("3. Monitor results and iterate strategy")

# ==================== 6. RUN ANALYSIS ====================
if __name__ == "__main__":
    main()

Key Features of This Analysis:

Top 3 Predictive Features Identified:

Based on correlation analysis, the script identifies:

1. Contract type - Month-to-month contracts have highest churn
2. Tech support - Lack of tech support increases churn risk
3. Monthly charges - High charges with low value perception drive churn

Four Visualizations Generated:

1. Correlation Matrix Heatmap - Shows relationships between features
2. Feature-Churn Correlation Bar Plot - Ranks feature importance
3. Top Feature Distribution - Shows how top feature relates to churn
4. Feature Relationship Scatter Plot - Visualizes interaction between top features

Actionable Recommendations Include:

Immediate Tactical Actions:

• Convert month-to-month to annual contracts with discounts
• Bundle tech support with premium plans
• Implement price-lock guarantees for high-paying customers

Strategic Initiatives:

1. Proactive Retention Campaign - Target at-risk customers
2. Customer Success Program - Improve customer experience
3. Win-back Campaign - Re-acquire recently churned customers

Implementation Timeline:

• Week 1-2: Data analysis and segmentation
• Week 3-4: Design retention offers
• Month 2: Pilot program
• Month 3: Full rollout
• Month 4-6: Optimization

To Use with Your Actual Data:

1. Replace create_sample_churn_data() with your data loading code:

df = pd.read_csv('your_churn_data.csv')

2. Ensure your dataset has a 'churn' column (binary: 'Yes'/'No' or 1/0)

3. The script automatically handles both numerical and categorical features

This analysis provides both technical insights (feature importance, correlations) and business-ready recommendations for reducing customer churn.