# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

# + [markdown] slideshow={"slide_type": "slide"}
# # Quantum Ensemble as Simple Averaging
#
# ### Theory and Implementation 
#
# ##### (Fixed $U_{(i,j)}$ for independent  quantum trajectories)

# + [markdown] slideshow={"slide_type": "skip"}
# $$\newcommand{\ket}[1]{\left|{#1}\right\rangle}$$
# $$\newcommand{\bra}[1]{\left\langle{#1}\right|}$$
# $$\newcommand{\braket}[1]{\left\langle{#1}\right\rangle}$$

# + [markdown] slideshow={"slide_type": "subslide"}
# ### (Step 1) State Preparation

# + [markdown] slideshow={"slide_type": "subslide"}
# # + Given a 2-qubits $control$ register $(d=2) \implies $  ensemble of $B=2^2$ classifiers
#
#
# # + $data$ register: *qubit encoding strategy* $\implies$  $N$ $2-$dimensional observations are encoded using $2 \times N$ qubits:
# $$
# \begin{align}
#     \text{data register: } \underset{features}{\big( \overset{4}{\underset{i=1}{\otimes}} \left| x_i \right\rangle \big)}\otimes  \underset{ labels}{\big( \overset{4}{\underset{i=1}{\otimes}} \left| y_i \right\rangle \big)}
# \end{align}
# $$
#
# where $x_i$ and $y_i$ are encoded as follows:
# $$
# \begin{align}
#     \left| x_i \right\rangle = x_{i,1}\left| 0 \right\rangle +  x_{i,2}\left| 1 \right\rangle
# \end{align}
# $$
#
# and
#
# if $\left| y_i \right\rangle=\left| 0 \right\rangle$ the $i$-th observation belongs to the class $0$. Otherwise if $\left| y_i \right\rangle=\left| 1 \right\rangle$ the $i$-th observation belongs to the class $1$. 

# + [markdown] slideshow={"slide_type": "subslide"}
# ### (Step 1) State Preparation

# + [markdown] slideshow={"slide_type": "fragment"}
# $$
# \begin{align*}
#     \left|\Phi_0\right\rangle &=
#     \big( H^{\otimes 2} \otimes S_{(x,y)} \big)\left|0\right\rangle \otimes \left|0\right\rangle \otimes \left|0\right\rangle \nonumber \\ 
#         & =
#    \left|c_1\right\rangle \otimes \left|c_2\right\rangle \otimes \left|x\right\rangle \left|y\right\rangle \nonumber\\
#     & =
#     \frac{1}{\sqrt{2}}\big(\left|0\right\rangle+\left|1\right\rangle\big) \otimes \frac{1}{\sqrt{2}}\big(\left|0\right\rangle+\left|1\right\rangle\big) \otimes \left|x_0,x_1,x_2,x_3\right\rangle \left|y_0,y_1,y_2,y_3\right\rangle \end{align*}
# $$

# + [markdown] slideshow={"slide_type": "fragment"}
# where $S_x$ is the routine which encodes in the amplitudes of a qubit a real vector $x$ and $H$ is the Hadamard transformation.

# + [markdown] slideshow={"slide_type": "subslide"}
# ### (Step 2) Sampling in Superposition
#
# The second step regards the generation of $2^d$ different transformations of the training set in superposition, each entangled with a state of the control register. To this end, $d$ steps are necessary, where each step consists in the entanglement of the $i$-th control qubit with two transformations of $\left|x,y\right\rangle$ based on two random unitaries, $U_{(i,1)}$ and $U_{(i,2)}$, for $i = 1,2$.
#
# As shown in the **Appendix A**, *Sampling in Superposition* step leads to the following quantum state:
# \begin{align}
# \ket{\Phi_{2}}
# =  \frac{1}{2}\Big[
#  \hspace{.2em} &\ket{00} U_{(2,1)}U_{(1,1)}\ket{x,y}    
# \nonumber \\ + & 
# \ket{01} U_{(2,1)}U_{(1,2)}\ket{x,y} 
# \nonumber \\ + & 
# \ket{10} U_{(2,2)}U_{(1,1)}\ket{x,y}
# \nonumber \\ + & 
# \ket{11} U_{(2,2)}U_{(1,2)}\ket{x,y} 
#          \Big]  \nonumber \\
#  & \hspace{-2.75em} = \frac{1}{\sqrt{4}} \sum_{b=1}^{4} \ket{b} V_b\ket{x,y}
# \end{align}

# + [markdown] slideshow={"slide_type": "subslide"}
# In order to obtain independend quantum trajectories we provide the following definition for $U_{(i,j)}$:
# $$U_{(1,1)} = \text{swap}(x_0,x_2) \times \text{swap}(y_0,y_2)$$
# $$U_{(1,2)} = \text{swap}(x_1,x_3) \times \text{swap}(y_1,y_3)$$
# $$U_{(2,1)} = \mathbf{I} $$
# $$U_{(2,2)} = \text{swap}(x_2,x_3) \times \text{swap}(y_2,y_3)$$
#
# where $ \mathbf{I}$ is the identity matrix. Thus, the step of *Sampling in Superposition* leads to:

# + [markdown] slideshow={"slide_type": "subslide"}
# \begin{align*}
#     \left|\Phi_{2}\right\rangle = \frac{1}{2}\Big[ 
#     & \left|11\right\rangle \left|x_0, x_3, x_1, x_2\right\rangle  \left|y_0, y_3, y_1, y_2\right\rangle  
#     \\ + &
#     \left|10\right\rangle \left|x_2, x_1, x_3, x_0\right\rangle  \left|y_2, y_1, y_3, y_0\right\rangle \nonumber\\
#     \hspace{.1em} 
#     + &
#     \left|01\right\rangle \left|x_0, x_3, x_2, x_1\right\rangle \left|y_0, y_3, y_2, y_1\right\rangle \\
#     + &
#     \left|00\right\rangle \left|x_2, x_1, x_0, x_3\right\rangle \left|y_2, y_1, y_0, y_3\right\rangle
#     \Big]
# \end{align*}

# + [markdown] slideshow={"slide_type": "subslide"}
# We can see that swap operations allows to entangle different observations (in terms of the indices of the qubits) to different state of the $control$ register. In particular, if considering the last qubit of the *features* and *labels* (sub-)registers, the above choice for $U_{(i,j)}$ guarantees that each quantum state of the control register is entangled with a different training observation. Using a compact representation:

# + [markdown] slideshow={"slide_type": "fragment"}
# \begin{align}
#      \left|\Phi_{2^{'}}\right\rangle = \frac{1}{2}\Big[
#     \left|11\right\rangle \left|x_2\right\rangle  \left|y_2\right\rangle  
#     + 
#     \left|10\right\rangle\left|x_0\right\rangle\left|y_0\right\rangle 
#     +
#     \left|01\right\rangle\left|x_1\right\rangle\left|y_1\right\rangle 
#     +
#     \left|00\right\rangle\left|x_3\right\rangle \left|y_3\right\rangle 
#     \Big] =
#     \frac{1}{\sqrt{4}}\sum_{i=0}^{3}\left|i\right\rangle\left|x_i,y_i\right\rangle
# \end{align}

# + [markdown] slideshow={"slide_type": "fragment"}
# Notice that, in this case the $i$-th basis state does not correspond to the integer representation of the binary state. 

# + [markdown] slideshow={"slide_type": "subslide"}
# ### (Step 3) Learning via interference

# + [markdown] slideshow={"slide_type": "fragment"}
# First, the $test$ register is initialised to encode the test set, $\tilde{x}$, considering also an additional register to store the final prediction:
#

# + [markdown] slideshow={"slide_type": "fragment"}
# \begin{align}
#     (S_{\tilde{x}} \otimes \mathbb{1}) \left|0\right\rangle \left|0\right\rangle =\left|x^{(test)}\right\rangle \left|0\right\rangle
# \end{align}

# + [markdown] slideshow={"slide_type": "fragment"}
# Then, the $data$ and $test$ registers interact via interference using the quantum version of the cosine classifier (gate $F$) to compute the estimates of the target variable:

# + [markdown] slideshow={"slide_type": "fragment"}
# \begin{align*}
#     \left|\Phi_{f}\right\rangle 
#                 = & \Big(\mathbb{1}^{\otimes 2} \otimes F \Big) \left|\Phi_{d}\right\rangle \nonumber \\ 
#                 = & (\mathbb{1}^{\otimes d} \otimes F )\Bigg[\frac{1}{\sqrt{2^d}}\sum_{b=1}^{2^d} \left|b\right\rangle \left|x_b, y_b\right\rangle\Bigg] \otimes \left|x^{(test)}\right\rangle \left|0\right\rangle  \nonumber \\ 
#                 = & \frac{1}{\sqrt{2^d}}\sum_{b=1}^{2^d} \left|b\right\rangle \left|x_b, y_b\right\rangle\left|x^{(test)}\right\rangle \left|\hat{f}_b\right\rangle
# \end{align*}

# + [markdown] slideshow={"slide_type": "fragment"}
# where $\hat{f}_b$ represents the $b$-th prediction for $\tilde{x}$ given the  $b$-th training set. 

# + [markdown] slideshow={"slide_type": "subslide"}
# ### (Step 4) Measurement

# + [markdown] slideshow={"slide_type": "fragment"}
# \begin{align*}
#     \left\langle M \right\rangle = &  
#     \frac{1}{2^d}\sum_{b=1}^{2^d} \left\langle\hat{f}_b|M|\hat{f}_b\right\rangle =
#     \frac{1}{2^d}\sum_{b=1}^{2^d}\left\langle M_b \right\rangle  \nonumber \\
#     = & \frac{1}{B} \sum_{b=1}^B \hat{f}_b = \hat{f}_{bag}(\tilde{x}|x,y)
# \end{align*}

# + [markdown] slideshow={"slide_type": "subslide"}
# ## Quantum Implementation

# + slideshow={"slide_type": "fragment"}
# Import pakages and functions
import sys
sys.path.insert(1, '../')

from Utils import *
from modeling import *

# + slideshow={"slide_type": "fragment"}
# load the toy dataset
X_data, Y_data, x_test = load_data_custom()

# Generate the quantum circuit
qc = ensemble_fixed_U(X_data, Y_data, x_test)

# + [markdown] slideshow={"slide_type": "fragment"}
# qc.draw(output='mpl', scale=.6,
#         style={'fontsize':15, 'dpi':200})

# + [markdown] slideshow={"slide_type": "slide"}
# # Quantum Ensemble as Simple Averaging -  Experiments
# -

# This notebook reproduces the results in Section 4.3 of the paper *Quantum Ensemble for Classification* where is shown that the quantum ensemble algorithm is able to compute the expectation value of multiple quantum trajectories in superposition with just one execution of the quantum cosine classifier. 
#
# For more details about the theoretical background see **Quantum Ensemble - Independent Trajectories.ipynb**

# +
# Import packages and functions
import sys
sys.path.insert(1, '../')

from Utils import *
from modeling import *

Create the toy dataset reported in Table 1 and execute the (classical) cosine classifiers.

# +
# Load data without normalisation
X_data, Y_data, x_test = load_data_custom(normalize = False)

# Create table as shown in the paper (Table 1)
data = pd.DataFrame(X_data, columns = [r'$X_1$', r'$X_2$'])

# Extract the value of the target variable as integer
y = [c[1] for c in Y_data]

# Compute the cosine distance between the training points and the test point
dist = [cosine_similarity([x], [x_test])[0][0] for x in X_data]

# Compute the value of the cosine distance classifier
# for the four training points from the test point
p = [cosine_classifier(x, x_test)[0][0] for x in X_data]

# Extract the probabilities for the test point to be classified in class 
# 1 according to the (classical) cosine classifies Equation (16)
probs = []
for i,j in zip(y,p):
    if i == 0:
        probs.append(1-j)
    else:
        probs.append(j)

# Create dataset as in paper (Table 1)
probs = np.array(probs)

# Rename columns
data[r'$y$'] = np.array(y)
data[r'$d($$\cdot$$, $ $x^{(test)})$'] = np.round(dist,2)
data[r'$P($$y^{(test)}$$=1$$|b$ $)$'] = probs

# Rename rows
data.index = [r'$x_1$', r'$x_2$', r'$x_3$', r'$x_4$',]

#Visualize dataset
data

# + slideshow={"slide_type": "subslide"}
# Load normalised data
X_data, Y_data, x_test = load_data_custom()

#Visualisation of quantum cosine classifier
quantum_cosine = quantum_cosine_classifier(X_data[0], x_test, Y_data[0] )
quantum_cosine.draw(output='mpl', scale=.7)

# + [markdown] slideshow={"slide_type": "subslide"}
# For each training point in *data* the quantum cosine classifier is executed to compute the prediction of the target variable for the test point $\tilde{x}$. Thus, given the measurements of the quantum circuts, the target probabilities are retrieved using the function *retrieve_proba*.

# + slideshow={"slide_type": "subslide"}
qc1 = quantum_cosine_classifier(X_data[0], x_test, Y_data[0] )
r1 = exec_simulator(qc1)
r1 = retrieve_proba(r1) 

qc2 = quantum_cosine_classifier(X_data[1], x_test, Y_data[1])
r2 = exec_simulator(qc2)
r2 = retrieve_proba(r2)

qc3 = quantum_cosine_classifier(X_data[2], x_test, Y_data[2])
r3 = exec_simulator(qc3)
r3 = retrieve_proba(r3)

qc4 = quantum_cosine_classifier(X_data[3], x_test, Y_data[3])
r4 = exec_simulator(qc4)
r4 = retrieve_proba(r4)

out = [r1, r2, r3, r4]

# + [markdown] slideshow={"slide_type": "subslide"}
# We compute the average of predictions provided by the four quantum cosine classifiers that corresponds to the classical ensemble prediction, using simple averaging as aggregation strategy.

# + slideshow={"slide_type": "subslide"}
p0 = [p[0] for p in out]
p1 = [p[1] for p in out]
r_avg = [np.mean(p0), np.mean(p1)]
print(np.mean(p0), np.mean(p1))

# + slideshow={"slide_type": "subslide"}
qc = ensemble_fixed_U(X_data, Y_data, x_test)
qc.draw(output='mpl', scale=.6, #filename='output/ensemble_circuit.png', 
        style={'fontsize':15, 'dpi':200})

# + slideshow={"slide_type": "subslide"}
r = exec_simulator(qc, n_shots=8192)
r_ens = retrieve_proba(r)
print(r_ens)

# +
# collect the results
output_simulator = [r1, r2, r3, r4, r_avg, r_ens]

data_pred = pd.DataFrame(output_simulator, 
                         columns=['p0', 'p1'], index=['qc1','qc2','qc3','qc4','AVG','Ensemble'])
data_pred
# data_pred.to_csv('output/sim_results.csv', index=False)

# + slideshow={"slide_type": "subslide"}
plot_cls(output_simulator, title= '')

# + [markdown] slideshow={"slide_type": "subslide"}
# The probability provided by the quantum cosine classifiers ($f_1$, $f_2$, $f_3$, $f_4$)  are pretty much the same to the classical cosine classifier (*data*). Furthermore, the average of the four classifiers is almost the same to the quantum ensemble prediction where it is required only one execution of the cosine classifier.

# + [markdown] slideshow={"slide_type": "subslide"}
# ### Multiple Experiments
#
# In order to show that the result of quantum ensemble holds regardless of the data reported in *Table 1*, the same experiment is performed on $20$ randomly generated datasets, and the average of the quantum cosine classifiers and the quantum ensemble prediction are compared

# + slideshow={"slide_type": "subslide"}
seed = 543

n_shots = 8192
N_runs = 20
y_labels =[[0,1], [1,0]]


p1_avg = []
p1_ens = []

np.random.seed(seed)
for run in np.arange(N_runs):
    # print(run)
    x1 = [np.random.randint(1, 9), np.random.randint(1, 9)]
    x2 = [np.random.randint(1, 9), np.random.randint(1, 9)]
    x3 = [np.random.randint(1, 9), np.random.randint(1, 9)]
    x4 = [np.random.randint(1, 9), np.random.randint(1, 9)]
    
    y1 = y_labels[np.random.randint(0, 2)]
    y2 = y_labels[np.random.randint(0, 2)]
    y3 = y_labels[np.random.randint(0, 2)]
    y4 = y_labels[np.random.randint(0, 2)]
    
    Y_data = [y1, y2, y3, y4]  
    X_data = [x1, x2, x3, x4]
    
    x_test = [np.random.randint(1, 9), np.random.randint(1, 9)]

    X_data, Y_data, x_test = load_data_custom(X_data, Y_data, x_test = x_test)

    qc1 = quantum_cosine_classifier(X_data[0], x_test, Y_data[0] )
    r1 = exec_simulator(qc1)
    r1 = retrieve_proba(r1)

    qc2 = quantum_cosine_classifier(X_data[1], x_test, Y_data[1])
    r2 = exec_simulator(qc2)
    r2 = retrieve_proba(r2)

    qc3 = quantum_cosine_classifier(X_data[2], x_test, Y_data[2])
    r3 = exec_simulator(qc3)
    r3 = retrieve_proba(r3)

    qc4 = quantum_cosine_classifier(X_data[3], x_test, Y_data[3])
    r4 = exec_simulator(qc4)
    r4 = retrieve_proba(r4)
     
    out = [r1, r2, r3, r4]

    p0 = [p[0] for p in out]
    p1 = [p[1] for p in out]
    
    r_avg = [np.mean(p0), np.mean(p1)]
    # print('AVG:', r_avg)

    
    qc = ensemble_fixed_U(X_data, Y_data, x_test)
    qc = transpile(qc, basis_gates = ['u1', 'u2', 'u3', 'cx'], optimization_level=3)

    r = exec_simulator(qc, n_shots=n_shots)
    r_ens = retrieve_proba(r)
    # print('Ensemble', r_ens)

    out = [r1, r2, r3, r4, r_avg, r_ens]
    
    p1_avg.append(r_avg[1])
    p1_ens.append(r_ens[1])

avg_vs_ensemble(p1_avg, p1_ens)

# +
# Execution on real device
IBMQ.load_account()
provider = IBMQ.get_provider(hub='ibm-q')
provider.backends()
backend_16 = provider.get_backend('ibmq_16_melbourne')
backend_5 = provider.get_backend('ibmq_rome')

def run_real_device(qc, backend, shots=8192):
    job = execute(qc, backend, shots=shots)
    results = job.result()
    r = results.get_counts(qc)
    return r

StarcoderdataJupyterScriptsDedupFiltered