Bottom-up, top-down and middle-out

Today let's understand what to do when forecasts made one series at a time do not add up.

Priya is a demand planner at a bakery chain. The chain has 9 stores in 3 regions: North, Central and South. She has to give finance the forecast of monthly units sold for 2025, for the whole chain, for each region and for each store. That's 13 series in all, and finance will add them up.

So she forecasts each of the 13 series on its own. But when finance adds up the store forecasts, the region forecasts and the total forecast, the three do not match. The same question, how many units will we sell, has three different answers.

There are three standard ways to fix this: bottom-up, top-down and middle-out. We will build all three on the same data and see what each one costs in accuracy.

The sales are simulated to look like real bakery data. The chart shows the chain's total units for each month from January 2019 to December 2024.

Switch between line and point. Every December is the highest month of its year, and sales trend upward over the years.

Why forecasts made series by series do not add up

Let's build Priya's data first, because every number from here on comes from it.

The 9 stores roll up into a hierarchy with three levels. At the bottom are the 9 stores. In the middle are the 3 regions, and each region is the sum of its 3 stores. At the top is the total, the sum of all 9 stores.

The sales are simulated for 84 months, January 2019 to December 2025. Each store has its own starting level, its own growth per month and its own noise, and all stores share one seasonal shape that peaks in December. N3 grows the fastest, at 1.5% a month.

RInteractive R
# Simulate 84 months of unit sales for 9 stores, three in each of three regions set.seed(2026) stores <- data.frame( store = c("N1", "N2", "N3", "C1", "C2", "C3", "S1", "S2", "S3"), region = rep(c("North", "Central", "South"), each = 3), start = c(620, 380, 100, 470, 290, 170, 460, 190, 100), growth = c(-0.001, 0.001, 0.015, 0.000, 0.002, 0.004, 0.001, 0.005, 0.008), noise = c(0.05, 0.08, 0.12, 0.06, 0.10, 0.14, 0.04, 0.09, 0.13) ) season <- c(0.86, 0.84, 0.92, 0.96, 1.00, 1.02, 0.98, 0.96, 1.00, 1.06, 1.12, 1.28) season <- season / mean(season) sales <- matrix(0, nrow = 84, ncol = 9, dimnames = list(NULL, stores$store)) for (i in 1:9) { trend <- stores$start[i] * (1 + stores$growth[i])^(0:83) sales[, i] <- round(trend * rep(season, 7) * exp(rnorm(84, 0, stores$noise[i]))) } head(sales, 3) #> N1 N2 N3 C1 C2 C3 S1 S2 S3 #> [1,] 547 325 79 395 249 167 395 160 79 #> [2,] 493 342 83 392 218 147 383 150 84 #> [3,] 573 350 90 474 284 193 426 171 82

  

Each column is a store and each row is a month, so these are the first three months of 2019.

Now we add the stores up to regions and a total. A region or a total is just the sum of the stores under it.

RInteractive R
# Add the stores up to 3 regions and 1 total, giving 13 series in all region_names <- c("North", "Central", "South") store_names <- stores$store region_of <- setNames(stores$region, stores$store) stores_in <- function(region) store_names[region_of == region] all_series <- cbind( Total = rowSums(sales), sapply(region_names, function(r) rowSums(sales[, stores_in(r)])), sales ) train_rows <- 1:72 test_rows <- 73:84 dim(all_series) #> [1] 84 13 all(all_series[, "Total"] == rowSums(all_series[, region_names])) #> [1] TRUE

  

The last line confirms that the actual sales add up in every one of the 84 months. Nothing forces the same to hold for forecasts.

Months 1 to 72, January 2019 to December 2024, are the training months. Months 73 to 84 are the year 2025. That year is the holdout: data that no model is fit on, kept back so we can measure each forecast against what was really sold.

Next, one model per series. We fit an ETS model to each of the 13 training series and forecast the 12 months of 2025. ETS is exponential smoothing with an error term, a trend and a seasonal term.

Setting model = "MAM" fixes every series to multiplicative error, additive trend and multiplicative seasonality, and damped = FALSE keeps the trend going in a straight line instead of flattening it out over the forecast. With the same model form for all 13 series, the only thing that will differ between the methods later is how the forecasts are combined.

A forecast made for one series with its own model, without any regard for the other series, is called a base forecast. To see whether the 13 base forecasts agree with each other, we add up the base forecasts of each series' children, the series directly below it, and compare the sum with the series' own base forecast.

RInteractive R
# Fit one ETS(M,A,M) model to each of the 13 training series and forecast 2025 library(forecast) base_fc <- sapply(colnames(all_series), function(name) { y <- ts(all_series[train_rows, name], start = c(2019, 1), frequency = 12) fit <- ets(y, model = "MAM", damped = FALSE) as.numeric(forecast(fit, h = 12)$mean) }) children <- list( Total = region_names, North = stores_in("North"), Central = stores_in("Central"), South = stores_in("South") ) january <- base_fc[1, ] sum_of_children <- sapply(children, function(kids) sum(january[kids])) round(data.frame( base_forecast = january[names(children)], sum_of_children = sum_of_children, gap = january[names(children)] - sum_of_children ), 1) #> base_forecast sum_of_children gap #> Total 2873.7 2856.0 17.7 #> North 1117.0 1126.1 -9.2 #> Central 922.8 921.2 1.6 #> South 816.3 814.4 1.9

  

Read the Total row first. For January 2025 the three region forecasts add up to 2,856.0, but the total's own forecast is 2,873.7. That is a gap of 17.7 units. Inside North, the three stores add up to 1,126.1 while North's forecast is 1,117.0, a gap of -9.2.

A set of forecasts is coherent when every region equals the sum of its stores and the total equals the sum of the regions. These base forecasts are not coherent. The 13 models were fit separately, so nothing made them agree.

The gaps do not go away over a year. Let's add the 12 forecast months at each level.

RInteractive R
# Add the 12 forecast months of 2025 at each level and compare with what was sold round(c( total_forecast = sum(base_fc[, "Total"]), sum_of_regions = sum(base_fc[, region_names]), sum_of_stores = sum(base_fc[, store_names]), sold_in_2025 = sum(all_series[test_rows, "Total"]) )) #> total_forecast sum_of_regions sum_of_stores sold_in_2025 #> 39612 39604 39720 39779

  

The total's model says 39,612 units for 2025. The region models add up to 39,604 and the store models add up to 39,720. What was actually sold is 39,779. So the same year has three different forecasts.

Finance cannot plan against three numbers. Priya needs one set of forecasts that is coherent. Turning base forecasts into a coherent set is called reconciliation, and there are three standard rules for it.