tutorials/Tutorial 1 - Training and Evaluation of Logistic Regression on Encrypted Data.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Tutorial 1: Training and Evaluation of Logistic Regression on Encrypted Data\n",
    "\n",
    "Welcome to this first use case tutorial, where we are going to show how to use TenSEAL for training and evaluating a logistic regression (LR) model on encrypted data (using homomorphic encryption) for heart disease prediction! If you haven't played with TenSEAL before, or need a quick overview of what homomorphic encryption is, I would suggest going through ['Tutorial 0 - Getting Started'](./Tutorial%200%20-%20Getting%20Started.ipynb) first.\n",
    "\n",
    "\n",
    "**Important note:** The goal of this tutorial isn't to show how efficient logistic regression is for this task, we will just go with whatever accuracy we get, but the training and evaluation on encrypted data should be comparable to when we use plain data.\n",
    "\n",
    "\n",
    "Authors:\n",
    "- Ayoub Benaissa - Twitter: [@y0uben11](https://twitter.com/y0uben11)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup\n",
    "\n",
    "All modules are imported here. Make sure everything is installed by running the cell below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import torch\n",
    "import tenseal as ts\n",
    "import pandas as pd\n",
    "import random\n",
    "from time import time\n",
    "\n",
    "# those are optional and are not necessary for training\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We now prepare the training and test data. The dataset was downloaded from Kaggle [here](https://www.kaggle.com/dileep070/heart-disease-prediction-using-logistic-regression). This dataset includes patients' information along with a 10-year risk of future coronary heart disease (CHD) as a label. The goal is to build a model that can predict this 10-year CHD risk based on patients' information. You can read more about the dataset in the link provided. \n",
    "\n",
    "Alternatively, we also provide the `random_data()` function below that generates random, linearly separable points. You can use it instead of the dataset from Kaggle, for those who just want to see how things work. The rest of the tutorial should work in the same way."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "############# Data summary #############\n",
      "x_train has shape: torch.Size([780, 9])\n",
      "y_train has shape: torch.Size([780, 1])\n",
      "x_test has shape: torch.Size([334, 9])\n",
      "y_test has shape: torch.Size([334, 1])\n",
      "#######################################\n"
     ]
    }
   ],
   "source": [
    "torch.random.manual_seed(73)\n",
    "random.seed(73)\n",
    "\n",
    "\n",
    "def split_train_test(x, y, test_ratio=0.3):\n",
    "    idxs = [i for i in range(len(x))]\n",
    "    random.shuffle(idxs)\n",
    "    # delimiter between test and train data\n",
    "    delim = int(len(x) * test_ratio)\n",
    "    test_idxs, train_idxs = idxs[:delim], idxs[delim:]\n",
    "    return x[train_idxs], y[train_idxs], x[test_idxs], y[test_idxs]\n",
    "\n",
    "\n",
    "def heart_disease_data():\n",
    "    data = pd.read_csv(\"./data/framingham.csv\")\n",
    "    # drop rows with missing values\n",
    "    data = data.dropna()\n",
    "    # drop some features\n",
    "    data = data.drop(columns=[\"education\", \"currentSmoker\", \"BPMeds\", \"diabetes\", \"diaBP\", \"BMI\"])\n",
    "    # balance data\n",
    "    grouped = data.groupby('TenYearCHD')\n",
    "    data = grouped.apply(lambda x: x.sample(grouped.size().min(), random_state=73).reset_index(drop=True))\n",
    "    # extract labels\n",
    "    y = torch.tensor(data[\"TenYearCHD\"].values).float().unsqueeze(1)\n",
    "    data = data.drop(\"TenYearCHD\", 'columns')\n",
    "    # standardize data\n",
    "    data = (data - data.mean()) / data.std()\n",
    "    x = torch.tensor(data.values).float()\n",
    "    return split_train_test(x, y)\n",
    "\n",
    "\n",
    "def random_data(m=1024, n=2):\n",
    "    # data separable by the line `y = x`\n",
    "    x_train = torch.randn(m, n)\n",
    "    x_test = torch.randn(m // 2, n)\n",
    "    y_train = (x_train[:, 0] >= x_train[:, 1]).float().unsqueeze(0).t()\n",
    "    y_test = (x_test[:, 0] >= x_test[:, 1]).float().unsqueeze(0).t()\n",
    "    return x_train, y_train, x_test, y_test\n",
    "\n",
    "\n",
    "# You can use whatever data you want without modification to the tutorial\n",
    "# x_train, y_train, x_test, y_test = random_data()\n",
    "x_train, y_train, x_test, y_test = heart_disease_data()\n",
    "\n",
    "print(\"############# Data summary #############\")\n",
    "print(f\"x_train has shape: {x_train.shape}\")\n",
    "print(f\"y_train has shape: {y_train.shape}\")\n",
    "print(f\"x_test has shape: {x_test.shape}\")\n",
    "print(f\"y_test has shape: {y_test.shape}\")\n",
    "print(\"#######################################\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Training a Logistic Regression Model\n",
    "\n",
    "We will start by training a logistic regression model (without any encryption), which can be viewed as a single layer neural network with a single node. We will be using this model as a means of comparison against encrypted training and evaluation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "class LR(torch.nn.Module):\n",
    "\n",
    "    def __init__(self, n_features):\n",
    "        super(LR, self).__init__()\n",
    "        self.lr = torch.nn.Linear(n_features, 1)\n",
    "        \n",
    "    def forward(self, x):\n",
    "        out = torch.sigmoid(self.lr(x))\n",
    "        return out"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "n_features = x_train.shape[1]\n",
    "model = LR(n_features)\n",
    "# use gradient descent with a learning_rate=1\n",
    "optim = torch.optim.SGD(model.parameters(), lr=1)\n",
    "# use Binary Cross Entropy Loss\n",
    "criterion = torch.nn.BCELoss()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Loss at epoch 1: 0.8504332900047302\n",
      "Loss at epoch 2: 0.6863385438919067\n",
      "Loss at epoch 3: 0.635811448097229\n",
      "Loss at epoch 4: 0.6193529367446899\n",
      "Loss at epoch 5: 0.6124349236488342\n"
     ]
    }
   ],
   "source": [
    "# define the number of epochs for both plain and encrypted training\n",
    "EPOCHS = 5\n",
    "\n",
    "def train(model, optim, criterion, x, y, epochs=EPOCHS):\n",
    "    for e in range(1, epochs + 1):\n",
    "        optim.zero_grad()\n",
    "        out = model(x)\n",
    "        loss = criterion(out, y)\n",
    "        loss.backward()\n",
    "        optim.step()\n",
    "        print(f\"Loss at epoch {e}: {loss.data}\")\n",
    "    return model\n",
    "\n",
    "model = train(model, optim, criterion, x_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy on plain test_set: 0.703592836856842\n"
     ]
    }
   ],
   "source": [
    "def accuracy(model, x, y):\n",
    "    out = model(x)\n",
    "    correct = torch.abs(y - out) < 0.5\n",
    "    return correct.float().mean()\n",
    "\n",
    "plain_accuracy = accuracy(model, x_test, y_test)\n",
    "print(f\"Accuracy on plain test_set: {plain_accuracy}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is worth to remember that a high accuracy isn't our goal. We just want to see that training on encrypted data doesn't affect the final result, so we will be comparing accuracies over encrypted data against the `plain_accuracy` we got here."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Encrypted Evaluation\n",
    "\n",
    "In this part, we will just focus on evaluating the logistic regression model with plain parameters (optionally encrypted parameters) on the encrypted test set. We first create a PyTorch-like LR model that can evaluate encrypted data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "class EncryptedLR:\n",
    "    \n",
    "    def __init__(self, torch_lr):\n",
    "        # TenSEAL processes lists and not torch tensors,\n",
    "        # so we take out the parameters from the PyTorch model\n",
    "        self.weight = torch_lr.lr.weight.data.tolist()[0]\n",
    "        self.bias = torch_lr.lr.bias.data.tolist()\n",
    "        \n",
    "    def forward(self, enc_x):\n",
    "        # We don't need to perform sigmoid as this model\n",
    "        # will only be used for evaluation, and the label\n",
    "        # can be deduced without applying sigmoid\n",
    "        enc_out = enc_x.dot(self.weight) + self.bias\n",
    "        return enc_out\n",
    "    \n",
    "    def __call__(self, *args, **kwargs):\n",
    "        return self.forward(*args, **kwargs)\n",
    "        \n",
    "    ################################################\n",
    "    ## You can use the functions below to perform ##\n",
    "    ## the evaluation with an encrypted model     ##\n",
    "    ################################################\n",
    "    \n",
    "    def encrypt(self, context):\n",
    "        self.weight = ts.ckks_vector(context, self.weight)\n",
    "        self.bias = ts.ckks_vector(context, self.bias)\n",
    "        \n",
    "    def decrypt(self, context):\n",
    "        self.weight = self.weight.decrypt()\n",
    "        self.bias = self.bias.decrypt()\n",
    "        \n",
    "\n",
    "eelr = EncryptedLR(model)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We now create a TenSEALContext for specifying the scheme and the parameters we are going to use. Here we choose small and secure parameters that allow us to make a single multiplication. That's enough for evaluating a logistic regression model, however, we will see that we need larger parameters when doing training on encrypted data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# parameters\n",
    "poly_mod_degree = 4096\n",
    "coeff_mod_bit_sizes = [40, 20, 40]\n",
    "# create TenSEALContext\n",
    "ctx_eval = ts.context(ts.SCHEME_TYPE.CKKS, poly_mod_degree, -1, coeff_mod_bit_sizes)\n",
    "# scale of ciphertext to use\n",
    "ctx_eval.global_scale = 2 ** 20\n",
    "# this key is needed for doing dot-product operations\n",
    "ctx_eval.generate_galois_keys()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will encrypt the whole test set before the evaluation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Encryption of the test-set took 1 seconds\n"
     ]
    }
   ],
   "source": [
    "t_start = time()\n",
    "enc_x_test = [ts.ckks_vector(ctx_eval, x.tolist()) for x in x_test]\n",
    "t_end = time()\n",
    "print(f\"Encryption of the test-set took {int(t_end - t_start)} seconds\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# (optional) encrypt the model's parameters\n",
    "# eelr.encrypt(ctx_eval)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you may have already noticed when we built the EncryptedLR class, we don't compute the sigmoid function on the encrypted output of the linear layer, simply because it's not needed, and computing sigmoid over encrypted data will increase the computation time and require larger encryption parameters. However, we will use sigmoid for the encrypted training part. We now proceed with the evaluation of the encrypted test set and compare the accuracy to the one on the plain test set."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Evaluated test_set of 334 entries in 1 seconds\n",
      "Accuracy: 225/334 = 0.6736526946107785\n",
      "Difference between plain and encrypted accuracies: 0.029940128326416016\n"
     ]
    }
   ],
   "source": [
    "def encrypted_evaluation(model, enc_x_test, y_test):\n",
    "    t_start = time()\n",
    "    \n",
    "    correct = 0\n",
    "    for enc_x, y in zip(enc_x_test, y_test):\n",
    "        # encrypted evaluation\n",
    "        enc_out = model(enc_x)\n",
    "        # plain comparison\n",
    "        out = enc_out.decrypt()\n",
    "        out = torch.tensor(out)\n",
    "        out = torch.sigmoid(out)\n",
    "        if torch.abs(out - y) < 0.5:\n",
    "            correct += 1\n",
    "    \n",
    "    t_end = time()\n",
    "    print(f\"Evaluated test_set of {len(x_test)} entries in {int(t_end - t_start)} seconds\")\n",
    "    print(f\"Accuracy: {correct}/{len(x_test)} = {correct / len(x_test)}\")\n",
    "    return correct / len(x_test)\n",
    "    \n",
    "\n",
    "encrypted_accuracy = encrypted_evaluation(eelr, enc_x_test, y_test)\n",
    "diff_accuracy = plain_accuracy - encrypted_accuracy\n",
    "print(f\"Difference between plain and encrypted accuracies: {diff_accuracy}\")\n",
    "if diff_accuracy < 0:\n",
    "    print(\"Oh! We got a better accuracy on the encrypted test-set! The noise was on our side...\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We saw that evaluating on the encrypted test set doesn't affect the accuracy that much. I've even seen examples where the encrypted evaluation performs better."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Training an Encrypted Logistic Regression Model on Encrypted Data\n",
    "\n",
    "In this part we will redefine a PyTorch-like model that can both forward encrypted data, as well as backpropagate to update the weights and thus train the encrypted logistic regression model on encrypted data. Below are more details about the training.\n",
    "\n",
    "#### Loss Function\n",
    "\n",
    "We are using the binary cross entropy loss function with regularization (more about the why of regularization will follow) where $y^{(i)}$ is the i'th expected label, $\\hat{y}^{(i)}$ is the i'th output of the logistic regression model and $\\theta$ is our n-sized weight vector.\n",
    "\n",
    "$$Loss(\\theta) = - \\frac{1}{m} \\sum_{i=1}^m [y^{(i)} log(\\hat{y}^{(i)}) + (1 - y^{(i)}) log (1 - \\hat{y}^{(i)})] + \\frac{\\lambda}{2m} \\sum_{j=1}^n \\theta_j^2$$\n",
    "\n",
    "#### Parameters Update\n",
    "\n",
    "For updating the parameter, the usual rule is as follows, where $x^{(i)}$ is the i'th input data:\n",
    "\n",
    "$$\\theta_j = \\theta_j - \\alpha \\; [ \\frac{1}{m} \\sum_{i=1}^m (\\hat{y}^{(i)} - y^{(i)}) x^{(i)} + \\frac{\\lambda}{m} \\theta_j]$$\n",
    "\n",
    "However, due to homomorphic encryption constraint, we preferred to use an $\\alpha = 1$ to reduce a multiplication and set $\\frac{\\lambda}{m} = 0.05$ which gets us to the following update rule:\n",
    "\n",
    "$$\\theta_j = \\theta_j - [ \\frac{1}{m} \\sum_{i=1}^m (\\hat{y}^{(i)} - y^{(i)}) x^{(i)} + 0.05 \\theta_j]$$\n",
    "\n",
    "#### Sigmoid Approximation\n",
    "\n",
    "Since we can't simply compute sigmoid on encrypted data, we need to approximate it using a low degree polynomial, the lower the degree the better, as we aim to perform as few multiplications as possible, to be able to use smaller parameters and thus optimize computation. This tutorial uses a degree 3 polynomial from https://eprint.iacr.org/2018/462.pdf, which approximates the sigmoid function in the range $[-5,5]$.\n",
    "\n",
    "$$\\sigma(x) = 0.5 + 0.197 x - 0.004 x^3$$\n",
    "\n",
    "#### Homomorphic Encryption Parameters\n",
    "\n",
    "From the input data to the parameter update, a ciphertext will need a multiplicative depth of 6, 1 for the dot product operation, 2 for the sigmoid approximation, and 3 for the backpropagation phase (one is actually hidden in the `self._delta_w += enc_x * out_minus_y` operation in the `backward()` function, which is multiplying a 1-sized vector with an n-sized one, which requires masking the first slot and replicating it n times in the first vector). With a scale of around 20 bits, we need 6 coefficients modulus with the same bit-size as the scale, plus the last coefficient, which needs more bits, we are already out of the 4096 polynomial modulus degree (which requires < 109 total bit count of the coefficients modulus, if we consider 128-bit security), so we will use 8192. This will allow us to batch up to 4096 values in a single ciphertext, but we are far away from this limitation, so we shouldn't even think about it.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "class EncryptedLR:\n",
    "    \n",
    "    def __init__(self, torch_lr):\n",
    "        self.weight = torch_lr.lr.weight.data.tolist()[0]\n",
    "        self.bias = torch_lr.lr.bias.data.tolist()\n",
    "        # we accumulate gradients and counts the number of iterations\n",
    "        self._delta_w = 0\n",
    "        self._delta_b = 0\n",
    "        self._count = 0\n",
    "        \n",
    "    def forward(self, enc_x):\n",
    "        enc_out = enc_x.dot(self.weight) + self.bias\n",
    "        enc_out = EncryptedLR.sigmoid(enc_out)\n",
    "        return enc_out\n",
    "    \n",
    "    def backward(self, enc_x, enc_out, enc_y):\n",
    "        out_minus_y = (enc_out - enc_y)\n",
    "        self._delta_w += enc_x * out_minus_y\n",
    "        self._delta_b += out_minus_y\n",
    "        self._count += 1\n",
    "        \n",
    "    def update_parameters(self):\n",
    "        if self._count == 0:\n",
    "            raise RuntimeError(\"You should at least run one forward iteration\")\n",
    "        # update weights\n",
    "        # We use a small regularization term to keep the output\n",
    "        # of the linear layer in the range of the sigmoid approximation\n",
    "        self.weight -= self._delta_w * (1 / self._count) + self.weight * 0.05\n",
    "        self.bias -= self._delta_b * (1 / self._count)\n",
    "        # reset gradient accumulators and iterations count\n",
    "        self._delta_w = 0\n",
    "        self._delta_b = 0\n",
    "        self._count = 0\n",
    "    \n",
    "    @staticmethod\n",
    "    def sigmoid(enc_x):\n",
    "        # We use the polynomial approximation of degree 3\n",
    "        # sigmoid(x) = 0.5 + 0.197 * x - 0.004 * x^3\n",
    "        # from https://eprint.iacr.org/2018/462.pdf\n",
    "        # which fits the function pretty well in the range [-5,5]\n",
    "        return enc_x.polyval([0.5, 0.197, 0, -0.004])\n",
    "    \n",
    "    def plain_accuracy(self, x_test, y_test):\n",
    "        # evaluate accuracy of the model on\n",
    "        # the plain (x_test, y_test) dataset\n",
    "        w = torch.tensor(self.weight)\n",
    "        b = torch.tensor(self.bias)\n",
    "        out = torch.sigmoid(x_test.matmul(w) + b).reshape(-1, 1)\n",
    "        correct = torch.abs(y_test - out) < 0.5\n",
    "        return correct.float().mean()    \n",
    "    \n",
    "    def encrypt(self, context):\n",
    "        self.weight = ts.ckks_vector(context, self.weight)\n",
    "        self.bias = ts.ckks_vector(context, self.bias)\n",
    "        \n",
    "    def decrypt(self):\n",
    "        self.weight = self.weight.decrypt()\n",
    "        self.bias = self.bias.decrypt()\n",
    "        \n",
    "    def __call__(self, *args, **kwargs):\n",
    "        return self.forward(*args, **kwargs)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "# parameters\n",
    "poly_mod_degree = 8192\n",
    "coeff_mod_bit_sizes = [40, 21, 21, 21, 21, 21, 21, 40]\n",
    "# create TenSEALContext\n",
    "ctx_training = ts.context(ts.SCHEME_TYPE.CKKS, poly_mod_degree, -1, coeff_mod_bit_sizes)\n",
    "ctx_training.global_scale = 2 ** 21\n",
    "ctx_training.generate_galois_keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Encryption of the training_set took 26 seconds\n"
     ]
    }
   ],
   "source": [
    "t_start = time()\n",
    "enc_x_train = [ts.ckks_vector(ctx_training, x.tolist()) for x in x_train]\n",
    "enc_y_train = [ts.ckks_vector(ctx_training, y.tolist()) for y in y_train]\n",
    "t_end = time()\n",
    "print(f\"Encryption of the training_set took {int(t_end - t_start)} seconds\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Below we study the distribution of `x.dot(weight) + bias` in both plain and encrypted domains. Making sure that it falls into the range $[-5,5]$, which is where our sigmoid approximation is good at, and we don't want to feed it data that is out of this range so that we don't get erroneous output, which can make our training behave unpredictably. But the weights will change during the training process, and we should try to keep them as small as possible while still learning. A technique often used with logistic regression, and we do exactly this (but serving another purpose which is *generalization*), is known as *regularization*, and you might already have spotted the additional term `self.weight * 0.05` in the `update_parameters()` function, which is the result of doing regularization.\n",
    "\n",
    "To recap, since our sigmoid approximation is only good in the range $[-5,5]$, we want to have all its inputs in that range. In order to do this, we need to keep our logistic regression parameters as small as possible, so we apply regularization.\n",
    "\n",
    "**Note:** Keeping the parameters small certainly reduces the magnitude of the output, but we can also get out of range if the data wasn't standardized. You may have spotted that we standardized the data with a mean of 0 and std of 1, this was both for better performance, as well as to keep the inputs to the sigmoid in the desired range."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Distribution on plain data:\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3daXAc553f8e8fAwxAHCQAAjzEWzJlmbKtC0t7vdbaqbVjSrtr2nt4JW/WTmyXSkmUxC82ZbpccZxV5YVjbzaVtbxc2VHZ2XJZScUX10uv5GzF97GkbOqgJEo0RYngCRAXcQ6Of15MDzgaDoAG0IPpbvw+VSjOdD+Y+bMB/PDg6X6eNndHRESSr6baBYiISDQU6CIiKaFAFxFJCQW6iEhKKNBFRFKitlpv3NHR4Tt37qzW24uIJNITTzzR6+6d5fZVLdB37tzJ0aNHq/X2IiKJZGYvz7VPQy4iIimhQBcRSQkFuohISijQRURSQoEuIpISCnQRkZRQoIuIpIQCXVatZ88N8Tc/Pc345HS1SxGJRNUmFolUU99Ijnu/8DMGxyZ58dIwf7b/9dUuSWTZ1EOXVelbx84yODbJG7eu438fPcPQ+GS1SxJZtlCBbmb7zOyEmZ00swNl9v97MzsWfDxjZtNm1h59uSLR+PZT57n5urV84u7XMT45w09/dbnaJYks24KBbmYZ4CHgLmAPcK+Z7Slu4+6fcfdb3f1W4OPA9929rxIFiyzXxNQ0T3cP8tbXdHDr9lYa6moU6JIKYXroe4GT7n7K3XPAo8D+edrfC3w1iuJEKuHZc0Pkpme4bXsr9bUZ7tjRxtGX1f+Q5AsT6FuAM0XPu4Nt1zCzRmAf8LU59t9nZkfN7GhPT89iaxWJxC9fGQDgtu1tANx83TpeuDjM1PRMNcsSWbYwgW5ltvkcbX8X+PFcwy3u/rC7d7l7V2dn2eV8RSru+QtDdDRn2bi2AYDXbW4hNzXDS70jVa5MZHnCBHo3sK3o+Vbg3Bxt70HDLRJzL/WOcH1H8+zzmzatBeDZ80PVKkkkEmEC/Qiw28x2mVmWfGgfKm1kZuuAtwHfirZEkWid6hnh+s6m2ee7Opowg9O9o1WsSmT5FpxY5O5TZvYA8BiQAR5x9+Nmdn+w/2DQ9L3A4+6uv1sltgZHJ7k8kmNXx9VAb6jLsHltAy9f1reuJFuomaLufhg4XLLtYMnzLwFfiqowkUp4KQjt4kAH2L6+kZf71EOXZNNMUVlVTvUMA3B9Z/Ortu9c36QeuiSeAl1WlVeCXvi29jWv2r5jfRO9wzmGJ6aqUZZIJBTosqqcHxino7me+trMq7Zvb28E4IyGXSTBFOiyqpwfGue61oZrtm9al992YXB8pUsSiYwCXVaV8wNjbF53baAXQv68Al0STIEuq8r5wXE2r1tzzfbO5npqDM4PjlWhKpFoKNBl1Rgan2R4YqrskEttpoYNLQ3qoUuiKdBl1Tg/kA/rTWV66ACbWxs0hi6JpkCXVaMwnHJdmTF0gM3rGjTkIommQJdVozCcsrm1fA9909o1nB8cx32uxURF4k2BLqvG+cFxzGBjS33Z/ZvXNTCam2ZoXJOLJJkU6LJq9FyZYH1TltpM+W/7ziDoe4cnVrIskcgo0GXV6B2eoKO5fO8cmN3Xe0WBLsmkQJdVY8FAb8kG7XIrVZJIpBTosmrkAz075/7OIOx7rujSRUkmBbqsGr1XcvP20Nsas2RqTD10SSwFuqwKIxNTjE1O0zHHFS4ANTVGe1NWJ0UlsRTosioUQnp909xDLpA/MapAl6RSoMuqUAjp+XrokL90sUdXuUhChQp0M9tnZifM7KSZHZijzdvN7JiZHTez70dbpsjy9FzJj4t3zjOGDtDRnNUYuiTWgjeJNrMM8BDwTqAbOGJmh9z92aI2rcDngX3u/oqZbahUwSJLMdtDXyDQO5vr6RmewN0xs5UoTSQyYXroe4GT7n7K3XPAo8D+kjbvB77u7q8AuPulaMsUWZ7ZMfR5LluEfODnpma4onuLSgKFCfQtwJmi593BtmI3Am1m9j0ze8LMPlDuhczsPjM7amZHe3p6llaxyBL0Dk/Q2lhH3RzT/gtmp/9rHF0SKEygl/u7s3Q5ulrgDuC3gXcB/8HMbrzmk9wfdvcud+/q7OxcdLEiS7XQNegFHbOTixTokjwLjqGT75FvK3q+FThXpk2vu48AI2b2A+AW4IVIqhRZpssj888SLSgMyfSP6sSoJE+YHvoRYLeZ7TKzLHAPcKikzbeAO82s1swagTcBz0VbqsjS9Q6H66G3NeYDvW9kstIliURuwR66u0+Z2QPAY0AGeMTdj5vZ/cH+g+7+nJn9PfAUMAN80d2fqWThIotxeXhiwUlFAK2NdYB66JJMYYZccPfDwOGSbQdLnn8G+Ex0pYlEY2p6hqHxKVobFw70hroMTdkMfSMKdEkezRSV1BsYyw+ftIfooQO0NWXpV6BLAinQJfUGguGTwnDKQtqbsvRpyEUSSIEuqdc/mu+ht4UYcim0Uw9dkkiBLqlXGA8PO+TS3pSd/SUgkiQKdEm9xQ65tDbWqYcuiaRAl9Rb7JBLe2OWKxNT5KZmKlmWSOQU6JJ6/SM5srU1NGYzodq3BUMzAzoxKgmjQJfU6x/N0dZYF3o53MJYu650kaRRoEvq9Y9Ohh5ugeLp/wp0SRYFuqRe/0huUYFe6KH3az0XSRgFuqRe/2iOtqZwV7gAtAVXw2jIRZJGgS6pN7DIIZfCmi+6dFGSRoEuqTYz48FJ0fCBnq2toaW+VisuSuIo0CXVroxPMePhJxUVaIEuSSIFuqRaoZcddtp/QVtTlj5N/5eEUaBLqhVObC5myCXfvo6+Ed1XVJJFgS6ptth1XAryKy6qhy7JokCXVCvcG3TRQy6NWU39l8QJFehmts/MTpjZSTM7UGb/281s0MyOBR+fjL5UkcW72kNf/JDLSG5aC3RJoix4T1EzywAPAe8EuoEjZnbI3Z8tafpDd/+dCtQosmT9ozkyNcbahlC3z53VWrRA14a1DZUoTSRyYXroe4GT7n7K3XPAo8D+ypYlEo2+kclFLcxVUJgtqhtdSJKECfQtwJmi593BtlK/bmZPmtl3zOzmci9kZveZ2VEzO9rT07OEckUWZ2A0t+jhFrh6VYwmF0mShAn0cl0bL3n+C2CHu98C/CXwzXIv5O4Pu3uXu3d1dnYurlKRJegfzdG+hEAvXBWjE6OSJGECvRvYVvR8K3CuuIG7D7n7cPD4MFBnZh2RVSmyRP0jk4u+ZBGKe+gacpHkCBPoR4DdZrbLzLLAPcCh4gZmtsmCQUoz2xu87uWoixVZrIGx3DIDXT10SY4FT/27+5SZPQA8BmSAR9z9uJndH+w/CPwB8C/NbAoYA+5x99JhGZEV5e6LvrlFwZpshvraGgbUQ5cECXUtVzCMcrhk28Gix58DPhdtaSLLMz45Q25qZkknRaEwW1Q9dEkOzRSV1BoYW9q0/4K2pqyGXCRRFOiSWoXhktY1Swz0xjqdFJVEUaBLahV61+uW2kNvVA9dkkWBLqk1GPSul3JSFPJDNTopKkmiQJfUGhgLhlyW0UMfGM0xM6MLtiQZFOiSWlfH0JfeQ5/x/G3sRJJAgS6pNTCaI1tbQ0Pd0r7NNblIkkaBLqk1MLq0lRYL2poKKy4q0CUZFOiSWgNjuSUPt8DVm2LoxKgkhQJdUqt/dHLJlyyChlwkeRToklqDwZDLUrVrxUVJGAW6pNZyh1xaGmqpMbSeiySGAl1Sa2B0aWuhF9TUGK2aLSoJokCXVBrLTTMxNbOsMXTQbFFJFgW6pFJhpcWlTvsv0HoukiQKdEml5a60WKAVFyVJFOiSSoVAX/6QS1Y3ipbEUKBLKhVCeDlXuUChh65Al2RQoEsqFVZaLEzfX6rWxizjkzOMT05HUZZIRYUKdDPbZ2YnzOykmR2Yp92vmdm0mf1BdCWKLN5yV1os0GxRSZIFA93MMsBDwF3AHuBeM9szR7tPA49FXaTIYg2MLW+lxYL2wgJdIzoxKvEX5rt9L3DS3U+5ew54FNhfpt2/Ab4GXIqwPpElGRiZpHXN0ldaLLi6QJd66BJ/YQJ9C3Cm6Hl3sG2WmW0B3gscnO+FzOw+MztqZkd7enoWW6tIaANjuWVfgw5Xh1z6FOiSAGECvVwXp/SeXP8N+Ji7z3vmyN0fdvcud+/q7OwMW6PIog0sc6XFgsLiXroWXZKgNkSbbmBb0fOtwLmSNl3Ao8Gftx3A3WY25e7fjKRKkUUaGJ1kx/rGZb/O7JCLFuiSBAgT6EeA3Wa2CzgL3AO8v7iBu+8qPDazLwHfVphLNQ2M5bilcd2yXydbW0NTNqMeuiTCgoHu7lNm9gD5q1cywCPuftzM7g/2zztuLlIN+dvPLX8MHTRbVJIjTA8ddz8MHC7ZVjbI3f2fL78skaUbn4xmpcWCtibNFpVk0ExRSZ3+iKb9F+RXXNSQi8SfAl1SpzBLdDm3nyumIRdJCgW6pE5UKy0WtGsJXUkIBbqkzuBYtEMurY1ZhsYnmZ4pnX4hEi8KdEmdQm96OfcTLdbWWIc7DI6ply7xpkCX1Lk6hh7RSdGmYPq/JhdJzCnQJXWiWmmxQAt0SVIo0CV1BkejWWmxQOu5SFIo0CV1+kdzkY2fg25yIcmhQJfUGRidnB0miULhl4OGXCTuFOiSOoNj+SGXqDTX11JbYxpykdhToEvqRD3kYmaaLSqJoECXVHF3+kcnZy81jEp7U53uKyqxp0CXVBnNTZObmqE9wjF0yF+6qJOiEncKdEmVwuSfqHvobY11sxOWROJKgS6pUuhFR91Db1MPXRJAgS6pUqkeemHIxV0LdEl8KdAlVWZ76BUYcpmcdkZy05G+rkiUQgW6me0zsxNmdtLMDpTZv9/MnjKzY2Z21MzeGn2pIgvrC65EqcSQC0C/FuiSGFsw0M0sAzwE3AXsAe41sz0lzf4BuMXdbwU+BHwx6kJFwugfyZGpMVoaQt0uN7Srs0V1YlTiK0wPfS9w0t1PuXsOeBTYX9zA3Yf96uBiE6CBRqmKvtEcbY111NREszBXQWFMXidGJc7CBPoW4EzR8+5g26uY2XvN7Hng78j30q9hZvcFQzJHe3p6llKvyLz6R3KRrYNe7OqKiwp0ia8wgV6uq3NND9zdv+HuNwHvAR4s90Lu/rC7d7l7V2dn5+IqFQmhbyQX+RUucHUMXUMuEmdhAr0b2Fb0fCtwbq7G7v4D4AYz61hmbSKL1j+ai/yEKMC6NeqhS/yFCfQjwG4z22VmWeAe4FBxAzN7jQV3EzCz24EscDnqYkUW0jcS/TouALWZGtY21KqHLrG24KUA7j5lZg8AjwEZ4BF3P25m9wf7DwK/D3zAzCaBMeCPXDMwZIXlF+bK0d4U3UqLxdqaNFtU4i3UtV3ufhg4XLLtYNHjTwOfjrY0kcUZGp9iesYrclIUCrNF1UOX+NJMUUmNwqSfqGeJFrQ11mlikcSaAl1So2+0Muu4FGiBLok7BbqkRt9wZVZaLGjVEroScwp0SY2+Ci3MVdDWmGV4Yorc1ExFXl9kuRTokhr9FVo6t6AwW3RgTMMuEk8KdEmNvtEc2UwNTdlMRV6/8ItCwy4SVwp0SY2+4RxtTXUEc9wipyV0Je4U6JIavcMTdDTXV+z1W2cX6FIPXeJJgS6pcXkkV9FAv7pAl3roEk8KdEmN3iuV7aHPDrmohy4xpUCXVHB3eodzdDRX5goXgDXZDPW1NeqhS2wp0CUVhsanyE3PVLSHDpotKvGmQJdUuDw8AUBHS+V66JA/MVq4EbVI3CjQJRV6g2n/le6hr2/O0jcyUdH3EFkqBbqkQm+hh17hQN/Q0kDPsAJd4kmBLqlQCPT1FTwpCtDZUs+loQl0/xaJIwW6pELvcA6zyq20WLChpZ6JqRmuTExV9H1ElkKBLqnQOzxBe2OW2kxlv6U7W/JDOpeGNOwi8aNAl1So9KSigkKg91xRoEv8hAp0M9tnZifM7KSZHSiz/4/N7Kng4ydmdkv0pYrMrXd4ouLj55AfcgG4dGW84u8lslgLBrqZZYCHgLuAPcC9ZranpNlLwNvc/Y3Ag8DDURcqMp9Kr+NS0NnSAKiHLvEUpoe+Fzjp7qfcPQc8CuwvbuDuP3H3/uDpz4Ct0ZYpMr+VGnJZ21BLtrZGgS6xFCbQtwBnip53B9vm8mHgO+V2mNl9ZnbUzI729PSEr1JkHsMTU4zkptm4tvKBbmZ0Ntcr0CWWwgR6ubsFlL0I18z+CflA/1i5/e7+sLt3uXtXZ2dn+CpF5nFxKD+evXFtw4q834a19VxSoEsMhQn0bmBb0fOtwLnSRmb2RuCLwH53vxxNeSILuzi4soGuHrrEVZhAPwLsNrNdZpYF7gEOFTcws+3A14E/cfcXoi9TZG4XrxQCvfJDLlDooesqF4mf2oUauPuUmT0APAZkgEfc/biZ3R/sPwh8ElgPfD64n+OUu3dVrmyRqy4M5nvLK9dDb6B/dJLc1AzZWk3lkPhYMNAB3P0wcLhk28Gixx8BPhJtaSLhXBwap6W+lqb6UN/Oy7Yh+Eugd3iC61rXrMh7ioSh7oUk3sWh8dmQXQmdzYXJRRpHl3hRoEviXRwaX7HhFrg6tHNhUOPoEi8KdEm8i0MTbFrBQL+uNf9e5wfHVuw9RcJQoEuizcw4l66Ms2EFA729KUt9bQ3n1UOXmFGgS6L1j+aYnPYVu2QR8rNFr2tdw9kB9dAlXhTokmiFXvJKDrkAbF7XwHkFusSMAl0SrdBL3tK2spcPbl63hnMDGnKReFGgS6Kd7Q8CfYWvB9/S2sClK+NMTc+s6PuKzEeBLol2dmCMhroa2psqf3OLYptb1zDjcFHXokuMKNAl0c72j7G1rZFgyYkVU5ghek7j6BIjCnRJtLMDYys+3AJw3br8SVgFusSJAl0S7ezA2IqfEIX8kAugE6MSKwp0SazR3BR9I7mq9NCb62tZ21CrHrrEigJdEqtwhcvWKvTQ8+/bSHf/aFXeW6QcBbokVvdAdQN9Z0cjL19WoEt8KNAlsbpnr0FvrMr771jfxJn+UV2LLrGhQJfEOt07wpq6zIqu41Js5/pGJqddi3RJbCjQJbFe6h1hZ0fTil+DXrBjfRMApy+PVOX9RUqFCnQz22dmJ8zspJkdKLP/JjP7qZlNmNmfRl+myLVO946wq6M6wy0AO2cDXePoEg8LBrqZZYCHgLuAPcC9ZranpFkf8G+Bz0ZeoUgZU9MzvNI3Ohuq1bChpZ6Guhpe7lUPXeIhTA99L3DS3U+5ew54FNhf3MDdL7n7EWCyAjWKXKO7f4ypGWdXR/UCvabG2NHepB66xEaYQN8CnCl63h1sWzQzu8/MjprZ0Z6enqW8hAgALwXj1tUMdIAd6xs1hi6xESbQy51x8qW8mbs/7O5d7t7V2dm5lJcQAeClnnyI7qxyoO/e2Mzp3hFyU7p0UaovTKB3A9uKnm8FzlWmHJFwXrh4hfamLOtXeNncUjdubGFqxjnVO1zVOkQgXKAfAXab2S4zywL3AIcqW5bI/J67cIXXbmyp2iWLBTdtWgvAiQtXqlqHCEDtQg3cfcrMHgAeAzLAI+5+3MzuD/YfNLNNwFFgLTBjZh8F9rj7UAVrl1VqZsZ58eIV3te1beHGFbaro4naGlOgSywsGOgA7n4YOFyy7WDR4wvkh2JEKu5M/yijuWlu2tRS7VLI1tZwQ2ezAl1iQTNFJXGeO58Pz5s2r61yJXk3bmrhxEUFulSfAl0S5/kLQ5jBjRubq10KADdtaqG7f4yhcU3DkOpSoEviPNU9yGs6m2nMhhoxrLg3bFkHwNPdg1WuRFY7Bbokirvzy1f6uW17a7VLmXXLtnwtv3ylv8qVyGqnQJdEefnyKP2jk9y2va3apcxat6aOGzqb+OUrA9UuRVY5Bbokyi/P5HvBceqhA9y2vY1jZwZwX9IkapFIKNAlUZ54uZ+mbIbdG6p/yWKx27a3cnkkp4W6pKoU6JIoPzl5mb272snUVHeGaKm33NABwI9O9la5ElnNFOiSGGcHxjjVO8Jbd8dvYbed6xvZ2raGH76gVUSlehTokhg/ejEflnfu7qhyJdcyM+7c3clPf3VZN42WqlGgS2I8fvwiW1rXsHtDPCYUlXrbjR1cmZjiH0/3VbsUWaUU6JIIg2OT/ODFHu56/aaqr7A4l7fduIHGbIa/ffJ8tUuRVUqBLonw2PELTE47v/3GzdUuZU5rshnedfMmDj99Xje8kKpQoEsifOXnr3BDZxO3bovX9eel3n3rdQyOTfL4sxeqXYqsQgp0ib1jZwZ48swAf/LmHbEdbin4zd2d7FzfyBd+cEqTjGTFKdAl9v7iuy+wbk0dv39H/Jfcz9QYH77zep7sHuTHJy9XuxxZZRToEmvfO3GJ77/Qw796+w20NNRVu5xQ/vCOrWxtW8N/+tvjTOoSRllBCnSJrcvDExz42tO8ZkMzH3zLzmqXE1pDXYZP/e7NvHhpmM8+fqLa5cgqokCXWOofyfGhLx+lfzTHX7zvVhrqMtUuaVHesWcjf/ym7fz190/x5Z+crnY5skqECnQz22dmJ8zspJkdKLPfzOy/B/ufMrPboy9VVosfvdjLez7/Y547N8Tn3n87b9i6rtolLcmn3n0z73jdBv7joeN87P88Rf9IrtolScoteMsXM8sADwHvBLqBI2Z2yN2fLWp2F7A7+HgT8FfBvyLzGstN0zs8wcuXRzl2pp+/P36BZ84Osb29ka/e9ybu2NFe7RKXrC5Tw8F/dgefffwFvvDDU3zrybPsu3kTv/GaDl63eS1bWtfQ2lgX+yt3JDnC3MNrL3DS3U8BmNmjwH6gOND3A//T89dp/czMWs1ss7tHPmXu+y/08OC3r751uUvDrtni8+8P8xqlTbykRbkr1Ba6aq30fcs1X+z7hnmN0lbhXmP+Wq95jRDHdGraGZucftW2129Zy5/tv5n3dW1L3DBLObWZGg7cdRO/d/sWvvyT0/zd0+f55rFzs/trDNbUZViTraWhrgYzMCz4N79GjAGUPpdE+6Nf28ZH7rw+8tcNE+hbgDNFz7u5tvddrs0W4FWBbmb3AfcBbN++fbG1AtBcX8trN5ashV3mO7x0U2kv6Nr9y3+N8nWUfI6V7g9TxwKvEaKQxb5vudhY+DUWjpriJhkz2pqydDbXc13rGt6wZR3rGpNxJcti3bixhf/83jfw4P7Xc7JnmFM9w3T3jzEwOsnY5DRjk9OMT06D53/xzbjjwWP34Ne4X/sLXZKpo7m+Iq8bJtDL/ZSWfleFaYO7Pww8DNDV1bWk78w7drRxx4743H5MZDFqaowbN7ZwY2mnRCQCYU6KdgPbip5vBc4toY2IiFRQmEA/Auw2s11mlgXuAQ6VtDkEfCC42uXNwGAlxs9FRGRuCw65uPuUmT0APAZkgEfc/biZ3R/sPwgcBu4GTgKjwL+oXMkiIlJOmDF03P0w+dAu3naw6LED/zra0kREZDE0U1REJCUU6CIiKaFAFxFJCQW6iEhKWLXuqmJmPcDLS/z0DqA3wnKiEte6IL61qa7FUV2Lk8a6drh7Z7kdVQv05TCzo+7eVe06SsW1LohvbaprcVTX4qy2ujTkIiKSEgp0EZGUSGqgP1ztAuYQ17ogvrWprsVRXYuzqupK5Bi6iIhcK6k9dBERKaFAFxFJidgGupn9oZkdN7MZM+sq2ffx4IbUJ8zsXXN8fruZfdfMXgz+jfyuGGb2v8zsWPBx2syOzdHutJk9HbQ7GnUdZd7vU2Z2tqi2u+doN+/NvytU22fM7PngZuLfMLPWOdpV/JjF8ebnZrbNzP6fmT0XfP//uzJt3m5mg0Vf309Wuq6i957361KlY/baomNxzMyGzOyjJW1W5JiZ2SNmdsnMninaFiqLIvl5dPdYfgCvA14LfA/oKtq+B3gSqAd2Ab8CMmU+/78AB4LHB4BPV7jePwc+Oce+00DHCh67TwF/ukCbTHDsrgeywTHdswK1/VOgNnj86bm+LpU+ZmH+/+SXhP4O+TtyvRn4+Qocn83A7cHjFuCFMnW9Hfj2Sn0/LebrUo1jVubreoH85JsVP2bAbwK3A88UbVswi6L6eYxtD93dn3P3E2V27QcedfcJd3+J/Brse+do9+Xg8ZeB91Sm0nyvBHgf8NVKvUcFzN78291zQOHm3xXl7o+7+1Tw9Gfk725VDWH+/7M3P3f3nwGtZra5kkW5+3l3/0Xw+ArwHPn78ybFih+zEr8F/MrdlzoLfVnc/QdAX8nmMFkUyc9jbAN9HnPdkLrURg/umhT8u6GCNd0JXHT3F+fY78DjZvZEcKPslfBA8CfvI3P8iRf2OFbSh8j35sqp9DEL8/+v6jEys53AbcDPy+z+dTN70sy+Y2Y3r1RNLPx1qfb31T3M3bGq1jELk0WRHLdQN7ioFDP7v8CmMrs+4e7fmuvTymyr2LWXIWu8l/l757/h7ufMbAPwXTN7PvhNXpG6gL8CHiR/XB4kPxz0odKXKPO5kRzHMMfMzD4BTAFfmeNlIj9mpWWW2bakm59Xgpk1A18DPuruQyW7f0F+SGE4OD/yTWD3StTFwl+Xah6zLPBu4ONldlfzmIURyXGraqC7+zuW8Glhb0h90cw2u/v54E++S5Wo0cxqgd8D7pjnNc4F/14ys2+Q//NqWeEU9tiZ2ReAb5fZVbEbe4c4Zh8Efgf4LQ8GEMu8RuTHrERsb35uZnXkw/wr7v710v3FAe/uh83s82bW4e4VX4QqxNelmjeMvwv4hbtfLN1RzWNGuCyK5LglccjlEHCPmdWb2S7yv2X/cY52HwwefxCYq8e/XO8Annf37nI7zazJzFoKj8mfFHymXNuolIxZvneO9wtz8+9K1LYP+BjwbncfnaPNShyzWN78PDgf8z+A59z9v87RZlPQDjPbS/7n+HIl6wreK8zXpZo3jJ/zL+VqHbNAmCyK5uex0md9l/pBPoi6gQngIvBY0b5PkD8jfAK4q2j7FwmuiAHWA/8AvBj8216hOr8E3F+y7TrgcPD4evJnrJ8EjpMfdqj0sfsb4GngqeCbYrTYGwwAAACgSURBVHNpXcHzu8lfRfGrlagreM+T5McKjwUfB6t1zMr9/4H7C19P8n8GPxTsf5qiq60qeHzeSv5P7aeKjtHdJXU9EByXJ8mfWH7LCn3tyn5dqn3MgvdtJB/Q64q2rfgxI/8L5TwwGeTXh+fKokr8PGrqv4hISiRxyEVERMpQoIuIpIQCXUQkJRToIiIpoUAXEUkJBbqISEoo0EVEUuL/A3Ar/hj8WM8fAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Distribution on encrypted data:\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3deXCc933f8fd3F8DiJkhcFG9KgiRTtmRLMG3HluP4qCUlY/qIGykZ243jqmqrup5OJlbGrSdTt5NR3LSZNHIYxdUkTTNVPPXFOrTko3GUjC2blKyLlChR1EGIJE6SuBdY4Ns/dhdaLRfEAtzd58DnNYPh7j7PPvvlA+CzP/z293t+5u6IiEj0JYIuQEREKkOBLiISEwp0EZGYUKCLiMSEAl1EJCbqgnrhrq4u37VrV1AvLyISSY8++uiIu3eX2hZYoO/atYvDhw8H9fIiIpFkZi8vt01dLiIiMaFAFxGJCQW6iEhMKNBFRGJCgS4iEhMKdBGRmFCgi4jEhAJd1p0Hnz7N1x8dQJeOlrgJbGKRSBC+9fNX+dzfPA7AyGSaf/GLVwRckUjllNVCN7ObzeyYmR03s7uX2ec9Zva4mR0xs7+vbJkil25x0fmjHzzHtVvaefdV3dz7d8eZnV8IuiyRilkx0M0sCdwL3ALsAW43sz1F+3QAXwE+5O7XAh+vQq0il+TnJ8/x0ug0n7lpN3fcdDnjsxl+dGw46LJEKqacFvpe4Li7n3D3OeABYF/RPr8OfMPdXwFw96HKlily6X7wzCB1CeO91/Ty9ss3sbG5nu8dPRN0WSIVU06gbwVOFtwfyD1W6Cpgo5n9yMweNbNPljqQmd1hZofN7PDwsFpGUls/eWGUG3ZsZENTPXXJBG+/vJNDL40FXZZIxZQT6FbiseLhAXXAjcAvAx8E/oOZXXXBk9zvc/d+d+/v7i559UeRqkhnFjh6apy37OxYeuzGnRs5OTbD4PhsgJWJVE45gT4AbC+4vw04VWKfB919yt1HgIeB6ytTosile/rVceYWFrlhx8alx/p3bQLg8EtngypLpKLKCfRDQJ+Z7TazBuA24EDRPt8GbjKzOjNrBt4GPFPZUkXW7smBcwC8eftrLfQ9l7VTlzCOnj4fVFkiFbXiOHR3z5jZXcBDQBK4392PmNmdue373f0ZM3sQeBJYBL7q7k9Xs3CR1XhucIKNzfX0tKWWHmuoS7C7q4VjZyYCrEykcsqaWOTuB4GDRY/tL7r/ZeDLlStNpHKeH5ykr7cNs9d/JHT15jYeP3kuoKpEKktT/yX23J3nBifo62m9YNs1m9sYODvDZDoTQGUilaVAl9gbnkgzPpspGeh9vW0AvDA0WeuyRCpOgS6x93wurK/KhXeh3V0tALw0OlXTmkSqQYEusXdiJBvWu7tbLti2Y1MzAC+PTte0JpFqUKBL7A2MTdNQl6C3rfGCbY31SbZsaFQLXWJBgS6xd/LsNNs6mkgkSk16hp2dLbw0okCX6FOgS+ydHJthW65rpZRdXc3qcpFYUKBL7J08O832jU3Lbt/Z2cLo1Bzjs/M1rEqk8hToEmvjs/Ocm55n+0Va6NtyYX/q3EytyhKpCgW6xNrJsWxXyvaNywf6lg4FusSDAl1i7eRYNqS3b1q+y2XrUqDrMroSbQp0ibWBsyu30LtbU9QnTS10iTwFusTawNkZWlN1dDTXL7tPImFs3tCoQJfIU6BLrJ05P8vmDY0XXGWx2JYNTepykchToEusnRmfZXP7hTNEi23taOJVtdAl4hToEmuD47P0lhHoWzqaODM+y8Ji8XK5ItGhQJfYWlh0hibSbN6QWnHfLR1NLCy6FoyWSFOgS2yNTKZZWPSyulzyoa9AlyhToEtsnTmfDedyulx6cldiHJpIV7UmkWpSoEtsncm1tjdvKCPQ27Mt9CG10CXCFOgSW4OrCPTOlhTJhDE4rha6RJcCXWLrzPlZ6hJGV8vKH4omE0Z3a4qhCbXQJboU6BJbZ8Zn6WlLLbuwRbGe9pRa6BJpZQW6md1sZsfM7LiZ3V1i+3vM7LyZPZ77+mLlSxVZncHxWXrL6G7J62lr1IeiEml1K+1gZkngXuADwABwyMwOuPvRol3/wd1/pQo1iqzJmfOzXNXbVvb+Pe0pfv7K2SpWJFJd5bTQ9wLH3f2Eu88BDwD7qluWyKUbGk+XNWQxr7etkdGpOeYyi1WsSqR6ygn0rcDJgvsDuceKvcPMnjCz75rZtaUOZGZ3mNlhMzs8PDy8hnJFyjM7v8BEOkN328ofiOblhy6OTKrbRaKpnEAv9YlS8QUvHgN2uvv1wH8HvlXqQO5+n7v3u3t/d3f36ioVWYWxqTkAOlsayn5Ob7tmi0q0lRPoA8D2gvvbgFOFO7j7uLtP5m4fBOrNrKtiVYqs0uhkLtBbV9FC12xRibhyAv0Q0Gdmu82sAbgNOFC4g5ltttwFp81sb+64o5UuVqRcI1PZUO5sLb+FrtmiEnUrjnJx94yZ3QU8BCSB+939iJndmdu+H/hV4F+aWQaYAW5zd12HVAKTb6GXM6kor7MlRcLQWHSJrBUDHZa6UQ4WPba/4PafAH9S2dJE1m50cvUt9GTC2NSSYnRKgS7RpJmiEkujU3M01idobkiu6nldrQ0MT8xVqSqR6lKgSyyNTKbpbEmtuJZose42tdAluhToEkujk3N0raK7Ja+zpUHj0CWyFOgSS6NT6VUNWczrak0xoi4XiSgFusTS6OTcqiYV5XW1pZiZX2B6LlOFqkSqS4EusePu2UBfQws9/yagVrpEkQJdYmd8NsPcwuKa+tC7ctd+GVY/ukSQAl1iJz8GvWstfei5iUijCnSJIAW6xM5o/sJca2qh57pcJtXlItGjQJfYWZoluopp/3mdaqFLhCnQJXbyreu19KE31CVob6zTWHSJJAW6xE7+wlwb1zBsEbIfjKrLRaJIgS6xMzqVpqO5nvrk2n68u1pTaqFLJCnQJXbWOqkor6tV0/8lmhToEjsjk2ub9p+XbaGry0WiR4EusTM6tbYLc+V1taY4PzPPXGaxglWJVJ8CXWJnNHfp3LXKj1/PLzQtEhUKdImVzMIiZ6fn1zSpKC8/w1T96BI1CnSJlbHp/CzRS+tDBwW6RI8CXWIlf5XErksc5QKa/i/Ro0CXWMkvH5e/auJa5Fv3Y1qKTiJGgS6xkp8leinj0FsakqTqEkvHEokKBbrESr7f+1L60M2MrtaUrokukVNWoJvZzWZ2zMyOm9ndF9nvrWa2YGa/WrkSRco3OjVHfdJob6y7pON0tTaohS6Rs2Kgm1kSuBe4BdgD3G5me5bZ7x7goUoXKVKu/Bh0M7uk43S2ppb640WiopwW+l7guLufcPc54AFgX4n9/g3wdWCogvWJrEp2LdG195/ndbaohS7RU06gbwVOFtwfyD22xMy2Ah8B9l/sQGZ2h5kdNrPDw8PDq61VZEUjU2tbHLpYZ2uK0ck53L0CVYnURjmBXupv1+Kf8j8CPu/uCxc7kLvf5+797t7f3d1dbo0iZRudTF/SGPS8rtYG5hYWmUhnKlCVSG2U88nRALC94P424FTRPv3AA7l+yy7gVjPLuPu3KlKlSJkq1uWSO8bo5BztjfWXfDyRWiinhX4I6DOz3WbWANwGHCjcwd13u/sud98F/B/gXynMpdam5zLMzC9UpstFa4tKBK3YQnf3jJndRXb0ShK4392PmNmdue0X7TcXqZVKTCrK69T0f4mgsgbruvtB4GDRYyWD3N3/2aWXJbJ6+UlFXRVooeePoaGLEiWaKSqxsdRCr0Af+sbm1/rQRaJCgS6xkW9NV6IPvaEuwYamevWhS6Qo0CU2RirYhw75xaLVQpfoUKBLbIxMpmlL1dFYn6zI8TpbU1rkQiJFgS6xUakx6HldrQ2Mal1RiRAFusTG6FS6Iv3neZ0tKfWhS6Qo0CU2RifnKtZ/DtnRMmen58ksLFbsmCLVpECX2BiZrMyFufKWlqKbVreLRIMCXWJhcdEZm0ovLfBcCfmLfGksukSFAl1i4dzMPIteuSGL8FoLXYEuUaFAl1gYrcBaosWWrrio6f8SEQp0iYWRCk77z+vKXXFRk4skKhToEgv5VnQlLsyV195UR13CNHRRIkOBLrFQyUvn5pkZna1aW1SiQ4EusTA6mSZh0NFcuUCHbItf0/8lKhToEgsjU3NsakmRTJRaAnftOltTjGj6v0SEAl1iYWSismPQ87paGtSHLpGhQJdYGJ2q7IW58tSHLlGiQJdYGJ1MLy3sXEmdrSlm5heYnstU/NgilaZAl1io9KVz8zo1/V8iRIEukTc7v8BEOlPRMeh5+WNqpItEgQJdIm8sNwplUwXHoOctTf9XC10iQIEukZdvPVejhb50gS5dz0UiQIEukZdvPVdj2GK+D13Xc5EoKCvQzexmMztmZsfN7O4S2/eZ2ZNm9riZHTazd1W+VJHShqvYQm+sT9KaqlOXi0RC3Uo7mFkSuBf4ADAAHDKzA+5+tGC3HwIH3N3N7Drga8A11ShYpFg1u1wgNxZdXS4SAeW00PcCx939hLvPAQ8A+wp3cPdJd/fc3RbAEamR0ck5WhqSNDUkq3L8zpYGjXKRSCgn0LcCJwvuD+Qeex0z+4iZPQv8LfDpUgcysztyXTKHh4eH11KvyAVGJtMVXdiiWFdrSl0uEgnlBHqpqx1d0AJ392+6+zXAh4EvlTqQu9/n7v3u3t/d3b26SkWWMTJZneu45HW2pvShqERCOYE+AGwvuL8NOLXczu7+MHCFmXVdYm0iZRmdnKta/zlkR8+MTaVZXFRPooRbOYF+COgzs91m1gDcBhwo3MHMrjQzy92+AWgARitdrEgp1e5y6WxpYNGzC1GLhNmKo1zcPWNmdwEPAUngfnc/YmZ35rbvBz4GfNLM5oEZ4NcKPiQVqZqFRWdsao7uKne5QPYCYNWYjSpSKSsGOoC7HwQOFj22v+D2PcA9lS1NZGVnp+dYdKrbQm99bXJRX2/VXkbkkmmmqERatcegFx5bY9El7BToEmnVnPafp0voSlQo0CXS8i30ana5dDQ3kDC0FJ2EngJdIm14Ihuy3VUM9GTC2NTSoMWiJfQU6BJpo1Nz1CeN9qayPt9fs86WFCMTaqFLuCnQJdJGJrJrieamQVRN9gJdaqFLuCnQJdJGJtN0tVV/bHj2ei5qoUu4KdAl0kan5uhsqV7/eV5na4NGuUjoKdAl0kYm0lUdg57X1ZpiIp1hdn6h6q8lslYKdIksd2dkaq4mXS75sehj6keXEFOgS2RNpDPMZRbpqkmXS/56Lgp0CS8FukRWfhhhTVro+eu5aPq/hJgCXSIrv+hELT4Uzf8VoBa6hJkCXSJraGIWgJ722oxyAU3/l3BToEtkDY1nw7W3rbHqr9XckKSxPqHJRRJqCnSJrKGJNA3JBB3N9VV/LTPLTv9XC11CTIEukTU0Pkt3W/Wn/ed1tTZosWgJNQW6RNbQRLom/ed5na26QJeEmwJdImtwfJaettoFek9biiEFuoSYAl0ia2giTU8NPhDN621vZHQqzfzCYs1eU2Q1FOgSSbPzC5yfmae3hl0uve2NuKMPRiW0FOgSSfmVimrbQs++eZw5P1uz1xRZDQW6RFItJxXl9bZn3zwGx9VCl3AqK9DN7GYzO2Zmx83s7hLbf8PMnsx9/djMrq98qSKvyU8qqnUfOrz2ZiISNisGupklgXuBW4A9wO1mtqdotxeBX3T364AvAfdVulCRQoPjtW+hd7Y0UJcwdblIaJXTQt8LHHf3E+4+BzwA7Cvcwd1/7O5nc3cfAbZVtkyR1xuaSFOXMDY1V/9Ki3mJhNHTllKXi4RWOYG+FThZcH8g99hyfgv4bqkNZnaHmR02s8PDw8PlVylSZHA8TXdbikSiNrNE83raG5f+OhAJm3ICvdRvjJfc0eyXyAb650ttd/f73L3f3fu7u7vLr1KkyNBEbScV5fW2pxToElrlBPoAsL3g/jbgVPFOZnYd8FVgn7uPVqY8kdKGJ9L0tNfuA9G8zWqhS4iVE+iHgD4z221mDcBtwIHCHcxsB/AN4BPu/lzlyxR5vVpP+8/raW9kfDbDzJwWi5bwqVtpB3fPmNldwENAErjf3Y+Y2Z257fuBLwKdwFdyV77LuHt/9cqW9WxmboGz0/Ns6Wiq+WtvXhqLPsuurpaav77IxawY6ADufhA4WPTY/oLbnwE+U9nSREo7fX4GgMs21L7LpVeBLiGmmaISOadz48Av21D7FvrS9H/1o0sIKdAlcl49l22hb+kIoIW+4bUWukjYKNAlck6fy4bp5gC6XNpSdbQ0JJf+ShAJEwW6RM7p8zN0tTaQqkvW/LXNjC0dTbx6dqbmry2yEgW6RM6p87OB9J/nbd3YxKnzCnQJHwW6RM7pczOB9J/nqYUuYaVAl8g5HXQLvaOJs9PzTM9lAqtBpBQFukTK+Ow8k+lMoC30bRuzbyanzqmVLuGiQJdIyY9wCbKFnp+hOqBuFwkZBbpESv7DyCBb6Fs78i10DV2UcFGgS6SEoYXe05YimTBePTcdWA0ipSjQJVJOnp2mLmFL11QJQl0yweb2RrXQJXQU6BIpr4xNs21jE8kar1RUbKuGLkoIKdAlUk6OTbN9U3PQZbB1Y9PSNWVEwkKBLpHyytg0O8IQ6B1NnBmfJbOwGHQpIksU6BIZ52fmOTc9H4pA37GpmYVFVz+6hIoCXSLj5Fh2VMnOzuADPV/DS6NTAVci8hoFukRGPtDD0Ie+O7dakQJdwkSBLpHxSogCvbstRVN9kpdGNBZdwkOBLpHxytg0G5vraW+sD7oUzIydnc28rBa6hIgCXSIjLCNc8nZ1tqjLRUJFgS6R8dLoFDs6W4IuY8murhZOjs2wsOhBlyICKNAlImbnFxg4O8MV3SEK9M5m5hYWOa3ViyQkygp0M7vZzI6Z2XEzu7vE9mvM7Cdmljaz3658mbLevTgyhTtc0d0adClLdub+Wnh5VB+MSjisGOhmlgTuBW4B9gC3m9meot3GgM8C/6XiFYoALwxPAnBlT3gCPT908USuNpGgldNC3wscd/cT7j4HPADsK9zB3Yfc/RAwX4UaRTg+NInZayEaBr3tKdpSdTw3qECXcCgn0LcCJwvuD+QeE6mZF4an2Laxicb6ZNClLDEzrtrcxnODE0GXIgKUF+ilrlO6po/1zewOMztsZoeHh4fXcghZp14YmgxV/3neVb2tPDc4gbtGukjwygn0AWB7wf1twKm1vJi73+fu/e7e393dvZZDyDq0uOicGJnkyhAGel9PG2en5xmZnAu6FJGyAv0Q0Gdmu82sAbgNOFDdskRe88rYNLPzi/T1hi/Qr97cBsDz6naREFgx0N09A9wFPAQ8A3zN3Y+Y2Z1mdieAmW02swHg3wH/3swGzKy9moXL+nH09DgA127ZEHAlF8q/yRxToEsI1JWzk7sfBA4WPba/4PYZsl0xIhV39NQ4dQkL1ZDFvO7WFBub6zXSRUJBM0Ul9I6eHufKntZQjXDJMzOu6m3j2TPjQZciokCX8Dt6apw3XBbeHrw3bt3A0VPjzGs5OgmYAl1CbXQyzZnxWfaEONCv27aBdGZR49ElcAp0CbUjp7JdGXu2hDfQr9/WAcBTA+cDrkTWOwW6hNpjr5zFLNsKDqudnc20N9bxhAJdAqZAl1B77JVzXN3bRlsIVilajplx3bYOnhw4F3Qpss4p0CW0Fhedn79ylrfs2Bh0KSu6fvsGjp2ZYCqdCboUWccU6BJax4cnmZjNcMOOjqBLWdHbdneSWXQeffls0KXIOqZAl9D66YtjAPTv2hRwJSvr37WRuoTxkxOjQZci65gCXULrH58fZmtHE7s6w7Mw9HKaG+q4fnsHP3lBgS7BUaBLKGUWFvnxC6O868ouzEpdwTl83nF5J0+9ep5J9aNLQBToEkpPvnqeidkM7+rrCrqUsv3ClZ0sLDr/+PxI0KXIOqVAl1D6f88MkTB455XRCfS37trEhqZ6vnfkTNClyDqlQJfQcXcOPnWad1zRyaaWhqDLKVt9MsH73tDDD58d0nVdJBAKdAmd5wYnOTEyxc1vvCzoUlbtg9du5vzMPD/LjdARqSUFuoTOtx9/lYTBB6/tDbqUVXt3XzetqTq+/thA0KXIOqRAl1CZyyzytcMDvPeaHnraGoMuZ9WaGpJ86M1bOPjUac7PzAddjqwzCnQJle8fHWRkMs1vvG1n0KWs2e1v3cHs/CLfVCtdakyBLqHh7nzlR8fZ2dnMu6/qDrqcNXvj1nb6d27kzx4+QTqzEHQ5so4o0CU0HjoyyJFT43z2vX0kE9GYTFSKmfHZ9/Vx+vwsXzusVrrUjgJdQmF6LsN/+tujXNHdwr43bwm6nEt2U18Xe3dv4g+/d4zRyXTQ5cg6oUCXUPj9g88ycHaG3//oddQlo/9jaWb85w+/kal0hi9882ncPeiSZB2I/m+ORN7/euRl/uqRl/nnN+1m7+7wX1mxXH29bfzOB6/hwSNn+G8/eD7ocmQdqAu6AFm/FhedP3v4BPc8+Cy/dHU3n7/5mqBLqrjP3LSbY4MT/PEPn2dmLsPv3HwN9TH4C0TCqayfLDO72cyOmdlxM7u7xHYzsz/ObX/SzG6ofKkSJ4++PMavf/UR7nnwWX75TZex/xM3xqKrpZiZ8Qcfu45PvH0nf/4PL/LRr/yYh58bVheMVMWKLXQzSwL3Ah8ABoBDZnbA3Y8W7HYL0Jf7ehvwp7l/ZZ2byyxybnqO4ck0zw9OcuTUeX747BAnhqfY1NLA73/0Tdz21u2RuUTuWiQSxpc+/EbecUUn//H/HuWT9/+MbRubeP8berlu2wb6etroaU/R2dIQyzc1qZ1yulz2Asfd/QSAmT0A7AMKA30f8D892+x4xMw6zOwydz9d6YL//rlhvvSd1166uKXjy9654O7rnnvhtuLn+vLbLtLYumh9Rc/1oq0X1nCx11zNc1dR0wr74stvW1h0ZuZfPw67Pmns3b2J33znbj52w1aaG9ZPr9+tb7qM972hhwOPn+LgU6d54NAr/MWPX38Rr4a6BI11CVL1SRqSCRIJMAwzyL/lmVn2du6xpfsSGb/21u185qbLK37ccn6btgInC+4PcGHru9Q+W4HXBbqZ3QHcAbBjx47V1gpAa6qOq3vbXv+gLX+3uOVX/INfuPnCbRd57gWvWbTvRY+7/HMv2HbBb2qVXueC5y4fEeUe14ANTfVsbGlgU0sDV3S3srurhYa69dsKTdUl+Xj/dj7ev53MwiIvjkzx4sgUQxNpRibTzMwvkJ5fJJ3J/utk31Tzb5buvP4xv7AhIOHX1ZqqynHLCfRSv9nFP0Hl7IO73wfcB9Df37+mn8Ibd27kxp3hXwVeZCV1yQR9vW30FTdQRNaonKbSALC94P424NQa9hERkSoqJ9APAX1mttvMGoDbgANF+xwAPpkb7fJ24Hw1+s9FRGR5K3a5uHvGzO4CHgKSwP3ufsTM7sxt3w8cBG4FjgPTwG9Wr2QRESmlrCEG7n6QbGgXPra/4LYD/7qypYmIyGqs3+EGIiIxo0AXEYkJBbqISEwo0EVEYsKCukiQmQ0DL6/x6V3ASAXLqZSw1gXhrU11rY7qWp041rXT3Uuu0RhYoF8KMzvs7v1B11EsrHVBeGtTXaujulZnvdWlLhcRkZhQoIuIxERUA/2+oAtYRljrgvDWprpWR3WtzrqqK5J96CIicqGottBFRKSIAl1EJCZCG+hm9nEzO2Jmi2bWX7Ttd3MLUh8zsw8u8/xNZvZ9M3s+92/FV8Uws78xs8dzXy+Z2ePL7PeSmT2V2+9wpeso8Xq/Z2avFtR26zL7XXTx7yrV9mUzeza3mPg3zaxjmf2qfs7CuPi5mW03s78zs2dyP///tsQ+7zGz8wXf3y9Wu66C177o9yWgc3Z1wbl43MzGzexzRfvU5JyZ2f1mNmRmTxc8VlYWVeT30d1D+QW8Abga+BHQX/D4HuAJIAXsBl4AkiWe/wfA3bnbdwP3VLnePwS+uMy2l4CuGp673wN+e4V9krlzdznQkDune2pQ2z8B6nK371nu+1Ltc1bO/5/sJaG/S3ZFrrcDP63B+bkMuCF3uw14rkRd7wG+U6ufp9V8X4I4ZyW+r2fITr6p+TkD3g3cADxd8NiKWVSp38fQttDd/Rl3P1Zi0z7gAXdPu/uLZK/BvneZ/f4yd/svgQ9Xp9JsqwT4p8D/rtZrVMHS4t/uPgfkF/+uKnf/nrtncncfIbu6VRDK+f8vLX7u7o8AHWZ2WTWLcvfT7v5Y7vYE8AzZ9XmjoubnrMj7gBfcfa2z0C+Juz8MjBU9XE4WVeT3MbSBfhHLLUhdrNdzqybl/u2pYk03AYPu/vwy2x34npk9mlsouxbuyv3Je/8yf+KVex6r6dNkW3OlVPuclfP/D/Qcmdku4C3AT0tsfoeZPWFm3zWza2tVEyt/X4L+ubqN5RtWQZ2zcrKoIuetrAUuqsXMfgBsLrHpC+7+7eWeVuKxqo29LLPG27l46/yd7n7KzHqA75vZs7l38qrUBfwp8CWy5+VLZLuDPl18iBLPrch5LOecmdkXgAzw18scpuLnrLjMEo+tafHzajCzVuDrwOfcfbxo82NkuxQmc5+PfAvoq0VdrPx9CfKcNQAfAn63xOYgz1k5KnLeAg10d3//Gp5W7oLUg2Z2mbufzv3JN1SNGs2sDvgocONFjnEq9++QmX2T7J9XlxRO5Z47M/tz4DslNlVtYe8yztmngF8B3ue5DsQSx6j4OSsS2sXPzayebJj/tbt/o3h7YcC7+0Ez+4qZdbl71S9CVcb3JcgF428BHnP3weINQZ4zysuiipy3KHa5HABuM7OUme0m+y77s2X2+1Tu9qeA5Vr8l+r9wLPuPlBqo5m1mFlb/jbZDwWfLrVvpRT1WX5kmdcrZ/HvatR2M/B54EPuPr3MPrU4Z6Fc/Dz3ecz/AJ5x9/+6zD6bc/thZnvJ/h6PVrOu3GuV830JcsH4Zf9SDuqc5ZSTRZX5faz2p75r/SIbRANAGhgEHirY9gWynwgfA24pePyr5EbEAJ3AD4Hnc/9uqlKdf4TQndsAAADBSURBVAHcWfTYFuBg7vblZD+xfgI4Qrbbodrn7q+Ap4Ancz8UlxXXlbt/K9lRFC/Uoq7cax4n21f4eO5rf1DnrNT/H7gz//0k+2fwvbntT1Ew2qqK5+ddZP/UfrLgHN1aVNddufPyBNkPln+hRt+7kt+XoM9Z7nWbyQb0hoLHan7OyL6hnAbmc/n1W8tlUTV+HzX1X0QkJqLY5SIiIiUo0EVEYkKBLiISEwp0EZGYUKCLiMSEAl1EJCYU6CIiMfH/Af+E0IWtdQ3qAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "normal_dist = lambda x, mean, var: np.exp(- np.square(x - mean) / (2 * var)) / np.sqrt(2 * np.pi * var)\n",
    "\n",
    "def plot_normal_dist(mean, var, rmin=-10, rmax=10):\n",
    "    x = np.arange(rmin, rmax, 0.01)\n",
    "    y = normal_dist(x, mean, var)\n",
    "    fig = plt.plot(x, y)\n",
    "    \n",
    "# plain distribution\n",
    "lr = LR(n_features)\n",
    "data = lr.lr(x_test)\n",
    "mean, var = map(float, [data.mean(), data.std() ** 2])\n",
    "plot_normal_dist(mean, var)\n",
    "print(\"Distribution on plain data:\")\n",
    "plt.show()\n",
    "\n",
    "# encrypted distribution\n",
    "def encrypted_out_distribution(eelr, enc_x_test):\n",
    "    w = eelr.weight\n",
    "    b = eelr.bias\n",
    "    data = []\n",
    "    for enc_x in enc_x_test:\n",
    "        enc_out = enc_x.dot(w) + b\n",
    "        data.append(enc_out.decrypt())\n",
    "    data = torch.tensor(data)\n",
    "    mean, var = map(float, [data.mean(), data.std() ** 2])\n",
    "    plot_normal_dist(mean, var)\n",
    "    print(\"Distribution on encrypted data:\")\n",
    "    plt.show()\n",
    "\n",
    "eelr = EncryptedLR(lr)\n",
    "eelr.encrypt(ctx_training)\n",
    "encrypted_out_distribution(eelr, enc_x_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Most of the data falls into $[-5,5]$, the sigmoid approximation should be good enough!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We finally reached the last part, which is about training an encrypted logistic regression model on encrypted data! You can see that we decrypt the weights and re-encrypt them again after every epoch, this is necessary since after updating the weights at the end of the epoch, we can no longer use them to perform enough multiplications, so we need to get them back to the initial ciphertext level. In a real scenario, this would translate to sending the weights back to the secret-key holder for decryption and re-encryption. In that case, it will result in just a few Kilobytes of communication per epoch."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Accuracy at epoch #0 is 0.523952066898346\n",
      "Accuracy at epoch #1 is 0.6796407103538513\n",
      "Accuracy at epoch #2 is 0.6796407103538513\n",
      "Accuracy at epoch #3 is 0.7005987763404846\n",
      "Accuracy at epoch #4 is 0.6916167736053467\n",
      "Accuracy at epoch #5 is 0.703592836856842\n",
      "\n",
      "Average time per epoch: 160 seconds\n",
      "Final accuracy is 0.703592836856842\n",
      "Difference between plain and encrypted accuracies: 0.0\n"
     ]
    }
   ],
   "source": [
    "eelr = EncryptedLR(LR(n_features))\n",
    "accuracy = eelr.plain_accuracy(x_test, y_test)\n",
    "print(f\"Accuracy at epoch #0 is {accuracy}\")\n",
    "\n",
    "times = []\n",
    "for epoch in range(EPOCHS):\n",
    "    eelr.encrypt(ctx_training)\n",
    "    \n",
    "    # if you want to keep an eye on the distribution to make sure\n",
    "    # the function approximation is still working fine\n",
    "    # WARNING: this operation is time consuming\n",
    "    # encrypted_out_distribution(eelr, enc_x_train)\n",
    "    \n",
    "    t_start = time()\n",
    "    for enc_x, enc_y in zip(enc_x_train, enc_y_train):\n",
    "        enc_out = eelr.forward(enc_x)\n",
    "        eelr.backward(enc_x, enc_out, enc_y)\n",
    "    eelr.update_parameters()\n",
    "    t_end = time()\n",
    "    times.append(t_end - t_start)\n",
    "    \n",
    "    eelr.decrypt()\n",
    "    accuracy = eelr.plain_accuracy(x_test, y_test)\n",
    "    print(f\"Accuracy at epoch #{epoch + 1} is {accuracy}\")\n",
    "\n",
    "\n",
    "print(f\"\\nAverage time per epoch: {int(sum(times) / len(times))} seconds\")\n",
    "print(f\"Final accuracy is {accuracy}\")\n",
    "\n",
    "diff_accuracy = plain_accuracy - accuracy\n",
    "print(f\"Difference between plain and encrypted accuracies: {diff_accuracy}\")\n",
    "if diff_accuracy < 0:\n",
    "    print(\"Oh! We got a better accuracy when training on encrypted data! The noise was on our side...\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Even after running this cell many times myself, I always feel the joy when I see it working on encrypted data, so I hope you are feeling this joy as well!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Congratulations!!! - Time to Join the Community!\n",
    "\n",
    "Congratulations on completing this notebook tutorial! If you enjoyed this and would like to join the movement towards privacy preserving, decentralized ownership of AI and the AI supply chain (data), you can do so in the following ways!\n",
    "\n",
    "### Star TenSEAL on GitHub\n",
    "\n",
    "The easiest way to help our community is just by starring the Repos! This helps raise awareness of the cool tools we're building.\n",
    "\n",
    "- [Star TenSEAL](https://github.com/OpenMined/TenSEAL)\n",
    "\n",
    "### Join our Slack!\n",
    "\n",
    "The best way to keep up to date on the latest advancements is to join our community! You can do so by filling out the form at [http://slack.openmined.org](http://slack.openmined.org). #lib_tenseal and #code_tenseal are the main channels for the TenSEAL project.\n",
    "\n",
    "### Join our Team!\n",
    "\n",
    "If you're excited about what we are working on TenSEAL, and if you're interested to work on homomorphic encryption related use cases, you should definitely join us!\n",
    "\n",
    "[Apply to the crypto team!](https://docs.google.com/forms/d/1T6MJ21V1lb7aEr4ilZOTYQXzxXP6KbpLumZVmTZMSuY/edit)\n",
    "\n",
    "\n",
    "### Donate\n",
    "\n",
    "If you don't have time to contribute to our codebase, but would still like to lend support, you can also become a Backer on our Open Collective. All donations go toward our web hosting and other community expenses such as hackathons and meetups!\n",
    "\n",
    "[OpenMined's Open Collective Page](https://opencollective.com/openmined)\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.9.6 (default, Sep 26 2022, 11:37:49) \n[Clang 14.0.0 (clang-1400.0.29.202)]"
  },
  "vscode": {
   "interpreter": {
    "hash": "31f2aee4e71d21fbe5cf8b01ff0e069b9275f58929596ceb00d14d90e3e16cd6"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}