Build a Stochastic Optimization Framework for Supply Chains

What You Will Build: A Working Stochastic Optimization Framework

You'll implement a functional stochastic optimization framework for automotive supply chain resilience. The system integrates Bayesian learning with robust optimization to dynamically adjust inventory policies under demand and supply uncertainty.

Recent research validates that integrated stochastic learning-optimization achieves 7.4% cost reduction in stable environments and 5.7% improvement during supply disruptions.

The working prototype processes historical demand data, updates probability distributions in real-time, and generates optimized procurement and production decisions. The framework follows the Bertsimas-Sim adjustable robust optimization approach, validated against realistic automotive operating conditions.

Prerequisites and Setup Checklist

Gather these requirements before starting. Missing components will cause failures in later steps.

• Python 3.9+ with NumPy, SciPy, and pandas installed
• Optimization solver: Gurobi (academic license free) or open-source CBC/GLPK
• Historical demand data: Minimum 12 months of daily or weekly records
• Supplier lead time data: Mean and variance for each supplier tier
• Computing resources: 8GB RAM minimum; 16GB recommended for simulation runs
• Access credentials: ERP or inventory management system API (read access)

Time estimate: 8-12 hours for initial implementation; 2-4 weeks for calibration and validation. Primary blockers include incomplete historical data and solver licensing issues.

Context: Why Stochastic Optimization for Automotive Supply Chain Resilience

Deterministic optimization assumes perfect forecasts. Reality delivers demand shocks, supplier failures, and logistics disruptions. Stochastic optimization explicitly models uncertainty, generating decisions that perform well across multiple scenarios.

The approach in this tutorial combines two proven methods. Bayesian inference continuously updates demand probability distributions as new data arrives. Robust optimization (Bertsimas-Sim) then generates decisions that remain feasible even when parameters deviate from expected values.

Alternative approaches exist. Scenario-based stochastic programming requires explicit scenario enumeration, which becomes computationally expensive. Pure machine learning approaches lack the mathematical guarantees that optimization provides. The integrated robust optimization and simulation framework from Aston University demonstrates that combining optimization with simulation improves supply chain resilience while supporting sustainable operational decision-making.

Step 1: Structure Your Data Pipeline

Action: Create a standardized data ingestion module that pulls historical demand, inventory levels, and supplier performance metrics.

Create a new Python file named data_pipeline.py. Your module must output three dataframes: demand history (indexed by date and SKU), inventory snapshots (daily positions by location), and supplier metrics (lead times, fill rates, capacity).

``` import pandas as pd from datetime import datetime, timedelta

class SupplyChainDataPipeline: def __init__(self, erp_connection): self.conn = erp_connection self.demand_df = None self.inventory_df = None self.supplier_df = None

def extract_demand_history(self, lookback_days=365): """Pull daily demand by SKU for Bayesian prior estimation""" query = f""" SELECT date, sku_id, quantity_demanded, quantity_fulfilled FROM demand_history WHERE date >= CURRENT_DATE - INTERVAL '{lookback_days} days' """ self.demand_df = pd.read_sql(query, self.conn) self.demand_df['date'] = pd.to_datetime(self.demand_df['date']) return self.demand_df ```

Expected result: Running pipeline.extract_demand_history() returns a dataframe with 365+ rows per active SKU. Checkpoint: Verify no null values in quantity columns and date coverage is continuous.

Common failure: Missing dates in demand history. Fix: Impute missing dates with zero demand or interpolated values before proceeding.

Step 2: Implement Bayesian Demand Estimation

Action: Build the Bayesian inference engine that updates demand distribution parameters as new observations arrive.

Create bayesian_estimator.py. Use conjugate priors for computational efficiency. For demand data following a normal distribution, use Normal-Inverse-Gamma conjugate pairs.

``` import numpy as np from scipy import stats

class BayesianDemandEstimator: def __init__(self, prior_mu=100, prior_kappa=1, prior_alpha=2, prior_beta=100): """Initialize with weakly informative priors""" self.mu_0 = prior_mu # Prior mean estimate self.kappa_0 = prior_kappa # Prior precision weight self.alpha_0 = prior_alpha # Prior shape for variance self.beta_0 = prior_beta # Prior scale for variance

def update_posterior(self, observations): """Update parameters given new demand observations""" n = len(observations) x_bar = np.mean(observations)

# Posterior parameters kappa_n = self.kappa_0 + n mu_n = (self.kappa_0 self.mu_0 + n x_bar) / kappa_n alpha_n = self.alpha_0 + n / 2 beta_n = (self.beta_0 + 0.5 np.sum((observations - x_bar)2) + (self.kappa_0 n (x_bar - self.mu_0)2) / (2 kappa_n))

return { 'mean': mu_n, 'variance': beta_n / (alpha_n - 1) if alpha_n > 1 else np.inf, 'confidence_interval': self._compute_ci(mu_n, kappa_n, alpha_n, beta_n) } ```

Expected result: Posterior mean converges toward sample mean as observations accumulate. With 100+ observations, posterior variance decreases significantly.

Checkpoint: Compare posterior mean to simple moving average. Difference should be less than 5% with sufficient data.

Common failure: Posterior variance explodes with outliers. Fix: Implement outlier detection (values beyond 3 standard deviations) and use robust estimation or winsorization.

Step 3: Build the Robust Optimization Model

Action: Implement the Bertsimas-Sim robust optimization formulation for inventory decisions under uncertainty.

Create robust_optimizer.py. The model minimizes worst-case costs while ensuring feasibility across the uncertainty set defined by your Bayesian estimates.

``` from gurobipy import Model, GRB, quicksum import numpy as np

class RobustInventoryOptimizer: def __init__(self, gamma=1.5): """ gamma: Bertsimas-Sim budget of uncertainty parameter Higher gamma = more conservative (robust) solutions """ self.gamma = gamma self.model = None

def build_model(self, demand_params, holding_cost, shortage_cost, order_cost, lead_time, planning_horizon): """Construct robust inventory optimization model""" self.model = Model("RobustInventory") T = planning_horizon

# Decision variables order_qty = self.model.addVars(T, name="order", lb=0) inventory = self.model.addVars(T, name="inventory", lb=0) shortage = self.model.addVars(T, name="shortage", lb=0)

# Auxiliary variables for robust counterpart z = self.model.addVars(T, name="z", lb=0) # Dual variables p = self.model.addVar(name="p", lb=0) # Budget constraint dual

# Objective: minimize expected cost plus robustness term self.model.setObjective( quicksum(holding_cost inventory[t] + shortage_cost shortage[t] + order_cost order_qty[t] for t in range(T)) + self.gamma p + quicksum(z[t] for t in range(T)), GRB.MINIMIZE )

# Inventory balance constraints with uncertainty for t in range(T): if t == 0: self.model.addConstr( inventory[t] + shortage[t] >= demand_params[t]['mean'] - order_qty[t] ) else: lt = min(lead_time, t) self.model.addConstr( inventory[t] + shortage[t] >= demand_params[t]['mean'] + inventory[t-1] - order_qty[t-lt] )

return self.model ```

Expected result: Model solves in under 60 seconds for 90-day planning horizon with 50 SKUs. Optimal objective value represents worst-case cost.

Checkpoint: Verify solution feasibility by substituting demand realizations at upper bound of uncertainty set. Inventory should remain non-negative.

Common failure: Model infeasible. Fix: Reduce gamma parameter or increase initial inventory bounds. Check that lead time does not exceed planning horizon.

Step 4: Integrate Simulation for Validation

Action: Build a discrete-event simulation module to test optimization decisions against realistic demand scenarios.

Create supply_chain_simulator.py. The simulator generates demand realizations from your Bayesian distributions and tracks inventory, orders, and costs over the planning horizon.

``` import numpy as np from collections import defaultdict

class SupplyChainSimulator: def __init__(self, optimizer, estimator, num_scenarios=1000): self.optimizer = optimizer self.estimator = estimator self.num_scenarios = num_scenarios self.results = defaultdict(list)

def run_simulation(self, initial_inventory, demand_params, holding_cost, shortage_cost, periods=365): """Monte Carlo simulation over planning horizon""" for scenario in range(self.num_scenarios): inventory = initial_inventory.copy() total_cost = 0 service_level_hits = 0

for t in range(periods): # Generate demand realization demand = np.random.normal( demand_params[t]['mean'], np.sqrt(demand_params[t]['variance']) ) demand = max(0, demand) # Non-negative demand

# Apply optimized order decision order = self.optimizer.get_order_decision(t, inventory)

# Update inventory inventory = inventory + order - demand

# Calculate costs if inventory >= 0: total_cost += holding_cost inventory service_level_hits += 1 else: total_cost += shortage_cost abs(inventory) inventory = 0 # Lost sales assumption

self.results['total_cost'].append(total_cost) self.results['service_level'].append(service_level_hits / periods)

return self._compute_statistics() ```

Expected result: After 1000 scenarios, mean service level should exceed 95% with well-calibrated parameters. Cost distribution provides confidence intervals for budgeting.

Checkpoint: Service level variance should decrease as gamma increases (more robust solutions). If variance increases, check demand parameter estimation.

Step 5: Implement Adaptive Learning Loop

Action: Connect the Bayesian estimator to real-time demand feeds for continuous parameter updates.

Create adaptive_controller.py. This module orchestrates the learning-optimization cycle: observe demand, update beliefs, re-optimize decisions.

``` import schedule import time from datetime import datetime

class AdaptiveSupplyChainController: def __init__(self, data_pipeline, estimator, optimizer, update_frequency='daily'): self.pipeline = data_pipeline self.estimator = estimator self.optimizer = optimizer self.update_frequency = update_frequency self.decision_log = []

def update_cycle(self): """Execute one learning-optimization cycle""" # Step 1: Pull latest demand observations new_demand = self.pipeline.get_recent_demand(lookback_days=7)

# Step 2: Update Bayesian posteriors for sku in new_demand['sku_id'].unique(): sku_demand = new_demand[new_demand['sku_id'] == sku]['quantity'] posterior = self.estimator.update_posterior(sku_demand.values) self.estimator.store_posterior(sku, posterior)

# Step 3: Re-optimize with updated parameters demand_params = self.estimator.get_all_posteriors() self.optimizer.build_model(demand_params, holding_cost=2.5, shortage_cost=15.0, order_cost=50.0, lead_time=3, planning_horizon=30) self.optimizer.solve()

# Step 4: Log decisions for audit decisions = self.optimizer.get_decisions() self.decision_log.append({ 'timestamp': datetime.now(), 'decisions': decisions, 'demand_params': demand_params })

return decisions ```

Expected result: Simulations over 365 periods show adaptive Bayesian learning with stochastic optimization improves service levels by 2.0% in stable scenarios.

Checkpoint: Decision log should show posterior means converging over time. Erratic swings indicate insufficient data or outlier contamination.

Step 6: Configure Multi-Echelon Extension

Action: Extend the single-location model to handle supplier tiers and distribution networks.

Modify robust_optimizer.py to include upstream supplier constraints and downstream distribution requirements. This matches the Aston University case study that applies adjustable robust optimization to multi-echelon automotive systems.

``` def build_multi_echelon_model(self, nodes, arcs, demand_params, supply_capacity, transport_cost): """Extend to multi-echelon network""" self.model = Model("MultiEchelonRobust")

# Index sets suppliers = [n for n in nodes if n['type'] == 'supplier'] plants = [n for n in nodes if n['type'] == 'plant'] warehouses = [n for n in nodes if n['type'] == 'warehouse']

# Flow variables for each arc flow = self.model.addVars( [(a['from'], a['to'], t) for a in arcs for t in range(T)], name="flow", lb=0 )

# Supplier capacity constraints with uncertainty for s in suppliers: for t in range(T): self.model.addConstr( quicksum(flow[s['id'], p['id'], t] for p in plants if (s['id'], p['id']) in arcs) ```

Expected result: Multi-echelon model captures supply-side uncertainty in addition to demand uncertainty. Solution balances inventory across network tiers.

Configuration and Customization Parameters

These parameters require adjustment for your specific supply chain context. Default values work for typical automotive components with moderate demand variability.

• gamma (robustness budget): Default 1.5. Increase to 2.0-2.5 for critical components where stockouts are unacceptable. Decrease to 1.0 for commodity items where cost minimization dominates.
• prior\_kappa (prior weight): Default 1. Set to 0.1 for new SKUs with little historical data (let data dominate quickly). Set to 5-10 for stable SKUs where historical patterns are reliable.
• update\_frequency: Default daily. Use weekly for slow-moving items. Use hourly for high-velocity components in JIT environments.
• planning\_horizon: Default 30 days. Extend to 90 days for long lead time suppliers. Reduce to 7 days for responsive supply networks.
• num\_scenarios: Default 1000. Increase to 5000 for high-stakes decisions. Reduce to 500 for rapid prototyping.

Must-change settings: Replace holding_cost, shortage_cost, and order_cost with your actual cost parameters. These directly impact optimization decisions.

Verification and Testing Protocol

Execute this test sequence to validate your implementation before production deployment.

Unit tests: Verify Bayesian estimator converges to known parameters when fed synthetic data with known mean and variance. Tolerance: 5% deviation after 100 observations.

Integration test: Run full pipeline with 6 months of historical data. Compare optimized decisions against naive (order-up-to) policy. Optimization should show measurable cost improvement.

Stress test: Inject demand shock (50% increase) at simulation day 180. Monitor recovery time. Framework should detect shift within 7-14 days and adjust parameters accordingly.

Edge cases to verify:

• Zero demand periods (holidays, shutdowns): Posterior should not collapse
• Supplier failure (capacity drops to zero): Model should remain feasible with alternative sourcing
• Demand spike beyond historical range: Check that robust solution maintains service level

Success definition: Simulation shows service level above 95% and cost within 10% of theoretical optimum across 1000 scenarios.

Common Errors and Fixes

Error: GurobiError: Model is infeasible

Cause: Uncertainty set too large relative to available capacity or initial inventory.

Fix: Reduce gamma parameter by 0.5 increments until feasible. Alternatively, increase initial inventory bounds or add slack variables with high penalty costs.

Error: Posterior variance increases over time instead of decreasing

Cause: Non-stationary demand pattern (trend or seasonality) violates model assumptions.

Fix: Detrend data before estimation. Implement seasonal decomposition or use time-varying prior parameters.

Error: MemoryError during simulation with large scenario counts

Cause: Storing full trajectory for each scenario exhausts RAM.

Fix: Compute summary statistics incrementally. Store only final metrics per scenario, not full time series.

Error: Optimization decisions oscillate wildly between periods

Cause: Insufficient regularization or overfitting to recent observations.

Fix: Increase prior\_kappa to weight historical patterns more heavily. Add smoothing constraints to limit period-over-period order changes.

Error: Service level drops during demand shocks despite robust formulation

Cause: Research confirms that Bayesian updating conservatism limits responsiveness during sudden demand shocks.

Fix: Implement change-point detection to trigger faster prior updates when distribution shifts are detected. Consider hybrid approach with scenario-based planning for extreme events.

Next Steps and Extensions

Your working framework provides the foundation for advanced capabilities. Consider these extensions based on your operational priorities.

Multi-objective optimization: Extend the model to balance cost, sustainability, and social factors. The Aston University framework demonstrates that robust optimization with adjustable parameters supports economic and social sustainability in stochastic automotive operations.

Real-time integration: Connect your adaptive controller to Supply Chain Disaster's hazard intelligence platform for early warning of supplier disruptions. Trigger immediate re-optimization when risk alerts arrive.

Machine learning enhancement: Replace Bayesian estimation with neural network forecasters for complex demand patterns. Maintain robust optimization layer for decision-making with uncertainty quantification.

For longitudinal studies on SCRES, track your framework's performance over 12-24 months. Document cost savings, service level improvements, and response times during actual disruptions. This empirical evidence validates the approach and supports continuous refinement.