Mixed Continuous-Discrete Bayesian Optimization¶

This tutorial demonstrates single-objective Bayesian optimization in Xopt for a mixed search space containing continuous and discrete variables.

We use a synthetic objective where a sinusoidal continuous landscape changes with a discrete phase value. This makes the optimizer balance:

• global exploration over discrete choices
• local optimization in a continuous coordinate

In [1]:

import math
import os
import random
import warnings

import matplotlib.pyplot as plt
import numpy as np
import torch

from xopt import Evaluator, Xopt
from xopt.generators.bayesian import ExpectedImprovementGenerator
from xopt.vocs import VOCS

warnings.filterwarnings("ignore")

SMOKE_TEST = os.environ.get("SMOKE_TEST")

# Keep smoke runs fast for docs validation while preserving behavior.
N_INITIAL = 3 if SMOKE_TEST else 5
N_STEPS = 2 if SMOKE_TEST else 25
N_RESTARTS = 1 if SMOKE_TEST else 8
N_MC_SAMPLES = 8 if SMOKE_TEST else 128
N_GRID = 30 if SMOKE_TEST else 100
SEED = 123

torch.manual_seed(SEED)
np.random.seed(SEED)
random.seed(SEED)
torch.set_num_threads(1)

import math import os import random import warnings import matplotlib.pyplot as plt import numpy as np import torch from xopt import Evaluator, Xopt from xopt.generators.bayesian import ExpectedImprovementGenerator from xopt.vocs import VOCS warnings.filterwarnings("ignore") SMOKE_TEST = os.environ.get("SMOKE_TEST") # Keep smoke runs fast for docs validation while preserving behavior. N_INITIAL = 3 if SMOKE_TEST else 5 N_STEPS = 2 if SMOKE_TEST else 25 N_RESTARTS = 1 if SMOKE_TEST else 8 N_MC_SAMPLES = 8 if SMOKE_TEST else 128 N_GRID = 30 if SMOKE_TEST else 100 SEED = 123 torch.manual_seed(SEED) np.random.seed(SEED) random.seed(SEED) torch.set_num_threads(1)

Define a mixed search space¶

The variable phase is discrete with values {0.0, 0.8, 1.6, 2.4}.

The objective is

$$ f(x_1, \phi) = (\cos(x_1 + \phi) - 0.25)^2 + 0.1\sin(3x_1 + \phi) + 0.05\phi, $$

where $\phi$ is optimized directly as a discrete phase offset.

This keeps the mixed optimization example simple while still showing joint continuous-discrete behavior.

In [2]:

PHASE_VALUES = [0.0, 0.8, 1.6, 2.4]


def evaluate_mixed_sinusoid(inputs):
    x1 = float(inputs["x1"])
    phase = float(inputs["phase"])

    value = (
        (math.cos(x1 + phase) - 0.25) ** 2
        + 0.1 * math.sin(3.0 * x1 + phase)
        + 0.05 * phase
    )
    return {"f": value}


vocs = VOCS(
    variables={
        "x1": [0.0, 2.0 * math.pi],
        "phase": set(PHASE_VALUES),
    },
    objectives={"f": "MINIMIZE"},
)

vocs

PHASE_VALUES = [0.0, 0.8, 1.6, 2.4] def evaluate_mixed_sinusoid(inputs): x1 = float(inputs["x1"]) phase = float(inputs["phase"]) value = ( (math.cos(x1 + phase) - 0.25) ** 2 + 0.1 * math.sin(3.0 * x1 + phase) + 0.05 * phase ) return {"f": value} vocs = VOCS( variables={ "x1": [0.0, 2.0 * math.pi], "phase": set(PHASE_VALUES), }, objectives={"f": "MINIMIZE"}, ) vocs

Out[2]:

VOCS(variables={'x1': ContinuousVariable(dtype=None, default_value=None, domain=[0.0, 6.283185307179586]), 'phase': DiscreteVariable(dtype=None, default_value=None, values={0.0, 0.8, 2.4, 1.6})}, objectives={'f': MinimizeObjective(dtype=None)}, constraints={}, constants={}, observables={})

Create Xopt objects¶

We use ExpectedImprovementGenerator and standard GP modeling options, with smoke-test-aware optimizer settings for quick docs checks.

In [3]:

evaluator = Evaluator(function=evaluate_mixed_sinusoid)
generator = ExpectedImprovementGenerator(vocs=vocs)
generator.gp_constructor.use_low_noise_prior = True
generator.numerical_optimizer.n_restarts = N_RESTARTS
generator.n_monte_carlo_samples = N_MC_SAMPLES

X = Xopt(evaluator=evaluator, generator=generator)

X.random_evaluate(N_INITIAL)
for _ in range(N_STEPS):
    X.step()

X.data

evaluator = Evaluator(function=evaluate_mixed_sinusoid) generator = ExpectedImprovementGenerator(vocs=vocs) generator.gp_constructor.use_low_noise_prior = True generator.numerical_optimizer.n_restarts = N_RESTARTS generator.n_monte_carlo_samples = N_MC_SAMPLES X = Xopt(evaluator=evaluator, generator=generator) X.random_evaluate(N_INITIAL) for _ in range(N_STEPS): X.step() X.data

Out[3]:

	x1	phase	f	xopt_runtime	xopt_error
0	1.907140	2.4	0.631500	0.000013	False
1	4.485319	0.8	0.224616	0.000003	False
2	4.857836	1.6	0.575022	0.000002	False
3	2.819172	0.8	1.351942	0.000002	False
4	1.762628	0.0	0.110287	0.000001	False
5	1.242293	0.0	-0.049969	0.000007	False
6	0.368800	0.0	0.555571	0.000006	False
7	1.155822	0.8	0.341090	0.000008	False
8	4.878955	0.0	0.094863	0.000007	False
9	4.318036	0.0	0.440033	0.000008	False
10	5.569488	0.0	0.171808	0.000008	False
11	5.232638	0.8	0.485608	0.000007	False
12	1.385346	0.0	-0.080614	0.000009	False
13	0.000000	2.4	1.162493	0.000006	False
14	3.896623	2.4	0.782253	0.000007	False
15	6.283185	0.0	0.562500	0.000007	False
16	6.283185	2.4	1.162493	0.000008	False
17	1.399218	0.0	-0.080760	0.000007	False
18	1.397922	0.0	-0.080769	0.000006	False
19	1.397266	0.0	-0.080771	0.000008	False
20	1.128814	2.4	1.455230	0.000008	False
21	2.702365	2.4	0.048639	0.000007	False
22	0.000000	0.8	0.311282	0.000007	False
23	3.701011	1.6	0.186780	0.000006	False
24	2.713436	1.6	0.456590	0.000006	False
25	1.396166	0.0	-0.080773	0.000007	False
26	1.395896	0.0	-0.080773	0.000007	False
27	1.395684	0.0	-0.080773	0.000008	False
28	1.395511	0.0	-0.080772	0.000007	False
29	1.395363	0.0	-0.080772	0.000010	False

Optimization diagnostics¶

We check three things:

1. best-so-far objective over evaluations
2. frequency of sampled phase values (rounded for display)
3. whether sampled phases exactly match the configured discrete set

In [4]:

data = X.data.copy()
data["eval_index"] = np.arange(len(data))
data["best_so_far"] = data["f"].cummin()

data["phase_raw"] = data["phase"].astype(float)
data["phase_display"] = data["phase_raw"].round(3)

allowed_phase = set(float(v) for v in PHASE_VALUES)
phase_is_configured = data["phase_raw"].apply(
    lambda value: any(np.isclose(value, v, atol=1e-8) for v in PHASE_VALUES)
)

fig, ax = plt.subplots(1, 2, figsize=(12, 4))
ax[0].plot(data["eval_index"], data["f"], "o-", alpha=0.5, label="objective")
ax[0].plot(data["eval_index"], data["best_so_far"], "k-", lw=2, label="best so far")
ax[0].set_xlabel("evaluation")
ax[0].set_ylabel("f")
ax[0].legend()

count_series = data["phase_display"].value_counts().sort_index()
ax[1].bar(count_series.index.astype(str), count_series.values, color="C1")
ax[1].set_xlabel("sampled phase (rounded)")
ax[1].set_ylabel("sample count")

fig.tight_layout()

print(f"Configured phase set: {sorted(allowed_phase)}")
print(
    f"Exact configured phase matches: {int(phase_is_configured.sum())}/{len(phase_is_configured)}"
)
if not phase_is_configured.all():
    print("Some sampled phase values are not exact configured values.")

data = X.data.copy() data["eval_index"] = np.arange(len(data)) data["best_so_far"] = data["f"].cummin() data["phase_raw"] = data["phase"].astype(float) data["phase_display"] = data["phase_raw"].round(3) allowed_phase = set(float(v) for v in PHASE_VALUES) phase_is_configured = data["phase_raw"].apply( lambda value: any(np.isclose(value, v, atol=1e-8) for v in PHASE_VALUES) ) fig, ax = plt.subplots(1, 2, figsize=(12, 4)) ax[0].plot(data["eval_index"], data["f"], "o-", alpha=0.5, label="objective") ax[0].plot(data["eval_index"], data["best_so_far"], "k-", lw=2, label="best so far") ax[0].set_xlabel("evaluation") ax[0].set_ylabel("f") ax[0].legend() count_series = data["phase_display"].value_counts().sort_index() ax[1].bar(count_series.index.astype(str), count_series.values, color="C1") ax[1].set_xlabel("sampled phase (rounded)") ax[1].set_ylabel("sample count") fig.tight_layout() print(f"Configured phase set: {sorted(allowed_phase)}") print( f"Exact configured phase matches: {int(phase_is_configured.sum())}/{len(phase_is_configured)}" ) if not phase_is_configured.all(): print("Some sampled phase values are not exact configured values.")

Configured phase set: [0.0, 0.8, 1.6, 2.4]
Exact configured phase matches: 30/30

Model and acquisition slices by discrete phase offset¶

To isolate the effect of the discrete offset, we visualize the model and acquisition over (x1, phase).

In [5]:

X.generator.visualize_model(variable_names=["x1", "phase"])

X.generator.visualize_model(variable_names=["x1", "phase"])

Out[5]:

(<Figure size 800x660 with 7 Axes>,
 array([[<Axes: title={'center': 'Posterior Mean [f]'}, xlabel='x1', ylabel='phase'>,
         <Axes: title={'center': 'Posterior SD [f]'}, xlabel='x1', ylabel='phase'>],
        [<Axes: title={'center': 'Acq. Function'}, xlabel='x1', ylabel='phase'>,
         <Axes: >]], dtype=object))

Best candidate and summary¶

The best observed point combines a continuous location x1 and a discrete phase offset.

In this example, Bayesian optimization jointly learns:

• which phase offset gives the best regime
• where to place continuous evaluations inside that regime

In [6]:

best_row = data.loc[data["f"].idxmin(), ["x1", "phase", "f"]]
summary = (
    data.groupby("phase_display", as_index=False)["f"]
    .min()
    .rename(columns={"f": "best_f_at_phase"})
)

print("Best observed candidate:")
display(best_row.to_frame().T)

print("Best objective reached at each sampled phase (rounded):")
display(summary.sort_values("best_f_at_phase"))

best_row = data.loc[data["f"].idxmin(), ["x1", "phase", "f"]] summary = ( data.groupby("phase_display", as_index=False)["f"] .min() .rename(columns={"f": "best_f_at_phase"}) ) print("Best observed candidate:") display(best_row.to_frame().T) print("Best objective reached at each sampled phase (rounded):") display(summary.sort_values("best_f_at_phase"))

Best observed candidate:

	x1	phase	f
25	1.396166	0.0	-0.080773

Best objective reached at each sampled phase (rounded):

	phase_display	best_f_at_phase
0	0.0	-0.080773
3	2.4	0.048639
2	1.6	0.186780
1	0.8	0.224616