Planning with the average can get expensive.

How it works

You could simply take average demand and plan with it. But that misleads you, because the costs don’t rise evenly: a shortfall weighs heavier than surplus stock. Anyone working only with the average overlooks these expensive outliers. Quicopt does the opposite: it evaluates every possible day individually with your actual costs, and averages only afterwards. That way the bad days are counted correctly, and you get the order quantity that is cheapest across all scenarios. In short: don’t average the demand — choose the order quantity whose expected cost across all possible days is lowest.

Eight products, one shelf

For a single product there is still a formula for the best order quantity. The moment several products share the same space, that ends: what product A gets, product B goes without. Here there are eight products, each with its own demand and its own shortage cost, sharing 1020 units of shelf space — less than their average demand adds up to.

Where the space goes

Product D is the striking one: it has by far the highest average demand and still gets the least — because running out of it only costs €2. Products C and F are stocked the other way, above their own average. No rule built on averages arrives at this split.

Extracting distribution parameters from spreadsheet data

In the model below, the average mu and the standard deviation sigma are fixed numbers. In reality they come from a spreadsheet: one row per sales day, one column per product.

import pandas as pd

sales    = pd.read_csv("sales_daily.csv", index_col=0)   # first column is the date
products = list(sales.columns)
mu       = sales.mean().tolist()            # average daily demand
sigma    = sales.std().tolist()             # how much it varies

That way mu and sigma come directly out of your sales history. If the history is an Excel file, export it as CSV once.

The complete model

That is all of it. The code runs against Quicopt’s API endpoint — no solver runs on your machine.

import pyomo.environ as pyo
from quicopt import Client
from quicopt.stochastic import Normal, set_distribution, set_scenarios, expectation, maximum

products = ["A", "B", "C", "D", "E", "F", "G", "H"]
mu       = [200.0, 150.0,  90.0, 300.0, 120.0,  60.0, 180.0, 100.0]  # average daily demand
sigma    = [ 20.0,  55.0,  30.0,  30.0,  45.0,  25.0,  18.0,  38.0]  # how much it varies
short    = [  4.0,   4.0,  12.0,   2.0,   9.0,  15.0,   3.0,   7.0]  # cost of one missing unit
hold     = 1.0                                                       # cost of one unit in stock
capacity = 1020.0                                                    # shared shelf space

m = pyo.ConcreteModel()
set_scenarios(m, 500, seed=7)                            # 500 possible days, drawn once

m.Q = pyo.Var(range(8), domain=pyo.NonNegativeIntegers)  # what we order, in whole units
m.d = pyo.Var(range(8))                                  # demand: one random variable per product
for i in range(8):
    set_distribution(m, m.d[i], Normal(mu[i], sigma[i]), name=f"demand_{products[i]}")

m.shelf = pyo.Constraint(expr=sum(m.Q[i] for i in range(8)) <= capacity)   # it has to fit

cost = sum(short[i] * maximum(m.d[i] - m.Q[i], 0.0)
           + hold * maximum(m.Q[i] - m.d[i], 0.0) for i in range(8))
m.cost = pyo.Objective(expr=expectation(cost), sense=pyo.minimize)   # expected cost, all 500 days

result = Client().solve(m)
print(result.display)

Output

▄▀▀▀▄        ▄                     █
█   █ █   █ ▄▄  ▄▀▀▀▄ ▄▀▀▀▄ █▀▀▀▄ ▀█▀▀
█ ▀▄▀ █   █  █  █   ▄ █   █ █▄▄▄▀  █  ▄
 ▀▀ ▀  ▀▀▀▀ ▀▀▀  ▀▀▀   ▀▀▀  █       ▀▀
├── status:     heuristic
├── feasible:   true
├── objective:  1014.523297602143
├── x:          x1=193, x2=140, x3=113, x4=73, x5=139, x6=82, …  (8 variables)
└── solve_time: 1.0007 s

planned with the average : 1347.0 EUR/day
optimized                : 1014.5 EUR/day

product   avg demand   by average   optimized
   A          200          170         193
   B          150          128         140
   C           90           76         113
   D          300          255          73
   E          120          102         139
   F           60           51          82
   G          180          153         169
   H          100           85         111

And now it gets interesting

Quicopt works purely primally: it evaluates your cost function scenario by scenario without presupposing any structure. So you can formulate costs and conditions that have no “smooth” mathematical formulation. The two examples below each add exactly one line to the model above, and nothing else changes:

m.cost = pyo.Objective(
    expr=expectation(cost) + 50.0 * sum(prob(m.d[i] - m.Q[i] >= 0.0) for i in range(8)),
    sense=pyo.minimize)

m.service = pyo.Constraint(
    expr=prob(sum(maximum(m.d[i] - m.Q[i], 0.0) for i in range(8)) <= 350.0) >= 0.90)

For context: the solver status reads “heuristic”: for this class of objective there is no proof of optimality, and the solver returns the best solution found. Several allocations come within less than a per mille of each other, so a later version of the solver may return a slightly different one.