Models

The models module contains the QBMBase abstract base class, which holds the visible and hidden units, weights, biases, and random number generator, and the BQRBM class which implements the bound-based quantum restricted Boltzmann machine on top of it. A BQRBM is trained with the train() method and sampled from with the sample() method, using either the simulation or the annealer backend depending on how it was instantiated. See the Getting Started page for a complete example, and the Annealer guide for setting up the annealer backend.

QBMBase

Bases: ABC

Abstract base class for Quantum Boltzmann Machines

Theory based on Quantum Boltzmann Machine by Amin et al. https://journals.aps.org/prx/abstract/10.1103/PhysRevX.8.021050

Source code in src/qbm/models/QBMBase.py

class QBMBase(ABC):
    """
    Abstract base class for Quantum Boltzmann Machines

    Theory based on Quantum Boltzmann Machine by Amin et al.
    https://journals.aps.org/prx/abstract/10.1103/PhysRevX.8.021050
    """

    b: np.ndarray
    W: np.ndarray

    def __init__(self, V_train: np.ndarray, n_hidden: int, seed: int | None) -> None:
        """
        Initializes the model.

        Args:
            V_train: Training data.
            n_hidden: Number of hidden units.
            seed: Seed for the random number generator.
        """
        self.V_train = V_train
        self.n_visible = V_train.shape[1]
        self.n_hidden = n_hidden
        self.n_qubits = self.n_visible + self.n_hidden
        self.seed = seed
        self.rng = get_rng(self.seed)
        self.grads: dict[str, np.ndarray] = {}

        self._initialize_weights_and_biases()

    def _apply_grads(self, learning_rate: float) -> None:
        """
        Applies the gradients from the positive and negative phases using the provided
        learning rate.

        Args:
            learning_rate: Learning rate to scale the gradients with.
        """
        self.b += learning_rate * (self.grads["b_pos"] - self.grads["b_neg"])
        self.W += learning_rate * (self.grads["W_pos"] - self.grads["W_neg"])

    def _binary_to_eigen(self, x: np.ndarray) -> np.ndarray:
        """
        Convert bit values {0, 1} to corresponding spin values {+1, -1}.

        Args:
            x: Input array of values {0, 1}.

        Returns:
            Output array of values {+1, -1}.
        """
        return (1 - 2 * x).astype(np.int8)

    def _eigen_to_binary(self, x: np.ndarray) -> np.ndarray:
        """
        Convert spin values {+1, -1} to corresponding bit values {0, 1}.

        Args:
            x: Input array of values {+1, -1}.

        Returns:
            Output array of values {0, 1}.
        """
        return ((1 - x) / 2).astype(np.int8)

    def _random_mini_batch_indices(self, mini_batch_size: int) -> list[np.ndarray]:
        """
        Generates random, non-intersecting sets of indices for creating mini-batches
        of the training data. The final mini-batch may be smaller than
        mini_batch_size if the training set size is not divisible by it.

        Args:
            mini_batch_size: Size of the mini-batches.

        Returns:
            List of numpy arrays, each array containing the indices corresponding to
            a mini-batch.
        """
        return np.split(
            self.rng.permutation(np.arange(self.V_train.shape[0])),
            np.arange(mini_batch_size, self.V_train.shape[0], mini_batch_size),
        )

    @abstractmethod
    def _compute_positive_grads(self, V_pos: np.ndarray) -> None:
        pass

    @abstractmethod
    def _compute_negative_grads(self, n_samples: int) -> None:
        pass

    @abstractmethod
    def _initialize_weights_and_biases(self, mu: float = 0, sigma: float = 0.1) -> None:
        pass

__init__(V_train, n_hidden, seed)

Initializes the model.

Parameters:

Name	Type	Description	Default
V_train	ndarray	Training data.	required
n_hidden	int	Number of hidden units.	required
seed	int | None	Seed for the random number generator.	required

Source code in src/qbm/models/QBMBase.py

def __init__(self, V_train: np.ndarray, n_hidden: int, seed: int | None) -> None:
    """
    Initializes the model.

    Args:
        V_train: Training data.
        n_hidden: Number of hidden units.
        seed: Seed for the random number generator.
    """
    self.V_train = V_train
    self.n_visible = V_train.shape[1]
    self.n_hidden = n_hidden
    self.n_qubits = self.n_visible + self.n_hidden
    self.seed = seed
    self.rng = get_rng(self.seed)
    self.grads: dict[str, np.ndarray] = {}

    self._initialize_weights_and_biases()

BQRBM

Bases: QBMBase

Bound-based Quantum Restricted Boltzmann Machine

Source code in src/qbm/models/BQRBM.py

class BQRBM(QBMBase):
    """
    Bound-based Quantum Restricted Boltzmann Machine
    """

    def __init__(
        self,
        V_train: np.ndarray,
        n_hidden: int,
        A_freeze: float,
        B_freeze: float,
        beta_initial: float = 1.0,
        beta_range: Sequence[float | int] = [0.1, 10],
        annealer_params: AnnealerParams | None = None,
        simulation_params: Mapping[str, Any] | None = None,
        seed: int | None = 0,
    ) -> None:
        """
        Initializes the model.

        Note:
            Exactly one of annealer_params or simulation_params must be provided.
            Whichever is provided determines whether the samples are generated by
            the annealer or the simulation, respectively.

            The simulation works by exact computation of ρ(s, T) = e^{-β * H(s)} / Z.
            H(s) is determined by A_freeze and B_freeze.

        Args:
            V_train: Training data set of visible vectors, shape (n_samples, n_visible).
                Values must be in {+1, -1}, or in {0, 1} which are converted to
                {+1, -1}.
            n_hidden: Number of hidden units.
            A_freeze: Value of A(s) at the freeze out point. Used for Γ = β * A_freeze.
                Units in GHz.
            B_freeze: Value of B(s) at the freeze out point. Used for
                b_i = -β * B_freeze * h_i and w_ij = -β * B_freeze * J_ij. Units in GHz.
            beta_initial: Initial value of the effective β. Units in 1/GHz.
            beta_range: Range of allowed β values, used for making sure β is not updated
                to an infeasible value (e.g. negative).
            annealer_params: Dictionary with keys:
                - "schedule": List of (t, s) tuples defining the anneal schedule.
                - "embedding": Dict mapping the logical to physical qubits.
                - "relative_chain_strength" [optional]: Relative chain strength value.
                - "qpu_params" [optional]: Parameters dict to unpack to DWaveSampler(),
                    e.g. {"region": "na-west-1", "solver": "Advantage_system4.1"}
            simulation_params: Dictionary with keys:
                - "beta": Effective β that the simulation generates samples at.
                - "h_range" [optional]: Allowed range of h values, defaults to
                    (-inf, inf).
                - "J_range" [optional]: Allowed range of J values, defaults to
                    (-inf, inf).
            seed: Seed for the random number generator. Used for random minibatches, as
                well as the exact sampler.

        Raises:
            ValueError: If neither or both of annealer_params and simulation_params are
                provided, if a required key is missing from the provided params dict,
                if beta_initial is not positive, if beta_range does not satisfy
                0 < min < max, or if the training data values are not in {-1, +1}.
        """
        self.qpu: DWaveSampler | None = None
        self.sampler: AnnealerSampler | None = None
        if beta_initial <= 0:
            raise ValueError(f"beta_initial must be positive (got {beta_initial})")
        if beta_range[0] <= 0 or beta_range[0] >= beta_range[1]:
            raise ValueError(
                f"beta_range must satisfy 0 < min < max (got {list(beta_range)})"
            )
        # Convert from binary to ±1 if necessary
        if set(np.unique(V_train)) == set([0, 1]):
            V_train = self._binary_to_eigen(V_train)
        if set(np.unique(V_train)) != set([-1, 1]):
            raise ValueError(
                "V_train values must be in {+1, -1} or {0, 1} "
                f"(got unique values {set(np.unique(V_train))})"
            )

        self.A_freeze = A_freeze
        self.B_freeze = B_freeze
        self.beta = beta_initial
        self.beta_range: Sequence[float | int] = beta_range
        self.beta_history: list[float] = [self.beta]
        super().__init__(V_train=V_train, n_hidden=n_hidden, seed=seed)

        # Check if requirements met to use annealer or simulation
        if (annealer_params is not None and simulation_params is not None) or (
            annealer_params is None and simulation_params is None
        ):
            raise ValueError(
                "You must pass one of either annealer_params or simulation_params, "
                "not both"
            )

        elif annealer_params is not None:
            annealer_params_keys = ["schedule", "embedding"]
            for k in annealer_params_keys:
                if k not in annealer_params:
                    raise ValueError(
                        "Missing key in annealer_params. "
                        f"Required keys are {annealer_params_keys}"
                    )

            self.annealer_params = annealer_params
            self._initialize_annealer()

        elif simulation_params is not None:
            simulation_params_keys = ["beta"]
            for k in simulation_params_keys:
                if k not in simulation_params:
                    raise ValueError(
                        "Missing key in simulation_params. "
                        f"Required keys are {simulation_params_keys}"
                    )

            self.simulation_params = simulation_params
            self._pauli_kron = get_pauli_kron(self.n_visible, self.n_hidden)
            self.h_range = np.array(
                self.simulation_params.get("h_range", [-np.inf, np.inf])
            )
            self.J_range = np.array(
                self.simulation_params.get("J_range", [-np.inf, np.inf])
            )

    def sample(
        self,
        n_samples: int,
        answer_mode: str = "raw",
        use_gauge: bool = True,
        binary: bool = False,
    ) -> SampleOutput:
        """
        Generate samples using the model, either exact or from the annealer.

        Args:
            n_samples: Number of samples to generate (num_reads param in sample_ising).
            answer_mode: "raw" or "histogram".
            use_gauge: If True will use a random gauge transformation (recommended for
                more robust sample generation).
            binary: If true will convert the state vector values from {+1, -1} to
                {0, 1}.

        Returns:
            Dictionary (exact) or Ocean SDK SampleSet object (annealer).
        """
        if hasattr(self, "simulation_params"):
            return self._sample_simulation(n_samples, binary=binary)
        else:
            return self._sample_annealer(
                n_samples, answer_mode=answer_mode, use_gauge=use_gauge, binary=binary
            )

    def train(
        self,
        n_epochs: int = 100,
        learning_rate: float | Sequence[float] | np.ndarray = 1e-1,
        learning_rate_beta: float | Sequence[float] | np.ndarray = 1e-1,
        mini_batch_size: int = 10,
        n_samples: int = 10_000,
        callback: Callback | None = None,
    ) -> None:
        """
        Fits the model to the training data.

        Args:
            n_epochs: Number of epochs to train for.
            learning_rate: Learning rate. If a list or array, then it will represent the
                learning rate over the epochs, must be of length n_epochs.
            learning_rate_beta: Learning rate for the effective temperature. If a list
                or array, then it will represent the learning rate over the epochs,
                must be of length n_epochs.
                Note: It might be useful to use a larger learning_rate_beta in the
                beginning to help the model find a good temperature, then drop it after
                a number of epochs.
            mini_batch_size: Size of the mini-batches.
            n_samples: Number of samples to generate after every epoch. Used for
                computing β gradient, as well as the callback.
            callback: A function called at the end of each epoch. It takes the arguments
                (model, samples), and returns a dictionary with required keys ["value",
                "print"], where the "print" value is a string to be printed at the end
                of each epoch.

        Raises:
            ValueError: If learning_rate or learning_rate_beta is a sequence whose
                length is not n_epochs.
        """
        if isinstance(learning_rate, float):
            learning_rates: list[float] = [learning_rate] * n_epochs
        else:
            learning_rates = np.asarray(learning_rate, dtype=np.float64).tolist()
        if len(learning_rates) != n_epochs:
            raise ValueError(
                f"learning_rate has length {len(learning_rates)}, "
                f"expected n_epochs = {n_epochs}"
            )

        if isinstance(learning_rate_beta, float):
            beta_learning_rates: list[float] = [learning_rate_beta] * n_epochs
        else:
            beta_learning_rates = np.asarray(
                learning_rate_beta, dtype=np.float64
            ).tolist()
        if len(beta_learning_rates) != n_epochs:
            raise ValueError(
                f"learning_rate_beta has length {len(beta_learning_rates)}, "
                f"expected n_epochs = {n_epochs}"
            )

        if not hasattr(self, "callback_history"):
            self.callback_history: list[Mapping[str, Any]] = []
        callback_output: Mapping[str, Any] = {}

        for epoch in range(1, n_epochs + 1):
            start_time = time()

            # Set the effective learning rates
            self.learning_rate: float = learning_rates[epoch - 1]
            self.learning_rate_beta: float = beta_learning_rates[epoch - 1]

            # Compute and apply gradient updates for each mini batch
            for mini_batch_indices in self._random_mini_batch_indices(mini_batch_size):
                V_pos = self.V_train[mini_batch_indices]
                self._compute_positive_grads(V_pos)
                self._compute_negative_grads(V_pos.shape[0])
                self._apply_grads(self.learning_rate)
                self._check_h_and_H_ranges()

            # Update β
            samples = self.sample(n_samples)
            self._update_beta(samples)
            self._check_h_and_H_ranges()

            # Callback function
            if callback is not None:
                callback_output = callback(self, self._get_state_vectors(samples))
                self.callback_history.append(callback_output)

            # Print diagnostics
            end_time = time()
            print(
                f"[{type(self).__name__}] epoch {epoch}:",
                f"β = {self.beta:.3f},",
                f"learning rate = {learning_rates[epoch - 1]:.2e},",
                f"β learning rate = {beta_learning_rates[epoch - 1]:.2e},",
                f"epoch duration = {timedelta(seconds=end_time - start_time)}",
            )
            if callback is not None and "print" in callback_output:
                print(callback_output["print"])

    def save(self, file_path: str | Path, reinitialize_annealer: bool = True) -> None:
        """
        Saves the BQRBM model at file_path. Necessary because of pickling issues
        with the qpu and sampler objects.

        Args:
            file_path: Path to save the model to. Must be a Path object or a string with
                ".pkl" file extension.
            reinitialize_annealer: If True and has attribute self.annealer_params, then
                will call self._initialize_annealer() after saving.
        """
        if hasattr(self, "annealer_params"):
            self.qpu = None
            self.sampler = None

        save_artifact(self, file_path)

        if hasattr(self, "annealer_params") and reinitialize_annealer:
            self._initialize_annealer()

    @staticmethod
    def load(file_path: str | Path, initialize_annealer: bool = True) -> "BQRBM":
        """
        Loads the BQRBM model at file_path. Necessary because of pickling issues
        with the qpu and sampler objects.

        Args:
            file_path: Path to the model to load. Must be a Path object or a string with
                ".pkl" file extension.
            initialize_annealer: If True and has attribute self.annealer_params, then
                will call self._initialize_annealer().

        Returns:
            BQRBM instance loaded from the file path.
        """
        model = load_artifact(file_path)
        if hasattr(model, "annealer_params") and initialize_annealer:
            model._initialize_annealer()

        return model

    @property
    def h(self) -> np.ndarray:
        """
        Ising h values. Correspond to b_i = -β * B_freeze * h_i
        """
        return -self.b / (self.beta * self.B_freeze)

    @property
    def J(self) -> np.ndarray:
        """
        Ising J values. Correspond to w_ij = -β * B_freeze * J_ij
        """
        J = np.zeros((self.n_qubits, self.n_qubits))
        J[: self.n_visible, self.n_visible :] = -self.W / (self.beta * self.B_freeze)
        return J

    def _check_h_and_H_ranges(self) -> None:
        """
        Raises an exception if h and J values do not fall within h_range and J_range.

        Raises:
            ValueError: If the learned h and J values are outside of the allowed range.
        """
        h_satisfied = np.logical_and(
            self.h > self.h_range.min(), self.h < self.h_range.max()
        ).all()
        J_satisfied = np.logical_and(
            self.J_range.min() < self.J, self.J_range.max() > self.J
        ).all()

        if not h_satisfied or not J_satisfied:
            raise ValueError("Learned h and J values outside of allowed range")

    def _compute_positive_grads(self, V_pos: np.ndarray) -> None:
        """
        Computes the gradients for the positive phase, i.e., the expectation values
        w.r.t. the clamped Hamiltonian.

        Args:
            V_pos: Training data set mini-batch of positive vectors, shape
                (mini_batch_size, n_visible).
        """
        b_hidden = self.b[self.n_visible :] + V_pos @ self.W
        D = np.sqrt((self.beta * self.A_freeze) ** 2 + b_hidden**2)
        H_pos = (b_hidden / D) * np.tanh(D)

        self.grads["b_pos"] = np.concatenate((V_pos.mean(axis=0), H_pos.mean(axis=0)))
        self.grads["W_pos"] = V_pos.T @ H_pos / V_pos.shape[0]

    def _compute_negative_grads(self, n_samples: int) -> None:
        """
        Computes the gradients for the negative phase, i.e., the expectation values
        w.r.t. the model distribution.

        Args:
            n_samples: Number of samples to use in the negative phase.
        """
        samples = self.sample(n_samples)
        state_vectors = self._get_state_vectors(samples)

        V_neg = state_vectors[:, : self.n_visible]
        b_hidden = self.b[self.n_visible :] + V_neg @ self.W
        D = np.sqrt((self.beta * self.A_freeze) ** 2 + b_hidden**2)
        H_neg = (b_hidden / D) * np.tanh(D)

        self.grads["b_neg"] = np.concatenate((V_neg.mean(axis=0), H_neg.mean(axis=0)))
        self.grads["W_neg"] = V_neg.T @ H_neg / V_neg.shape[0]

    def _get_state_vectors(self, samples: SampleOutput) -> np.ndarray:
        """
        Get the state vectors from the samples (depending on exact or annealer
        generated).

        Args:
            samples: Return value out of BQRBM.sample().

        Returns:
            Array of state vectors, shape (n_samples, n_qubits).
        """
        if isinstance(samples, SampleSet):
            return samples.record.sample
        return samples["state_vectors"]

    def _initialize_annealer(self) -> None:
        """
        Initializes the D-Wave sampler using the fixed embedding provided to the object
        instantiation.
        """
        self.qpu = DWaveSampler(**self.annealer_params.get("qpu_params", {}))
        self.qpu.validate_anneal_schedule(self.annealer_params["schedule"])
        self.sampler = FixedEmbeddingComposite(
            self.qpu, self.annealer_params["embedding"]
        )
        self.h_range = np.array(self.qpu.properties["h_range"])
        self.J_range = np.array(self.qpu.properties["j_range"])

    def _initialize_weights_and_biases(self, mu: float = 0, sigma: float = 0.1) -> None:
        """
        Initializes the weights and biases. The biases are initialized to zero,
        and the weights are drawn from a normal distribution.

        Args:
            mu: Mean of the normal distribution of the weights.
            sigma: Standard deviation of the normal distribution of the weights.
        """
        self.b = np.zeros(self.n_qubits)
        self.W = self.rng.normal(mu, sigma, (self.n_visible, self.n_hidden))

    def _mean_classical_energy(
        self, V: np.ndarray, H: np.ndarray, VW: np.ndarray
    ) -> float:
        """
        Computes the mean classical energy w.r.t. the weights and biases over the
        provided visible and hidden unit state vectors.

        Args:
            V: Numpy array where the rows are visible units.
            H: Numpy array where the rows are hidden units.
            VW: V @ W (used to avoid double computation).

        Returns:
            Mean energy.
        """
        return (
            -(V @ self.b[: self.n_visible]).sum()
            - (H @ self.b[self.n_visible :]).sum()
            - np.einsum("kj,kj", VW, H)
        ) / V.shape[0]

    def _sample_annealer(
        self,
        n_samples: int,
        answer_mode: str = "raw",
        use_gauge: bool = True,
        binary: bool = False,
    ) -> SampleSet:
        """
        Obtain a sample set using the annealer.

        Args:
            n_samples: Number of samples to generate (num_reads param in sample_ising).
            answer_mode: "raw" or "histogram".
            use_gauge: If True will use a random gauge transformation (recommended for
                more robust sample generation).
            binary: If true will convert the state vector values from {-1, +1} to
                {0, 1}.

        Returns:
            Ocean SDK SampleSet object.
        """
        # Compute the h's and J's
        h = self.h
        J = self.J

        # Apply a random gauge
        gauge: np.ndarray | None = None
        if use_gauge:
            gauge = self.rng.choice([-1, 1], self.n_qubits)
            if gauge is None:
                raise RuntimeError("Failed to generate a random gauge")
            h *= gauge
            J *= np.outer(gauge, gauge)

        # Compute the chain strength
        chain_strength = self.annealer_params.get("relative_chain_strength")
        if chain_strength is not None:
            chain_strength *= max(np.abs(h).max(), np.abs(J).max())
            chain_strength = min(chain_strength, self.J_range.max())

        # Get samples from the annealer
        if self.sampler is None:
            raise RuntimeError("Annealer sampler is not initialized")
        samples = self.sampler.sample_ising(
            h,
            J,
            num_reads=n_samples,
            anneal_schedule=self.annealer_params.get("schedule"),
            chain_strength=chain_strength,
            answer_mode=answer_mode,
            auto_scale=False,
        )

        # Undo the gauge
        if gauge is not None:
            samples.record.sample *= gauge

        # Convert to binary if specified
        if binary:
            samples.record.sample = self._eigen_to_binary(samples.record.sample)

        return samples

    def _sample_simulation(
        self, n_samples: int, binary: bool = False
    ) -> SimulationSamples:
        """
        Sample using the exact computed probabilities.

        Args:
            n_samples: Number of samples to generate.
            binary: If true will convert the state vector values from {-1, +1} to
                {0, 1}.

        Returns:
            Dict with keys:
                - "E": Energies, i.e., the diagonal of the Hamiltonian.
                - "p": Probabilities, i.e., the diagonal of the density matrix.
                - "states": Integer state numbers of sampled states.
                - "state_vectors": Array of sampled state vectors, shape
                    (n_samples, n_qubits).
        """
        # Compute the h's and J's
        h = self.h
        J = self.J

        # Compute the Hamiltonian and density matrix
        H = compute_H(
            h, J, self.A_freeze, self.B_freeze, self.n_qubits, self._pauli_kron
        )
        rho = compute_rho(
            H, self.simulation_params["beta"], diagonal=(self.A_freeze == 0)
        )

        # Sample using the probabilities on the diagonal of rho
        probabilities = np.diag(rho).copy()
        states = self.rng.choice(
            range(2**self.n_qubits), size=n_samples, p=probabilities
        )
        state_vectors = self._binary_to_eigen(
            np.vstack([Discretizer.int_to_bit_vector(x, self.n_qubits) for x in states])
        )
        samples: SimulationSamples = {
            "E": np.diag(H).copy(),
            "p": probabilities,
            "states": states,
            "state_vectors": state_vectors,
        }

        # Convert to binary if specified
        if binary:
            samples["state_vectors"] = self._eigen_to_binary(samples["state_vectors"])

        return samples

    def _update_beta(self, samples: SampleOutput) -> None:
        """
        Updates the effective β = 1 / kT estimator. Used for scaling the coefficients
        sent to the annealer.

        Note:
            In its current form, this method only works when s_freeze = 1, i.e., when
            training and sampling from classical Boltzmann distributions.

        Args:
            samples: Samples to use for computing the mean energies w.r.t. the model.
        """
        # Compute the train energy
        VW_train = self.V_train @ self.W
        b_eff = self.b[self.n_visible :] + VW_train
        D = np.sqrt((self.beta * self.A_freeze) ** 2 + b_eff**2)
        H_train = (b_eff / D) * np.tanh(D)
        E_train = self._mean_classical_energy(self.V_train, H_train, VW_train)

        # Compute the model energy
        state_vectors = self._get_state_vectors(samples)
        V_model = state_vectors[:, : self.n_visible]
        VW_model = V_model @ self.W
        H_model = state_vectors[:, self.n_visible :]
        E_model = self._mean_classical_energy(V_model, H_model, VW_model)

        # Update the params
        self.beta = np.clip(
            self.beta + self.learning_rate_beta * (E_train - E_model),
            self.beta_range[0],
            self.beta_range[1],
        )
        self.beta_history.append(self.beta)

J property

Ising J values. Correspond to w_ij = -β * B_freeze * J_ij

h property

Ising h values. Correspond to b_i = -β * B_freeze * h_i

__init__(V_train, n_hidden, A_freeze, B_freeze, beta_initial=1.0, beta_range=[0.1, 10], annealer_params=None, simulation_params=None, seed=0)

Initializes the model.

Note

Exactly one of annealer_params or simulation_params must be provided. Whichever is provided determines whether the samples are generated by the annealer or the simulation, respectively.

The simulation works by exact computation of ρ(s, T) = e^{-β * H(s)} / Z. H(s) is determined by A_freeze and B_freeze.

Parameters:

Name	Type	Description	Default
V_train	ndarray	Training data set of visible vectors, shape (n_samples, n_visible). Values must be in {+1, -1}, or in {0, 1} which are converted to {+1, -1}.	required
n_hidden	int	Number of hidden units.	required
A_freeze	float	Value of A(s) at the freeze out point. Used for Γ = β * A_freeze. Units in GHz.	required
B_freeze	float	Value of B(s) at the freeze out point. Used for b_i = -β * B_freeze * h_i and w_ij = -β * B_freeze * J_ij. Units in GHz.	required
beta_initial	float	Initial value of the effective β. Units in 1/GHz.	1.0
beta_range	Sequence[float | int]	Range of allowed β values, used for making sure β is not updated to an infeasible value (e.g. negative).	[0.1, 10]
annealer_params	AnnealerParams | None	Dictionary with keys: - "schedule": List of (t, s) tuples defining the anneal schedule. - "embedding": Dict mapping the logical to physical qubits. - "relative_chain_strength" [optional]: Relative chain strength value. - "qpu_params" [optional]: Parameters dict to unpack to DWaveSampler(), e.g. {"region": "na-west-1", "solver": "Advantage_system4.1"}	None
simulation_params	Mapping[str, Any] | None	Dictionary with keys: - "beta": Effective β that the simulation generates samples at. - "h_range" [optional]: Allowed range of h values, defaults to (-inf, inf). - "J_range" [optional]: Allowed range of J values, defaults to (-inf, inf).	None
seed	int | None	Seed for the random number generator. Used for random minibatches, as well as the exact sampler.	0

Raises:

Type	Description
ValueError	If neither or both of annealer_params and simulation_params are provided, if a required key is missing from the provided params dict, if beta_initial is not positive, if beta_range does not satisfy 0 < min < max, or if the training data values are not in {-1, +1}.

Source code in src/qbm/models/BQRBM.py

def __init__(
    self,
    V_train: np.ndarray,
    n_hidden: int,
    A_freeze: float,
    B_freeze: float,
    beta_initial: float = 1.0,
    beta_range: Sequence[float | int] = [0.1, 10],
    annealer_params: AnnealerParams | None = None,
    simulation_params: Mapping[str, Any] | None = None,
    seed: int | None = 0,
) -> None:
    """
    Initializes the model.

    Note:
        Exactly one of annealer_params or simulation_params must be provided.
        Whichever is provided determines whether the samples are generated by
        the annealer or the simulation, respectively.

        The simulation works by exact computation of ρ(s, T) = e^{-β * H(s)} / Z.
        H(s) is determined by A_freeze and B_freeze.

    Args:
        V_train: Training data set of visible vectors, shape (n_samples, n_visible).
            Values must be in {+1, -1}, or in {0, 1} which are converted to
            {+1, -1}.
        n_hidden: Number of hidden units.
        A_freeze: Value of A(s) at the freeze out point. Used for Γ = β * A_freeze.
            Units in GHz.
        B_freeze: Value of B(s) at the freeze out point. Used for
            b_i = -β * B_freeze * h_i and w_ij = -β * B_freeze * J_ij. Units in GHz.
        beta_initial: Initial value of the effective β. Units in 1/GHz.
        beta_range: Range of allowed β values, used for making sure β is not updated
            to an infeasible value (e.g. negative).
        annealer_params: Dictionary with keys:
            - "schedule": List of (t, s) tuples defining the anneal schedule.
            - "embedding": Dict mapping the logical to physical qubits.
            - "relative_chain_strength" [optional]: Relative chain strength value.
            - "qpu_params" [optional]: Parameters dict to unpack to DWaveSampler(),
                e.g. {"region": "na-west-1", "solver": "Advantage_system4.1"}
        simulation_params: Dictionary with keys:
            - "beta": Effective β that the simulation generates samples at.
            - "h_range" [optional]: Allowed range of h values, defaults to
                (-inf, inf).
            - "J_range" [optional]: Allowed range of J values, defaults to
                (-inf, inf).
        seed: Seed for the random number generator. Used for random minibatches, as
            well as the exact sampler.

    Raises:
        ValueError: If neither or both of annealer_params and simulation_params are
            provided, if a required key is missing from the provided params dict,
            if beta_initial is not positive, if beta_range does not satisfy
            0 < min < max, or if the training data values are not in {-1, +1}.
    """
    self.qpu: DWaveSampler | None = None
    self.sampler: AnnealerSampler | None = None
    if beta_initial <= 0:
        raise ValueError(f"beta_initial must be positive (got {beta_initial})")
    if beta_range[0] <= 0 or beta_range[0] >= beta_range[1]:
        raise ValueError(
            f"beta_range must satisfy 0 < min < max (got {list(beta_range)})"
        )
    # Convert from binary to ±1 if necessary
    if set(np.unique(V_train)) == set([0, 1]):
        V_train = self._binary_to_eigen(V_train)
    if set(np.unique(V_train)) != set([-1, 1]):
        raise ValueError(
            "V_train values must be in {+1, -1} or {0, 1} "
            f"(got unique values {set(np.unique(V_train))})"
        )

    self.A_freeze = A_freeze
    self.B_freeze = B_freeze
    self.beta = beta_initial
    self.beta_range: Sequence[float | int] = beta_range
    self.beta_history: list[float] = [self.beta]
    super().__init__(V_train=V_train, n_hidden=n_hidden, seed=seed)

    # Check if requirements met to use annealer or simulation
    if (annealer_params is not None and simulation_params is not None) or (
        annealer_params is None and simulation_params is None
    ):
        raise ValueError(
            "You must pass one of either annealer_params or simulation_params, "
            "not both"
        )

    elif annealer_params is not None:
        annealer_params_keys = ["schedule", "embedding"]
        for k in annealer_params_keys:
            if k not in annealer_params:
                raise ValueError(
                    "Missing key in annealer_params. "
                    f"Required keys are {annealer_params_keys}"
                )

        self.annealer_params = annealer_params
        self._initialize_annealer()

    elif simulation_params is not None:
        simulation_params_keys = ["beta"]
        for k in simulation_params_keys:
            if k not in simulation_params:
                raise ValueError(
                    "Missing key in simulation_params. "
                    f"Required keys are {simulation_params_keys}"
                )

        self.simulation_params = simulation_params
        self._pauli_kron = get_pauli_kron(self.n_visible, self.n_hidden)
        self.h_range = np.array(
            self.simulation_params.get("h_range", [-np.inf, np.inf])
        )
        self.J_range = np.array(
            self.simulation_params.get("J_range", [-np.inf, np.inf])
        )

load(file_path, initialize_annealer=True) staticmethod

Loads the BQRBM model at file_path. Necessary because of pickling issues with the qpu and sampler objects.

Parameters:

Name	Type	Description	Default
file_path	str | Path	Path to the model to load. Must be a Path object or a string with ".pkl" file extension.	required
initialize_annealer	bool	If True and has attribute self.annealer_params, then will call self._initialize_annealer().	True

Returns:

Type	Description
BQRBM	BQRBM instance loaded from the file path.

Source code in src/qbm/models/BQRBM.py

@staticmethod
def load(file_path: str | Path, initialize_annealer: bool = True) -> "BQRBM":
    """
    Loads the BQRBM model at file_path. Necessary because of pickling issues
    with the qpu and sampler objects.

    Args:
        file_path: Path to the model to load. Must be a Path object or a string with
            ".pkl" file extension.
        initialize_annealer: If True and has attribute self.annealer_params, then
            will call self._initialize_annealer().

    Returns:
        BQRBM instance loaded from the file path.
    """
    model = load_artifact(file_path)
    if hasattr(model, "annealer_params") and initialize_annealer:
        model._initialize_annealer()

    return model

sample(n_samples, answer_mode='raw', use_gauge=True, binary=False)

Generate samples using the model, either exact or from the annealer.

Parameters:

Name	Type	Description	Default
n_samples	int	Number of samples to generate (num_reads param in sample_ising).	required
answer_mode	str	"raw" or "histogram".	'raw'
use_gauge	bool	If True will use a random gauge transformation (recommended for more robust sample generation).	True
binary	bool	If true will convert the state vector values from {+1, -1} to {0, 1}.	False

Returns:

Type	Description
SampleOutput	Dictionary (exact) or Ocean SDK SampleSet object (annealer).

Source code in src/qbm/models/BQRBM.py

def sample(
    self,
    n_samples: int,
    answer_mode: str = "raw",
    use_gauge: bool = True,
    binary: bool = False,
) -> SampleOutput:
    """
    Generate samples using the model, either exact or from the annealer.

    Args:
        n_samples: Number of samples to generate (num_reads param in sample_ising).
        answer_mode: "raw" or "histogram".
        use_gauge: If True will use a random gauge transformation (recommended for
            more robust sample generation).
        binary: If true will convert the state vector values from {+1, -1} to
            {0, 1}.

    Returns:
        Dictionary (exact) or Ocean SDK SampleSet object (annealer).
    """
    if hasattr(self, "simulation_params"):
        return self._sample_simulation(n_samples, binary=binary)
    else:
        return self._sample_annealer(
            n_samples, answer_mode=answer_mode, use_gauge=use_gauge, binary=binary
        )

save(file_path, reinitialize_annealer=True)

Saves the BQRBM model at file_path. Necessary because of pickling issues with the qpu and sampler objects.

Parameters:

Name	Type	Description	Default
file_path	str | Path	Path to save the model to. Must be a Path object or a string with ".pkl" file extension.	required
reinitialize_annealer	bool	If True and has attribute self.annealer_params, then will call self._initialize_annealer() after saving.	True

Source code in src/qbm/models/BQRBM.py

def save(self, file_path: str | Path, reinitialize_annealer: bool = True) -> None:
    """
    Saves the BQRBM model at file_path. Necessary because of pickling issues
    with the qpu and sampler objects.

    Args:
        file_path: Path to save the model to. Must be a Path object or a string with
            ".pkl" file extension.
        reinitialize_annealer: If True and has attribute self.annealer_params, then
            will call self._initialize_annealer() after saving.
    """
    if hasattr(self, "annealer_params"):
        self.qpu = None
        self.sampler = None

    save_artifact(self, file_path)

    if hasattr(self, "annealer_params") and reinitialize_annealer:
        self._initialize_annealer()

train(n_epochs=100, learning_rate=0.1, learning_rate_beta=0.1, mini_batch_size=10, n_samples=10000, callback=None)

Fits the model to the training data.

Parameters:

Name	Type	Description	Default
n_epochs	int	Number of epochs to train for.	100
learning_rate	float | Sequence[float] | ndarray	Learning rate. If a list or array, then it will represent the learning rate over the epochs, must be of length n_epochs.	0.1
learning_rate_beta	float | Sequence[float] | ndarray	Learning rate for the effective temperature. If a list or array, then it will represent the learning rate over the epochs, must be of length n_epochs. Note: It might be useful to use a larger learning_rate_beta in the beginning to help the model find a good temperature, then drop it after a number of epochs.	0.1
mini_batch_size	int	Size of the mini-batches.	10
n_samples	int	Number of samples to generate after every epoch. Used for computing β gradient, as well as the callback.	10000
callback	Callback | None	A function called at the end of each epoch. It takes the arguments (model, samples), and returns a dictionary with required keys ["value", "print"], where the "print" value is a string to be printed at the end of each epoch.	None

Raises:

Type	Description
ValueError	If learning_rate or learning_rate_beta is a sequence whose length is not n_epochs.

Source code in src/qbm/models/BQRBM.py

def train(
    self,
    n_epochs: int = 100,
    learning_rate: float | Sequence[float] | np.ndarray = 1e-1,
    learning_rate_beta: float | Sequence[float] | np.ndarray = 1e-1,
    mini_batch_size: int = 10,
    n_samples: int = 10_000,
    callback: Callback | None = None,
) -> None:
    """
    Fits the model to the training data.

    Args:
        n_epochs: Number of epochs to train for.
        learning_rate: Learning rate. If a list or array, then it will represent the
            learning rate over the epochs, must be of length n_epochs.
        learning_rate_beta: Learning rate for the effective temperature. If a list
            or array, then it will represent the learning rate over the epochs,
            must be of length n_epochs.
            Note: It might be useful to use a larger learning_rate_beta in the
            beginning to help the model find a good temperature, then drop it after
            a number of epochs.
        mini_batch_size: Size of the mini-batches.
        n_samples: Number of samples to generate after every epoch. Used for
            computing β gradient, as well as the callback.
        callback: A function called at the end of each epoch. It takes the arguments
            (model, samples), and returns a dictionary with required keys ["value",
            "print"], where the "print" value is a string to be printed at the end
            of each epoch.

    Raises:
        ValueError: If learning_rate or learning_rate_beta is a sequence whose
            length is not n_epochs.
    """
    if isinstance(learning_rate, float):
        learning_rates: list[float] = [learning_rate] * n_epochs
    else:
        learning_rates = np.asarray(learning_rate, dtype=np.float64).tolist()
    if len(learning_rates) != n_epochs:
        raise ValueError(
            f"learning_rate has length {len(learning_rates)}, "
            f"expected n_epochs = {n_epochs}"
        )

    if isinstance(learning_rate_beta, float):
        beta_learning_rates: list[float] = [learning_rate_beta] * n_epochs
    else:
        beta_learning_rates = np.asarray(
            learning_rate_beta, dtype=np.float64
        ).tolist()
    if len(beta_learning_rates) != n_epochs:
        raise ValueError(
            f"learning_rate_beta has length {len(beta_learning_rates)}, "
            f"expected n_epochs = {n_epochs}"
        )

    if not hasattr(self, "callback_history"):
        self.callback_history: list[Mapping[str, Any]] = []
    callback_output: Mapping[str, Any] = {}

    for epoch in range(1, n_epochs + 1):
        start_time = time()

        # Set the effective learning rates
        self.learning_rate: float = learning_rates[epoch - 1]
        self.learning_rate_beta: float = beta_learning_rates[epoch - 1]

        # Compute and apply gradient updates for each mini batch
        for mini_batch_indices in self._random_mini_batch_indices(mini_batch_size):
            V_pos = self.V_train[mini_batch_indices]
            self._compute_positive_grads(V_pos)
            self._compute_negative_grads(V_pos.shape[0])
            self._apply_grads(self.learning_rate)
            self._check_h_and_H_ranges()

        # Update β
        samples = self.sample(n_samples)
        self._update_beta(samples)
        self._check_h_and_H_ranges()

        # Callback function
        if callback is not None:
            callback_output = callback(self, self._get_state_vectors(samples))
            self.callback_history.append(callback_output)

        # Print diagnostics
        end_time = time()
        print(
            f"[{type(self).__name__}] epoch {epoch}:",
            f"β = {self.beta:.3f},",
            f"learning rate = {learning_rates[epoch - 1]:.2e},",
            f"β learning rate = {beta_learning_rates[epoch - 1]:.2e},",
            f"epoch duration = {timedelta(seconds=end_time - start_time)}",
        )
        if callback is not None and "print" in callback_output:
            print(callback_output["print"])