Task2_VQA.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Variational classifiers can be used to classify and categorise data into their respective groups.\n",
    "\n",
    "Classical variational classifiers are used commonly today. In this notebook we implement a quantum variational classifier, using quantum circuits and training them to categorise new data. These quantum variational classifiers can be used for supervised machine learning.\n",
    "\n",
    "We will work on one dataset, learning the parity function, which will also show how to use basis encoding.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Parity"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import necessary libraries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pennylane as qml\n",
    "from pennylane import numpy as np\n",
    "from pennylane.optimize import AdamOptimizer"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Initiate our device:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "dev = qml.device(\"default.qubit\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We first encode our data. This problem needs a basis encoding state preparation (since it is basic qubit states with 1 and 0 e.g. |1100>).\n",
    "\n",
    "Pennylane's `BasisState()` method can help us encode our data with basis encoding:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def state_prep(x):\n",
    "    qml.BasisState(x, wires=range(4))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And now we have to define our quantum ansatze. \n",
    "\n",
    "An ansatze in general, is kind of like an assumption about the form of an unknown function. Here, it is a set of parametrised quantum gates used to train a quantum machine learning algorithm, or exclusively here, a variational quantum algorithm (VQA).\n",
    "\n",
    "Our parity function is a simple problem, therefore a set of rotation gates and entangling such as R or RZ-RY-RZ (as the tutorial implements) are enough; but based on [this paper](http://arxiv.org/pdf/2111.13730), we will implement an RX-CX-RX ansatz and a more efficient RX-RZ-CX ansatz with alternating entanglement. This latter ansatz has more effective parameters to train, and is shown to reduce error by ~95% compared to RX-CX-RX.\n",
    "\n",
    "We will implement all; please run each cell separately to try them out:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# R ansatz (tutorial implementation)\n",
    "def layer(layer_weights):\n",
    "    for wire in range(4):\n",
    "        qml.Rot(*layer_weights[wire], wires=wire)\n",
    "\n",
    "    for wires in ([0, 1], [1, 2], [2, 3], [3, 0]):\n",
    "        qml.CNOT(wires)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# RX-CX-RX ansatz\n",
    "def layer(layer_weights):\n",
    "    for wire in range(4):\n",
    "        qml.RX(layer_weights[wire][0], wires=wire)\n",
    "\n",
    "    for wires in ([0, 1], [2, 3]):\n",
    "        qml.CNOT(wires)\n",
    "\n",
    "    for wire in range(4):\n",
    "        qml.RX(layer_weights[wire][1], wires=wire)\n",
    "\n",
    "    for wires in ([1, 2], [3, 0]):\n",
    "        qml.CNOT(wires)\n",
    "\n",
    "    for wire in range(4):\n",
    "        qml.RX(layer_weights[wire][2], wires=wire)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# RX-RZ-CX ansatz\n",
    "def layer(layer_weights):\n",
    "    for wire in range(4):\n",
    "        qml.RX(layer_weights[wire][0], wires=wire)\n",
    "        qml.RZ(layer_weights[wire][0], wires=wire)\n",
    "\n",
    "    for wires in ([0, 1], [2, 3]):\n",
    "        qml.CNOT(wires)\n",
    "\n",
    "    for wire in range(4):\n",
    "        qml.RX(layer_weights[wire][1], wires=wire)\n",
    "        qml.RZ(layer_weights[wire][1], wires=wire)\n",
    "\n",
    "    for wires in ([1, 2], [3, 0]):\n",
    "        qml.CNOT(wires)\n",
    "\n",
    "    for wire in range(4):\n",
    "        qml.RX(layer_weights[wire][2], wires=wire)\n",
    "        qml.RZ(layer_weights[wire][2], wires=wire)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We then define our quantum circuit. In pennylane we have to make our circuit a \"QNode\" so that it differs from a normal python function.\n",
    "\n",
    "Our circuit will prepare the state (basis encoding) then for the number of layers, add our ansatz layer, and at last measure the expected value:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "@qml.qnode(dev)\n",
    "def circuit(weights, x):\n",
    "    state_prep(x)\n",
    "\n",
    "    for layer_weights in weights:\n",
    "        layer(layer_weights)\n",
    "\n",
    "    return qml.expval(qml.PauliZ(0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Our full model will be our circuit from the previous cell, plus a \"classical bias\" term (trainable) as a form of post-processing for the circuit:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def variational_classifier(weights, bias, x):\n",
    "    return circuit(weights, x) + bias"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we have to define our loss and metric for the model to be evaluated and trained upon. We will define MSE and Accuracy:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def square_loss(labels, predictions):\n",
    "    return np.mean((labels - qml.math.stack(predictions)) ** 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def accuracy(labels, predictions):\n",
    "    acc = sum(abs(label - pred) < 1e-4 for label, pred in zip(labels, predictions))\n",
    "    acc = acc / len(labels)\n",
    "    return acc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And now our cost function (cost function simply is the loss function applied to all of the predictions, while loss function is applied to each prediction separately):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cost(weights, bias, X, Y):\n",
    "    predictions = [variational_classifier(weights, bias, x) for x in X]\n",
    "    return square_loss(Y, predictions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will load our data and split them into X (input) and Y (label). We also do some preprocessing (shifting labels from [0, 1] to [-1, 1]) so that label distinguishing can be more significant for the model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x = [0 0 0 1], y = 1\n",
      "x = [0 0 1 0], y = 1\n",
      "x = [0 1 0 0], y = 1\n",
      "x = [0 1 0 1], y = -1\n",
      "x = [0 1 1 0], y = -1\n",
      "x = [0 1 1 1], y = 1\n",
      "x = [1 0 0 0], y = 1\n",
      "x = [1 0 0 1], y = -1\n",
      "x = [1 0 1 1], y = 1\n",
      "x = [1 1 1 1], y = -1\n"
     ]
    }
   ],
   "source": [
    "data = np.loadtxt(\"variational_classifier/data/parity_train.txt\", dtype=int)\n",
    "X = np.array(data[:, :-1])\n",
    "Y = np.array(data[:, -1])\n",
    "Y = Y * 2 - 1\n",
    "\n",
    "for x, y in zip(X, Y):\n",
    "    print(f\"x = {x}, y = {y}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now to initialise some random weights and our bias. We put `require_grad` to `True` to let the optimiser know these are the parameters we want to train.\n",
    "\n",
    "We can change the hyperparameters (num_layers, weights_init, bias_init) as we like:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Weights: [[[ 0.00552532 -0.03709256  0.01247802]\n",
      "  [ 0.02290623  0.02074381  0.03773278]\n",
      "  [-0.00223397 -0.00724203  0.0029735 ]\n",
      "  [-0.00875566  0.04342514  0.02304621]]\n",
      "\n",
      " [[-0.03637625 -0.00276099  0.01079679]\n",
      "  [-0.03550565  0.02629753 -0.00946896]\n",
      "  [-0.0218446  -0.00500055 -0.01964589]\n",
      "  [ 0.02062538  0.00982668 -0.00893293]]\n",
      "\n",
      " [[-0.0161272   0.00262536 -0.0242512 ]\n",
      "  [ 0.00319982 -0.01510446  0.00699792]\n",
      "  [ 0.01955084 -0.0027717   0.00207713]\n",
      "  [ 0.00601182  0.01936411  0.01739248]]]\n",
      "Bias:  0.7\n"
     ]
    }
   ],
   "source": [
    "np.random.seed(17)\n",
    "num_qubits = 4\n",
    "num_layers = 3\n",
    "weights_init = 0.02 * np.random.randn(num_layers, num_qubits, 3, requires_grad=True)\n",
    "bias_init = np.array(0.7, requires_grad=True)\n",
    "\n",
    "print(\"Weights:\", weights_init)\n",
    "print(\"Bias: \", bias_init)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We initiate our optimiser, which is an Adam optimiser, and our batch size to feed data to our model cost:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "opt = AdamOptimizer(0.3)\n",
    "batch_size = 4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And now we train our model. The optimiser will update the weights and bias to acquire the best possible parameters, and then we make predictions and compute the accuracy each iteration/epoch. We can increase or decrease the number of iterations to change the learnability of the model. \n",
    "\n",
    "Note that increasing the iterations too much might lead to overfitting, or too few iterations causes underfitting:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Iter:    1 | Cost: 1.4913523 | Accuracy: 0.5000000\n",
      "Iter:    2 | Cost: 1.2636099 | Accuracy: 0.4000000\n",
      "Iter:    3 | Cost: 1.2582003 | Accuracy: 0.4000000\n",
      "Iter:    4 | Cost: 0.9442554 | Accuracy: 0.6000000\n",
      "Iter:    5 | Cost: 0.9396397 | Accuracy: 0.7000000\n",
      "Iter:    6 | Cost: 0.8886298 | Accuracy: 0.7000000\n",
      "Iter:    7 | Cost: 0.9455230 | Accuracy: 0.7000000\n",
      "Iter:    8 | Cost: 0.8868692 | Accuracy: 0.7000000\n",
      "Iter:    9 | Cost: 0.8266910 | Accuracy: 0.8000000\n",
      "Iter:   10 | Cost: 0.8261738 | Accuracy: 0.6000000\n",
      "Iter:   11 | Cost: 0.8971381 | Accuracy: 0.6000000\n",
      "Iter:   12 | Cost: 0.9153667 | Accuracy: 0.8000000\n",
      "Iter:   13 | Cost: 0.8212975 | Accuracy: 0.7000000\n",
      "Iter:   14 | Cost: 0.7611873 | Accuracy: 0.6000000\n",
      "Iter:   15 | Cost: 0.8117961 | Accuracy: 0.7000000\n",
      "Iter:   16 | Cost: 1.0210193 | Accuracy: 0.6000000\n",
      "Iter:   17 | Cost: 1.1843587 | Accuracy: 0.6000000\n",
      "Iter:   18 | Cost: 1.1922951 | Accuracy: 0.6000000\n",
      "Iter:   19 | Cost: 1.0539110 | Accuracy: 0.6000000\n",
      "Iter:   20 | Cost: 1.0469186 | Accuracy: 0.6000000\n",
      "Iter:   21 | Cost: 1.0559071 | Accuracy: 0.6000000\n",
      "Iter:   22 | Cost: 1.0682091 | Accuracy: 0.6000000\n",
      "Iter:   23 | Cost: 1.0049225 | Accuracy: 0.6000000\n",
      "Iter:   24 | Cost: 0.9551770 | Accuracy: 0.6000000\n",
      "Iter:   25 | Cost: 0.8805330 | Accuracy: 0.7000000\n",
      "Iter:   26 | Cost: 0.8124167 | Accuracy: 0.8000000\n",
      "Iter:   27 | Cost: 0.7520682 | Accuracy: 0.8000000\n",
      "Iter:   28 | Cost: 0.7007398 | Accuracy: 0.8000000\n",
      "Iter:   29 | Cost: 0.6346478 | Accuracy: 0.8000000\n",
      "Iter:   30 | Cost: 0.5743598 | Accuracy: 0.9000000\n",
      "Iter:   31 | Cost: 0.5796448 | Accuracy: 0.8000000\n",
      "Iter:   32 | Cost: 0.6280195 | Accuracy: 0.9000000\n",
      "Iter:   33 | Cost: 0.6807400 | Accuracy: 0.8000000\n",
      "Iter:   34 | Cost: 0.6426360 | Accuracy: 1.0000000\n",
      "Iter:   35 | Cost: 0.5835949 | Accuracy: 0.9000000\n",
      "Iter:   36 | Cost: 0.6486932 | Accuracy: 0.9000000\n",
      "Iter:   37 | Cost: 0.6121696 | Accuracy: 0.9000000\n",
      "Iter:   38 | Cost: 0.7297439 | Accuracy: 0.7000000\n",
      "Iter:   39 | Cost: 0.8720292 | Accuracy: 0.7000000\n",
      "Iter:   40 | Cost: 0.9451530 | Accuracy: 0.7000000\n",
      "Iter:   41 | Cost: 1.0202595 | Accuracy: 0.5000000\n",
      "Iter:   42 | Cost: 1.0515379 | Accuracy: 0.5000000\n",
      "Iter:   43 | Cost: 1.0774587 | Accuracy: 0.6000000\n",
      "Iter:   44 | Cost: 0.9633594 | Accuracy: 0.7000000\n",
      "Iter:   45 | Cost: 0.8242365 | Accuracy: 0.6000000\n",
      "Iter:   46 | Cost: 0.7740563 | Accuracy: 0.7000000\n",
      "Iter:   47 | Cost: 0.7953407 | Accuracy: 0.6000000\n",
      "Iter:   48 | Cost: 0.7742674 | Accuracy: 0.5000000\n",
      "Iter:   49 | Cost: 0.7104017 | Accuracy: 0.6000000\n",
      "Iter:   50 | Cost: 0.5490838 | Accuracy: 0.7000000\n"
     ]
    }
   ],
   "source": [
    "weights = weights_init\n",
    "bias = bias_init\n",
    "for iteration in range(50):\n",
    "\n",
    "    batch_index = np.random.randint(0, len(X), (batch_size,))  # random batch of data based on batch size\n",
    "    X_batch = X[batch_index]\n",
    "    Y_batch = Y[batch_index]\n",
    "\n",
    "    weights, bias = opt.step(cost, weights, bias, X=X_batch, Y=Y_batch)  # updating weights and bias\n",
    "\n",
    "    train_preds = [np.sign(variational_classifier(weights, bias, x)) for x in X]  # predictions based on the new weights and bias\n",
    "\n",
    "    # compute new cost and accuracy\n",
    "    le_cost = cost(weights, bias, X, Y)\n",
    "    train_acc = accuracy(Y, train_preds)\n",
    "\n",
    "    print(f\"Iter: {iteration+1:4d} | Cost: {le_cost:0.7f} | Accuracy: {train_acc:0.7f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we can see, the model already achieved perfect accuracy and low cost early in ~15-25 iterations onwards (with R ansatz). This shows our choice of ansatz and optimiser were pretty good.\n",
    "\n",
    "Results of different ansatzes:\n",
    "- R ansatz:         Final accuracy 100%, Final cost 0.0088925   (faster model. pretty amazing ansatz!)\n",
    "- RX-CX-RX ansatz:  Final accuracy 50%, Final cost 1.0372201    (alternating accuracy between 50% & 60%, cost between [0.9, 1.2]; not a good ansatz!)\n",
    "- RX-RZ-CX ansatz:  Final accuracy 100%, Final cost 0.1980340   (mixed results. sometimes better than all, sometimes worse that R; but always better than RX-CX-RX.)\n",
    "\n",
    "In summary, it seems if we want a good model, better to go with the R ansatz. But we'll test the models first then conclude.\n",
    "\n",
    "Do note that the hyperparameters we choose also affect the results."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And finally we test our model on some test data to see how it does on unseen data. Here we also do the same preprocessing as in train data, and then run our model and calculate predictions and accuracy on unseen data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x = [0 0 0 0], y = -1,   pred = 1.0\n",
      "x = [0 0 1 1], y = -1,   pred = 1.0\n",
      "x = [1 0 1 0], y = -1,   pred = 1.0\n",
      "x = [1 1 1 0], y = 1,   pred = 1.0\n",
      "x = [1 1 0 0], y = -1,   pred = -1.0\n",
      "x = [1 1 0 1], y = 1,   pred = 1.0\n",
      "x = [1 0 0 0], y = 1,   pred = 1.0\n",
      "x = [1 1 1 1], y = -1,   pred = -1.0\n",
      "x = [0 1 0 1], y = -1,   pred = 1.0\n",
      "x = [0 0 1 0], y = 1,   pred = 1.0\n",
      "x = [0 1 0 0], y = 1,   pred = -1.0\n",
      "Accuracy on unseen data: 0.5454545454545454\n"
     ]
    }
   ],
   "source": [
    "data = np.loadtxt(\"variational_classifier/data/parity_test.txt\", dtype=int)\n",
    "X_test = np.array(data[:, :-1])\n",
    "Y_test = np.array(data[:, -1])\n",
    "Y_test = Y_test * 2 - 1\n",
    "\n",
    "test_preds = [np.sign(variational_classifier(weights, bias, x)) for x in X_test]\n",
    "\n",
    "for x, y, pred in zip(X_test, Y_test, test_preds):\n",
    "    print(f\"x = {x}, y = {y},   pred = {pred}\")\n",
    "\n",
    "test_acc = accuracy(Y_test, test_preds)\n",
    "print(\"Accuracy on unseen data:\", test_acc)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- R ansatz: 100% accuracy\n",
    "- RX-CX-RX ansatz: 45% accuracy\n",
    "- RX-RZ-CX ansatz: mixed accuracy\n",
    "\n",
    "Therefore to conclude, our best choice of ansatz would be R, but we can also use RX-RZ-CX as its accuracy and cost are similar to R, but take into account the speed, as R was the fastest of all.\n",
    "\n",
    "RX-CX-RX ansatz was the worst, so it's better not to use it."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "nickel",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}