DESCRIPTIVE DIAGNOSTIC PREDICTIVE PRESCRIPTIVE
"What happened?" "Why did it happen?" "What will happen?" "What should we do?"
 Sales dropped Ad campaign failed Next quarter sales Increase marketing
 last quarter. in Target segment. will drop 12%. in Tier-2 cities.

 ◄─────────────── Looking backward ────────────────
 ──────────────── Looking forward ───────────────►

 Structured Data Semi-Structured Unstructured
 CRM databases JSON API responses Customer reviews
 Transaction logs IoT sensor telemetry Social media posts
 Sales spreadsheets XML feeds Call center recordings
 ERP system exports Web server logs Email text

import pandas as pd
import numpy as np

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

# ── Temporal features ───────────────────────────────────────────
df['day_of_week'] = pd.to_datetime(df['date']).dt.dayofweek # 0=Mon, 6=Sun
df['is_weekend'] = df['day_of_week'].isin([5, 6]).astype(int) # Binary flag
df['month'] = pd.to_datetime(df['date']).dt.month
df['quarter'] = pd.to_datetime(df['date']).dt.quarter

# ── Lag features (what happened last week/month?) ───────────────
df['sales_lag_7'] = df['units_sold'].shift(7) # Sales from 7 days ago
df['sales_lag_30'] = df['units_sold'].shift(30) # Sales from 30 days ago

# ── Rolling window aggregates ───────────────────────────────────
df['rolling_7d_avg'] = df['units_sold'].rolling(7).mean()
df['rolling_30d_std'] = df['units_sold'].rolling(30).std() # Volatility!

# ── Polynomial terms (to capture curvature) ─────────────────────
df['price_squared'] = df['price'] ** 2
df['log_price'] = np.log(df['price'] + 1) # +1 to avoid log(0)!

# ── One-hot encoding for categorical variables ──────────────────
df = pd.get_dummies(df, columns=['product_category', 'region'], drop_first=True)

Instead of one fixed train/val split, rotate through K splits:

Fold 1: [VAL][TRAIN][TRAIN][TRAIN][TRAIN]
Fold 2: [TRAIN][VAL][TRAIN][TRAIN][TRAIN]
Fold 3: [TRAIN][TRAIN][VAL][TRAIN][TRAIN]
Fold 4: [TRAIN][TRAIN][TRAIN][VAL][TRAIN]
Fold 5: [TRAIN][TRAIN][TRAIN][TRAIN][VAL]

Average performance across 5 folds → more reliable estimate of real-world accuracy!

# Simple data drift detection
from scipy.stats import ks_2samp

# training data distribution
train_prices = training_df['price'].values
# current production data distribution 
live_prices = live_df['price'].values

stat, p_value = ks_2samp(train_prices, live_prices)
if p_value < 0.05:
 print(" DATA DRIFT DETECTED! Input distribution has changed significantly.")
 print("Consider retraining the model!")

PROBLEM: A national retailer wants to reduce:
 - Out-of-stock events (lose sales when popular items run out)
 - Overstock (waste money storing things nobody buys)

SOLUTION BUILT:
 Data: 5 years of sales + weather data + advertising budgets + holidays
 Model: Gradient Boosted Trees (LightGBM)
 Training frequency: Weekly retraining with latest data

RESULTS:
 ↓ 30% reduction in out-of-stock events during holiday season
 ↓ 12% reduction in inventory holding costs
 = Millions saved + happier customers!

DEPLOYMENT: Predictions feed directly into automated restocking system.
When model says "Product X will sell 500 units next week in Store Y",
the procurement system automatically places the order!

Why not minimize sum of |errors|?
 Because squaring has beautiful math properties!
 - Squaring makes all errors positive (no positive and negative canceling)
 - Squaring penalizes BIG errors more than small ones (outlier sensitivity)
 - Squared function is differentiable → can find minimum with calculus
 - Results in the BLUE (Best Linear Unbiased Estimator) under Gauss-Markov

R² = 0.0 → model explains NOTHING (worse than just predicting the mean!)
R² = 0.5 → model explains 50% of variation in y
R² = 0.95 → model explains 95% of variation in y 
R² = 1.0 → perfect fit (almost certainly overfitting or a trivial relationship!)

t-statistic interpretation:
 t = 0 → slope estimate = 0, no evidence of relationship
 t = 2.0 → slope is 2 standard errors from zero → likely significant (p ≈ 0.05)
 t = 4.0 → very strong evidence of real relationship (p < 0.001)
 t = -3.5 → strong negative relationship evidence

p-value:
 p < 0.05 → Significant at 5% level (conventional threshold)
 p < 0.01 → Significant at 1% level (stronger evidence)
 p < 0.001 → Very strong evidence

import pandas as pd
import numpy as np
import statsmodels.api as sm
import seaborn as sns
import matplotlib.pyplot as plt

# ── Load data ────────────────────────────────────────────────────
data = pd.read_csv('car_data.csv')
# x = engine_size (litres), y = fuel_consumption (L/100km)

# ── Fit the model ─────────────────────────────────────────────────
X = sm.add_constant(data['engine_size']) # Add β₀ intercept column
model = sm.OLS(data['fuel_consumption'], X).fit()

# ── Detailed results ─────────────────────────────────────────────
print(model.summary())
# Output: β₀=2.3, β₁=0.07, R²=0.68, p-value < 0.001

# ── Scatter plot with fitted line ────────────────────────────────
fig, axes = plt.subplots(1, 3, figsize=(15, 5))

# Plot 1: Data + regression line
axes[0].scatter(data['engine_size'], data['fuel_consumption'], alpha=0.5)
axes[0].plot(data['engine_size'], model.fittedvalues, 'r-', linewidth=2)
axes[0].set_xlabel('Engine Size (L)')
axes[0].set_ylabel('Fuel Consumption (L/100km)')
axes[0].set_title('Simple Linear Regression')
# Model: ŷ = 2.3 + 0.07x

# Plot 2: Residuals vs Fitted (check homoscedasticity)
axes[1].scatter(model.fittedvalues, model.resid, alpha=0.5)
axes[1].axhline(0, color='black', linestyle='--')
axes[1].set_xlabel('Fitted Values')
axes[1].set_ylabel('Residuals')
axes[1].set_title('Residuals vs Fitted\n(Should look like a random cloud!)')

# Plot 3: Q-Q plot (check normality of residuals)
sm.qqplot(model.resid, line='45', ax=axes[2])
axes[2].set_title('Q-Q Plot\n(Points should follow the 45° line)')

plt.tight_layout()
plt.savefig('regression_diagnostics.png', dpi=150)
plt.show()

β₀ = 2.3: When engine_size = 0L (theoretical baseline),
 expected fuel consumption is 2.3 L/100km.
 (Not physically meaningful here, just math!)

β₁ = 0.07: Each additional 1-litre increase in engine size
 is associated with 0.07 more L/100km consumption.
 A 3.0L engine vs 1.0L engine: 0.07 × 2 = 0.14 L/100km more.

R² = 0.68: Engine size explains 68% of the variability in fuel consumption.
 Remaining 32% is due to: driver behavior, transmission type,
 aerodynamics, load weight, tire pressure, etc.

Interpretation (each coefficient "holds all others fixed"):

β₁ = +120: Each extra m² adds ₹1,20,000 to price
 (holding bedrooms, age, distance constant)

β₂ = +10,000: Each extra bedroom adds ₹10,00,000
 (holding size, age, distance constant)

β₃ = -8,000: Each extra year of age REDUCES price by ₹8,00,000
 (building wear & tear effect)

β₄ = -2,500: Each extra km from city center REDUCES price by ₹25,000
 (location premium)

Intercept = 50,000: A 0 m², 0 bedroom, brand new house at distance 0...
 Not meaningful! Just the mathematical anchor.

Example:
 Perfect scenario: size and bedrooms are independent
 Problem scenario: larger houses ALWAYS have more bedrooms

 Consequence: model might say:
 β₁ (size) = +200 and β₂ (bedrooms) = -50
 (Rooms look NEGATIVE even though more rooms = higher price!)
 The model is confused because the two predictors are intertwined.

from statsmodels.stats.outliers_influence import variance_inflation_factor

X = housing[['size', 'bedrooms', 'age', 'dist_center']]
for i, col in enumerate(X.columns):
 vif = variance_inflation_factor(X.values, i)
 print(f"VIF for {col}: {vif:.2f}")

# Output:
# VIF for size: 1.45 → OK! No multicollinearity issue
# VIF for bedrooms: 1.82 → OK!
# VIF for age: 1.21 → OK!
# VIF for dist_center: 9.7 → CONCERN! dist_center correlates with other vars

from sklearn.linear_model import Ridge, Lasso, ElasticNet
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import cross_val_score
import numpy as np

X = housing[['size', 'bedrooms', 'age', 'dist_center']].values
y = housing['price'].values

# Always scale before regularization!
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

# Try Ridge with different λ values
for alpha in [0.01, 0.1, 1.0, 10.0, 100.0]:
 model = Ridge(alpha=alpha)
 scores = cross_val_score(model, X_scaled, y, cv=5, scoring='r2')
 print(f"Ridge α={alpha:6.2f} | CV R²: {scores.mean():.3f} ± {scores.std():.3f}")

# Lasso for feature selection
lasso = Lasso(alpha=0.1)
lasso.fit(X_scaled, y)
feature_names = ['size', 'bedrooms', 'age', 'dist_center']
for name, coef in zip(feature_names, lasso.coef_):
 print(f"{name:15s}: {coef:.2f} {'← ELIMINATED' if coef == 0 else ''}")

β₁₂ = +500

Interpretation:
 In non-luxury areas: effect of size = β₁ + β₁₂ × 0 = β₁ alone
 In luxury areas: effect of size = β₁ + β₁₂ × 1 = β₁ + 500

 → Every extra m² is worth ₹500 more in luxury areas than non-luxury!
 Size matters MORE in luxury neighborhoods.

import statsmodels.formula.api as smf

# Interaction via formula syntax (much cleaner than manual calculation!)
model = smf.ols(
 'price ~ size + bedrooms + age + dist_center'
 ' + size:is_luxury_area' # interaction term
 ' + np.power(size, 2)', # quadratic term (diminishing returns for size)
 data=housing
).fit()

print(model.summary())

Raw coefficients:
 β(size) = 120 (₹ per m²)
 β(bedrooms) = 10,000 (₹ per bedroom)
 β(age) = -8,000 (₹ per year)

Who matters more? You can't compare 120 and 10,000 directly
because they use completely different units!

Solution: STANDARDIZE the predictors first!
 z = (x - mean) / std

Standardized coefficients:
 β*(size) = 0.65 → 1-std increase in size → 0.65 std increase in price
 β*(bedrooms) = 0.42 → 1-std increase in bedrooms → 0.42 std increase in price
 β*(age) = -0.38 → 1-std increase in age → 0.38 std DECREASE in price

Now comparison is meaningful!
SIZE is the strongest driver of price (0.65 > 0.42 > 0.38)

from sklearn.preprocessing import StandardScaler
import pandas as pd

scaler = StandardScaler()
X_scaled = pd.DataFrame(
 scaler.fit_transform(X),
 columns=X.columns
)
model_std = sm.OLS(y, sm.add_constant(X_scaled)).fit()
print("\nStandardized (Beta) Coefficients:")
print(model_std.params.drop('const').sort_values(key=abs, ascending=False))

Size: t=24.0, p<0.001 → Extremely significant! CI=[110,130] excludes 0 
Distance: t=-0.83, p=0.406 → NOT significant! CI=[-8500,3500] INCLUDES 0 
 → Distance to city center doesn't improve this model
 → Consider removing it (AIC/BIC will tell you for sure)

Original: city ∈ {Delhi, Mumbai, Bangalore}
Create 2 dummies (NOT 3 — always k-1 to avoid perfect multicollinearity!):

 D1_Mumbai D2_Bangalore
Delhi: 0 0 ← reference category
Mumbai: 1 0
Bangalore: 0 1

Model: price = β₀ + β₁(size) + β₂(D1_Mumbai) + β₃(D2_Bangalore) + ...

Interpretation:
 β₂ = +15,00,000 → Mumbai properties cost ₹15L MORE than Delhi, same size/age/etc.
 β₃ = +8,00,000 → Bangalore properties cost ₹8L MORE than Delhi, same features.

 Want Mumbai vs Bangalore comparison? β₂ - β₃ = ₹7,00,000 difference.

Famous spurious correlations:
 Ice cream sales → Drowning deaths (r = 0.85)
 Nicolas Cage films → Pool drownings (r = 0.92!)
 iPhone units sold → Suicide rates (r = 0.88)

Cause? NONE! All driven by confounders or randomness.
Ice cream and drowning are both caused by HOT WEATHER.

Regression gives you ASSOCIATION.
Causation requires:
 1. Temporal precedence (X happened before Y)
 2. Mechanism (how does X cause Y?)
 3. Ruling out confounders (RCT or natural experiment)

Dataset 1: Linear relationship (regression is appropriate!)
Dataset 2: Curved relationship (need polynomial term!)
Dataset 3: Perfect linear with one massive outlier (outlier drives the line!)
Dataset 4: All x-values the same except one (terrible design!)

LESSON: Always visualize your data before and after fitting!

import pandas as pd
import numpy as np
import statsmodels.api as sm
import seaborn as sns
import matplotlib.pyplot as plt
from statsmodels.graphics.gofplots import ProbPlot

# Fit model
data = sns.load_dataset('mpg').dropna()
X = sm.add_constant(data[['horsepower', 'weight']])
model = sm.OLS(data['mpg'], X).fit()

fitted = model.fittedvalues
residuals = model.resid
std_resid = (residuals - residuals.mean()) / residuals.std()

fig, axes = plt.subplots(2, 3, figsize=(18, 12))
fig.suptitle('Complete Regression Diagnostic Dashboard', fontsize=16, y=1.02)

# ── Plot 1: Scatter with fit line (for simple case) ──────────────
axes[0,0].scatter(data['horsepower'], data['mpg'], alpha=0.5, color='steelblue')
axes[0,0].set_xlabel('Horsepower')
axes[0,0].set_ylabel('MPG')
axes[0,0].set_title('Data + Relationship')

# ── Plot 2: Residuals vs Fitted ──────────────────────────────────
axes[0,1].scatter(fitted, residuals, alpha=0.5, color='orange')
axes[0,1].axhline(0, color='black', linestyle='--', linewidth=1)
axes[0,1].set_xlabel('Fitted Values')
axes[0,1].set_ylabel('Residuals')
axes[0,1].set_title('Residuals vs Fitted\n(Want: no pattern!)')
# Add a smoother to see trends
import statsmodels.nonparametric.smoothers_lowess as lowess
smooth = lowess.lowess(residuals, fitted, frac=0.5)
axes[0,1].plot(smooth[:,0], smooth[:,1], 'r-', linewidth=2)

# ── Plot 3: Q-Q Plot ─────────────────────────────────────────────
pp = ProbPlot(residuals)
pp.qqplot(line='45', ax=axes[0,2], alpha=0.5)
axes[0,2].set_title('Q-Q Plot of Residuals\n(Points should follow 45° line)')

# ── Plot 4: Scale-Location (check homoscedasticity) ─────────────
axes[1,0].scatter(fitted, np.sqrt(np.abs(std_resid)), alpha=0.5, color='green')
axes[1,0].set_xlabel('Fitted Values')
axes[1,0].set_ylabel('√|Standardized Residuals|')
axes[1,0].set_title('Scale-Location\n(Horizontal band = OK)')

# ── Plot 5: Cook's Distance ──────────────────────────────────────
influence = model.get_influence()
(cook_d, _) = influence.cooks_distance
axes[1,1].stem(cook_d, markerfmt=',', linefmt='grey', basefmt='black')
n = len(cook_d)
threshold = 4 / n # Common rule of thumb: 4/n
axes[1,1].axhline(threshold, color='red', linestyle='--', label=f'Threshold (4/n = {threshold:.3f})')
axes[1,1].set_xlabel('Observation Index')
axes[1,1].set_ylabel("Cook's Distance")
axes[1,1].set_title("Cook's Distance\n(Points above red line = influential!)")
axes[1,1].legend()

# ── Plot 6: Coefficient Forest Plot ─────────────────────────────
params = model.params.drop('const')
conf = model.conf_int().drop('const')
axes[1,2].errorbar(
 params.values,
 range(len(params)),
 xerr=[params.values - conf[0].values,
 conf[1].values - params.values],
 fmt='o', color='purple', capsize=5, capthick=2
)
axes[1,2].axvline(0, color='black', linestyle='--', alpha=0.5)
axes[1,2].set_yticks(range(len(params)))
axes[1,2].set_yticklabels(params.index)
axes[1,2].set_xlabel('Coefficient Value')
axes[1,2].set_title('Coefficient Plot\n(CI crossing 0 = not significant)')

plt.tight_layout()
plt.savefig('full_diagnostic_dashboard.png', dpi=150, bbox_inches='tight')
plt.show()

import seaborn as sns
import matplotlib.pyplot as plt

fig, axes = plt.subplots(1, 2, figsize=(12, 4))

# Histogram with KDE overlay
sns.histplot(data['price'], bins=30, kde=True, ax=axes[0], color='steelblue')
axes[0].set_title('Price Distribution')
axes[0].set_xlabel('Price (₹)')

# Log-transformed (often better for right-skewed financial data)
sns.histplot(np.log1p(data['price']), bins=30, kde=True, ax=axes[1], color='orange')
axes[1].set_title('Log(Price) Distribution\n(More symmetric!)')
axes[1].set_xlabel('log(1 + Price)')
plt.tight_layout()

Normal (bell curve) → Good for linear regression directly
Right-skewed (long tail →) → Consider log transformation
Bimodal (two humps) → Two distinct populations! Separate them!
Uniform (flat) → Unusual, check data quality

# Compare price across product categories
sns.boxplot(x='product_category', y='price', data=df,
 order=df.groupby('product_category')['price'].median().sort_values().index)
plt.xticks(rotation=45)
plt.title('Price Distribution by Category')

# Violin plot: more information than box plot!
sns.violinplot(x='region', y='units_sold', data=df, inner='box')
plt.title('Sales Distribution by Region\n(Wider = more data at that value)')

# Check ALL pairwise relationships in one shot!
numeric_cols = ['size', 'bedrooms', 'age', 'dist_center', 'price']
pair_grid = sns.pairplot(
 housing[numeric_cols],
 diag_kind='kde', # Show distribution on diagonal
 plot_kws={'alpha': 0.5}, # Transparency for overlapping points
 corner=True # Show only lower triangle
)
pair_grid.fig.suptitle('Pairwise Relationships — Housing Data', y=1.02)

corr_matrix = housing[numeric_cols].corr()

plt.figure(figsize=(8, 6))
mask = np.triu(np.ones_like(corr_matrix, dtype=bool)) # Hide upper triangle
sns.heatmap(
 corr_matrix,
 mask=mask,
 annot=True, # Show correlation values in cells
 fmt='.2f', # 2 decimal places
 cmap='RdYlGn', # Red=negative, Green=positive
 center=0, # Center colorscale at 0
 vmin=-1, vmax=1,
 square=True,
 linewidths=0.5
)
plt.title('Correlation Matrix\n(Red=negative, Green=positive correlation)')

(size, bedrooms) = 0.82 → High correlation! Multicollinearity risk!
(age, price) = -0.65 → Older buildings tend to be cheaper
(dist_center, size) = -0.15 → Weak relationship, no concern

# Monthly sales with trend and seasonality
plt.figure(figsize=(14, 6))

# Raw values
plt.plot(df['date'], df['monthly_sales'], alpha=0.4, label='Monthly Sales', color='steelblue')

# 12-month rolling average (smooth out seasonality)
rolling_avg = df.set_index('date')['monthly_sales'].rolling('365D').mean()
plt.plot(rolling_avg, label='12-Month Rolling Avg', color='red', linewidth=2)

# Mark significant events
for event_date, event_name in events.items():
 plt.axvline(pd.to_datetime(event_date), color='green', linestyle='--', alpha=0.7)
 plt.text(pd.to_datetime(event_date), plt.ylim()[1]*0.95, event_name,
 rotation=45, fontsize=8)

plt.legend()
plt.title('Monthly Sales — Trend + Seasonality Analysis')
plt.xlabel('Date')
plt.ylabel('Units Sold')

BEFORE modeling, ask yourself:

 Have I checked all variable distributions? (histograms)
 Are there obvious outliers? (box plots, scatter plots)
 Are there missing values? (df.isnull().sum())
 Are predictors correlated with each other? (correlation heatmap)
 Are predictors correlated with the target? (pair plot)
 Is the target variable normally distributed or skewed?
 For time series: is there seasonality or trend?
 Are there categorical variables? How many unique values?
 Are there data quality issues (impossible values, duplicates)?
 Have I documented my findings to share with the team?

Interaction Type What It Does Tool
─────────────────────────────────────────────────────────────────
ZOOM & PAN Inspect dense data at detail Plotly, D3.js, Bokeh
BRUSHING Select subset in one chart → Vega-Altair, D3.js
 highlights appear in linked charts
FILTERING Dropdowns/sliders to show Dash, Shiny, Tableau
 subset of data
TOOLTIPS Hover for exact values + metadata Every major library
PARAMETER SLIDERS Change model params → live update Dash, Shiny, ipywidgets
RANGE SELECTOR Zoom into a time window Plotly, Bokeh
DRILL-DOWN Click a bar → see component details Power BI, Tableau
CROSS-FILTER Select in one chart → filters all Tableau, Dash

from dash import Dash, dcc, html, Input, Output, callback
import plotly.express as px
import pandas as pd

app = Dash(__name__)
df = px.data.gapminder() # World development dataset

app.layout = html.Div([
 html.H1(" Global Development Explorer", style={'textAlign': 'center'}),

 # ── Controls row ─────────────────────────────────────────────
 html.Div([
 html.Label("Select Year:"),
 dcc.Slider(
 id='year-slider',
 min=df['year'].min(), max=df['year'].max(),
 step=5,
 value=2007,
 marks={str(y): str(y) for y in df['year'].unique()}
 ),
 html.Label("Color by:"),
 dcc.Dropdown(
 id='color-dropdown',
 options=[{'label': c, 'value': c} for c in ['continent', 'gdpPercap']],
 value='continent',
 clearable=False
 ),
 ], style={'padding': '20px'}),

 # ── Charts row ────────────────────────────────────────────────
 html.Div([
 dcc.Graph(id='bubble-chart', style={'width': '50%', 'display': 'inline-block'}),
 dcc.Graph(id='bar-chart', style={'width': '50%', 'display': 'inline-block'}),
 ]),

 # ── The scatter shows what was clicked ────────────────────────
 dcc.Graph(id='time-series'),
])

@callback(
 Output('bubble-chart', 'figure'),
 Output('bar-chart', 'figure'),
 Input('year-slider', 'value'),
 Input('color-dropdown', 'value')
)
def update_charts(selected_year, color_col):
 filtered = df[df['year'] == selected_year]

 bubble = px.scatter(
 filtered,
 x='gdpPercap', y='lifeExp',
 size='pop', color=color_col,
 hover_name='country',
 log_x=True,
 size_max=60,
 title=f'GDP vs Life Expectancy ({selected_year})'
 )

 bar = px.bar(
 filtered.groupby('continent')['pop'].sum().reset_index(),
 x='continent', y='pop',
 color='continent',
 title=f'Population by Continent ({selected_year})'
 )

 return bubble, bar

@callback(
 Output('time-series', 'figure'),
 Input('bubble-chart', 'clickData')
)
def show_country_history(click_data):
 if click_data is None:
 country = 'India'
 else:
 country = click_data['points'][0]['hovertext']

 country_data = df[df['country'] == country]
 fig = px.line(country_data, x='year', y='lifeExp',
 title=f'Life Expectancy Over Time: {country}',
 markers=True)
 return fig

if __name__ == '__main__':
 app.run(debug=True, port=8050)

import altair as alt
import pandas as pd

cars = pd.read_csv('cars.csv')

# Brush selection (interactive) — selecting in one chart highlights in both!
brush = alt.selection_interval()

# Chart 1: Horsepower vs MPG
scatter1 = alt.Chart(cars).mark_circle(size=60).encode(
 x='Horsepower:Q',
 y='Miles_per_Gallon:Q',
 color=alt.condition(brush, 'Origin:N', alt.value('lightgray')),
 tooltip=['Name:N', 'Origin:N', 'Horsepower:Q', 'Miles_per_Gallon:Q']
).add_params(brush).properties(
 title='Horsepower vs MPG',
 width=300, height=300
)

# Chart 2: Weight vs Acceleration
scatter2 = alt.Chart(cars).mark_circle(size=60).encode(
 x='Weight_in_lbs:Q',
 y='Acceleration:Q',
 color=alt.condition(brush, 'Origin:N', alt.value('lightgray')),
 tooltip=['Name:N', 'Origin:N', 'Weight_in_lbs:Q', 'Acceleration:Q']
).properties(
 title='Weight vs Acceleration\n(Highlights match your selection!)',
 width=300, height=300
)

# Side-by-side linked charts
(scatter1 | scatter2).save('linked_brushing.html')
# Open linked_brushing.html → brush a region in left chart → right chart updates!

A regional sales manager needs to:
 1. See overall company trends
 2. Compare regions
 3. Drill down into specific stores
 4. Monitor KPIs against targets

Without interaction (static):
 Need 10+ different reports emailed every Monday.
 By the time you've read them, the data is stale.

With interaction (Tableau dashboard):
 One dashboard. All live data.

 User story:
 1. Sees national revenue map → notices Northeast is RED (below target)
 2. Brushes over Northeast stores → all other charts filter to NE only
 3. Bar chart now shows NE product categories → Electronics is dragging
 4. Clicks Electronics bar → time series shows Electronics sales dropped sharply in Week 3
 5. Adds tooltip → Week 3 note: "Competitor sale event"

 This insight took 2 minutes. It would have taken 3 hours with static reports!

Types of data sources in real enterprise systems:

Internal:
 CRM (Salesforce/HubSpot) → customer behavior, sales pipeline
 ERP (SAP/Oracle) → inventory, orders, financials
 App databases (MySQL) → user actions, transactions
 Data warehouse snapshots → historical aggregates

External:
 REST APIs → Twitter/X sentiment, weather data, competitor prices
 IoT sensors → factory sensors, vehicle GPS, smart meters
 Web scraping → market research, product reviews
 Data vendors → demographic data, financial feeds (Bloomberg)

Real-time pathway (sub-second latency):
 Events → Kafka topics → Spark Streaming / Flink → Feature Store

Batch pathway (minutes-hours latency):
 Database dumps → Airflow DAG (schedule) → S3/HDFS

CDC (Change Data Capture):
 Debezium monitors database logs → streams INSERT/UPDATE/DELETE events to Kafka
 Your data warehouse is always in sync with the operational database!

# Simple Kafka producer example (sending IoT sensor data)
from kafka import KafkaProducer
import json, time, random

producer = KafkaProducer(
 bootstrap_servers=['kafka1:9092', 'kafka2:9092'],
 value_serializer=lambda v: json.dumps(v).encode('utf-8')
)

while True:
 event = {
 'device_id': 'sensor_042',
 'temperature': 25.0 + random.gauss(0, 2),
 'humidity': 65.0 + random.gauss(0, 5),
 'timestamp': time.time()
 }
 producer.send('iot-readings', event)
 time.sleep(1) # One reading per second

THE "MEDALLION ARCHITECTURE" (modern data lake pattern):

Bronze Layer (Raw): Exact copy of source data. Never delete. Never modify.
 HDFS: /data/bronze/salesforce/2026/02/26/*.json

Silver Layer (Clean): Standardized, deduplicated, parsed, type-corrected.
 HDFS: /data/silver/customers/dt=2026-02-26/*.parquet

Gold Layer (Enriched): Business-ready aggregates and features for models.
 Snowflake table: gold.customer_features

Feature Store: Serves pre-computed features to models at low latency.
 Training: reads historical features from S3
 Serving: reads real-time features from Redis

# Example: Training a churn prediction model with Spark MLlib on a cluster
from pyspark.ml import Pipeline
from pyspark.ml.feature import VectorAssembler, StandardScaler
from pyspark.ml.classification import GBTClassifier
from pyspark.ml.evaluation import BinaryClassificationEvaluator
from pyspark.sql import SparkSession

spark = SparkSession.builder.appName("ChurnModel").getOrCreate()

# Load 50 million customer records from data warehouse
df = spark.read.parquet("s3://company-datalake/gold/customer_features/")

feature_cols = ['usage_30d', 'support_tickets', 'days_since_last_login',
 'monthly_spend', 'plan_tier_encoded', 'num_products']

assembler = VectorAssembler(inputCols=feature_cols, outputCol='features_raw')
scaler = StandardScaler(inputCol='features_raw', outputCol='features')
gbt = GBTClassifier(labelCol='churned', featuresCol='features',
 maxIter=100, maxDepth=5)

pipeline = Pipeline(stages=[assembler, scaler, gbt])

train_df, test_df = df.randomSplit([0.8, 0.2], seed=42)
model = pipeline.fit(train_df)

# Evaluate
evaluator = BinaryClassificationEvaluator(labelCol='churned')
auc = evaluator.evaluate(model.transform(test_df))
print(f"Test AUC: {auc:.4f}") # e.g., 0.8923

# Save model
model.write().overwrite().save("s3://company-models/churn/v3/")

import mlflow
import mlflow.sklearn
from sklearn.ensemble import GradientBoostingClassifier
from sklearn.metrics import roc_auc_score, precision_score, recall_score

mlflow.set_experiment("churn-prediction")

with mlflow.start_run(run_name="GBT-v3-hypertuned"):
 # Log hyperparameters
 mlflow.log_param("n_estimators", 200)
 mlflow.log_param("max_depth", 5)
 mlflow.log_param("learning_rate", 0.05)
 mlflow.log_param("feature_set", "v3-with-support-tickets")

 # Train
 model = GradientBoostingClassifier(n_estimators=200, max_depth=5, learning_rate=0.05)
 model.fit(X_train, y_train)

 # Log metrics
 y_pred_proba = model.predict_proba(X_test)[:, 1]
 y_pred = model.predict(X_test)
 mlflow.log_metric("roc_auc", roc_auc_score(y_test, y_pred_proba))
 mlflow.log_metric("precision", precision_score(y_test, y_pred))
 mlflow.log_metric("recall", recall_score(y_test, y_pred))

 # Log the model artifact itself
 mlflow.sklearn.log_model(model, "churn_model")

 # Log feature importance plot
 plt.figure(figsize=(10, 6))
 feature_imp = pd.Series(model.feature_importances_, index=feature_cols)
 feature_imp.sort_values().plot(kind='barh', color='steelblue')
 plt.title('Feature Importances')
 plt.tight_layout()
 mlflow.log_figure(plt.gcf(), "feature_importance.png")

print(" Experiment logged! View at: http://localhost:5000")

from fastapi import FastAPI
from pydantic import BaseModel
import joblib
import numpy as np
import uvicorn

app = FastAPI(title="Churn Prediction API", version="3.0")

# Load model at startup (once, not per request!)
model = joblib.load("churn_model_v3.pkl")
scaler = joblib.load("scaler_v3.pkl")

class CustomerFeatures(BaseModel):
 usage_30d: float # Total usage in last 30 days
 support_tickets: int # Open support tickets
 days_since_last_login: int # Recency
 monthly_spend: float # Monthly revenue in ₹
 plan_tier_encoded: int # 0=Free, 1=Basic, 2=Pro, 3=Enterprise
 num_products: int # Number of products subscribed to

class PredictionResponse(BaseModel):
 churn_probability: float
 risk_segment: str # Low / Medium / High
 recommended_action: str

@app.post("/predict", response_model=PredictionResponse)
async def predict_churn(features: CustomerFeatures):
 # Convert to numpy array
 X = np.array([[
 features.usage_30d, features.support_tickets,
 features.days_since_last_login, features.monthly_spend,
 features.plan_tier_encoded, features.num_products
 ]])
 X_scaled = scaler.transform(X)

 prob = model.predict_proba(X_scaled)[0, 1]

 # Business rules on top of model output
 if prob < 0.3:
 segment, action = "Low", "No action needed"
 elif prob < 0.7:
 segment, action = "Medium", "Send retention email campaign"
 else:
 segment, action = "High", "Assign to retention specialist — call within 24h"

 return PredictionResponse(
 churn_probability=round(prob, 4),
 risk_segment=segment,
 recommended_action=action
 )

@app.get("/health")
async def health_check():
 return {"status": "healthy", "model_version": "3.0"}

# uvicorn main:app --host 0.0.0.0 --port 8000

PROBLEM:
 Customer complaints about inaccurate delivery time estimates
 "Your website said 2 days but it took 6!"

DATA SOURCES:
 GPS logs → actual transit times (historical)
 Weather APIs → real-time and forecast weather along routes
 Traffic APIs → current and historical congestion data
 Order system → distance, carrier, package weight

PIPELINE:
 Kafka → ingests GPS and weather streams in real-time
 Snowflake → stores all historical transit + features
 Spark → daily model training (5 million deliveries/week of data)
 GBT model → features: distance, day_of_week, carrier, weather_score, traffic_index
 FastAPI + Kubernetes → prediction endpoint (50ms latency)

MODEL FEATURES:
 distance_km (obvious)
 day_of_week (weekends are slower)
 is_holiday (holidays massively slow delivery)
 carrier_id (each carrier has different SLAs)
 live_traffic_speed_on_route (real-time!)
 weather_severity_forecast (rain/snow along route)
 origin_region_load (is the warehouse currently busy?)

RESULTS:
 ↓ 15% reduction in late deliveries
 ↓ 40% reduction in "where is my package?" customer contacts
 Routing team uses same pipeline to A/B test new routes

MONITORING:
 Alert fires if model's MAE (mean absolute error) rises above 4 hours
 Weekly automated retraining keeps model fresh as traffic patterns evolve

GDPR (Europe): Must get consent, can request deletion, must explain predictions
HIPAA (US Healthcare): PHI (Protected Health Info) cannot leave approved systems
PDP Bill (India): Similar consent and data localization requirements

Practical steps:
 Anonymize PII before sending to model training
 Audit trails for who accessed which data
 Right-to-explanation: can you explain WHY credit was denied?
 Model cards documenting fairness, limitations, use cases

SCENARIO: Fraud detection at payment gateway

Real-time (< 100ms budget):
 Can't run a full neural network! Too slow!
 Use: Lightweight logistic regression or decision tree
 Sacrifice some accuracy for speed

Batch nightly:
 Run full deep learning on all day's transactions
 Flag suspicious patterns for human review

HYBRID approach (best of both worlds):
 Real-time: fast model catches 80% of obvious fraud
 Batch: deep model re-scores everything, catches remaining 15%
 Human review: handles the remaining edge cases

# Instead of predicting price, predict log(price)
model = sm.OLS(np.log(y), X).fit()
# To get predictions on original scale: np.exp(model.predict(X_new))

# If variance is proportional to fitted values
weights = 1.0 / model.fittedvalues
wls_model = sm.WLS(y, X, weights=weights).fit()

model_robust = sm.OLS(y, X).fit(cov_type='HC3')

from scipy.stats import ks_2samp
stat, p = ks_2samp(train_data['age'], live_data['age'])
if p < 0.05: print("Data drift detected in 'age' feature!")

# Re-run evaluation on a fresh labeled dataset
# Make sure you're computing precision correctly
# Check for bugs in the evaluation code first!

# Compare training data distribution vs recent production data
for col in feature_columns:
 stat, p = ks_2samp(train_df[col], live_df[col])
 if p < 0.05:
 print(f" Drift in feature: {col}")

# Look at False Positives (predicted churn but they stayed)
fp_customers = test_df[(predictions == 1) & (actuals == 0)]
print("False Positives profile:")
print(fp_customers[features].describe())
# Compare to True Positives — what's different?