Integration with qiskit-symb

doc/requirements.txt

@@ -23,3 +23,4 @@ moabb==1.1.0
 pyriemann==0.6
 docplex>=2.21.207
 firebase_admin==6.5.0
+qiskit-symb

pyriemann_qiskit/classification.py

@@ -81,7 +81,7 @@ class QuanticClassifierBase(BaseEstimator, ClassifierMixin):
         If true, will output all intermediate results and logs.
     shots : int, default=1024
         Number of repetitions of each circuit, for sampling.
-    gen_feature_map : Callable[int, QuantumCircuit | FeatureMap], \
+    gen_feature_map : Callable[[int, str], QuantumCircuit | FeatureMap], \
                       default=Callable[int, ZZFeatureMap]
         Function generating a feature map to encode data into a quantum state.
     seed : int | None, default=None
@@ -357,7 +357,7 @@ class QuanticSVM(QuanticClassifierBase):
         If true, will output all intermediate results and logs.
     shots : int, default=1024
         Number of repetitions of each circuit, for sampling.
-    gen_feature_map : Callable[int, QuantumCircuit | FeatureMap], \
+    gen_feature_map : Callable[[int, str], QuantumCircuit | FeatureMap], \
                       default=Callable[int, ZZFeatureMap]
         Function generating a feature map to encode data into a quantum state.
     seed : int | None, default=None
@@ -420,6 +420,7 @@ def _init_algo(self, n_features):
         if self.quantum:
             quantum_kernel = get_quantum_kernel(
                 self._feature_map,
+                self.gen_feature_map,
                 self._quantum_instance,
                 self.use_fidelity_state_vector_kernel,
             )
@@ -496,7 +497,7 @@ class QuanticVQC(QuanticClassifierBase):
         If true, will output all intermediate results and logs
     shots : int, default=1024
         Number of repetitions of each circuit, for sampling
-    gen_feature_map : Callable[int, QuantumCircuit | FeatureMap], \
+    gen_feature_map : Callable[[int, str], QuantumCircuit | FeatureMap], \
                       default=Callable[int, ZZFeatureMap]
         Function generating a feature map to encode data into a quantum state.
     seed : int | None, default=None

pyriemann_qiskit/pipelines.py

@@ -458,7 +458,7 @@ class QuantumMDMVotingClassifier(BasePipeline):
         If true, will output all intermediate results and logs.
     shots : int (default:1024)
         Number of repetitions of each circuit, for sampling.
-    gen_feature_map : Callable[int, QuantumCircuit | FeatureMap] \
+    gen_feature_map : Callable[[int, str], QuantumCircuit | FeatureMap] \
                       (default : Callable[int, ZZFeatureMap])
         Function generating a feature map to encode data into a quantum state.
 

pyriemann_qiskit/utils/hyper_params_factory.py

@@ -25,8 +25,11 @@ def gen_x_feature_map(reps=2):
 
     Returns
     -------
-    ret : XFeatureMap
-        An instance of XFeatureMap.
+    ret : Callable[[int, str], XFeatureMap]
+        A callable that takes into arguments:
+        - the number of features
+        - the prefix string for the parameters
+        And returns an instance of XFeatureMap.
 
     Raises
     ------
@@ -41,16 +44,18 @@ def gen_x_feature_map(reps=2):
     Notes
     -----
     .. versionadded:: 0.2.0
+    .. versionchanged:: 0.4.0
+        Possibility to specify parameter prefix.
     """
     if reps < 1:
         raise ValueError(f"Parameter reps must be superior or equal to 1 (Got {reps})")
 
-    return lambda n_features: PauliFeatureMap(
+    return lambda n_features, param_prefix="x": PauliFeatureMap(
         feature_dimension=n_features,
         paulis=["X"],
         reps=reps,
         data_map_func=None,
-        parameter_prefix="x",
+        parameter_prefix=param_prefix,
         insert_barriers=False,
         name="XFeatureMap",
     )
@@ -70,8 +75,11 @@ def gen_z_feature_map(reps=2):
 
     Returns
     -------
-    ret : ZFeatureMap
-        An instance of ZFeatureMap.
+    ret : Callable[[int, str], ZFeatureMap]
+        A callable that takes into arguments:
+        - the number of features
+        - the prefix string for the parameters
+        And returns an instance of ZFeatureMap.
 
     Raises
     ------
@@ -86,11 +94,17 @@ def gen_z_feature_map(reps=2):
     Notes
     -----
     .. versionadded:: 0.2.0
+    .. versionchanged:: 0.4.0
+        Possibility to specify parameter prefix.
     """
     if reps < 1:
         raise ValueError(f"Parameter reps must be superior or equal to 1 (Got {reps})")
 
-    return lambda n_features: ZFeatureMap(feature_dimension=n_features, reps=reps)
+    return lambda n_features, param_prefix="x": ZFeatureMap(
+        feature_dimension=n_features,
+        reps=reps,
+        parameter_prefix=param_prefix,
+    )
 
 
 def gen_zz_feature_map(reps=2, entanglement="linear"):
@@ -114,8 +128,11 @@ def gen_zz_feature_map(reps=2, entanglement="linear"):
 
     Returns
     -------
-    ret : ZZFeatureMap
-        An instance of ZZFeatureMap.
+    ret : Callable[[int, str], ZZFeatureMap]
+        A callable that takes into arguments:
+        - the number of features
+        - the prefix string for the parameters
+        And returns an instance of ZZFeatureMap.
 
     Raises
     ------
@@ -132,12 +149,17 @@ def gen_zz_feature_map(reps=2, entanglement="linear"):
     Notes
     -----
     .. versionadded:: 0.0.1
+    .. versionchanged:: 0.4.0
+        Possibility to specify parameter prefix.
     """
     if reps < 1:
         raise ValueError(f"Parameter reps must be superior or equal to 1 (Got {reps})")
 
-    return lambda n_features: ZZFeatureMap(
-        feature_dimension=n_features, reps=reps, entanglement=entanglement
+    return lambda n_features, param_prefix="x": ZZFeatureMap(
+        feature_dimension=n_features,
+        reps=reps,
+        entanglement=entanglement,
+        parameter_prefix=param_prefix,
     )
 
 

pyriemann_qiskit/utils/quantum_provider.py

@@ -1,7 +1,11 @@
 """Module containing helpers for IBM quantum backends
    providers and simulators."""
 
+import joblib
 import logging
+import numpy as np
+import os
+import pickle
 
 from qiskit_aer import AerSimulator
 from qiskit_aer.quantum_info import AerStatevector
@@ -11,6 +15,117 @@
     FidelityStatevectorKernel,
     FidelityQuantumKernel,
 )
+from qiskit_symb.quantum_info import Statevector
+
+
+class SymbFidelityStatevectorKernel:
+
+    """Symbolic Statevector kernel
+
+    An implementation of the quantum kernel for classically simulated
+    state vectors [1]_ using qiskit-symb for symbolic representation
+    of statevectors [2]_.
+
+    Here, the kernel function is defined as the overlap of two simulated quantum
+    statevectors produced by a parametrized quantum circuit (called feature map) [1]_.
+
+    Notes
+    -----
+    .. versionadded:: 0.4.0
+
+    Parameters
+    ----------
+    feature_map: QuantumCircuit | FeatureMap
+        An instance of a feature map.
+    gen_feature_map: Callable[[int, str], QuantumCircuit | FeatureMap], \
+                      default=Callable[int, ZZFeatureMap]
+        Function generating a feature map to encode data into a quantum state.
+    n_jobs: int
+        The number of job for fidelity evaluation.
+        If null or negative, the number of jobs is set to 1
+        If set to 1, evaluation will run on the main thread.
+
+    References
+    ----------
+    .. [1] \
+    https://github.com/qiskit-community/qiskit-machine-learning/blob/30dad803e9457f955464220eddc1e55a65452bbc/qiskit_machine_learning/kernels/fidelity_statevector_kernel.py#L31
+    .. [2] https://github.com/SimoneGasperini/qiskit-symb/issues/6
+
+    """
+
+    def __init__(self, feature_map, gen_feature_map, n_jobs=1):
+        self.n_jobs = n_jobs if n_jobs >= 1 else 1
+        cached_file = os.path.join(
+            "symb_statevectors", f"{feature_map.name}-{feature_map.reps}"
+        )
+
+        if os.path.isfile(cached_file):
+            print("Loading symbolic Statevector from cache")
+            file = open(cached_file, "rb")
+            sv = pickle.load(file)
+        else:
+            print("Computing symbolic Statevector")
+            fm2 = gen_feature_map(feature_map.num_qubits, "b")
+            self.circuit = feature_map.compose(fm2.inverse()).decompose()
+            sv = Statevector(self.circuit)
+            print(f"Dumping to {cached_file}")
+            file = open(cached_file, "wb")
+            pickle.dump(sv, file)
+
+        self.function = sv.to_lambda()
+
+    """Evaluate the quantum kernel.
+
+    Returns
+    -------
+    kernel : ndarray, shape (len(x_vec), len(y_vec))
+        The kernel matrix.
+
+    Notes
+    -----
+    .. versionadded:: 0.4.0
+    """
+
+    def evaluate(self, x_vec, y_vec=None):
+        if y_vec is None:
+            y_vec = x_vec
+
+        x_vec_len = len(x_vec)
+        y_vec_len = len(y_vec)
+
+        is_sim = x_vec_len == y_vec_len and (x_vec == y_vec).all()
+
+        kernel_matrix = np.zeros((x_vec_len, y_vec_len))
+
+        chunck = x_vec_len // self.n_jobs
+
+        def compute_fidelity_partial_matrix(i_thread):
+            for i in range(i_thread * chunck, (i_thread + 1) * chunck):
+                x = x_vec[i]
+                for j in range(i if is_sim else y_vec_len):
+                    y = y_vec[j]
+                    if isinstance(x, np.float64):
+                        # Pegagos implementation
+                        fidelity = abs(self.function(x, y)[0, 0]) ** 2
+                    else:
+                        fidelity = abs(self.function(*x, *y)[0, 0]) ** 2
+
+                    kernel_matrix[i, j] = fidelity
+                    if is_sim:
+                        kernel_matrix[j, i] = fidelity
+            return kernel_matrix
+
+        if self.n_jobs == 1:
+            return compute_fidelity_partial_matrix(0)
+        else:
+            print("n_jobs greater than 1, parallelizing")
+            results = joblib.Parallel(n_jobs=self.n_jobs)(
+                joblib.delayed(compute_fidelity_partial_matrix)(i_thread)
+                for i_thread in range(self.n_jobs)
+            )
+            for result in results:
+                kernel_matrix += result
+            return kernel_matrix
 
 
 def get_provider():
@@ -93,12 +208,16 @@ def get_device(provider, min_qubits):
     )
 
 
-def get_quantum_kernel(feature_map, quantum_instance, use_fidelity_state_vector_kernel):
+def get_quantum_kernel(
+    feature_map, gen_feature_map, quantum_instance, use_fidelity_state_vector_kernel
+):
     """Get a quantum kernel
 
     Return an instance of FidelityQuantumKernel or
     FidelityStatevectorKernel (in the case of a simulation).
 
+    For simulation with a small number of qubits (< 9), qiskit-symb is used.
+
     Parameters
     ----------
     feature_map: QuantumCircuit | FeatureMap
@@ -116,26 +235,40 @@ def get_quantum_kernel(feature_map, quantum_instance, use_fidelity_state_vector_
     Notes
     -----
     .. versionadded:: 0.3.0
+    .. versionchanged:: 0.4.0
+        Add support for qiskit-symb
     """
     if use_fidelity_state_vector_kernel and isinstance(
         quantum_instance._backend, AerSimulator
     ):
-        logging.log(
-            logging.WARN,
-            """FidelityQuantumKernel skipped because of time.
+        # For simulation:
+        if feature_map.num_qubits <= 9:
+            # With a small number of qubits, let's use qiskit-symb
+            # See:
+            # https://medium.com/qiskit/qiskit-symb-a-qiskit-ecosystem-package-for-symbolic-quantum-computation-b6b4407fa705
+            kernel = SymbFidelityStatevectorKernel(
+                feature_map, gen_feature_map, n_jobs=4
+            )
+            logging.log(
+                logging.WARN,
+                """Using SymbFidelityStatevectorKernel""",
+            )
+        else:
+            # For a larger number of qubits,
+            # we will not use FidelityQuantumKernel as it is slow. See
+            # https://github.com/qiskit-community/qiskit-machine-learning/issues/547#issuecomment-1486527297
+            kernel = FidelityStatevectorKernel(
+                feature_map=feature_map,
+                statevector_type=AerStatevector,
+                shots=quantum_instance.options["shots"],
+            )
+            logging.log(
+                logging.WARN,
+                """FidelityQuantumKernel skipped because of time.
                     Using FidelityStatevectorKernel with AerStatevector.
                     Seed cannot be set with FidelityStatevectorKernel.
                     Increase the number of shots to diminish the noise.""",
-        )
-
-        # if this is a simulation,
-        # we will not use FidelityQuantumKernel as it is slow. See
-        # https://github.com/qiskit-community/qiskit-machine-learning/issues/547#issuecomment-1486527297
-        kernel = FidelityStatevectorKernel(
-            feature_map=feature_map,
-            statevector_type=AerStatevector,
-            shots=quantum_instance.options["shots"],
-        )
+            )
     else:
         kernel = FidelityQuantumKernel(
             feature_map=feature_map, fidelity=ComputeUncompute(quantum_instance)

requirements.txt

@@ -14,3 +14,4 @@ docplex>=2.21.207
 firebase_admin==6.5.0
 tqdm
 pandas
+qiskit-symb

setup.py

@@ -49,7 +49,8 @@
                         'firebase_admin==6.5.0',
                         'scikit-learn==1.5.1',
                         'tqdm',
-                        'pandas'
+                        'pandas',
+                        'qiskit-symb'
                         ],
       extras_require={'docs': [
                                 'sphinx-gallery',