# ====== Base class imports ======
import os
import time
import warnings
import numpy as np
import pandas as pd
from sklearn.ensemble import RandomForestClassifier
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score, auc, classification_report, confusion_matrix, f1_score
from sklearn.model_selection import GridSearchCV, RandomizedSearchCV
from sklearn.neural_network import MLPClassifier
from sklearn.svm import SVC
try:
from xgboost import XGBClassifier
XGBOOST_AVAILABLE = True
except Exception:
XGBOOST_AVAILABLE = False
XGBClassifier = None # type: ignore
# from qiskit.primitives import Sampler
from functools import reduce
# ====== Qiskit imports ======
from qiskit import QuantumCircuit
from qiskit.quantum_info import Pauli
from qiskit_ibm_runtime.exceptions import IBMRuntimeError, RuntimeJobFailureError
from sklearn import svm
import qbiocode.utils.qutils as qutils
# ====== Additional local imports ======
from qbiocode.evaluation.model_evaluation import modeleval
[docs]
def compute_pqk(
X_train,
X_test,
y_train,
y_test,
args,
model="PQK",
data_key="",
verbose=False,
encoding="Z",
primitive="estimator",
entanglement="linear",
reps=2,
):
"""
This function generates quantum circuits, computes projections of the data onto these circuits,
and evaluates the performance of classical machine learning models on the projected data.
It uses a feature map to encode the data into quantum states and then measures the expectation values
of Pauli operators to obtain the features. The classical models are trained on the projected training data and
evaluated on the projected test data. The function returns evaluation metrics and model parameters.
This function requires a quantum backend (simulator or real quantum hardware) for execution.
It supports various configurations such as encoding methods, entanglement strategies, and repetitions
of the feature map. The results are saved to files for training and test projections, which are reused
if they already exist to avoid redundant computations.
This function is part of the main quantum machine learning pipeline (QProfiler.py) and is intended for use in supervised learning tasks.
It leverages quantum computing to enhance feature extraction and classification performance on complex datasets.
The function returns the performance results, including accuracy, F1-score, AUC, runtime, as well as model parameters, and other relevant metrics.
Args:
X_train (np.ndarray): Training data features.
X_test (np.ndarray): Test data features.
y_train (np.ndarray): Training data labels.
y_test (np.ndarray): Test data labels.
args (dict): Arguments containing backend and other configurations.
model (str): Model type, default is 'PQK'.
data_key (str): Key for the dataset, default is ''.
verbose (bool): If True, print additional information, default is False.
encoding (str): Encoding method for the quantum circuit, default is 'Z'.
primitive (str): Primitive type to use, default is 'estimator'.
entanglement (str): Entanglement strategy, default is 'linear'.
reps (int): Number of repetitions for the feature map, default is 2.
Returns:
modeleval (pd.DataFrame): A DataFrame containing evaluation metrics and model parameters for all models.
"""
classical_models = ["svc"]
beg_time = time.time()
feat_dimension = X_train.shape[1]
projection_dir = os.path.expanduser(args.get("pqk_projection_dir", "pqk_projections"))
if not os.path.exists(projection_dir):
os.makedirs(projection_dir)
file_projection_train = os.path.join(
projection_dir, "pqk_projection_" + data_key + "_train.npy"
)
file_projection_test = os.path.join(
projection_dir, "pqk_projection_" + data_key + "_test.npy"
)
checkpoint_dir = os.path.join(projection_dir, "checkpoints")
os.makedirs(checkpoint_dir, exist_ok=True)
checkpoint_every = int(args.get("pqk_checkpoint_every", 1))
session_chunk_size = int(args.get("pqk_session_chunk_size", 100))
max_runtime_retries = int(args.get("pqk_runtime_max_retries", 3))
def _checkpoint_path(final_path):
base = os.path.basename(final_path)
return os.path.join(checkpoint_dir, base.replace(".npy", ".partial.npy"))
def _save_projection_array(path, projections):
tmp_path = path + ".tmp.npy"
np.save(tmp_path, np.asarray(projections))
os.replace(tmp_path, path)
def _load_checkpoint(path, expected_len):
if not os.path.exists(path):
return []
projections = np.load(path, allow_pickle=False)
if len(projections) > expected_len:
raise ValueError(
f"Checkpoint {path} has {len(projections)} rows, "
f"but the dataset only has {expected_len} rows."
)
return list(projections)
def _validate_projection_file(path, expected_len):
if not os.path.exists(path):
return
projections = np.load(path, allow_pickle=False)
if len(projections) != expected_len:
raise ValueError(
f"Projection file {path} has {len(projections)} rows, "
f"but the current dataset expects {expected_len} rows. "
"Remove this projection file or use a different pqk_projection_dir."
)
def _is_closed_session_error(exc):
return isinstance(exc, IBMRuntimeError) and (
"Session has been closed" in str(exc) or '"code":1217' in str(exc)
)
def _is_retryable_runtime_error(exc):
return isinstance(exc, RuntimeJobFailureError) and (
"Temporary Internal Error" in str(exc) or "Error code 9707" in str(exc)
)
def _close_session(session):
if not isinstance(session, type(None)):
session.close()
def _refresh_runtime(session):
_close_session(session)
_, new_session, new_prim = qutils.get_backend_session(
args, "estimator", num_qubits=num_qubits
)
return new_session, new_prim
# This function ensures that all multiplicative factors of data features inside single qubit gates are 1.0
def data_map_func(x: np.ndarray):
"""
Define a function map from R^n to R.
Args:
x: data
Returns:
the mapped value (float or Parameter expression)
"""
coeff = x[0] / 2 if len(x) == 1 else reduce(lambda m, n: (m * n) / 2, x)
# Check if coeff is a numeric type before converting to float
# If it's a Parameter expression, return it as-is for Qiskit to handle
try:
return float(coeff)
except (TypeError, ValueError):
# If conversion fails, it's likely a Parameter expression
return coeff
# choose a method for mapping your features onto the circuit
feature_map, _ = qutils.get_feature_map(
feature_map=encoding,
feat_dimension=X_train.shape[1],
reps=reps,
entanglement=entanglement,
data_map_func=data_map_func,
)
# Build quantum circuit
circuit = QuantumCircuit(feature_map.num_qubits)
circuit.compose(feature_map, inplace=True)
num_qubits = circuit.num_qubits
_validate_projection_file(file_projection_train, len(X_train))
_validate_projection_file(file_projection_test, len(X_test))
if (not os.path.exists(file_projection_train)) | (not os.path.exists(file_projection_test)):
# Generate the backend, session and primitive
backend, session, prim = qutils.get_backend_session(
args, "estimator", num_qubits=num_qubits
)
# Transpile
if args["backend"] != "simulator":
circuit = qutils.transpile_circuit(
circuit, opt_level=3, backend=backend, PT=True, initial_layout=None
)
# Set the global phase to 0 to avoid header size issues
circuit.global_phase = 0
for f_tr, dat in [
(file_projection_train, X_train.copy()),
(file_projection_test, X_test.copy()),
]:
if not os.path.exists(f_tr):
projections = []
# Identity operator on all qubits
id = "I" * feat_dimension
# We group all commuting observables
# These groups are the Pauli X, Y and Z operators on individual qubits
# Apply the circuit layout to the observable if mapped to device
if args["backend"] != "simulator":
observables_x = []
observables_y = []
observables_z = []
for i in range(feat_dimension):
observables_x.append(
Pauli(id[:i] + "X" + id[(i + 1) :]).apply_layout(
circuit.layout, num_qubits=backend.num_qubits
)
)
observables_y.append(
Pauli(id[:i] + "Y" + id[(i + 1) :]).apply_layout(
circuit.layout, num_qubits=backend.num_qubits
)
)
observables_z.append(
Pauli(id[:i] + "Z" + id[(i + 1) :]).apply_layout(
circuit.layout, num_qubits=backend.num_qubits
)
)
else:
observables_x = [
Pauli(id[:i] + "X" + id[(i + 1) :]) for i in range(feat_dimension)
]
observables_y = [
Pauli(id[:i] + "Y" + id[(i + 1) :]) for i in range(feat_dimension)
]
observables_z = [
Pauli(id[:i] + "Z" + id[(i + 1) :]) for i in range(feat_dimension)
]
checkpoint_file = _checkpoint_path(f_tr)
projections = _load_checkpoint(checkpoint_file, len(dat))
if projections:
print(
f"Resuming {os.path.basename(f_tr)} from "
f"datapoint {len(projections)}"
)
datapoints_in_session = 0
for i in range(len(projections), len(dat)):
if i % 100 == 0:
print(f"at datapoint {str(i)}")
if (
session is not None
and session_chunk_size > 0
and datapoints_in_session >= session_chunk_size
):
session, prim = _refresh_runtime(session)
datapoints_in_session = 0
# Get training sample
parameters = dat[i]
# We define the primitive unified blocs (PUBs) consisting of the embedding circuit,
# set of observables and the circuit parameters
pub_x = (circuit, observables_x, parameters)
pub_y = (circuit, observables_y, parameters)
pub_z = (circuit, observables_z, parameters)
retry_count = 0
while True:
try:
job = prim.run([pub_x, pub_y, pub_z])
job_result = job.result()
job_result_x = job_result[0].data.evs
job_result_y = job_result[1].data.evs
job_result_z = job_result[2].data.evs
break
except Exception as exc:
_save_projection_array(checkpoint_file, projections)
if session is not None and _is_closed_session_error(exc):
session, prim = _refresh_runtime(session)
datapoints_in_session = 0
continue
if (
session is not None
and _is_retryable_runtime_error(exc)
and retry_count < max_runtime_retries
):
retry_count += 1
print(
f"Retrying datapoint {i} after temporary runtime "
f"failure ({retry_count}/{max_runtime_retries})"
)
session, prim = _refresh_runtime(session)
datapoints_in_session = 0
continue
raise
# Record <X>, <Y> and <Z> on all qubits for the current datapoint
projections.append([job_result_x, job_result_y, job_result_z])
datapoints_in_session += 1
if checkpoint_every > 0 and len(projections) % checkpoint_every == 0:
_save_projection_array(checkpoint_file, projections)
_save_projection_array(f_tr, projections)
if os.path.exists(checkpoint_file):
os.remove(checkpoint_file)
if not isinstance(session, type(None)):
session.close()
# Load computed projections
projections_train = np.load(file_projection_train)
projections_train = np.array(projections_train).reshape(len(projections_train), -1)
projections_test = np.load(file_projection_test)
projections_test = np.array(projections_test).reshape(len(projections_test), -1)
model = create_svc_model(args["seed"])
method_pqk = "pqk"
model.fit(projections_train, y_train)
y_predicted = model.predict(projections_test)
hyperparameters = {
"feature_map": feature_map.__class__.__name__,
"feature_map_reps": reps,
"entanglement": entanglement,
"best_params": model.best_params_,
# Add other hyperparameters as needed
}
model_params = hyperparameters
return modeleval(
y_test, y_predicted, beg_time, params=model_params, args=args, model=method_pqk, verbose=verbose
)
[docs]
def create_svc_model(seed):
svc_param_distributions = {
"C": [0.1, 1, 10, 100],
"gamma": [0.001, 0.01, 0.1, 1],
"kernel": ["linear", "rbf", "poly", "sigmoid"],
}
# Initialize the SVC
svc = SVC(random_state=seed)
# Initialize RandomizedSearchCV
svc_model = RandomizedSearchCV(
estimator=svc,
param_distributions=svc_param_distributions,
n_iter=40,
cv=5,
random_state=seed,
n_jobs=-1,
)
return svc_model
