Optuna Hyperparameter Tuning

Optuna Hyperparameter Tuning for Optimization Algorithms

You are an expert in automated algorithm configuration with Optuna. This skill covers tuning the parameters of metaheuristics and MIP solvers: defining search spaces, choosing samplers (TPE, CMA-ES, NSGA-II), building multi-instance objectives, pruning bad configurations early, persisting and parallelizing studies, and validating tuned configurations on held-out instances. Use the pattern catalog below to set up tuning runs that produce parameters which generalize, not parameters that memorize the training instances.

Initial Assessment

Establish these facts before writing any tuning code:

• What is being tuned? A metaheuristic (SA, GA, ALNS), a MIP solver (Gurobi, CP-SAT), or a hybrid? The target determines the objective metric and the cost per evaluation.
• Parameter inventory. List every parameter with its type (continuous, integer, categorical, conditional), a plausible range, and the current default. Fewer than 10 parameters is the normal case; more than 15 suggests the algorithm design should be simplified first.
• Objective metric. Solution quality at a fixed budget? Time to proven optimality? Gap at a time limit? Anytime behavior? The metric must match how the algorithm will be used and reported later.
• Cost per trial. One trial = (instances per trial) x (seeds per instance) x (single-run budget). Compute the total wall-clock cost of the tuning run before starting it.
• Instance set. How many training instances exist? Are they representative of the instances used in the final experiments? Is a train/test split possible (it should be)?
• Stochasticity. Is the tuned algorithm randomized? If yes, plan multiple seeds per instance inside every trial; a single seed makes the objective so noisy that TPE chases luck.
• Tuning budget. How many trials are affordable? Under 30 trials, model-based sampling barely beats random search; plan the budget before choosing the sampler.
• Hardware. Single machine or cluster? Parallel trials need shared storage (SQLite is fine for a handful of workers, JournalStorage or an RDBMS beyond that).
• Reproducibility requirements. Will the tuning protocol be described in a paper? Then fix and record: sampler seed, instance list, seed lists, per-run budgets, and the Optuna version.
• Reporting obligation. Fair comparisons require that every compared algorithm receives the same tuning effort (same budget, same protocol). Decide this now, not after the experiments.

The Algorithm Configuration Problem

Offline parameter tuning is itself an optimization problem. Given an algorithm $A$ with parameter space $\Theta$, a distribution $\mathcal{D}$ over problem instances, and a cost function $c(\theta, \pi, s)$ (the result of running $A(\theta)$ on instance $\pi$ with random seed $s$), the goal is

$$ \theta^{*} \in \arg\min_{\theta \in \Theta} ; \mathbb{E}_{\pi \sim \mathcal{D}, ; s} \big[ , c(\theta, \pi, s) , \big]. $$

Three properties make this problem awkward and drive every pattern in this skill:

1. The expectation is approximated by a finite sample of training instances and seeds. The empirical minimizer overfits this sample exactly like a machine-learning model overfits a training set. Generalization is checked with held-out instances, never assumed.
2. Evaluating one $\theta$ is expensive (a full algorithm run per instance per seed), so the tuner must be sample-efficient: model-based sampling (TPE) and early stopping of bad configurations (pruning, racing) both attack this.
3. The evaluation is noisy for randomized algorithms, so single observations of $c(\theta, \pi, s)$ must never be trusted; aggregate over seeds and instances.

Families of tuning approaches

How TPE works, in two sentences

The Tree-structured Parzen Estimator (Bergstra, Bardenet, Bengio & Kegl 2011, "Algorithms for Hyper-Parameter Optimization") splits observed trials into a good set (cost below a quantile $\gamma$) and a bad set, fits kernel-density estimates $l(\theta)$ over the good and $g(\theta)$ over the bad configurations, and proposes the candidate maximizing $l(\theta)/g(\theta)$ — the configuration most characteristic of good trials and least characteristic of bad ones. It handles continuous, integer, log-scaled, and categorical dimensions natively, which is exactly the mix that metaheuristic parameter spaces have.

Sampler selection

Budgeting rule of thumb

Total cost = trials x instances/trial x seeds/instance x run budget. As a floor, give the tuner 10-20 trials per tuned parameter, after the n_startup_trials random warm-up (default 10). Worked example: 5 parameters, 6 training instances, 2 seeds, 30 s per run, 100 trials = 100 x 6 x 2 x 30 s = 10 hours sequential. If that is too much, cut in this order: (1) per-run budget (tune on a reduced budget, validate on the real one), (2) seeds via pruning (Pattern 4), (3) instance count — never below 4-5, or the tuned configuration will be instance-specific.

Pattern Catalog: Search Spaces and Objectives

Each pattern below is a self-contained, runnable block: motivation, code, pitfall. Helper functions are repeated where needed so every block stands alone.

Pattern 1 — Minimal tuning loop for a metaheuristic

The base pattern: an objective(trial) function that (a) draws a configuration from the search space, (b) runs the algorithm under that configuration on training data over several seeds, (c) returns an aggregate score. Here: simulated annealing for the quadratic assignment problem (QAP), with an O(n) swap-delta evaluation so the tuning run itself is not needlessly slow.

import numpy as np
import optuna


def random_qap(n: int, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Random symmetric QAP instance (flow, distance) with zero diagonals."""
    rng = np.random.default_rng(seed)
    flow = np.triu(rng.integers(0, 100, (n, n)), k=1)
    dist = np.triu(rng.integers(1, 50, (n, n)), k=1)
    return flow + flow.T, dist + dist.T


def qap_cost(perm: np.ndarray, flow: np.ndarray, dist: np.ndarray) -> float:
    """Cost of assigning facility i to location perm[i]."""
    return float(np.sum(flow * dist[np.ix_(perm, perm)]))


def swap_delta(perm: np.ndarray, flow: np.ndarray, dist: np.ndarray,
               r: int, s: int) -> float:
    """O(n) cost change for swapping positions r and s (symmetric QAP)."""
    mask = np.ones(len(perm), dtype=bool)
    mask[[r, s]] = False
    df = flow[r, mask] - flow[s, mask]
    dd = dist[perm[s], perm[mask]] - dist[perm[r], perm[mask]]
    return 2.0 * float(df @ dd)


def sa_qap(flow: np.ndarray, dist: np.ndarray, t0: float, alpha: float,
           moves_per_temp: int, budget: int, seed: int) -> float:
    """Simulated annealing for QAP: swap moves, geometric cooling."""
    rng = np.random.default_rng(seed)
    perm = rng.permutation(len(flow))
    cost = qap_cost(perm, flow, dist)
    best, temp, evals = cost, t0, 0
    while evals < budget:
        for _ in range(min(moves_per_temp, budget - evals)):
            r, s = rng.choice(len(perm), size=2, replace=False)
            delta = swap_delta(perm, flow, dist, r, s)
            evals += 1
            if delta < 0 or rng.random() < np.exp(-delta / temp):
                perm[[r, s]] = perm[[s, r]]
                cost += delta
                best = min(best, cost)
        temp = max(alpha * temp, 1e-9)
    return best


def objective(trial: optuna.Trial) -> float:
    """Mean best cost of SA over 3 seeds on one training instance."""
    t0 = trial.suggest_float("t0", 1.0, 1e4, log=True)
    alpha = trial.suggest_float("alpha", 0.80, 0.999)
    moves = trial.suggest_int("moves_per_temp", 10, 1000, log=True)
    flow, dist = random_qap(25, seed=0)
    costs = [sa_qap(flow, dist, t0, alpha, moves, budget=30_000, seed=s)
             for s in (1, 2, 3)]
    return float(np.mean(costs))


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(direction="minimize",
                            sampler=optuna.samplers.TPESampler(seed=42))
study.optimize(objective, n_trials=50)
print(study.best_params)
print(round(study.best_value, 1))
# Expected: a (t0, alpha, moves_per_temp) combination whose mean cost over the
# three seeds beats the default (t0=100, alpha=0.95, moves=100) by a clear margin.

Search-space conventions that matter: use log=True for scale-type parameters (temperatures, penalty coefficients, mutation rates, population sizes) so TPE searches orders of magnitude uniformly; use step= on integers when only multiples make sense; keep ranges generous on the first run, then narrow once study.trials_dataframe() shows where good trials cluster.

Pitfall: This block tunes on a single instance, which is acceptable only as a teaching device. Parameters tuned on one instance encode that instance's quirks (size, density, cost scale) and routinely lose to defaults elsewhere; real tuning uses Pattern 3 (multi-instance objective) plus Pattern 6 (held-out validation). Equally common: returning the cost of a single seed — TPE's good/bad split is then contaminated by luck, and the "best trial" is often just the luckiest seed.

Pattern 2 — Conditional search spaces

Metaheuristic spaces are rarely flat: choosing the crossover operator activates that operator's own parameters; choosing penalty-based constraint handling activates a penalty coefficient. Optuna supports this directly — call suggest_* inside if branches. Sample a parameter only when the configuration actually uses it, and give branch-specific parameters distinct names.

import numpy as np
import optuna


def ga_knapsack(values: np.ndarray, weights: np.ndarray, capacity: float,
                pop_size: int, crossover: str, swap_prob: float, mode: str,
                penalty_coef: float, mutation_rate: float,
                generations: int, seed: int) -> float:
    """Vectorized binary GA for 0-1 knapsack; returns best feasible value."""
    rng = np.random.default_rng(seed)
    n = len(values)
    drop_order = np.argsort(values / weights)  # worst ratio first

    def repair(pop: np.ndarray) -> np.ndarray:
        for ind in pop:  # rows are views; edits write through
            load = float(ind @ weights)
            for j in drop_order:
                if load <= capacity:
                    break
                if ind[j]:
                    ind[j] = False
                    load -= weights[j]
        return pop

    pop = rng.random((pop_size, n)) < 0.5
    if mode == "repair":
        pop = repair(pop)
    best = 0.0
    for _ in range(generations):
        load = pop @ weights
        value = pop @ values
        excess = np.maximum(0.0, load - capacity)
        fit = value - (penalty_coef * excess if mode == "penalty" else 0.0)
        feasible = load <= capacity
        if feasible.any():
            best = max(best, float(value[feasible].max()))
        a = rng.integers(0, pop_size, size=pop_size)  # binary tournament
        b = rng.integers(0, pop_size, size=pop_size)
        parents = pop[np.where(fit[a] >= fit[b], a, b)]
        half = pop_size // 2
        p1, p2 = parents[:half], parents[half:2 * half]
        if crossover == "uniform":
            mask = rng.random((half, n)) < swap_prob
        else:  # one_point
            cut = rng.integers(1, n, size=(half, 1))
            mask = np.arange(n) >= cut
        pop = np.vstack([np.where(mask, p2, p1), np.where(mask, p1, p2)])
        pop = pop ^ (rng.random(pop.shape) < mutation_rate)
        if mode == "repair":
            pop = repair(pop)
    return best


rng = np.random.default_rng(0)
VALUES = rng.integers(10, 100, 50).astype(float)
WEIGHTS = rng.integers(5, 50, 50).astype(float)
CAPACITY = 0.4 * float(WEIGHTS.sum())


def objective(trial: optuna.Trial) -> float:
    """Conditional space: operator choice activates operator parameters."""
    pop_size = trial.suggest_int("pop_size", 20, 200, step=20)
    mutation_rate = trial.suggest_float("mutation_rate", 1e-3, 0.1, log=True)
    crossover = trial.suggest_categorical("crossover", ["uniform", "one_point"])
    swap_prob = (trial.suggest_float("uniform_swap_prob", 0.1, 0.5)
                 if crossover == "uniform" else 0.5)
    mode = trial.suggest_categorical("mode", ["penalty", "repair"])
    coef = (trial.suggest_float("penalty_coef", 1.0, 100.0, log=True)
            if mode == "penalty" else 0.0)
    vals = [ga_knapsack(VALUES, WEIGHTS, CAPACITY, pop_size, crossover,
                        swap_prob, mode, coef, mutation_rate,
                        generations=100, seed=s) for s in (1, 2)]
    return float(np.mean(vals))


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(direction="maximize",
                            sampler=optuna.samplers.TPESampler(seed=13))
study.optimize(objective, n_trials=40)
print(study.best_params, round(study.best_value, 1))
# Expected: on this loosely constrained instance, repair-based handling with a
# moderate mutation rate typically wins; the conditional penalty_coef appears
# in best_params only when mode == "penalty" was selected.

Pitfall: Never reuse one parameter name across branches (e.g., a generic "prob" sampled under both crossover choices). Optuna treats parameters by name; merging two semantically different quantities into one distribution corrupts TPE's density model. Also know that with the default univariate TPE, each conditional parameter is modeled only from the trials where it was sampled — fine in principle, but with many branches each branch gets few observations. For heavily conditional spaces, set TPESampler(multivariate=True, group=True) (see Advanced Techniques) or consider irace, which was designed for such spaces.

Pattern 3 — Multi-instance objective with gap normalization

A configuration must work across the instance distribution, so each trial evaluates a set of training instances. Raw objective values must not be averaged across instances of different size or cost scale — the largest instance would dominate the mean. Normalize each run to a relative gap against a per-instance reference cost (best known value, or the best of a few default-configuration runs), then aggregate.

import numpy as np
import optuna


def random_qap(n: int, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Random symmetric QAP instance (flow, distance) with zero diagonals."""
    rng = np.random.default_rng(seed)
    flow = np.triu(rng.integers(0, 100, (n, n)), k=1)
    dist = np.triu(rng.integers(1, 50, (n, n)), k=1)
    return flow + flow.T, dist + dist.T


def sa_qap(flow: np.ndarray, dist: np.ndarray, t0: float, alpha: float,
           steps: int, seed: int) -> float:
    """Compact SA for QAP: swap moves, per-move geometric cooling."""
    rng = np.random.default_rng(seed)
    perm = rng.permutation(len(flow))
    cost = float(np.sum(flow * dist[np.ix_(perm, perm)]))
    best, temp = cost, t0
    for _ in range(steps):
        r, s = rng.choice(len(perm), size=2, replace=False)
        cand = perm.copy()
        cand[[r, s]] = cand[[s, r]]
        new = float(np.sum(flow * dist[np.ix_(cand, cand)]))
        if new < cost or rng.random() < np.exp((cost - new) / temp):
            perm, cost = cand, new
            best = min(best, cost)
        temp = max(alpha * temp, 1e-12)
    return best


TRAIN = [random_qap(n, seed) for n, seed in [(15, 1), (20, 2), (25, 3), (30, 4)]]
SEEDS = (101, 102)
STEPS = 4000

# Per-instance reference: best of five default-configuration runs.
REF = [min(sa_qap(f, d, t0=200.0, alpha=0.999, steps=STEPS, seed=s)
           for s in range(5)) for f, d in TRAIN]


def objective(trial: optuna.Trial) -> float:
    """Mean relative gap to reference over all training instances and seeds."""
    t0 = trial.suggest_float("t0", 1.0, 1e4, log=True)
    alpha = trial.suggest_float("alpha", 0.990, 0.99999)
    gaps = []
    for (flow, dist), ref in zip(TRAIN, REF):
        for s in SEEDS:
            cost = sa_qap(flow, dist, t0, alpha, STEPS, s)
            gaps.append((cost - ref) / ref)
    return float(np.mean(gaps))


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(direction="minimize",
                            sampler=optuna.samplers.TPESampler(seed=21))
study.optimize(objective, n_trials=40)
print(study.best_params, round(study.best_value, 4))
# Expected: best mean gap at or below 0.0 -- the tuned configuration matches or
# beats the default-configuration reference across all four training instances.

Aggregation choice: the mean gap rewards configurations that are good everywhere but is sensitive to one catastrophic instance; the median is robust to outliers but can hide a configuration that fails on a minority of instances. For solver runtimes use the shifted geometric mean (Pattern 7). When instance hardness varies wildly, report both mean and worst-case gap, and consider tuning on the mean while monitoring the max.

Pitfall: Averaging raw costs instead of gaps is the classic error: a configuration that shaves 2% off the n=30 instance moves the mean more than one that shaves 20% off the n=15 instance, so the tuner silently specializes to the biggest instance. The second error is a reference set that is too weak (e.g., a single random-restart cost): gaps then compress toward large negative values and differences between good configurations become invisible. Spend real compute on the references once; they are reused by every trial.

Pattern Catalog: Pruning, Storage, and Parallel Studies

Pattern 4 — Instance-wise pruning (racing on a budget)

When each trial evaluates several instances sequentially, a configuration that is clearly bad after two instances should not waste budget on the remaining ones. Report the running aggregate after each instance with trial.report(value, step) and stop when trial.should_prune() fires. This is a pragmatic, model-based cousin of racing methods like F-Race (Birattari, Stuetzle, Paquete & Varrentrapp 2002, "A Racing Algorithm for Configuring Metaheuristics").

import numpy as np
import optuna


def random_qap(n: int, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Random symmetric QAP instance (flow, distance) with zero diagonals."""
    rng = np.random.default_rng(seed)
    flow = np.triu(rng.integers(0, 100, (n, n)), k=1)
    dist = np.triu(rng.integers(1, 50, (n, n)), k=1)
    return flow + flow.T, dist + dist.T


def sa_qap(flow: np.ndarray, dist: np.ndarray, t0: float, alpha: float,
           steps: int, seed: int) -> float:
    """Compact SA for QAP: swap moves, per-move geometric cooling."""
    rng = np.random.default_rng(seed)
    perm = rng.permutation(len(flow))
    cost = float(np.sum(flow * dist[np.ix_(perm, perm)]))
    best, temp = cost, t0
    for _ in range(steps):
        r, s = rng.choice(len(perm), size=2, replace=False)
        cand = perm.copy()
        cand[[r, s]] = cand[[s, r]]
        new = float(np.sum(flow * dist[np.ix_(cand, cand)]))
        if new < cost or rng.random() < np.exp((cost - new) / temp):
            perm, cost = cand, new
            best = min(best, cost)
        temp = max(alpha * temp, 1e-12)
    return best


TRAIN = [random_qap(n, seed) for n, seed in [(15, 1), (20, 2), (25, 3), (30, 4)]]
SEEDS = (101, 102)
STEPS = 4000
REF = [min(sa_qap(f, d, t0=200.0, alpha=0.999, steps=STEPS, seed=s)
           for s in range(5)) for f, d in TRAIN]


def objective(trial: optuna.Trial) -> float:
    """Evaluate instances in a FIXED order; prune on the running mean gap."""
    t0 = trial.suggest_float("t0", 1.0, 1e4, log=True)
    alpha = trial.suggest_float("alpha", 0.990, 0.99999)
    gaps: list[float] = []
    for step, ((flow, dist), ref) in enumerate(zip(TRAIN, REF)):
        for s in SEEDS:
            cost = sa_qap(flow, dist, t0, alpha, STEPS, s)
            gaps.append((cost - ref) / ref)
        trial.report(float(np.mean(gaps)), step)
        if trial.should_prune():
            raise optuna.TrialPruned()
    return float(np.mean(gaps))


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(
    direction="minimize",
    sampler=optuna.samplers.TPESampler(seed=11),
    pruner=optuna.pruners.MedianPruner(n_startup_trials=10, n_warmup_steps=1),
)
study.optimize(objective, n_trials=60)
n_pruned = sum(t.state == optuna.trial.TrialState.PRUNED for t in study.trials)
print(round(study.best_value, 4), n_pruned)
# Expected: roughly a third to half of the 60 trials are pruned after 2-3 of
# the 4 instances, cutting total tuning time substantially with little or no
# loss in best_value compared to the unpruned study of Pattern 3.

Pruner choice: MedianPruner (prune if worse than the median of completed trials at the same step) is the safe default for instance-wise evaluation. SuccessiveHalvingPruner and HyperbandPruner (Jamieson & Talwalkar 2016, successive halving; Li et al. 2018, Hyperband) allocate budget more aggressively and pay off when trials report many steps (e.g., one step per instance-seed pair, or per restart of an anytime algorithm).

Pitfall: Intermediate values are only comparable across trials if every trial sees the instances in the same fixed order — never shuffle the instance list per trial. Order the instances from cheapest to most expensive so pruning saves the most time, but be aware this biases early pruning decisions toward performance on small instances; n_warmup_steps>=1 and a conservative n_startup_trials protect against killing configurations that only shine on the larger instances. Pruning also interacts with noise: with one seed per instance, the running mean is so noisy that good configurations get pruned on bad luck.

Pattern 5 — Persistent storage, resuming, and parallel workers

Tuning runs of hours must survive crashes and should use idle cores. Back the study with persistent storage; any number of worker processes can then attach to the same study by name and call optimize concurrently. Each worker process runs this same script.

import numpy as np
import optuna


def one_plus_one_es(sigma0: float, adapt: float, seed: int) -> float:
    """(1+1)-ES on a 10-D sphere; returns best objective after 500 evaluations."""
    rng = np.random.default_rng(seed)
    x = rng.normal(0.0, 1.0, 10)
    fx = float(x @ x)
    sigma = sigma0
    for _ in range(500):
        y = x + sigma * rng.normal(0.0, 1.0, 10)
        fy = float(y @ y)
        if fy <= fx:
            x, fx, sigma = y, fy, sigma * adapt
        else:
            sigma /= adapt ** 0.25
    return fx


def objective(trial: optuna.Trial) -> float:
    """Tune the ES step-size schedule; mean final quality over 3 seeds."""
    sigma0 = trial.suggest_float("sigma0", 1e-3, 10.0, log=True)
    adapt = trial.suggest_float("adapt", 1.01, 2.0)
    return float(np.mean([one_plus_one_es(sigma0, adapt, s) for s in (1, 2, 3)]))


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(
    study_name="es_sphere_v1",
    storage="sqlite:///es_sphere.db",   # file-backed; resumable
    direction="minimize",
    sampler=optuna.samplers.TPESampler(seed=7, multivariate=True),
    load_if_exists=True,                # attach instead of failing if it exists
)
study.optimize(objective, n_trials=25)

df = study.trials_dataframe(attrs=("number", "value", "params", "state"))
df.to_csv("es_sphere_trials.csv", index=False)
print(len(study.trials), study.best_params)
# Expected: re-running this script continues the same study (trial count grows
# by 25 each run); best sigma0 is moderate (around 0.5-3) with adapt near
# 1.2-1.6, the classic 1/5th-success-rule regime.

Pitfall: SQLite serializes writes through file locks; with more than roughly 4-8 concurrent workers (or a network filesystem) you will see database is locked errors. Switch to optuna.storages.JournalStorage with a journal file backend, or a PostgreSQL/MySQL storage URL, for larger worker counts. Second trap: giving every worker the same TPESampler(seed=7) makes workers propose near-identical configurations in parallel. Either give each worker a different sampler seed (and accept that exact reproducibility now requires replaying the stored study) or set constant_liar=True in TPESampler so concurrent workers are steered apart. Finally, never run two different objective definitions against the same study_name — the parameter records will be inconsistent and the surrogate model meaningless.

Pattern Catalog: Validation and Solver Tuning

Pattern 6 — Train/test instance split and the overtuning check

Tuning is learning, so it is evaluated like learning: split the instance set, tune on the training part, and report performance of the single chosen configuration on the held-out part with fresh seeds. The gap between training score and test score is the measured overtuning bias. Hooker (1995, "Testing Heuristics: We Have It All Wrong") made the underlying point early: results demonstrated on the development set alone are not evidence.

import numpy as np
import optuna


def random_qap(n: int, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Random symmetric QAP instance (flow, distance) with zero diagonals."""
    rng = np.random.default_rng(seed)
    flow = np.triu(rng.integers(0, 100, (n, n)), k=1)
    dist = np.triu(rng.integers(1, 50, (n, n)), k=1)
    return flow + flow.T, dist + dist.T


def sa_qap(flow: np.ndarray, dist: np.ndarray, t0: float, alpha: float,
           steps: int, seed: int) -> float:
    """Compact SA for QAP: swap moves, per-move geometric cooling."""
    rng = np.random.default_rng(seed)
    perm = rng.permutation(len(flow))
    cost = float(np.sum(flow * dist[np.ix_(perm, perm)]))
    best, temp = cost, t0
    for _ in range(steps):
        r, s = rng.choice(len(perm), size=2, replace=False)
        cand = perm.copy()
        cand[[r, s]] = cand[[s, r]]
        new = float(np.sum(flow * dist[np.ix_(cand, cand)]))
        if new < cost or rng.random() < np.exp((cost - new) / temp):
            perm, cost = cand, new
            best = min(best, cost)
        temp = max(alpha * temp, 1e-12)
    return best


SIZES_SEEDS = [(15, 1), (20, 2), (25, 3), (30, 4), (18, 5), (22, 6), (26, 7), (28, 8)]
ALL = [random_qap(n, s) for n, s in SIZES_SEEDS]
TRAIN, TEST = ALL[:5], ALL[5:]
STEPS = 4000


def references(instances: list[tuple[np.ndarray, np.ndarray]]) -> list[float]:
    """Best of five default-configuration runs per instance."""
    return [min(sa_qap(f, d, 200.0, 0.999, STEPS, s) for s in range(5))
            for f, d in instances]


REF_TRAIN, REF_TEST = references(TRAIN), references(TEST)


def evaluate(config: dict, instances: list, refs: list[float],
             seeds: tuple[int, ...]) -> float:
    """Mean relative gap of a configuration over instances x seeds."""
    gaps = [(sa_qap(f, d, config["t0"], config["alpha"], STEPS, s) - ref) / ref
            for (f, d), ref in zip(instances, refs) for s in seeds]
    return float(np.mean(gaps))


def objective(trial: optuna.Trial) -> float:
    config = {"t0": trial.suggest_float("t0", 1.0, 1e4, log=True),
              "alpha": trial.suggest_float("alpha", 0.990, 0.99999)}
    return evaluate(config, TRAIN, REF_TRAIN, seeds=(101, 102))


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(direction="minimize",
                            sampler=optuna.samplers.TPESampler(seed=3))
study.optimize(objective, n_trials=60)

DEFAULT = {"t0": 200.0, "alpha": 0.999}
TEST_SEEDS = (501, 502, 503)  # disjoint from the tuning seeds
print("train gap, tuned:  ", round(study.best_value, 4))
print("test  gap, tuned:  ", round(evaluate(study.best_params, TEST, REF_TEST, TEST_SEEDS), 4))
print("test  gap, default:", round(evaluate(DEFAULT, TEST, REF_TEST, TEST_SEEDS), 4))
# Expected: the tuned configuration beats the default on the test set, but its
# test gap is worse than its train gap; that difference is the overtuning bias
# and should be reported, not hidden.

Pitfall: The number printed by study.best_value is an optimistically biased estimate — it is the minimum of 60 noisy observations (the winner's curse). Never report it as the performance of the tuned algorithm. Re-run the chosen configuration on test instances with seeds that were never used during tuning. Also keep the split honest: if test instances come from a different size range or generator than training instances, you are measuring extrapolation, not generalization — sometimes that is exactly what you want to measure, but say so. For the final tuned-vs-default comparison across instances, use paired statistical tests as described in algorithm-benchmarking-statistics.

Pattern 7 — Tuning Gurobi parameters for time-to-optimality

Solver parameter tuning has a different objective shape: minimize runtime to proven optimality at a time limit, with a penalty for timeouts (PAR2 = twice the limit, a convention from the algorithm-configuration literature, e.g., Hutter, Hoos & Leyton-Brown 2011, "Sequential Model-Based Optimization for General Algorithm Configuration"). Aggregate runtimes with the shifted geometric mean, the standard in MIP benchmarking (Achterberg 2007, "Constraint Integer Programming").

import gurobipy as gp
import numpy as np
import optuna
from gurobipy import GRB


def random_mkp(n: int, m: int, seed: int) -> tuple[np.ndarray, np.ndarray, np.ndarray]:
    """Correlated multidimensional knapsack instance (Chu-Beasley style)."""
    rng = np.random.default_rng(seed)
    w = rng.integers(1, 1000, size=(m, n))
    cap = (0.25 * w.sum(axis=1)).astype(np.int64)
    profit = (w.mean(axis=0) / 2 + rng.integers(1, 500, size=n)).astype(np.int64)
    return profit, w, cap


def add_capacity_constraints(model: gp.Model, x: gp.tupledict,
                             w: np.ndarray, cap: np.ndarray) -> None:
    """One knapsack constraint per resource dimension."""
    m, n = w.shape
    for k in range(m):
        model.addConstr(
            gp.quicksum(int(w[k, j]) * x[j] for j in range(n)) <= int(cap[k]),
            name=f"capacity[{k}]",
        )


def par2_runtime(profit: np.ndarray, w: np.ndarray, cap: np.ndarray,
                 params: dict, time_limit: float) -> float:
    """Runtime to proven optimality; PAR2 penalty (2x limit) on timeout."""
    model = gp.Model("mkp")
    model.Params.OutputFlag = 0
    model.Params.Threads = 1          # deterministic-ish, fair timing
    model.Params.TimeLimit = time_limit
    for name, value in params.items():
        model.setParam(name, value)
    n = len(profit)
    x = model.addVars(n, vtype=GRB.BINARY, name="x")
    add_capacity_constraints(model, x, w, cap)
    model.setObjective(gp.quicksum(int(profit[j]) * x[j] for j in range(n)),
                       GRB.MAXIMIZE)
    model.optimize()
    if model.Status == GRB.OPTIMAL:
        return model.Runtime
    return 2.0 * time_limit


def shifted_geomean(times: list[float], shift: float = 1.0) -> float:
    """Shifted geometric mean, the standard MIP runtime aggregate."""
    arr = np.asarray(times, dtype=float) + shift
    return float(np.exp(np.log(arr).mean()) - shift)


TRAIN = [random_mkp(n=120, m=10, seed=s) for s in (1, 2, 3, 4)]


def objective(trial: optuna.Trial) -> float:
    params = {
        "MIPFocus": trial.suggest_categorical("MIPFocus", [0, 1, 2, 3]),
        "Heuristics": trial.suggest_float("Heuristics", 0.0, 0.5),
        "Cuts": trial.suggest_categorical("Cuts", [-1, 0, 1, 2]),
        "Presolve": trial.suggest_categorical("Presolve", [-1, 0, 1, 2]),
        "VarBranch": trial.suggest_categorical("VarBranch", [-1, 0, 1, 2, 3]),
    }
    times = [par2_runtime(p, w, c, params, time_limit=30.0) for p, w, c in TRAIN]
    return shifted_geomean(times)


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(direction="minimize",
                            sampler=optuna.samplers.TPESampler(seed=9))
study.optimize(objective, n_trials=40)
print(study.best_params, round(study.best_value, 2))
# Expected: on these small instances most configurations finish well inside the
# limit, so differences are modest; on instances sized to take 10-60 s at
# default settings, tuned parameters often give 1.5-3x geomean speedups.

Choose training instances that take roughly 10-60 seconds at default settings: instances that solve in milliseconds make all configurations look identical (runtime differences are noise), and instances that always time out make the objective a constant PAR2 plateau with no gradient for the tuner. If proving optimality is out of reach, switch the metric to the MIP gap at the time limit, read via model.MIPGap after checking model.SolCount > 0.

Pitfall: Wall-clock runtimes are machine- and load-dependent: a tuning run shared with other jobs learns the scheduler, not the solver. Pin Threads=1, run on an idle machine, and prefer Gurobi's deterministic work units (model.Params.WorkLimit and the model.Work attribute) over seconds when comparability across machines matters. Also remember Gurobi ships its own tuning tool (grbtune / model.tune()); use Optuna instead when the metric is custom (gap at limit, primal integral, multi-instance aggregates) or when solver parameters are tuned jointly with parameters of a surrounding matheuristic.

Advanced Techniques

Multivariate and grouped TPE

The default TPE is univariate: it models each parameter independently, which wastes information when parameters interact (cooling rate and moves-per-temperature in SA jointly determine total cooling; population size and mutation rate in a GA trade off against each other). TPESampler(multivariate=True) fits a joint kernel density over all parameters; adding group=True partitions the space into sub-spaces that co-occur in trials, which is the correct treatment for the conditional spaces of Pattern 2. For parallel studies add constant_liar=True so simultaneously running trials repel each other. These three flags together are the recommended non-default setup for metaheuristic tuning.

Warm-starting the study with known configurations

The published default configuration and configurations from earlier tuning runs are valuable evidence; force them to be evaluated first with enqueue_trial. This guarantees the baseline is in the study (so "tuned beats default" is checked inside the same protocol) and gives TPE a head start in a sensible region.

import optuna


def objective(trial: optuna.Trial) -> float:
    """Search space of Pattern 5; any multi-instance objective plugs in here."""
    sigma0 = trial.suggest_float("sigma0", 1e-3, 10.0, log=True)
    adapt = trial.suggest_float("adapt", 1.01, 2.0)
    return (sigma0 - 1.0) ** 2 + (adapt - 1.3) ** 2  # synthetic score for demo


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(direction="minimize",
                            sampler=optuna.samplers.TPESampler(seed=0))
study.enqueue_trial({"sigma0": 1.0, "adapt": 1.5})  # literature default
study.enqueue_trial({"sigma0": 0.3, "adapt": 1.2})  # previous tuning result
study.optimize(objective, n_trials=30)
print(study.trials[0].params, round(study.trials[0].value, 4))
# Expected: trial 0 is exactly the enqueued default configuration, so the
# tuned-vs-default comparison is available from the study itself.

Multi-objective tuning: quality versus budget

Solution quality and computational effort usually trade off, and collapsing them into one scalar hides the trade-off from the person who has to choose. Declare two objectives and let NSGA-II approximate the Pareto front of configurations; then pick the knee point, or the cheapest configuration within tolerance of the best quality.

import numpy as np
import optuna


def es_run(sigma0: float, adapt: float, evals: int, seed: int) -> float:
    """(1+1)-ES on a 10-D sphere with a configurable evaluation budget."""
    rng = np.random.default_rng(seed)
    x = rng.normal(0.0, 1.0, 10)
    fx = float(x @ x)
    sigma = sigma0
    for _ in range(evals):
        y = x + sigma * rng.normal(0.0, 1.0, 10)
        fy = float(y @ y)
        if fy <= fx:
            x, fx, sigma = y, fy, sigma * adapt
        else:
            sigma /= adapt ** 0.25
    return fx


def objective(trial: optuna.Trial) -> tuple[float, float]:
    """Objective 1: mean final quality. Objective 2: evaluation budget."""
    sigma0 = trial.suggest_float("sigma0", 1e-3, 10.0, log=True)
    adapt = trial.suggest_float("adapt", 1.01, 2.0)
    evals = trial.suggest_int("evals", 100, 2000, log=True)
    quality = float(np.mean([es_run(sigma0, adapt, evals, s) for s in (1, 2, 3)]))
    return quality, float(evals)


optuna.logging.set_verbosity(optuna.logging.WARNING)
study = optuna.create_study(directions=["minimize", "minimize"],
                            sampler=optuna.samplers.NSGAIISampler(seed=5))
study.optimize(objective, n_trials=80)
front = sorted((t.values[1], t.values[0], t.params) for t in study.best_trials)
for budget, quality, params in front[:5]:
    print(int(budget), round(quality, 6), params)
# Expected: study.best_trials holds the Pareto set; quality improves
# monotonically along the front as the evaluation budget grows.

CMA-ES as the sampler for continuous spaces

When every tuned parameter is numeric (temperatures, rates, weights — no categoricals, no conditionals), optuna.samplers.CmaEsSampler typically outperforms TPE after the first few dozen trials, because CMA-ES learns the covariance of good configurations and follows correlated valleys in the parameter space. The restrictions are real: categorical and conditional parameters fall back to random sampling inside CmaEsSampler, so a mixed space silently loses most of the benefit — split the tuning (categoricals first, then continuous refinement) or stay with multivariate TPE. The mechanics of CMA-ES itself, and using it directly without Optuna, are covered in evolution-strategies.

irace, racing, and when not to use Optuna

irace (Lopez-Ibanez, Dubois-Lacoste, Perez Caceres, Birattari & Stuetzle 2016, "The irace package: Iterated racing for automatic algorithm configuration") implements iterated F-Race: candidate configurations race over a growing in

…(truncated)