examples/04_model_select_and_optimize/tuning_spark_als.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "source": [
    "<i>Copyright (c) Microsoft Corporation. All rights reserved.</i>\n",
    "\n",
    "<i>Licensed under the MIT License.</i>"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "# Hyperparameter tuning (Spark based recommender)"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Hyperparameter tuning for Spark based recommender algorithm is important to select a model with the optimal performance. This notebook introduces good practices in performing hyperparameter tuning for building recommender models with the utility functions provided in the [Microsoft/Recommenders](https://github.com/Microsoft/Recommenders.git) repository.\n",
    "\n",
    "Three different approaches are introduced and comparatively studied.\n",
    "* Spark native/custom constructs (`ParamGridBuilder`, `TrainValidationSplit`).\n",
    "* `hyperopt` package with Tree of Parzen Estimator algorithm. \n",
    "* Brute-force random search of parameter values sampled with pre-defined space. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 0 Global settings and import"
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "source": [
    "# set the environment path to find Recommenders\n",
    "%matplotlib notebook\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import sys\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "import pyspark\n",
    "import pyspark.sql.functions as F\n",
    "from pyspark.ml.recommendation import ALS\n",
    "from pyspark.ml.tuning import ParamGridBuilder, TrainValidationSplit\n",
    "from pyspark.ml.evaluation import Evaluator, RegressionEvaluator\n",
    "from pyspark.ml.pipeline import Estimator, Model\n",
    "from pyspark import keyword_only  \n",
    "from pyspark.ml.param.shared import *\n",
    "from pyspark.ml.util import *\n",
    "from pyspark.mllib.evaluation import RankingMetrics\n",
    "from pyspark.sql.types import ArrayType, IntegerType\n",
    "\n",
    "from hyperopt import fmin, tpe, hp, STATUS_OK, Trials\n",
    "from hyperopt.pyll.stochastic import sample\n",
    "\n",
    "from recommenders.utils.timer import Timer\n",
    "from recommenders.utils.spark_utils import start_or_get_spark\n",
    "from recommenders.evaluation.spark_evaluation import SparkRankingEvaluation, SparkRatingEvaluation\n",
    "from recommenders.datasets.movielens import load_spark_df\n",
    "from recommenders.datasets.spark_splitters import spark_random_split\n",
    "\n",
    "print(\"System version: {}\".format(sys.version))\n",
    "print(\"Pandas version: {}\".format(pd.__version__))\n",
    "print(\"PySpark version: {}\".format(pyspark.__version__))"
   ],
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "System version: 3.5.5 |Anaconda custom (64-bit)| (default, May 13 2018, 21:12:35) \n",
      "[GCC 7.2.0]\n",
      "Pandas version: 0.23.0\n",
      "PySpark version: 2.3.1\n"
     ]
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "source": [
    "MOVIELENS_DATA_SIZE = \"100k\"\n",
    "\n",
    "NUMBER_CORES = 1\n",
    "NUMBER_ITERATIONS = 25\n",
    "\n",
    "COL_USER = \"userID\"\n",
    "COL_ITEM = \"itemID\"\n",
    "COL_TIMESTAMP = \"timestamp\"\n",
    "COL_RATING = \"rating\"\n",
    "COL_PREDICTION = \"prediction\"\n",
    "\n",
    "HEADER = {\n",
    "    \"col_user\": COL_USER,\n",
    "    \"col_item\": COL_ITEM,\n",
    "    \"col_rating\": COL_RATING,\n",
    "    \"col_prediction\": COL_PREDICTION,\n",
    "}\n",
    "\n",
    "HEADER_ALS = {\n",
    "    \"userCol\": COL_USER,\n",
    "    \"itemCol\": COL_ITEM,\n",
    "    \"ratingCol\": COL_RATING\n",
    "}\n",
    "\n",
    "SUBSET_RATIO = 0.5\n",
    "\n",
    "RANK = [10, 15, 20, 30, 40]\n",
    "REG = [ 0.1, 0.01, 0.001, 0.0001, 0.00001]"
   ],
   "outputs": [],
   "metadata": {
    "tags": [
     "parameters"
    ]
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 1 Data preparation"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "A Spark session is created. Note in this case, to study the running time for different approaches, the Spark session in local mode uses only one core for running. This eliminates the impact of parallelization of parameter tuning. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "source": [
    "spark = start_or_get_spark(url=\"local[{}]\".format(NUMBER_CORES))"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "MovieLens 100k dataset is used for running the demonstration."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "source": [
    "data = load_spark_df(spark, size=MOVIELENS_DATA_SIZE, header=(COL_USER, COL_ITEM, COL_RATING))"
   ],
   "outputs": [
    {
     "output_type": "stream",
     "name": "stderr",
     "text": [
      "100%|██████████| 4.81k/4.81k [00:01<00:00, 2.47kKB/s]\n"
     ]
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "To reduce time spent on the comparitive study, 50% of the data is used for the experimentation below."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "source": [
    "data, _ = spark_random_split(data, ratio=SUBSET_RATIO)"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "The dataset is split into 3 subsets randomly with a given split ratio. The hyperparameter tuning is performed on the training and the validating data, and then the optimal recommender selected is evaluated on the testing dataset."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "source": [
    "train, valid, test = spark_random_split(data, ratio=[3, 1, 1])"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 2 Hyper parameter tuning with Azure Machine Learning Services"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "The `hyperdrive` module in the [Azure Machine Learning Services](https://azure.microsoft.com/en-us/services/machine-learning-service/) runs [hyperparameter tuning and optimizing for machine learning model selection](https://docs.microsoft.com/en-us/azure/machine-learning/service/how-to-tune-hyperparameters). At the moment, the service supports running hyperparameter tuning on heterogenous computing targets such as cluster of commodity compute nodes with or without GPU devices (see detailed documentation [here](https://docs.microsoft.com/en-us/azure/machine-learning/service/how-to-set-up-training-targets)). It is feasible to run parameter tuning on a cluster of VM nodes. In this case, the service containerizes individual and independent Spark session on each node of the cluster to run the parameter tuning job in parallel, instead of inside a single Spark session where the training is executed in a distributed manner.  \n",
    "\n",
    "Detailed instructions of tuning hyperparameter of non-Spark workloads by using Azure Machine Learning Services can be found in [this](./hypertune_aml_wide_and_deep_quickstart.ipynb) notebook. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 3 Hyper parameter tuning with Spark ML constructs"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "### 3.1 Spark native construct"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Spark ML lib implements modules such as `CrossValidator` and `TrainValidationSplit` for tuning hyperparameters (see [here](https://spark.apache.org/docs/2.2.0/ml-tuning.html)). However, by default, it does not support custom machine learning algorithms, data splitting methods, and evaluation metrics, like what are offered as utility functions in the Recommenders repository. \n",
    "\n",
    "For example, the Spark native constuct can be used for tuning a recommender against the `rmse` metric which is one of the available regression metrics in Spark."
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Firstly, a Spark ALS object needs to be created. In this case, for illustration purpose, it is an ALS model object."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "source": [
    "# NOTE the parameters of interest, rank and regParam, are left unset, \n",
    "# because their values will be assigned in the parameter grid builder.\n",
    "als = ALS(\n",
    "    maxIter=15,\n",
    "    implicitPrefs=False,\n",
    "    alpha=0.1,\n",
    "    coldStartStrategy='drop',\n",
    "    nonnegative=False,\n",
    "    **HEADER_ALS\n",
    ")"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Then, a parameter grid can be defined as follows. Without loss of generity, only `rank` and `regParam` are considered."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "source": [
    "paramGrid = ParamGridBuilder() \\\n",
    "    .addGrid(als.rank, RANK) \\\n",
    "    .addGrid(als.regParam, REG) \\\n",
    "    .build()"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Given the settings above, a `TrainValidationSplit` constructor can be created for fitting the best model in the given parameter range. In this case, the `RegressionEvaluator` is using `RMSE`, by default, as an evaluation metric. \n",
    "\n",
    "Since the data splitter is embedded in the `TrainValidationSplit` object, to make sure the splitting ratio is consistent across different approaches, the split ratio is set to be 0.75 and in the model training the training dataset and validating dataset are combined. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "source": [
    "tvs = TrainValidationSplit(\n",
    "    estimator=als,\n",
    "    estimatorParamMaps=paramGrid,\n",
    "    # A regression evaluation method is used. \n",
    "    evaluator=RegressionEvaluator(labelCol='rating'),\n",
    "    # 75% of the data will be used for training, 25% for validation.\n",
    "    # NOTE here the splitting is random. The Spark splitting utilities (e.g. chrono splitter)\n",
    "    # are therefore not available here. \n",
    "    trainRatio=0.75\n",
    ")"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "source": [
    "with Timer() as time_spark:\n",
    "    # Run TrainValidationSplit, and choose the best set of parameters.\n",
    "    # NOTE train and valid is union because in Spark TrainValidationSplit does splitting by itself.\n",
    "    model = tvs.fit(train.union(valid))\n",
    "\n"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "The model parameters in the grid and the best metrics can be then returned. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "source": [
    "for idx, item in enumerate(model.getEstimatorParamMaps()):\n",
    "    print('Run {}:'.format(idx))\n",
    "    print('\\tValidation Metric: {}'.format(model.validationMetrics[idx]))\n",
    "    for key, value in item.items():\n",
    "        print('\\t{0}: {1}'.format(repr(key), value))"
   ],
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "Run 0:\n",
      "\tValidation Metric: 1.0505385750367227\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 10\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.1\n",
      "Run 1:\n",
      "\tValidation Metric: 1.0444319735752456\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 15\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.1\n",
      "Run 2:\n",
      "\tValidation Metric: 1.040060458376737\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 20\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.1\n",
      "Run 3:\n",
      "\tValidation Metric: 1.0293843505140208\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 30\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.1\n",
      "Run 4:\n",
      "\tValidation Metric: 1.0216137585741758\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 40\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.1\n",
      "Run 5:\n",
      "\tValidation Metric: 1.4359689750708238\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 10\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.01\n",
      "Run 6:\n",
      "\tValidation Metric: 1.4607579006632527\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 15\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.01\n",
      "Run 7:\n",
      "\tValidation Metric: 1.462040851503185\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 20\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.01\n",
      "Run 8:\n",
      "\tValidation Metric: 1.3991557293601262\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 30\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.01\n",
      "Run 9:\n",
      "\tValidation Metric: 1.3410599805798034\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 40\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.01\n",
      "Run 10:\n",
      "\tValidation Metric: 1.988790687632913\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 10\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.001\n",
      "Run 11:\n",
      "\tValidation Metric: 1.87653183445857\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 15\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.001\n",
      "Run 12:\n",
      "\tValidation Metric: 1.9156975302076813\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 20\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.001\n",
      "Run 13:\n",
      "\tValidation Metric: 1.8551322988886354\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 30\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.001\n",
      "Run 14:\n",
      "\tValidation Metric: 1.8640288497082245\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 40\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.001\n",
      "Run 15:\n",
      "\tValidation Metric: 3.0577123716478947\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 10\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.0001\n",
      "Run 16:\n",
      "\tValidation Metric: 2.849817587112537\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 15\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.0001\n",
      "Run 17:\n",
      "\tValidation Metric: 2.5637113568725733\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 20\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.0001\n",
      "Run 18:\n",
      "\tValidation Metric: 2.549950528259629\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 30\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.0001\n",
      "Run 19:\n",
      "\tValidation Metric: 2.9097893609292096\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 40\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 0.0001\n",
      "Run 20:\n",
      "\tValidation Metric: 4.9701569583584115\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 10\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 1e-05\n",
      "Run 21:\n",
      "\tValidation Metric: 4.230109404473068\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 15\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 1e-05\n",
      "Run 22:\n",
      "\tValidation Metric: 3.5723074921735707\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 20\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 1e-05\n",
      "Run 23:\n",
      "\tValidation Metric: 3.469775656947356\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 30\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 1e-05\n",
      "Run 24:\n",
      "\tValidation Metric: 4.426604995574413\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='rank', doc='rank of the factorization'): 40\n",
      "\tParam(parent='ALS_496eab0f4b1e8da03092', name='regParam', doc='regularization parameter (>= 0).'): 1e-05\n"
     ]
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "source": [
    "model.validationMetrics"
   ],
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "[1.0505385750367227,\n",
       " 1.0444319735752456,\n",
       " 1.040060458376737,\n",
       " 1.0293843505140208,\n",
       " 1.0216137585741758,\n",
       " 1.4359689750708238,\n",
       " 1.4607579006632527,\n",
       " 1.462040851503185,\n",
       " 1.3991557293601262,\n",
       " 1.3410599805798034,\n",
       " 1.988790687632913,\n",
       " 1.87653183445857,\n",
       " 1.9156975302076813,\n",
       " 1.8551322988886354,\n",
       " 1.8640288497082245,\n",
       " 3.0577123716478947,\n",
       " 2.849817587112537,\n",
       " 2.5637113568725733,\n",
       " 2.549950528259629,\n",
       " 2.9097893609292096,\n",
       " 4.9701569583584115,\n",
       " 4.230109404473068,\n",
       " 3.5723074921735707,\n",
       " 3.469775656947356,\n",
       " 4.426604995574413]"
      ]
     },
     "metadata": {},
     "execution_count": 12
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "To get the best model, just do"
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "source": [
    "model_best_spark = model.bestModel"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "### 3.2 Custom `Estimator`, `Transformer`, and `Evaluator` for Spark ALS\n",
    "\n",
    "One can also customize Spark modules to allow tuning hyperparameters for a desired model and evaluation metric, given that the native Spark ALS does not allow tuning hyperparameters for ranking metrics such as precision@k, recall@k, etc. This can be done by creating custom `Estimator`, `Transformer` and `Evaluator`. The benefit is that, after the customization, the tuning process can make use of `trainValidSplit` directly, which distributes the tuning in a Spark session."
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "#### Customized `Estimator` and `Transformer` for top k recommender based on Spark ALS\n",
    "\n",
    "The following shows how to implement a PySpark `Estimator` and `Transfomer` for recommending top k items from ALS model. The latter generates top k recommendations from the model object. Both of the two are designed by following the protocol of Spark APIs, to make sure that they can be run with the hyperparameter tuning constructs in Spark."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "source": [
    "class ALSTopK(\n",
    "    ALS,\n",
    "    Estimator,\n",
    "    HasInputCol,\n",
    "    HasPredictionCol\n",
    "):    \n",
    "    rank = Param(Params._dummy(), \"rank\", \"rank of the factorization\",\n",
    "                 typeConverter=TypeConverters.toInt)\n",
    "    numUserBlocks = Param(Params._dummy(), \"numUserBlocks\", \"number of user blocks\",\n",
    "                          typeConverter=TypeConverters.toInt)\n",
    "    numItemBlocks = Param(Params._dummy(), \"numItemBlocks\", \"number of item blocks\",\n",
    "                          typeConverter=TypeConverters.toInt)\n",
    "    implicitPrefs = Param(Params._dummy(), \"implicitPrefs\", \"whether to use implicit preference\",\n",
    "                          typeConverter=TypeConverters.toBoolean)\n",
    "    alpha = Param(Params._dummy(), \"alpha\", \"alpha for implicit preference\",\n",
    "                  typeConverter=TypeConverters.toFloat)\n",
    "    userCol = Param(Params._dummy(), \"userCol\", \"column name for user ids. Ids must be within \" +\n",
    "                    \"the integer value range.\", typeConverter=TypeConverters.toString)\n",
    "    itemCol = Param(Params._dummy(), \"itemCol\", \"column name for item ids. Ids must be within \" +\n",
    "                    \"the integer value range.\", typeConverter=TypeConverters.toString)\n",
    "    ratingCol = Param(Params._dummy(), \"ratingCol\", \"column name for ratings\",\n",
    "                      typeConverter=TypeConverters.toString)\n",
    "    nonnegative = Param(Params._dummy(), \"nonnegative\",\n",
    "                        \"whether to use nonnegative constraint for least squares\",\n",
    "                        typeConverter=TypeConverters.toBoolean)\n",
    "    intermediateStorageLevel = Param(Params._dummy(), \"intermediateStorageLevel\",\n",
    "                                     \"StorageLevel for intermediate datasets. Cannot be 'NONE'.\",\n",
    "                                     typeConverter=TypeConverters.toString)\n",
    "    finalStorageLevel = Param(Params._dummy(), \"finalStorageLevel\",\n",
    "                              \"StorageLevel for ALS model factors.\",\n",
    "                              typeConverter=TypeConverters.toString)\n",
    "    coldStartStrategy = Param(Params._dummy(), \"coldStartStrategy\", \"strategy for dealing with \" +\n",
    "                              \"unknown or new users/items at prediction time. This may be useful \" +\n",
    "                              \"in cross-validation or production scenarios, for handling \" +\n",
    "                              \"user/item ids the model has not seen in the training data. \" +\n",
    "                              \"Supported values: 'nan', 'drop'.\",\n",
    "                              typeConverter=TypeConverters.toString)\n",
    "\n",
    "    @keyword_only\n",
    "    def __init__(\n",
    "        self,\n",
    "        rank=10, maxIter=10, regParam=0.1, numUserBlocks=10, numItemBlocks=10,\n",
    "        implicitPrefs=False, alpha=1.0, userCol=\"user\", itemCol=\"item\", seed=None, k=10,\n",
    "        ratingCol=\"rating\", nonnegative=False, checkpointInterval=10,\n",
    "        intermediateStorageLevel=\"MEMORY_AND_DISK\",\n",
    "        finalStorageLevel=\"MEMORY_AND_DISK\", coldStartStrategy=\"nan\"\n",
    "    ):\n",
    "        super(ALS, self).__init__()\n",
    "        self._java_obj = self._new_java_obj(\"org.apache.spark.ml.recommendation.ALS\", self.uid)\n",
    "        self._setDefault(rank=10, maxIter=10, regParam=0.1, numUserBlocks=10, numItemBlocks=10,\n",
    "                         implicitPrefs=False, alpha=1.0, userCol=\"user\", itemCol=\"item\",\n",
    "                         ratingCol=\"rating\", nonnegative=False, checkpointInterval=10,\n",
    "                         intermediateStorageLevel=\"MEMORY_AND_DISK\",\n",
    "                         finalStorageLevel=\"MEMORY_AND_DISK\", coldStartStrategy=\"nan\")\n",
    "\n",
    "        kwargs = self._input_kwargs \n",
    "        kwargs = {x: kwargs[x] for x in kwargs if x not in {'k'}}\n",
    "        self.setParams(**kwargs)\n",
    "        \n",
    "        # The manually added parameter is not present in ALS Java implementation. \n",
    "        self.k = k\n",
    "        \n",
    "    def setRank(self, value):\n",
    "        \"\"\"\n",
    "        Sets the value of :py:attr:`rank`.\n",
    "        \"\"\"\n",
    "        return self._set(rank=value)\n",
    "\n",
    "    def getRank(self):\n",
    "        \"\"\"\n",
    "        Gets the value of rank or its default value.\n",
    "        \"\"\"\n",
    "        return self.getOrDefault(self.rank)\n",
    "    \n",
    "    def setParams(self, rank=10, maxIter=10, regParam=0.1, numUserBlocks=10, numItemBlocks=10,\n",
    "                  implicitPrefs=False, alpha=1.0, userCol=\"user\", itemCol=\"item\", seed=None,\n",
    "                  ratingCol=\"rating\", nonnegative=False, checkpointInterval=10,\n",
    "                  intermediateStorageLevel=\"MEMORY_AND_DISK\",\n",
    "                  finalStorageLevel=\"MEMORY_AND_DISK\", coldStartStrategy=\"nan\"):\n",
    "        \"\"\"\n",
    "        setParams(self, rank=10, maxIter=10, regParam=0.1, numUserBlocks=10, numItemBlocks=10, \\\n",
    "                 implicitPrefs=False, alpha=1.0, userCol=\"user\", itemCol=\"item\", seed=None, \\\n",
    "                 ratingCol=\"rating\", nonnegative=False, checkpointInterval=10, \\\n",
    "                 intermediateStorageLevel=\"MEMORY_AND_DISK\", \\\n",
    "                 finalStorageLevel=\"MEMORY_AND_DISK\", coldStartStrategy=\"nan\")\n",
    "        Sets params for ALS.\n",
    "        \"\"\"\n",
    "        kwargs = self._input_kwargs\n",
    "        kwargs = {x: kwargs[x] for x in kwargs if x not in {'k'}}\n",
    "        return self._set(**kwargs)\n",
    "        \n",
    "    def _fit(self, dataset):\n",
    "        kwargs = self._input_kwargs    \n",
    "        # Exclude k as it is not a parameter for ALS.\n",
    "        kwargs = {x: kwargs[x] for x in kwargs if x not in {'k'}}\n",
    "        kwargs['rank'] = self.getRank()\n",
    "        kwargs['regParam'] = self.getOrDefault(self.regParam)\n",
    "        als = ALS(\n",
    "            **kwargs\n",
    "        )\n",
    "        als_model = als.fit(dataset)\n",
    "        \n",
    "        user_col = kwargs['userCol']\n",
    "        item_col = kwargs['itemCol']\n",
    "        \n",
    "        k = self.k\n",
    "                     \n",
    "        topk_model = ALSTopKModel()\n",
    "        topk_model.setParams(\n",
    "            als_model,\n",
    "            user_col, \n",
    "            item_col, \n",
    "            k\n",
    "        )\n",
    "        \n",
    "        return topk_model\n",
    "    \n",
    "    \n",
    "class ALSTopKModel(\n",
    "    Model,\n",
    "    HasInputCol,\n",
    "    HasPredictionCol,\n",
    "    HasLabelCol\n",
    "):    \n",
    "    def setParams(self, model, userCol, itemCol, k):\n",
    "        self.model = model\n",
    "        self.userCol = userCol\n",
    "        self.itemCol = itemCol\n",
    "        self.k = k\n",
    "    \n",
    "    def _transform(self, dataset):\n",
    "        predictionCol = self.getPredictionCol()\n",
    "        labelCol = self.getLabelCol()\n",
    "        \n",
    "        users = dataset.select(self.userCol).distinct()\n",
    "        topk_recommendation = self.model.recommendForUserSubset(users, self.k)  \n",
    "        \n",
    "        extract_value = F.udf((lambda x: [y[0] for y in x]), ArrayType(IntegerType()))\n",
    "        topk_recommendation = topk_recommendation.withColumn(predictionCol, extract_value(F.col(\"recommendations\")))        \n",
    "        \n",
    "        dataset = (\n",
    "            dataset\n",
    "            .groupBy(self.userCol)\n",
    "            .agg(F.collect_list(F.col(self.itemCol)).alias(labelCol))\n",
    "        )\n",
    "            \n",
    "        topk_recommendation_all = dataset.join(\n",
    "            topk_recommendation, \n",
    "            on=self.userCol,\n",
    "            how=\"outer\"\n",
    "        )\n",
    "        \n",
    "        return topk_recommendation_all.select(self.userCol, labelCol, predictionCol)"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "#### Customized precision@k evaluation metric\n",
    "\n",
    "In addition to the custom `Estimator` and `Transformer`, it may also be desired to customize an `Evaluator` to allow \"beyond-rating\" metrics. The codes as following illustrates a precision@k evaluator. Other types of evaluators can be developed in a similar way."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "source": [
    "# Define a custom Evaulator. Here precision@k is used.\n",
    "class PrecisionAtKEvaluator(Evaluator):\n",
    "\n",
    "    def __init__(self, predictionCol=\"prediction\", labelCol=\"label\", k=10):\n",
    "        self.predictionCol = predictionCol\n",
    "        self.labelCol = labelCol\n",
    "        self.k = k\n",
    "\n",
    "    def _evaluate(self, dataset):\n",
    "        \"\"\"\n",
    "        Returns a random number. \n",
    "        Implement here the true metric\n",
    "        \"\"\"      \n",
    "        # Drop Nulls.\n",
    "        dataset = dataset.na.drop()\n",
    "        metrics = RankingMetrics(dataset.select(self.predictionCol, self.labelCol).rdd)\n",
    "        return metrics.precisionAt(self.k)\n",
    "\n",
    "    def isLargerBetter(self):\n",
    "        return True"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Then a new ALS top-k recommender can be created, and the Spark native construct, `TrainValidationSplit` module, can be used to find the optimal model w.r.t the precision@k metric."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "source": [
    "alstopk = ALSTopK(\n",
    "    userCol=COL_USER,\n",
    "    itemCol=COL_ITEM,\n",
    "    ratingCol=COL_RATING,\n",
    "    k=10\n",
    ")\n",
    "\n",
    "# Here for illustration purpose, a small grid is used.\n",
    "paramGrid = ParamGridBuilder() \\\n",
    "    .addGrid(alstopk.rank, [10, 20]) \\\n",
    "    .addGrid(alstopk.regParam, [0.1, 0.01]) \\\n",
    "    .build()\n",
    "\n",
    "tvs = TrainValidationSplit(\n",
    "    estimator=alstopk,\n",
    "    estimatorParamMaps=paramGrid,\n",
    "    # A regression evaluation method is used. \n",
    "    evaluator=PrecisionAtKEvaluator(),\n",
    "    # 75% of the data will be used for training, 25% for validation.\n",
    "    # NOTE here the splitting is random. The Spark splitting utilities (e.g. chrono splitter)\n",
    "    # are therefore not available here. \n",
    "    trainRatio=0.75\n",
    ")"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "source": [
    "# Run TrainValidationSplit, and choose the best set of parameters.\n",
    "# NOTE train and valid is union because in Spark TrainValidationSplit does splitting by itself.\n",
    "model_precision = tvs.fit(train.union(valid))\n",
    "\n",
    "model_precision.getEstimatorParamMaps()"
   ],
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "[{Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 10,\n",
       "  Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.1},\n",
       " {Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 10,\n",
       "  Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.01},\n",
       " {Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 20,\n",
       "  Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.1},\n",
       " {Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 20,\n",
       "  Param(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.01}]"
      ]
     },
     "metadata": {},
     "execution_count": 17
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "source": [
    "def best_param(model, is_larger_better=True):\n",
    "    if is_larger_better:\n",
    "        best_metric = max(model.validationMetrics)\n",
    "    else:\n",
    "        best_metric = min(model.validationMetrics)\n",
    "        \n",
    "    parameters = model.getEstimatorParamMaps()[model.validationMetrics.index(best_metric)]\n",
    "     \n",
    "    return list(parameters.values())"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "source": [
    "params = best_param(model_precision)"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "source": [
    "model_precision.bestModel.transform(valid).limit(5).show()\n",
    "\n",
    "for idx, item in enumerate(model_precision.getEstimatorParamMaps()):\n",
    "    print('Run {}:'.format(idx))\n",
    "    print('\\tValidation Metric: {}'.format(model_precision.validationMetrics[idx]))\n",
    "    for key, value in item.items():\n",
    "        print('\\t{0}: {1}'.format(repr(key), value))"
   ],
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "+------+--------------------+--------------------+\n",
      "|userID|               label|          prediction|\n",
      "+------+--------------------+--------------------+\n",
      "|   148|     [116, 135, 189]|[126, 735, 222, 7...|\n",
      "|   463|[15, 100, 103, 12...|[87, 922, 343, 52...|\n",
      "|   471|      [95, 102, 946]|[663, 114, 488, 9...|\n",
      "|   496|[94, 143, 196, 28...|[475, 645, 69, 18...|\n",
      "|   833|[50, 68, 79, 92, ...|[428, 20, 135, 18...|\n",
      "+------+--------------------+--------------------+\n",
      "\n",
      "Run 0:\n",
      "\tValidation Metric: 0.01411637931034483\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 10\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.1\n",
      "Run 1:\n",
      "\tValidation Metric: 0.007866379310344828\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 10\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.01\n",
      "Run 2:\n",
      "\tValidation Metric: 0.01799568965517241\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 20\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.1\n",
      "Run 3:\n",
      "\tValidation Metric: 0.018103448275862084\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='rank', doc='rank of the factorization'): 20\n",
      "\tParam(parent='ALSTopK_4f48b7cc6cf2badfcea7', name='regParam', doc='regularization parameter (>= 0).'): 0.01\n"
     ]
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 4 Hyperparameter tuning with `hyperopt`"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "`hyperopt` is an open source Python package that is designed for tuning parameters for generic function with any pre-defined loss. More information about `hyperopt` can be found [here](https://github.com/hyperopt/hyperopt). `hyperopt` supports parallelization on MongoDB but not Spark. In our case, the tuning is performed in a sequential mode on a local computer.\n",
    "\n",
    "In `hyperopt`, an *objective* function is defined for optimizing the hyper parameters. In this case, the objective is similar to that in the Spark native construct situation, which is *to the RMSE metric for an ALS recommender*. Parameters of `rank` and `regParam` are used as hyperparameters. \n",
    "\n",
    "The objective function shown below demonstrates a RMSE loss for an ALS recommender. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "source": [
    "# Customize an objective function\n",
    "def objective(params):\n",
    "    with Timer() as time_run_start:\n",
    "    \n",
    "        rank = params['rank']\n",
    "        reg = params['reg']\n",
    "        train = params['train'] \n",
    "        valid = params['valid'] \n",
    "        col_user = params['col_user'] \n",
    "        col_item = params['col_item']\n",
    "        col_rating = params['col_rating'] \n",
    "        col_prediction = params['col_prediction'] \n",
    "        k = params['k']\n",
    "        relevancy_method = params['relevancy_method']\n",
    "\n",
    "        als = ALS(\n",
    "            rank=rank,\n",
    "            maxIter=15,\n",
    "            implicitPrefs=False,\n",
    "            alpha=0.1,\n",
    "            regParam=reg,\n",
    "            coldStartStrategy='drop',\n",
    "            nonnegative=False,\n",
    "            seed=42,\n",
    "            **HEADER_ALS\n",
    "        )\n",
    "\n",
    "        model = als.fit(train) \n",
    "        prediction = model.transform(valid)\n",
    "\n",
    "        rating_eval = SparkRatingEvaluation(\n",
    "            valid, \n",
    "            prediction, \n",
    "            **HEADER\n",
    "        )\n",
    "\n",
    "        rmse = rating_eval.rmse()\n",
    "    \n",
    "    # Return the objective function result.\n",
    "    return {\n",
    "        'loss': rmse,\n",
    "        'status': STATUS_OK,\n",
    "        'eval_time': time_run_start.interval\n",
    "    }"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "A search space is usually defined for hyperparameter exploration. Design of search space is empirical, and depends on the understanding of how distribution of parameter of interest affects the model performance measured by the loss function. \n",
    "\n",
    "In the ALS algorithm, the two hyper parameters, rank and reg, affect model performance in a way that\n",
    "* The higher the rank, the better the model performance but also the higher risk of overfitting.\n",
    "* The reg parameter prevents overfitting in certain way. \n",
    "\n",
    "Therefore, in this case, a uniform distribution and a lognormal distribution sampling spaces are used for rank and reg, respectively. A narrow search space is used for illustration purpose, that is, the range of rank is from 10 to 20, while that of reg is from $e^{-5}$ to $e^{-1}$. Together with the randomly sampled hyper parameters, other parameters use for building / evaluating the recommender, like `k`, column names, data, etc., are kept as constants."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "source": [
    "# define a search space\n",
    "space = {\n",
    "    'rank': hp.quniform('rank', 10, 40, 5),\n",
    "    'reg': hp.loguniform('reg', -5, -1),\n",
    "    'train': train, \n",
    "    'valid': valid, \n",
    "    'col_user': COL_USER, \n",
    "    'col_item': COL_ITEM, \n",
    "    'col_rating': COL_RATING, \n",
    "    'col_prediction': \"prediction\", \n",
    "    'k': 10,\n",
    "    'relevancy_method': \"top_k\"\n",
    "}"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "### 4.1 Hyperparameter tuning with TPE"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "`fmin` of `hyperopt` is used for running the trials for searching optimal hyper parameters. In `hyperopt`, there are different strategies for intelligently optimize hyper parameters. For example, `hyperopt` avails [Tree of Parzen Estimators (TPE) method](https://papers.nips.cc/paper/4443-algorithms-for-hyper-parameter-optimization.pdf) for searching optimal parameters. \n",
    "\n",
    "The TPE method models a surface response of $p(x|y)$ by transforming a generative process, replacing the distributions of the configuration prior with non-parametric densities, where $p$ is the probability of configuration space $x$ given the loss $y$. For different configuration space, the TPE method does different replacements. That is, uniform $\\to$ truncated Gaussian mixture, log-uniform $\\to$ exponentiated truncated Gaussian mixture, categorical $\\to$ re-weighted categorical, etc. Using different observations ${x(1), ..., x(k)}$ in the non-parametric densities, these substitutions represent a learning algorithm that can produce a variety of densities over the configuration space $X$. By maintaining sorted lists of observed variables in $H$, the runtime of each iteration of the TPE algorithm can scale linearly in $|H|$ and linearly in the number of variables (dimensions) being optimized. In a nutshell, the algorithm recognizes the irrelevant variables in the configuration space, and thus reduces iterations in searching for the optimal ones. Details of the TPE algorithm can be found in the reference paper.\n",
    "\n",
    "The following runs the trials with the pre-defined objective function and search space. TPE is used as the optimization method. Totally there will be 10 evaluations run for searching the best parameters."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "source": [
    "with Timer() as time_hyperopt:\n",
    "    # Trials for recording each iteration of the hyperparameter searching.\n",
    "    trials = Trials()\n",
    "\n",
    "    best = fmin(\n",
    "        fn=objective,\n",
    "        space=space,\n",
    "        algo=tpe.suggest,\n",
    "        trials=trials,\n",
    "        max_evals=NUMBER_ITERATIONS\n",
    "    )\n",
    "                  \n"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "source": [
    "trials.best_trial"
   ],
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "{'book_time': datetime.datetime(2019, 7, 17, 12, 28, 19, 108000),\n",
       " 'exp_key': None,\n",
       " 'misc': {'cmd': ('domain_attachment', 'FMinIter_Domain'),\n",
       "  'idxs': {'rank': [21], 'reg': [21]},\n",
       "  'tid': 21,\n",
       "  'vals': {'rank': [35.0], 'reg': [0.20807325457673764]},\n",
       "  'workdir': None},\n",
       " 'owner': None,\n",
       " 'refresh_time': datetime.datetime(2019, 7, 17, 12, 28, 28, 873000),\n",
       " 'result': {'eval_time': 9.763921976089478,\n",
       "  'loss': 0.9948670591255364,\n",
       "  'status': 'ok'},\n",
       " 'spec': None,\n",
       " 'state': 2,\n",
       " 'tid': 21,\n",
       " 'version': 0}"
      ]
     },
     "metadata": {},
     "execution_count": 24
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "source": [
    "parameters = ['rank', 'reg']\n",
    "cols = len(parameters)\n",
    "f, axes = plt.subplots(nrows=1, ncols=cols, figsize=(15,5))\n",
    "cmap = plt.cm.jet\n",
    "for i, val in enumerate(parameters):\n",
    "    xs = np.array([t['misc']['vals'][val] for t in trials.trials]).ravel()\n",
    "    ys = [t['result']['loss'] for t in trials.trials]\n",
    "    xs, ys = zip(*sorted(zip(xs, ys)))\n",
    "    ys = np.array(ys)\n",
    "    axes[i].scatter(xs, ys, s=20, linewidth=0.01, alpha=0.75, c=cmap(float(i)/len(parameters)))\n",
    "    axes[i].set_title(val)"
   ],
   "outputs": [
    {
     "output_type": "display_data",
     "data": {
      "application/javascript": "/* Put everything inside the global mpl namespace */\nwindow.mpl = {};\n\n\nmpl.get_websocket_type = function() {\n    if (typeof(WebSocket) !== 'undefined') {\n        return WebSocket;\n    } else if (typeof(MozWebSocket) !== 'undefined') {\n        return MozWebSocket;\n    } else {\n        alert('Your browser does not have WebSocket support.' +\n              'Please try Chrome, Safari or Firefox ≥ 6. ' +\n              'Firefox 4 and 5 are also supported but you ' +\n              'have to enable WebSockets in about:config.');\n    };\n}\n\nmpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n    this.id = figure_id;\n\n    this.ws = websocket;\n\n    this.supports_binary = (this.ws.binaryType != undefined);\n\n    if (!this.supports_binary) {\n        var warnings = document.getElementById(\"mpl-warnings\");\n        if (warnings) {\n            warnings.style.display = 'block';\n            warnings.textContent = (\n                \"This browser does not support binary websocket messages. \" +\n                    \"Performance may be slow.\");\n        }\n    }\n\n    this.imageObj = new Image();\n\n    this.context = undefined;\n    this.message = undefined;\n    this.canvas = undefined;\n    this.rubberband_canvas = undefined;\n    this.rubberband_context = undefined;\n    this.format_dropdown = undefined;\n\n    this.image_mode = 'full';\n\n    this.root = $('<div/>');\n    this._root_extra_style(this.root)\n    this.root.attr('style', 'display: inline-block');\n\n    $(parent_element).append(this.root);\n\n    this._init_header(this);\n    this._init_canvas(this);\n    this._init_toolbar(this);\n\n    var fig = this;\n\n    this.waiting = false;\n\n    this.ws.onopen =  function () {\n            fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n            fig.send_message(\"send_image_mode\", {});\n            if (mpl.ratio != 1) {\n                fig.send_message(\"set_dpi_ratio\", {'dpi_ratio': mpl.ratio});\n            }\n            fig.send_message(\"refresh\", {});\n        }\n\n    this.imageObj.onload = function() {\n            if (fig.image_mode == 'full') {\n                // Full images could contain transparency (where diff images\n                // almost always do), so we need to clear the canvas so that\n                // there is no ghosting.\n                fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n            }\n            fig.context.drawImage(fig.imageObj, 0, 0);\n        };\n\n    this.imageObj.onunload = function() {\n        fig.ws.close();\n    }\n\n    this.ws.onmessage = this._make_on_message_function(this);\n\n    this.ondownload = ondownload;\n}\n\nmpl.figure.prototype._init_header = function() {\n    var titlebar = $(\n        '<div class=\"ui-dialog-titlebar ui-widget-header ui-corner-all ' +\n        'ui-helper-clearfix\"/>');\n    var titletext = $(\n        '<div class=\"ui-dialog-title\" style=\"width: 100%; ' +\n        'text-align: center; padding: 3px;\"/>');\n    titlebar.append(titletext)\n    this.root.append(titlebar);\n    this.header = titletext[0];\n}\n\n\n\nmpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n\n}\n\n\nmpl.figure.prototype._root_extra_style = function(canvas_div) {\n\n}\n\nmpl.figure.prototype._init_canvas = function() {\n    var fig = this;\n\n    var canvas_div = $('<div/>');\n\n    canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n\n    function canvas_keyboard_event(event) {\n        return fig.key_event(event, event['data']);\n    }\n\n    canvas_div.keydown('key_press', canvas_keyboard_event);\n    canvas_div.keyup('key_release', canvas_keyboard_event);\n    this.canvas_div = canvas_div\n    this._canvas_extra_style(canvas_div)\n    this.root.append(canvas_div);\n\n    var canvas = $('<canvas/>');\n    canvas.addClass('mpl-canvas');\n    canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n\n    this.canvas = canvas[0];\n    this.context = canvas[0].getContext(\"2d\");\n\n    var backingStore = this.context.backingStorePixelRatio ||\n\tthis.context.webkitBackingStorePixelRatio ||\n\tthis.context.mozBackingStorePixelRatio ||\n\tthis.context.msBackingStorePixelRatio ||\n\tthis.context.oBackingStorePixelRatio ||\n\tthis.context.backingStorePixelRatio || 1;\n\n    mpl.ratio = (window.devicePixelRatio || 1) / backingStore;\n\n    var rubberband = $('<canvas/>');\n    rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n\n    var pass_mouse_events = true;\n\n    canvas_div.resizable({\n        start: function(event, ui) {\n            pass_mouse_events = false;\n        },\n        resize: function(event, ui) {\n            fig.request_resize(ui.size.width, ui.size.height);\n        },\n        stop: function(event, ui) {\n            pass_mouse_events = true;\n            fig.request_resize(ui.size.width, ui.size.height);\n        },\n    });\n\n    function mouse_event_fn(event) {\n        if (pass_mouse_events)\n            return fig.mouse_event(event, event['data']);\n    }\n\n    rubberband.mousedown('button_press', mouse_event_fn);\n    rubberband.mouseup('button_release', mouse_event_fn);\n    // Throttle sequential mouse events to 1 every 20ms.\n    rubberband.mousemove('motion_notify', mouse_event_fn);\n\n    rubberband.mouseenter('figure_enter', mouse_event_fn);\n    rubberband.mouseleave('figure_leave', mouse_event_fn);\n\n    canvas_div.on(\"wheel\", function (event) {\n        event = event.originalEvent;\n        event['data'] = 'scroll'\n        if (event.deltaY < 0) {\n            event.step = 1;\n        } else {\n            event.step = -1;\n        }\n        mouse_event_fn(event);\n    });\n\n    canvas_div.append(canvas);\n    canvas_div.append(rubberband);\n\n    this.rubberband = rubberband;\n    this.rubberband_canvas = rubberband[0];\n    this.rubberband_context = rubberband[0].getContext(\"2d\");\n    this.rubberband_context.strokeStyle = \"#000000\";\n\n    this._resize_canvas = function(width, height) {\n        // Keep the size of the canvas, canvas container, and rubber band\n        // canvas in synch.\n        canvas_div.css('width', width)\n        canvas_div.css('height', height)\n\n        canvas.attr('width', width * mpl.ratio);\n        canvas.attr('height', height * mpl.ratio);\n        canvas.attr('style', 'width: ' + width + 'px; height: ' + height + 'px;');\n\n        rubberband.attr('width', width);\n        rubberband.attr('height', height);\n    }\n\n    // Set the figure to an initial 600x600px, this will subsequently be updated\n    // upon first draw.\n    this._resize_canvas(600, 600);\n\n    // Disable right mouse context menu.\n    $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n        return false;\n    });\n\n    function set_focus () {\n        canvas.focus();\n        canvas_div.focus();\n    }\n\n    window.setTimeout(set_focus, 100);\n}\n\nmpl.figure.prototype._init_toolbar = function() {\n    var fig = this;\n\n    var nav_element = $('<div/>')\n    nav_element.attr('style', 'width: 100%');\n    this.root.append(nav_element);\n\n    // Define a callback function for later on.\n    function toolbar_event(event) {\n        return fig.toolbar_button_onclick(event['data']);\n    }\n    function toolbar_mouse_event(event) {\n        return fig.toolbar_button_onmouseover(event['data']);\n    }\n\n    for(var toolbar_ind in mpl.toolbar_items) {\n        var name = mpl.toolbar_items[toolbar_ind][0];\n        var tooltip = mpl.toolbar_items[toolbar_ind][1];\n        var image = mpl.toolbar_items[toolbar_ind][2];\n        var method_name = mpl.toolbar_items[toolbar_ind][3];\n\n        if (!name) {\n            // put a spacer in here.\n            continue;\n        }\n        var button = $('<button/>');\n        button.addClass('ui-button ui-widget ui-state-default ui-corner-all ' +\n                        'ui-button-icon-only');\n        button.attr('role', 'button');\n        button.attr('aria-disabled', 'false');\n        button.click(method_name, toolbar_event);\n        button.mouseover(tooltip, toolbar_mouse_event);\n\n        var icon_img = $('<span/>');\n        icon_img.addClass('ui-button-icon-primary ui-icon');\n        icon_img.addClass(image);\n        icon_img.addClass('ui-corner-all');\n\n        var tooltip_span = $('<span/>');\n        tooltip_span.addClass('ui-button-text');\n        tooltip_span.html(tooltip);\n\n        button.append(icon_img);\n        button.append(tooltip_span);\n\n        nav_element.append(button);\n    }\n\n    var fmt_picker_span = $('<span/>');\n\n    var fmt_picker = $('<select/>');\n    fmt_picker.addClass('mpl-toolbar-option ui-widget ui-widget-content');\n    fmt_picker_span.append(fmt_picker);\n    nav_element.append(fmt_picker_span);\n    this.format_dropdown = fmt_picker[0];\n\n    for (var ind in mpl.extensions) {\n        var fmt = mpl.extensions[ind];\n        var option = $(\n            '<option/>', {selected: fmt === mpl.default_extension}).html(fmt);\n        fmt_picker.append(option)\n    }\n\n    // Add hover states to the ui-buttons\n    $( \".ui-button\" ).hover(\n        function() { $(this).addClass(\"ui-state-hover\");},\n        function() { $(this).removeClass(\"ui-state-hover\");}\n    );\n\n    var status_bar = $('<span class=\"mpl-message\"/>');\n    nav_element.append(status_bar);\n    this.message = status_bar[0];\n}\n\nmpl.figure.prototype.request_resize = function(x_pixels, y_pixels) {\n    // Request matplotlib to resize the figure. Matplotlib will then trigger a resize in the client,\n    // which will in turn request a refresh of the image.\n    this.send_message('resize', {'width': x_pixels, 'height': y_pixels});\n}\n\nmpl.figure.prototype.send_message = function(type, properties) {\n    properties['type'] = type;\n    properties['figure_id'] = this.id;\n    this.ws.send(JSON.stringify(properties));\n}\n\nmpl.figure.prototype.send_draw_message = function() {\n    if (!this.waiting) {\n        this.waiting = true;\n        this.ws.send(JSON.stringify({type: \"draw\", figure_id: this.id}));\n    }\n}\n\n\nmpl.figure.prototype.handle_save = function(fig, msg) {\n    var format_dropdown = fig.format_dropdown;\n    var format = format_dropdown.options[format_dropdown.selectedIndex].value;\n    fig.ondownload(fig, format);\n}\n\n\nmpl.figure.prototype.handle_resize = function(fig, msg) {\n    var size = msg['size'];\n    if (size[0] != fig.canvas.width || size[1] != fig.canvas.height) {\n        fig._resize_canvas(size[0], size[1]);\n        fig.send_message(\"refresh\", {});\n    };\n}\n\nmpl.figure.prototype.handle_rubberband = function(fig, msg) {\n    var x0 = msg['x0'] / mpl.ratio;\n    var y0 = (fig.canvas.height - msg['y0']) / mpl.ratio;\n    var x1 = msg['x1'] / mpl.ratio;\n    var y1 = (fig.canvas.height - msg['y1']) / mpl.ratio;\n    x0 = Math.floor(x0) + 0.5;\n    y0 = Math.floor(y0) + 0.5;\n    x1 = Math.floor(x1) + 0.5;\n    y1 = Math.floor(y1) + 0.5;\n    var min_x = Math.min(x0, x1);\n    var min_y = Math.min(y0, y1);\n    var width = Math.abs(x1 - x0);\n    var height = Math.abs(y1 - y0);\n\n    fig.rubberband_context.clearRect(\n        0, 0, fig.canvas.width, fig.canvas.height);\n\n    fig.rubberband_context.strokeRect(min_x, min_y, width, height);\n}\n\nmpl.figure.prototype.handle_figure_label = function(fig, msg) {\n    // Updates the figure title.\n    fig.header.textContent = msg['label'];\n}\n\nmpl.figure.prototype.handle_cursor = function(fig, msg) {\n    var cursor = msg['cursor'];\n    switch(cursor)\n    {\n    case 0:\n        cursor = 'pointer';\n        break;\n    case 1:\n        cursor = 'default';\n        break;\n    case 2:\n        cursor = 'crosshair';\n        break;\n    case 3:\n        cursor = 'move';\n        break;\n    }\n    fig.rubberband_canvas.style.cursor = cursor;\n}\n\nmpl.figure.prototype.handle_message = function(fig, msg) {\n    fig.message.textContent = msg['message'];\n}\n\nmpl.figure.prototype.handle_draw = function(fig, msg) {\n    // Request the server to send over a new figure.\n    fig.send_draw_message();\n}\n\nmpl.figure.prototype.handle_image_mode = function(fig, msg) {\n    fig.image_mode = msg['mode'];\n}\n\nmpl.figure.prototype.updated_canvas_event = function() {\n    // Called whenever the canvas gets updated.\n    this.send_message(\"ack\", {});\n}\n\n// A function to construct a web socket function for onmessage handling.\n// Called in the figure constructor.\nmpl.figure.prototype._make_on_message_function = function(fig) {\n    return function socket_on_message(evt) {\n        if (evt.data instanceof Blob) {\n            /* FIXME: We get \"Resource interpreted as Image but\n             * transferred with MIME type text/plain:\" errors on\n             * Chrome.  But how to set the MIME type?  It doesn't seem\n             * to be part of the websocket stream */\n            evt.data.type = \"image/png\";\n\n            /* Free the memory for the previous frames */\n            if (fig.imageObj.src) {\n                (window.URL || window.webkitURL).revokeObjectURL(\n                    fig.imageObj.src);\n            }\n\n            fig.imageObj.src = (window.URL || window.webkitURL).createObjectURL(\n                evt.data);\n            fig.updated_canvas_event();\n            fig.waiting = false;\n            return;\n        }\n        else if (typeof evt.data === 'string' && evt.data.slice(0, 21) == \"data:image/png;base64\") {\n            fig.imageObj.src = evt.data;\n            fig.updated_canvas_event();\n            fig.waiting = false;\n            return;\n        }\n\n        var msg = JSON.parse(evt.data);\n        var msg_type = msg['type'];\n\n        // Call the  \"handle_{type}\" callback, which takes\n        // the figure and JSON message as its only arguments.\n        try {\n            var callback = fig[\"handle_\" + msg_type];\n        } catch (e) {\n            console.log(\"No handler for the '\" + msg_type + \"' message type: \", msg);\n            return;\n        }\n\n        if (callback) {\n            try {\n                // console.log(\"Handling '\" + msg_type + \"' message: \", msg);\n                callback(fig, msg);\n            } catch (e) {\n                console.log(\"Exception inside the 'handler_\" + msg_type + \"' callback:\", e, e.stack, msg);\n            }\n        }\n    };\n}\n\n// from http://stackoverflow.com/questions/1114465/getting-mouse-location-in-canvas\nmpl.findpos = function(e) {\n    //this section is from http://www.quirksmode.org/js/events_properties.html\n    var targ;\n    if (!e)\n        e = window.event;\n    if (e.target)\n        targ = e.target;\n    else if (e.srcElement)\n        targ = e.srcElement;\n    if (targ.nodeType == 3) // defeat Safari bug\n        targ = targ.parentNode;\n\n    // jQuery normalizes the pageX and pageY\n    // pageX,Y are the mouse positions relative to the document\n    // offset() returns the position of the element relative to the document\n    var x = e.pageX - $(targ).offset().left;\n    var y = e.pageY - $(targ).offset().top;\n\n    return {\"x\": x, \"y\": y};\n};\n\n/*\n * return a copy of an object with only non-object keys\n * we need this to avoid circular references\n * http://stackoverflow.com/a/24161582/3208463\n */\nfunction simpleKeys (original) {\n  return Object.keys(original).reduce(function (obj, key) {\n    if (typeof original[key] !== 'object')\n        obj[key] = original[key]\n    return obj;\n  }, {});\n}\n\nmpl.figure.prototype.mouse_event = function(event, name) {\n    var canvas_pos = mpl.findpos(event)\n\n    if (name === 'button_press')\n    {\n        this.canvas.focus();\n        this.canvas_div.focus();\n    }\n\n    var x = canvas_pos.x * mpl.ratio;\n    var y = canvas_pos.y * mpl.ratio;\n\n    this.send_message(name, {x: x, y: y, button: event.button,\n                             step: event.step,\n                             guiEvent: simpleKeys(event)});\n\n    /* This prevents the web browser from automatically changing to\n     * the text insertion cursor when the button is pressed.  We want\n     * to control all of the cursor setting manually through the\n     * 'cursor' event from matplotlib */\n    event.preventDefault();\n    return false;\n}\n\nmpl.figure.prototype._key_event_extra = function(event, name) {\n    // Handle any extra behaviour associated with a key event\n}\n\nmpl.figure.prototype.key_event = function(event, name) {\n\n    // Prevent repeat events\n    if (name == 'key_press')\n    {\n        if (event.which === this._key)\n            return;\n        else\n            this._key = event.which;\n    }\n    if (name == 'key_release')\n        this._key = null;\n\n    var value = '';\n    if (event.ctrlKey && event.which != 17)\n        value += \"ctrl+\";\n    if (event.altKey && event.which != 18)\n        value += \"alt+\";\n    if (event.shiftKey && event.which != 16)\n        value += \"shift+\";\n\n    value += 'k';\n    value += event.which.toString();\n\n    this._key_event_extra(event, name);\n\n    this.send_message(name, {key: value,\n                             guiEvent: simpleKeys(event)});\n    return false;\n}\n\nmpl.figure.prototype.toolbar_button_onclick = function(name) {\n    if (name == 'download') {\n        this.handle_save(this, null);\n    } else {\n        this.send_message(\"toolbar_button\", {name: name});\n    }\n};\n\nmpl.figure.prototype.toolbar_button_onmouseover = function(tooltip) {\n    this.message.textContent = tooltip;\n};\nmpl.toolbar_items = [[\"Home\", \"Reset original view\", \"fa fa-home icon-home\", \"home\"], [\"Back\", \"Back to  previous view\", \"fa fa-arrow-left icon-arrow-left\", \"back\"], [\"Forward\", \"Forward to next view\", \"fa fa-arrow-right icon-arrow-right\", \"forward\"], [\"\", \"\", \"\", \"\"], [\"Pan\", \"Pan axes with left mouse, zoom with right\", \"fa fa-arrows icon-move\", \"pan\"], [\"Zoom\", \"Zoom to rectangle\", \"fa fa-square-o icon-check-empty\", \"zoom\"], [\"\", \"\", \"\", \"\"], [\"Download\", \"Download plot\", \"fa fa-floppy-o icon-save\", \"download\"]];\n\nmpl.extensions = [\"eps\", \"jpeg\", \"pdf\", \"png\", \"ps\", \"raw\", \"svg\", \"tif\"];\n\nmpl.default_extension = \"png\";var comm_websocket_adapter = function(comm) {\n    // Create a \"websocket\"-like object which calls the given IPython comm\n    // object with the appropriate methods. Currently this is a non binary\n    // socket, so there is still some room for performance tuning.\n    var ws = {};\n\n    ws.close = function() {\n        comm.close()\n    };\n    ws.send = function(m) {\n        //console.log('sending', m);\n        comm.send(m);\n    };\n    // Register the callback with on_msg.\n    comm.on_msg(function(msg) {\n        //console.log('receiving', msg['content']['data'], msg);\n        // Pass the mpl event to the overridden (by mpl) onmessage function.\n        ws.onmessage(msg['content']['data'])\n    });\n    return ws;\n}\n\nmpl.mpl_figure_comm = function(comm, msg) {\n    // This is the function which gets called when the mpl process\n    // starts-up an IPython Comm through the \"matplotlib\" channel.\n\n    var id = msg.content.data.id;\n    // Get hold of the div created by the display call when the Comm\n    // socket was opened in Python.\n    var element = $(\"#\" + id);\n    var ws_proxy = comm_websocket_adapter(comm)\n\n    function ondownload(figure, format) {\n        window.open(figure.imageObj.src);\n    }\n\n    var fig = new mpl.figure(id, ws_proxy,\n                           ondownload,\n                           element.get(0));\n\n    // Call onopen now - mpl needs it, as it is assuming we've passed it a real\n    // web socket which is closed, not our websocket->open comm proxy.\n    ws_proxy.onopen();\n\n    fig.parent_element = element.get(0);\n    fig.cell_info = mpl.find_output_cell(\"<div id='\" + id + \"'></div>\");\n    if (!fig.cell_info) {\n        console.error(\"Failed to find cell for figure\", id, fig);\n        return;\n    }\n\n    var output_index = fig.cell_info[2]\n    var cell = fig.cell_info[0];\n\n};\n\nmpl.figure.prototype.handle_close = function(fig, msg) {\n    var width = fig.canvas.width/mpl.ratio\n    fig.root.unbind('remove')\n\n    // Update the output cell to use the data from the current canvas.\n    fig.push_to_output();\n    var dataURL = fig.canvas.toDataURL();\n    // Re-enable the keyboard manager in IPython - without this line, in FF,\n    // the notebook keyboard shortcuts fail.\n    IPython.keyboard_manager.enable()\n    $(fig.parent_element).html('<img src=\"' + dataURL + '\" width=\"' + width + '\">');\n    fig.close_ws(fig, msg);\n}\n\nmpl.figure.prototype.close_ws = function(fig, msg){\n    fig.send_message('closing', msg);\n    // fig.ws.close()\n}\n\nmpl.figure.prototype.push_to_output = function(remove_interactive) {\n    // Turn the data on the canvas into data in the output cell.\n    var width = this.canvas.width/mpl.ratio\n    var dataURL = this.canvas.toDataURL();\n    this.cell_info[1]['text/html'] = '<img src=\"' + dataURL + '\" width=\"' + width + '\">';\n}\n\nmpl.figure.prototype.updated_canvas_event = function() {\n    // Tell IPython that the notebook contents must change.\n    IPython.notebook.set_dirty(true);\n    this.send_message(\"ack\", {});\n    var fig = this;\n    // Wait a second, then push the new image to the DOM so\n    // that it is saved nicely (might be nice to debounce this).\n    setTimeout(function () { fig.push_to_output() }, 1000);\n}\n\nmpl.figure.prototype._init_toolbar = function() {\n    var fig = this;\n\n    var nav_element = $('<div/>')\n    nav_element.attr('style', 'width: 100%');\n    this.root.append(nav_element);\n\n    // Define a callback function for later on.\n    function toolbar_event(event) {\n        return fig.toolbar_button_onclick(event['data']);\n    }\n    function toolbar_mouse_event(event) {\n        return fig.toolbar_button_onmouseover(event['data']);\n    }\n\n    for(var toolbar_ind in mpl.toolbar_items){\n        var name = mpl.toolbar_items[toolbar_ind][0];\n        var tooltip = mpl.toolbar_items[toolbar_ind][1];\n        var image = mpl.toolbar_items[toolbar_ind][2];\n        var method_name = mpl.toolbar_items[toolbar_ind][3];\n\n        if (!name) { continue; };\n\n        var button = $('<button class=\"btn btn-default\" href=\"#\" title=\"' + name + '\"><i class=\"fa ' + image + ' fa-lg\"></i></button>');\n        button.click(method_name, toolbar_event);\n        button.mouseover(tooltip, toolbar_mouse_event);\n        nav_element.append(button);\n    }\n\n    // Add the status bar.\n    var status_bar = $('<span class=\"mpl-message\" style=\"text-align:right; float: right;\"/>');\n    nav_element.append(status_bar);\n    this.message = status_bar[0];\n\n    // Add the close button to the window.\n    var buttongrp = $('<div class=\"btn-group inline pull-right\"></div>');\n    var button = $('<button class=\"btn btn-mini btn-primary\" href=\"#\" title=\"Stop Interaction\"><i class=\"fa fa-power-off icon-remove icon-large\"></i></button>');\n    button.click(function (evt) { fig.handle_close(fig, {}); } );\n    button.mouseover('Stop Interaction', toolbar_mouse_event);\n    buttongrp.append(button);\n    var titlebar = this.root.find($('.ui-dialog-titlebar'));\n    titlebar.prepend(buttongrp);\n}\n\nmpl.figure.prototype._root_extra_style = function(el){\n    var fig = this\n    el.on(\"remove\", function(){\n\tfig.close_ws(fig, {});\n    });\n}\n\nmpl.figure.prototype._canvas_extra_style = function(el){\n    // this is important to make the div 'focusable\n    el.attr('tabindex', 0)\n    // reach out to IPython and tell the keyboard manager to turn it's self\n    // off when our div gets focus\n\n    // location in version 3\n    if (IPython.notebook.keyboard_manager) {\n        IPython.notebook.keyboard_manager.register_events(el);\n    }\n    else {\n        // location in version 2\n        IPython.keyboard_manager.register_events(el);\n    }\n\n}\n\nmpl.figure.prototype._key_event_extra = function(event, name) {\n    var manager = IPython.notebook.keyboard_manager;\n    if (!manager)\n        manager = IPython.keyboard_manager;\n\n    // Check for shift+enter\n    if (event.shiftKey && event.which == 13) {\n        this.canvas_div.blur();\n        event.shiftKey = false;\n        // Send a \"J\" for go to next cell\n        event.which = 74;\n        event.keyCode = 74;\n        manager.command_mode();\n        manager.handle_keydown(event);\n    }\n}\n\nmpl.figure.prototype.handle_save = function(fig, msg) {\n    fig.ondownload(fig, null);\n}\n\n\nmpl.find_output_cell = function(html_output) {\n    // Return the cell and output element which can be found *uniquely* in the notebook.\n    // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n    // IPython event is triggered only after the cells have been serialised, which for\n    // our purposes (turning an active figure into a static one), is too late.\n    var cells = IPython.notebook.get_cells();\n    var ncells = cells.length;\n    for (var i=0; i<ncells; i++) {\n        var cell = cells[i];\n        if (cell.cell_type === 'code'){\n            for (var j=0; j<cell.output_area.outputs.length; j++) {\n                var data = cell.output_area.outputs[j];\n                if (data.data) {\n                    // IPython >= 3 moved mimebundle to data attribute of output\n                    data = data.data;\n                }\n                if (data['text/html'] == html_output) {\n                    return [cell, data, j];\n                }\n            }\n        }\n    }\n}\n\n// Register the function which deals with the matplotlib target/channel.\n// The kernel may be null if the page has been refreshed.\nif (IPython.notebook.kernel != null) {\n    IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n}\n",
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {}
    },
    {
     "output_type": "display_data",
     "data": {
      "text/html": [
       "<img src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABdwAAAH0CAYAAAAnhe8sAAAgAElEQVR4XuzdDdBdZXkv/CtKnKCvSGJrlQrahoSA1DGUDvLRIq3ThNI0HGmt7aslflDQSTWmVURHC9W2ehBsbPxqx0rp26FVKyihSTqtBUqAUqhoG76SDOkoULUnAaUQThDeWbsmKk+yedZ+1t73x/o9M87px9rrvu7ff527177Y7D3r8ccffzz8ESBAgAABAgQIECBAgAABAgQIECBAgAABAjMSmGXgPiM/LyZAgAABAgQIECBAgAABAgQIECBAgAABAgMBA3cPAgECBAgQIECAAAECBAgQIECAAAECBAgQ6EDAwL0DRLcgQIAAAQIECBAgQIAAAQIECBAgQIAAAQIG7p4BAgQIECBAgAABAgQIECBAgAABAgQIECDQgYCBeweIbkGAAAECBAgQIECAAAECBAgQIECAAAECBAzcPQMECBAgQIAAAQIECBAgQIAAAQIECBAgQKADAQP3DhDdggABAgQIECBAgAABAgQIECBAgAABAgQIGLh7BggQIECAAAECBAgQIECAAAECBAgQIECAQAcCBu4dILoFAQIECBAgQIAAAQIECBAgQIAAAQIECBAwcPcMECBAgAABAgQIECBAgAABAgQIECBAgACBDgQM3DtAdAsCBAgQIECAAAECBAgQIECAAAECBAgQIGDg7hkgQIAAAQIECBAgQIAAAQIECBAgQIAAAQIdCBi4d4DoFgQIECBAgAABAgQIECBAgAABAgQIECBAwMDdM0CAAAECBAgQIECAAAECBAgQIECAAAECBDoQMHDvANEtCBAgQIAAAQIECBAgQIAAAQIECBAgQICAgbtngAABAgQIECBAgAABAgQIECBAgAABAgQIdCBg4N4BolsQIECAAAECBAgQIECAAAECBAgQIECAAAEDd88AAQIECBAgQIAAAQIECBAgQIAAAQIECBDoQMDAvQNEtyBAgAABAgQIECBAgAABAgQIECBAgAABAgbungECBAgQIECAAAECBAgQIECAAAECBAgQINCBgIF7B4huQYAAAQIECBAgQIAAAQIECBAgQIAAAQIEDNw9AwQIECBAgAABAgQIECBAgAABAgQIECBAoAMBA/cOEN2CAAECBAgQIECAAAECBAgQIECAAAECBAgYuHsGCBAgQIAAAQIECBAgQIAAAQIECBAgQIBABwIG7h0gugUBAgQIECBAgAABAgQIECBAgAABAgQIEDBw9wwQIECAAAECBAgQIECAAAECBAgQIECAAIEOBAzcO0B0CwIECBAgQIAAAQIECBAgQIAAAQIECBAgYODuGSBAgAABAgQIECBAgAABAgQIECBAgAABAh0IGLh3gOgWBAgQIECAAAECBAgQIECAAAECBAgQIEDAwN0zQIAAAQIECBAgQIAAAQIECBAgQIAAAQIEOhAwcO8A0S0IECBAgAABAgQIECBAgAABAgQIECBAgICBu2eAAAECBAgQIECAAAECBAgQIECAAAECBAh0IGDg3gGiWxAgQIAAAQIECBAgQIAAAQIECBAgQIAAAQN3zwABAgQIECBAgAABAgQIECBAgAABAgQIEOhAwMC9A0S3IECAAAECBAgQIECAAAECBAgQIECAAAECBu6eAQIECBAgQIAAAQIECBAgQIAAAQIECBAg0IGAgXsHiG5BgAABAgQIECBAgAABAgQIECBAgAABAgQM3D0DBAgQIECAAAECBAgQIECAAAECBAgQIECgAwED9w4Q3YIAAQIECBAgQIAAAQIECBAgQIAAAQIECBi4ewYIECBAgEALgRe+8IVx9NFHx7p161q8yqUECBAgQIAAAQIECBAgQIBAHwQM3PuQsj0SIECAQGcCBu6dUboRAQIECBAgQIAAAQIECBCoTsDAvbpIbYgAAQIEGoGHHnoonv70p3eOYeDeOakbEiBAgAABAgRmJDCuvm9GRXkxAQIECPRWwMC9t9HbOAECBOoROP/88+OCCy6IW265Jf7gD/4g/uEf/iHmzJkTV155ZXzwgx+MG2+8Mb7+9a/Hj/zIj8Txxx8f73//++MFL3jBXoBLLrkkXvva18YXv/jF+PSnPx2f+cxn4vHHH4+TTz451q5dG4cccsjea/c1cP/oRz8ab37zm+Nd73rXoA5/BAgQIECAAAEC4xHYX9933333xc033xy/93u/F9ddd93gwxdHHnlknHfeefHKV77yB4pp/vdve9vb4ktf+lI8+9nPjte97nWD3vCss86Ku+++O5p+zx8BAgQIEBhVwMB9VDmvI0CAAIFsBPa88WreKL3qVa+Kl7/85fHf//3fsXv37vi3f/u3WLx4ccybNy+aN2If+9jH4rbbbhv854d+6IcGe9gzcP/xH//xOO200+IXf/EX46tf/ergjdhLXvKSwSB+z9/3D9yboXxzzYc//OH4kz/5k1ixYkU2JgohQIAAAQIECNQosL++76CDDoqlS5fGcccdFytXroxnPetZ8Vd/9VeDPu9Tn/rU3j7tK1/5yuCahQsXxjvf+c7BvxH58Y9/fNAbbt++3cC9xofGnggQIDBhAQP3CYNbjgABAgS6F9jzxus973nP0E+Yf+c734ldu3YNPunefBK++VT69w/c3/SmN8VHPvKRvQVeeOGF8fa3v30wqH/uc587+J/vGbg3n4J/zWteE3//938ff/M3fxM/93M/1/3G3JEAAQIECBAgQOAHBPbX9zWfZj/wwAPjpptuigMOOGDva5YtWzb4tyC/9rWvxVOe8pTBp93/9m//djBc3/Phi8ceeyx+4id+YjB09wl3DxwBAgQIzFTAwH2mgl5PgAABAskF9rzx+vKXvxwvfvGL99bz4IMPxnvf+97BQLx5U9UM3Pf8nXPOOYNPuzd/ez7hvmHDhliyZMneazZu3Dj4pFTzlTTNJ6Gav2bg/rznPW/wX99zzz2DN2xHH310cgMFECBAgAABAgT6ILCvvm/r1q2xYMGCwVcJvuUtb/kBhj/90z+N5kMVzTC9Gco3H7xo+rovfOELP3Bd87WAzb0N3PvwFNkjAQIExitg4D5eX3cnQIAAgQkI7Hnj9Y1vfCN++Id/eO+Kv/RLvzT4Pvd3v/vd8VM/9VPR/KvGs2bNil/4hV8Y/KcZtDd/ewbu//Iv/xLHHnvs3tdfffXVccopp8Q//uM/xste9rLB/7wZuH/729+OnTt3Dr7n8xOf+MQEdmgJAgQIECBAgACBRmBffd+mTZvipJNOGgp07bXXxk//9E8PPv3e/HZPM4j//r/ma2Xe+MY3Grh7zAgQIEBgxgIG7jMmdAMCBAgQSC2w543XN7/5zb3/avADDzwQc+fOjd/93d8d/GfP3yOPPBLPeMYz4tWvfvXIA/fmE+2/8iu/Eq9//esH/2k+Kd/8K8r+CBAgQIAAAQIExiuwr77vzjvvjEWLFg1+IPUVr3jFPgs44ogj4pnPfOZ+P+G+574+4T7e/NydAAECfRAwcO9DyvZIgACBygX29cbrW9/61uDHsv7wD/8w3vGOd+wVaL6jvfkhrTPPPHNGA/d169bFZz/72fj1X//1wXeB/vmf/3k89alPrVza9ggQIECAAAECaQX21fc1FTU/gtp8rcxVV101tEDf4Z42P6sTIECgDwIG7n1I2R4JECBQucD+3nidfPLJ8e///u/xgQ98YPBVMNdcc0188pOfjIcffjiWL18+44F7w9p8h/sv//Ivx6mnnhqXXXZZPO1pT6tc2/YIECBAgAABAukE9tf3NV8B2PRjTf+3YsWK+NEf/dHYsWNH3H777fGv//qv0fzgffPX/ObPS1/60sGA/l3vetfgh1abr5PZvHlz/Md//MfgP4cddli6DVqZAAECBIoXMHAvPkIbIECAAIH9vfFqftS0+eGsL37xi/Hoo4/GiSeeOPgxrdNOO23wneyjfod785UyzSfc9/w13/W+bNmywXeHfu5znxu8cfNHgAABAgQIECDQvcD++r5mpa985Svx+7//+9H0Zs3v7Tz72c+Oo446avBvI5599tl7i7nuuuvid37nd+LWW28dfAXha17zmsHXEp577rlx//33D/4tSX8ECBAgQGBUAQP3UeW8jgABAgQIECBAgAABAgQIEKhC4Od//udj+/btcdddd1WxH5sgQIAAgXQCBu7p7K1MgAABAgQIECBAgAABAgQITFhg9erVsXjx4jj00EMHXzvzl3/5l4N/S7H56sHXve51E67GcgQIECBQm4CBe22J2g8BAgQIECBAgAABAgQIECCwX4HmKwe/8IUvxH/+53/GrFmzBl87s2rVqnj1q19NjQABAgQIzFjAwH3GhG5AgAABAgQIECBAgAABAgQIECBAgAABAgQiDNw9BQQIECBAgAABAgQIECBAgAABAgQIECBAoAMBA/cOEN2CAAECBAgQIECAAAECBAgQIECAAAECBAgYuHsGCBAgQIAAAQIECBAgQIAAAQIECBAgQIBABwIG7h0gukUZAo899ljce++98cxnPnPwwzj+CBAgQIAAAQJPJvD444/Ht7/97TjkkEPiKU95ypNd7n+fSECflwjesgQIECBAoGABfV7B4WVeuoF75gEprzuBr33ta3HooYd2d0N3IkCAAAECBHoj8NWvfjWe//zn92a/pW1Un1daYuolQIAAAQL5COjz8smilkoM3GtJ0j6eVOCBBx6Igw8+OJqD9KCDDnrS611AgAABAgQIEPjWt741+Af2999/fzzrWc8CkqmAPi/TYJRFgAABAgQyFtDnZRxO4aUZuBceoPKnL9AcpM0b5eYNmYH79N1cSYAAAQIE+iygfygjfTmVkZMqCRAgQIBATgL6h5zSqKsWA/e68rSbIQIOUo8HAQIECBAg0FZA/9BWLM31ckrjblUCBAgQIFCygP6h5PTyrt3APe98VNehgIO0Q0y3IkCAAAECPRHQP5QRtJzKyEmVBAgQIEAgJwH9Q05p1FWLgXtdedrNEAEHqceDAAECBAgQaCugf2grluZ6OaVxtyoBAgQIEChZQP9Qcnp5127gnnc+qutQwEHaIaZbESBAgACBngjoH8oIWk5l5KRKAgQIECCQk4D+Iac06qrFwL2uPO1miICD1ONBgAABAgQItBXQP7QVS3O9nNK4W5UAAQIECJQsoH8oOb28azdwzzsf1XUo4CDtENOtCBAgQIBATwT0D2UELacyclIlAQIECBDISUD/kFMaddVi4F5XnnYzRMBB6vEgQIAAAQIE2groH9qKpbleTmncrUqAAAECBEoW0D+UnF7etRu4552P6joUcJB2iOlWBAgQIECgJwL6hzKCllMZOamSAAECBAjkJKB/yCmNumoxcK8rT7sZIuAg9XgQIECAAAECbQX0D23F0lwvpzTuViVAgAABAiUL6B9KTi/v2g3c885HdR0KOEg7xHQrAgQIECDQEwH9QxlBy6mMnFRJgAABAgRyEtA/5JRGXbUYuNeVp90MEXCQejwIECBAgACBtgL6h7Ziaa6XUxp3qxIgQIAAgZIF9A8lp5d37Qbueeejug4FHKQdYroVAQIECBDoiYD+oYyg5VRGTqokQIAAAQI5CegfckqjrloM3OvK026GCDhIPR4ECBAgQIBAWwH9Q1uxNNfLKY27VQkQIECAQMkC+oeS08u7dgP3vPNRXYcCDtIOMd2KAAECBAj0RED/UEbQciojJ1USIECAAIGcBPQPOaVRVy0G7nXlWcRurr322rjwwgvjlltuifvuuy8uv/zyOP3006dV+6ZNm+Lkk0+Oo48+Om699dZpvWbPRQ7SVlwuJkCAAAECBCJC/9DuMdDntfNyNQECBAgQIJBOQJ+Xzr72lQ3ca084w/2tX78+msH5McccE2eccca0B+4PPPDA4DWHH354fP3rXzdwzzBbJREgQIAAgdoEvBFrl6g+r52XqwkQIECAAIF0Avq8dPa1r2zgXnvCme9v1qxZ0x64v+pVr4oFCxbEU5/61LjiiisM3DPPVnn9Edi1a3ds2LA1tm3bGfPnz42lSw+POXNm9wfATgkQqFrAG7HR4629z9sdu2NL3BU7Y2fMjbmxIBbG7PB//0Z/YrySAAECBAhMVkCfN1nvPq1m4N6ntDPc63TfiH3qU5+Kj370o3HDDTfE+973vmkN3B955JFo/rPnrzlIDz300Gg+KX/QQQdlqKEkAuUJNMP2lSvXD4bte/6aofvatacaupcXp4oJENiHgDdioz8WNfd5zbD9qrgydsSOvUDzYl6cFssM3Ud/ZLySAAECBAhMVECfN1HuXi1m4N6ruPPb7HTeiG3ZsiVOOumk+Kd/+qdYuHBhnH/++dMauDfXXXDBBVM2beCe33OgonIFrrji9liz5qYpG1i16rhYvnxRuRtTOQECBL4r4I3Y6I9CzX3ebbE5bowbpuAcHyfEkXHU6GheSYAAAQIECExMQJ83MereLWTg3rvI89rwk70R+853vhMvfelL4/Wvf32cc845g+KnO3D3Cfe8slZNnQIXXXR9rFu3Zcrmli1bGKtXH1/npu2KAIFeCXgjNnrcNfd5m+K6uDPumIJzRCyKE+Ok0dG8kgABAgQIEJiYgD5vYtS9W8jAvXeR57XhJ3sjdv/998fcuXMH39u+5++xxx6Lxx9/fPA/+7u/+7v42Z/92WltykE6LSYXEWgl4BPurbhcTIBAgQL6h9FDq7nP8wn30Z8LryRAgAABArkI6PNySaK+Ogzc68u0qB092RuxZrh+2223/cCemu9y/+IXvxif/exn48d+7MfiGc94xrT27CCdFpOLCLQS8B3urbhcTIBAgQL6h9FDq7nP8x3uoz8XXkmAAAECBHIR0OflkkR9dRi415dp9jt68MEHY+vWrYM6Fy9eHBdffHGccsopMW/evDjssMPivPPOi3vuuScuvfTSfe5lul8p88QXO0izfzQUWKhAM3TfuHHb4IdTmx9MXbJkvh9MLTRLZRMgMFVA/9DuqehTn9cM3bfGlsEPpzY/mHp4LPCDqe0eF1cTIECAAIGkAvq8pPxVL27gXnW8eW7u6quvHgzYn/h35plnxiWXXBIrVqyI7du3R3Pdvv4M3PPMVVUECBAgQKBGAW/E2qWqz2vn5WoCBAgQIEAgnYA+L5197SsbuNeesP3tFXCQehgIECBAgACBtgL6h7Ziaa6XUxp3qxIgQIAAgZIF9A8lp5d37Qbueeejug4FHKQdYroVAQIECBDoiYD+oYyg5VRGTqokQIAAAQI5CegfckqjrloM3OvK026GCDhIPR4ECBAgQIBAWwH9Q1uxNNfLKY27VQkQIECAQMkC+oeS08u7dgP3vPNRXYcCDtIOMd2KAAECBAj0RED/UEbQciojJ1USIECAAIGcBPQPOaVRVy0G7nXlaTdDBBykHg8CBAgQIECgrYD+oa1YmuvllMbdqgQIECBAoGQB/UPJ6eVdu4F73vmorkMBB2mHmG5FgAABAgR6IqB/KCNoOZWRkyoJECBAgEBOAvqHnNKoqxYD97rytJshAg5SjwcBAgQIECDQVkD/0FYszfVySuNuVQIECBAgULKA/qHk9PKu3cA973xU16GAg7RDTLciQIAAAQI9EdA/lBG0nMrISZUECBAgQCAnAf1DTmnUVYuBe1152s0QAQepx4MAAQIECBBoK6B/aCuW5no5pXG3KgECBAgQKFlA/1ByennXbuCedz6q61DAQdohplsRIECAAIGeCOgfyghaTmXkpEoCBAgQIJCTgP4hpzTqqsXAva487WaIgIPU40GAAAECBAi0FdA/tBVLc72c0rhblQABAgQIlCygfyg5vbxrN3DPOx/VdSjgIO0Q060IECBAgEBPBPQPZQQtpzJyUiUBAgQIEMhJQP+QUxp11WLgXleedjNEwEHq8SBAoEaBXbt2x4YNW2Pbtp0xf/7cWLr08JgzZ3aNW7UnAkkE9A9J2FsvmiKn3bE7tsRdsTN2xtyYGwtiYcwO52/r8LyAAAECBAgkEkjRPyTaqmUnLGDgPmFwy6UTcJCms7cyAQLjEWiG7StXrh8M2/f8NUP3tWtPNXQfD7m79lBA/1BG6JPOqRm2XxVXxo7YsRdoXsyL02KZoXsZj4wqCRAgQIBATLp/QN4fAQP3/mTd+506SHv/CAAgUJ3AFVfcHmvW3DRlX6tWHRfLly+qbr82RCCFgP4hhXr7NSed022xOW6MG6YUenycEEfGUe034BUECBAgQIDAxAUm3T9MfIMWTCZg4J6M3sKTFnCQTlrcegQIjFvgoouuj3XrtkxZZtmyhbF69fHjXt79CfRCQP9QRsyTzmlTXBd3xh1TcI6IRXFinFQGmioJECBAgEDPBSbdP/Scu1fbN3DvVdz93qyDtN/52z2BGgV8wr3GVO0pNwH9Q26J7LueSefkE+5lPBeqJECAAAECwwQm3T9Ioz8CBu79ybr3O3WQ9v4RAECgOgHf4V5dpDaUoYD+IcNQ9lHSpHPyHe5lPBeqJECAAAECBu6egRQCBu4p1K2ZRGDSb8SSbNKiBAj0TqAZum/cuG3ww6nND6YuWTLfD6b27imw4XEK6B/GqdvdvVPk1Azdt8aWwQ+nNj+Yengs8IOp3UXqTgQIECBAYOwCKfqHsW/KAlkIGLhnEYMiJiHgIJ2EsjUIECBAgEBdAvqHMvKUUxk5qZIAAQIECOQkoH/IKY26ajFwrytPuxki4CD1eBAgQIAAAQJtBfQPbcXSXC+nNO5WJUCAAAECJQvoH0pOL+/aDdzzzkd1HQo4SDvEdCsCBAgQINATAf1DGUHLqYycVEmAAAECBHIS0D/klEZdtRi415Wn3QwRcJB6PAgQIECAAIG2AvqHtmJprpdTGnerEiBAgACBkgX0DyWnl3ftBu5556O6DgUcpB1iuhUBAgQIEOiJgP6hjKDlVEZOqiRAgAABAjkJ6B9ySqOuWgzc68rTboYIOEg9HgQIECBAgEBbAf1DW7E018spjbtVCRAgQIBAyQL6h5LTy7t2A/e881FdhwIO0g4x3YoAAQIECPREQP9QRtByKiMnVRIgQIAAgZwE9A85pVFXLQbudeVpN0MEHKQeDwIECBAgQKCtgP6hrVia6+WUxt2qBAgQIECgZAH9Q8np5V27gXve+aiuQwEHaYeYbkWAAAECBHoioH8oI2g5lZGTKgkQIECAQE4C+oec0qirFgP3uvK0myECDlKPBwECBAgQINBWQP/QVizN9XJK425VAgQIECBQsoD+oeT08q7dwD3vfFTXoYCDtENMtyJAgAABAj0R0D+UEbScyshJlQQIECBAICcB/UNOadRVi4F7XXnazRABB6nHgwABAgQIEGgroH9oK5bmejmlcbcqAQIECBAoWUD/UHJ6eddu4J53PqrrUMBB2iGmWxEgQIAAgZ4I6B/KCFpOZeSkSgIECBAgkJOA/iGnNOqqxcC9rjztZoiAg9TjQYAAAQIECLQV0D+0FUtzvZzSuFuVAAECBAiULKB/KDm9vGs3cM87H9V1KOAg7RDTrQgQIECAQE8E9A9lBC2nMnJSJQECBAgQyElA/5BTGnXVYuBeV552M0TAQerxIECAAAECBNoK6B/aiqW5Xk5p3K1KgAABAgRKFtA/lJxe3rUbuOedj+o6FHCQdojpVgQIECBAoCcC+ocygpZTGTmpkgABAgQI5CSgf8gpjbpqMXCvK0+7GSLgIPV4ECBAgAABAm0F9A9txdJcL6c07lYlQIAAAQIlC+gfSk4v79oN3PPOR3UdCjhIO8R0KwIECBAg0BMB/UMZQcupjJxUSYAAAQIEchLQP+SURl21GLjXlafdDBFwkHo8CBAgQIAAgbYC+oe2Ymmul1Mad6sSIECAAIGSBfQPJaeXd+0G7nnno7oOBRykHWK6FQECBAgQ6ImA/qGMoOVURk6qJECAAAECOQnoH3JKo65aDNzrytNuhgg4SD0eBAgQIECAQFsB/UNbsTTXyymNu1UJECBAgEDJAvqHktPLu3YD97zzUV2HAg7SDjHdigABAgQI9ERA/1BG0HIqIydVEiBAgACBnAT0DzmlUVctBu515Wk3QwQcpB4PAgQIECBAoK2A/qGtWJrr5ZTG3aoECBAgQKBkAf1DyenlXbuBe975qK5DAQdph5huRYAAAQIEeiKgfygjaDmVkZMqCRAgQIBATgL6h5zSqKsWA/e68rSbIQIOUo8HAQIECBAg0FZA/9BWLM31ckrjblUCBAgQIFCygP6h5PTyrt3APe98VNehgIO0Q0y3IkCAAAECPRHQP5QRtJzKyEmVBAgQIEAgJwH9Q05p1FWLgXtdedrNEAEHqceDAAECBAgQaCugf2grluZ6OaVxtyoBAgQIEChZQP9Qcnp5127gnnc+qutQwEHaIaZbESBAgACBngjoH8oIWk5l5KRKAgQIECCQk4D+Iac06qrFwL2uPO1miICD1ONBgAABAgQItBXQP7QVS3O9nNK4W5UAAQIECJQsoH8oOb28azdwzzsf1XUo4CDtENOtCBAgQIBATwT0D2UELacyclIlAQIECBDISUD/kFMaddVi4F5XnnYzRMBB6vEgQIAAAQIE2groH9qKpbleTmncrUqAAAECBEoW0D+UnF7etRu4552P6joUcJB2iOlWBAgQIECgJwL6hzKCllMZOamSAAECBAjkJKB/yCmNumoxcK8rT7sZIuAg9XgQIECAAAECbQX0D23F0lwvpzTuViVAgAABAiUL6B9KTi/v2g3c885HdR0KOEg7xHQrAgQIECDQEwH9QxlBy6mMnFRJgAABAgRyEtA/5JRGXbUYuNeVp90MEXCQejwIECBAgACBtgL6h7Ziaa6XUxp3qxIgQIAAgZIF9A8lp5d37Qbueeejug4FHKQdYroVAQIECBDoiYD+oYyg5VRGTqokQIAAAQI5CegfckqjrloM3OvK026GCDhIPR4ECBAgQIBAWwH9Q1uxNNfLKY27VQkQIDWkVQUAACAASURBVECAQMkC+oeS08u7dgP3vPNRXYcCDtIOMd2KAAECBAj0RED/UEbQciojJ1USIECAAIGcBPQPOaVRVy0G7nXlaTdDBBykHg8CBAgQIECgrYD+oa1YmuvllMbdqgQIECBAoGQB/UPJ6eVdu4F73vmorkMBB2mHmG5FgAABAgR6IqB/KCNoOZWRkyoJECBAgEBOAvqHnNKoqxYD97rytJshAg5SjwcBAgQIECDQVkD/0FYszfVySuNuVQIECBAgULKA/qHk9PKu3cA973xU16GAg7RDTLciQIAAAQI9EdA/lBG0nMrISZUECBAgQCAnAf1DTmnUVYuBe1152s0QAQepx4MAAQIECBBoK6B/aCuW5no5pXG3KgECBAgQKFlA/1ByennXbuCedz6q61DAQdohplsRIECAAIGeCOgfyghaTmXkpEoCBAgQIJCTgP4hpzTqqsXAva487WaIgIPU40GAAAECBAi0FdA/tBVLc72c0rhblQABAgQIlCygfyg5vbxrN3DPOx/VdSjgIO0Q060IECBAgEBPBPQPZQQtpzJyUiUBAgQIEMhJQP+QUxp11WLgXleedjNEwEHq8SBAgAABAgTaCugf2oqluV5OadytSoAAAQIEShbQP5ScXt61G7jnnY/qOhRwkHaI6VYECBAgQKAnAvqHMoKWUxk5qZIAAQIECOQkoH/IKY26ajFwrytPuxki4CD1eBAgQIAAAQJtBfQPbcXSXC+nNO5WJUCAAAECJQvoH0pOL+/aDdzzzkd1HQo4SDvEdCsCBAgQINATAf1DGUHLqYycVEmAAAECBHIS0D/klEZdtRi415Wn3QwRcJB6PAgQIECAAIG2AvqHtmJprpdTGnerEiBAgACBkgX0DyWnl3ftBu5556O6DgUcpB1iuhUBAgQIEOiJgP6hjKDlVEZOqiRAgAABAjkJ6B9ySqOuWgzc68rTboYIOEg9HgQIECBAgEBbAf1DW7E018spjbtVCRAgQIBAyQL6h5LTy7t2A/e881FdhwIO0g4x3YoAAQIECPREQP9QRtByKiMnVRIgQIAAgZwE9A85pVFXLQbudeVZxG6uvfbauPDCC+OWW26J++67Ly6//PI4/fTT91v7ddddF+eee27ccccd8dBDD8ULXvCCOPvss+Otb31rq/06SFtxuZgAAQIECBCICP1Du8dAn9fOy9UECBAgQIBAOgF9Xjr72lc2cK894Qz3t379+ti0aVMcc8wxccYZZzzpwP1LX/rSYNj+4he/OJ7xjGdEM4BvBu4f+tCH4jd/8zenvUMH6bSpXEiAAAECBAh8V0D/0O5R0Oe183I1AQIECBAgkE5An5fOvvaVDdxrTzjz/c2aNetJB+772sIrXvGKwfD9L/7iL6a9QwfptKlcSIAAAQIECBi4z/gZ0OfNmNANCBAgQIAAgTEKmBONEbfntzZw7/kDkHr7o7wRaz7xfuqpp8b73ve+eMMb3jDtLThIp03lQgIECBAgQMDAfcbPgD5vxoRuQIAAAQIECIxRwJxojLg9v7WBe88fgNTbb/NG7PnPf35885vfjEcffTTOP//8ePe73z20/EceeSSa/+z5aw7SQw89NB544IE46KCDUm/d+gQIECBAgEABAt6IjR6SPm90O68kQIAAAQIExi+gzxu/cV9XMHDva/KZ7LvNG7G77747HnzwwbjxxhvjHe94R6xduzZ+7dd+bb87aYbyF1xwwZT/vYF7JuErgwABAgQIFCDgjdjoIenzRrfzSgIECBAgQGD8Avq88Rv3dQUD974mn8m+27wR+/6Sm6+Tab6//c4779zvTnzCPZOQlUGAAAECBAoW8EZs9PD0eaPbeSUBAgQIECAwfgF93viN+7qCgXtfk89k36O+EXvve98bn/zkJ2P79u3T3omDdNpULiRAgAABAgS+K6B/GP1R0OeNbueVBAgQIECAwPgF9HnjN+7rCgbufU0+4b6br4XZunXroILFixfHxRdfHKecckrMmzcvDjvssDjvvPPinnvuiUsvvXRwzUc+8pHB/3zRokWD//66666LVatWxW/91m8Nfjh1un8O0ulKuY4AAQIECBDYI6B/aPcs6PPaebmaAAECBAgQSCegz0tnX/vKBu61J5zh/q6++urBgP2Jf2eeeWZccsklsWLFisEn15vrmr8//uM/jk984hPRfIf7AQccEPPnz4+zzjorzj777HjKU54y7R06SKdN5UICBAgQIEDguwL6h3aPgj6vnZerCRAgQIAAgXQC+rx09rWvbOBee8L2t1fAQephIECAAAECBNoK6B/aiqW5Xk5p3K1KgAABAgRKFtA/lJxe3rUbuOedj+o6FHCQdojpVgQIECBAoCcC+ocygpZTGTmpkgABAgQI5CSgf8gpjbpqMXCvK0+7GSLgIPV4ECBAgAABAm0F9A9txdJcL6c07lYlQIAAAQIlC+gfSk4v79oN3PPOR3UdCjhIO8R0KwIECBAg0BMB/UMZQcupjJxUSYAAAQIEchLQP+SURl21GLjXlafdDBFwkHo8CBAgQIAAgbYC+oe2Ymmul1Mad6sSIECAAIGSBfQPJaeXd+0G7nnno7oOBRykHWK6FQECBAgQ6ImA/qGMoOVURk6qJECAAAECOQnoH3JKo65aDNzrytNuhgg4SD0eBAgQIECAQFsB/UNbsTTXyymNu1UJECBAgEDJAvqHktPLu3YD97zzUV2HAg7SDjHdigABAgQI9ERA/1BG0HIqIydVEiBAgACBnAT0DzmlUVctBu515Wk3QwQcpB4PAgQIECBAoK2A/qGtWJrr5ZTG3aoECBAgQKBkAf1DyenlXbuBe975qK5DAQdph5huRYAAAQIEeiKgfygjaDmVkZMqCRAgQIBATgL6h5zSqKsWA/e68rSbIQIOUo8HAQIECBAg0FZA/9BWLM31ckrjblUCBAgQIFCygP6h5PTyrt3APe98VNehgIO0Q0y3IkCAAAECPRHQP5QRtJzKyEmVBAgQIEAgJwH9Q05p1FWLgXtdedrNEAEHqceDAAECBAgQaCugf2grluZ6OaVxtyoBAgQIEChZQP9Qcnp5127gnnc+qutQwEHaIWZmt9q1a3ds2LA1tm3bGfPnz42lSw+POXNmZ1alcggQIECgRAH9QxmplZbT7tgdW+Ku2Bk7Y27MjQWxMGaH3qWMp02VBAgQIFCLQGn9Qy3ufdiHgXsfUrbHgYCDtM4HoRm2r1y5fjBs3/PXDN3Xrj3V0L3OyO2KAAECExXQP0yUe+TFSsqpGbZfFVfGjtixd7/zYl6cFssM3Ud+AryQAAECBAi0Fyipf2i/O69IKWDgnlLf2hMVcJBOlHtii11xxe2xZs1NU9Zbteq4WL580cTqsBABAgQI1Cmgfygj15Jyui02x41xwxTY4+OEODKOKgNclQQIECBAoAKBkvqHCrh7tQUD917F3e/NOkjrzP+ii66Pdeu2TNncsmULY/Xq4+vctF0RIECAwMQE9A8To57RQiXltCmuizvjjin7PSIWxYlx0owcvJgAAQIECBCYvkBJ/cP0d+XKHAQM3HNIQQ0TEXCQToR54ov4hPvEyS1IgACBXgnoH8qIu6ScfMK9jGdKlQQIECBQv0BJ/UP9adS1QwP3uvK0myECDtI6H4/mO9zf+Mar4uab742HH340DjzwgDj22EPiYx87zXe41xm5XREgQGCiAvqHiXKPvFhJOe35Dvf/iv+Kh+PhaP775odTfzV+LZ4eTx/ZwAsJECBAgACBdgIl9Q/tdubq1AIG7qkTsP7EBBykE6Oe6ELNwP2cc9bFzTffFw8/vDsOPHB2HHvs8+LjH/9FA/eJJmExAgQI1Cmgfygj19Jyeigeik/HX8XO2Dn4odQD48B4djzbD6eW8bipkgABAgQqESitf6iEvRfbMHDvRcw22Qg4SOt8DnylTJ252hUBAgRyEdA/5JLE8DpKy8nXypTxXKmSAAECBOoWKK1/qDuNunZn4F5XnnYzRMBBWufj4UdT68zVrggQIJCLgP4hlyTqGrj74dQynitVEiBAgEDdAvq8uvNNuTsD95T61p6ogIN0otwTW8wn3CdGbSECBAj0UkD/UEbspeXkE+5lPFeqJECAAIG6BUrrH+pOo67dGbjXlafdDBFwkNb5eDTf4b5y5frYtm3n3g3Onz831q491Xe41xm5XREgQGCiAvqHiXKPvFhpOe354dQdsWPvnufFPN/hPvIT4IUECBAgQKC9QGn9Q/sdekUqAQP3VPLWnbiAg3Ti5BNbsBm6b9y4bTB0b4btS5bMN2yfmL6FCBAgULeA/qGMfEvMqRm6b40t0Qzdm2H74bFg8AOq/ggQIECAAIHJCJTYP0xGxiozFTBwn6mg1xcj4CAtJiqFEiBAgACBbAT0D9lEMbQQOZWRkyoJECBAgEBOAvqHnNKoqxYD97rytJshAg5SjwcBAgQIECDQVkD/0FYszfVySuNuVQIECBAgULKA/qHk9PKu3cA973xU16GAg7RDTLciQIAAAQI9EdA/lBG0nMrISZUECBAgQCAnAf1DTmnUVYuBe1152s0QAQepx4MAAQIECBBoK6B/aCuW5no5pXG3KgECBAgQKFlA/1ByennXbuCedz6q61DAQdohplsRIECAAIGeCOgfyghaTmXkpEoCBAgQIJCTgP4hpzTqqsXAva487WaIgIPU40GAAAECBAi0FdA/tBVLc72c0rhblQABAgQIlCygfyg5vbxrN3DPOx/VdSjgIO0Q060IECBAgEBPBPQPZQQtpzJyUiUBAgQIEMhJQP+QUxp11WLgXleedjNEwEHq8SBAgAABAgTaCugf2oqluV5OadytSoAAAQIEShbQP5ScXt61G7jnnY/qOhRwkHaI6VYECBAgQKAnAvqHMoKWUxk5qZIAAQIECOQkoH/IKY26ajFwrytPuxki4CD1eBAgQIAAAQJtBfQPbcXSXC+nNO5WJUCAAAECJQvoH0pOL+/aDdzzzkd1HQo4SDvEdCsCBAgQINATAf1DGUHLqYycVEmAAAECBHIS0D/klEZdtRi415Wn3QwRcJB6PAgQIECAAIG2AvqHtmJprpdTGnerEiBAgACBkgX0DyWnl3ftBu5556O6DgUcpB1iuhUBAgQIEOiJgP6hjKDlVEZOqiRAgAABAjkJ6B9ySqOuWgzc68rTboYIOEg9HgQIECBAgEBbAf1DW7E018spjbtVCRAgQIBAyQL6h5LTy7t2A/e881FdhwIO0g4x3YoAAQIECPREQP9QRtByKiMnVRIgQIAAgZwE9A85pVFXLQbudeVpN0MEHKQeDwIECBAgQKCtgP6hrVia6+WUxt2qBAgQIECgZAH9Q8np5V27gXve+aiuQwEHaYeYbkWAAAECBHoioH8oI2g5lZGTKgkQIECAQE4C+oec0qirFgP3uvK0myECDlKPBwECBAgQINBWQP/QVizN9XJK425VAgQIECBQsoD+oeT08q7dwD3vfFTXoYCDtENMtyJAgAABAj0R0D+UEbScyshJlQQIECBAICcB/UNOadRVi4F7XXnazRABB6nHgwABAgQIEGgroH9oK5bmejmlcbcqAQIECBAoWUD/UHJ6eddu4J53PqrrUMBB2iGmWxEgQIAAgZ4I6B/KCFpOZeSkSgIECBAgkJOA/iGnNOqqxcC9rjztZoiAg9TjQYAAAQIECLQV0D+0FUtzvZzSuFuVAAECBAiULKB/KDm9vGs3cM87H9V1KOAg7RDTrQgQIECAQE8E9A9lBC2nMnJSJQECBAgQyElA/5BTGnXVYuBeV552M0TAQerxIECAAAECBNoK6B/aiqW5Xk5p3K1KgAABAgRKFtA/lJxe3rUbuOedj+o6FHCQdojpVgQIECBAoCcC+ocygpZTGTmpkgABAgQI5CSgf8gpjbpqMXCvK0+7GSLgIPV4ECBAgAABAm0F9A9txdJcL6c07lYlQIAAAQIlC+gfSk4v79oN3PPOR3UdCjhIO8R0KwIECBAg0BMB/UMZQcupjJxUSYAAAQIEchLQP+SURl21GLjXlafdDBFwkHo8CBAgQIAAgbYC+oe2Ymmul1Mad6sSIECAAIGSBfQPJaeXd+0G7nnno7oOBRykHWK6FQECBAgQ6ImA/qGMoOVURk6qJECAAAECOQnoH3JKo65aDNzrytNuhgg4SD0eBAgQIECAQFsB/UNbsTTXyymNu1UJECBAgEDJAvqHktPLu3YD97zzUV2HAg7SDjHdigABAgQI9ERA/1BG0HIqIydVEiBAgACBnAT0DzmlUVctBu515Wk3QwQcpB4PAgQIECBAoK2A/qGtWJrr5ZTG3aoECBAgQKBkAf1DyenlXbuBe975qK5DAQdph5huRYAAAQIEeiKgfygjaDmVkZMqCRAgQIBATgL6h5zSqKsWA/e68rSbIQIOUo8HAQIECBAg0FZA/9BWLM31ckrjblUCBAgQIFCygP6h5PTyrt3APe98VNehgIO0Q0y3IkCAAAECPRHQP5QRtJzKyEmVBAgQIEAgJwH9Q05p1FWLgXtdedrNEAEHqceDAAECBAgQaCugf2grluZ6OaVxtyoBAgQIEChZQP9Qcnp5127gnnc+qutQwEHaIaZbESBAgACBngjoH8oIWk5l5KRKAgQIECCQk4D+Iac06qrFwL2uPO1miICD1ONBgAABAgQItBXQP7QVS3O9nNK4W5UAAQIECJQsoH8oOb28azdwzzsf1XUo4CDtENOtCBAgQIBATwT0D2UELacyclIlAQIECBDISUD/kFMaddVi4F5XnnYzRMBB6vEgQIAAAQIE2groH9qKpbleTmncrUqAAAECBEoW0D+UnF7etRu4552P6joUcJB2iOlWBAgQIECgJwL6hzKCllMZOamSAAECBAjkJKB/yCmNumoxcK8rT7sZIuAg9XgQIECAAAECbQX0D23F0lwvpzTuViVAgAABAiUL6B9KTi/v2g3c885HdR0KOEg7xHQrAgQIECDQEwH9QxlBy6mMnFRJgAABAgRyEtA/5JRGXbUYuNeVp90MEXCQejwIECBAgACBtgL6h7Ziaa6XUxp3qxIgQIAAgZIF9A8lp5d37Qbueeejug4FHKQdYroVAQIECBDoiYD+oYyg5VRGTqokQIAAAQI5CegfckqjrloM3OvK026GCDhIPR4ECBAgQIBAWwH9Q1uxNNfLKY27VQkQIECAQMkC+oeS08u7dgP3vPNRXYcCDtIOMd2KAAECBAj0RED/UEbQciojJ1USIECAAIGcBPQPOaVRVy0G7nXlaTdDBBykHg8CBAgQIECgrYD+oa1YmuvllMbdqgQIECBAoGQB/UPJ6eVdu4F73vmorkMBB2mHmG5FgAABAgR6IqB/KCNoOZWRkyoJECBAgEBOAvqHnNKoqxYD97ryLGI31157bVx44YVxyy23xH333ReXX355nH766fut/XOf+1x87GMfi1tvvTUeeeSReNGLXhTnn39+LFmypNV+HaStuFxMgAABAgQIRIT+od1joM9r5+VqAgQIECBAIJ2APi+dfe0rG7jXnnCG+1u/fn1s2rQpjjnmmDjjjDOedOC+atWqOOSQQ+KUU06Jgw8+OD71qU/FBz/4wfjnf/7nWLx48bR36CCdNpULCRAgQIAAge8K6B/aPQr6vHZeriZAgAABAgTSCejz0tnXvrKBe+0JZ76/WbNmPenAfV9baD7l/qu/+qvxnve8Z9o7dJBOm8qFBAgQIECAgIH7jJ8Bfd6MCd2AAAECBAgQGKOAOdEYcXt+awP3nj8Aqbc/yhuxxx57LF74whfG29/+9li5cuW0t+AgnTaVCwkQIECAAAED9xk/A/q8GRO6AQECBAgQIDBGAXOiMeL2/NYG7j1/AFJvf5Q3Ys33v7///e+P22+/PZ7znOfsdwvN9703/9nz1xykhx56aDzwwANx0EEHpd669QkQIECAAIECBLwRGz0kfd7odl5JgAABAgQIjF9Anzd+476uYODe1+Qz2XfbN2KXXXZZvOENb4jPf/7z8fKXv3zoLpofVr3gggumXGPgnkn4yiBAgAABAgUIeCM2ekj6vNHtvJIAAQIECBAYv4A+b/zGfV3BwL2vyWey7zZvxP76r/86Xvva18ZnPvOZOO200550Bz7h/qRELiBAgAABAgSeRMAbsdEfEX3e6HZeSYAAAQIECIxfQJ83fuO+rmDg3tfkM9n3dN+INZ9sf93rXhfN/3v66aePVL2DdCQ2LyJAgAABAr0W0D+MHr8+b3Q7ryRAgAABAgTGL6DPG79xX1cwcO9r8gn3/eCDD8bWrVsHFSxevDguvvjiOOWUU2LevHlx2GGHxXnnnRf33HNPXHrppYNrmiH7b/zGb8SaNWviFa94xd7KDzzwwHjWs5417Z04SKdN5UICBAgQIEDguwL6h3aPgj6vnZerCRAgQIAAgXQC+rx09rWvbOBee8IZ7u/qq68eDNif+HfmmWfGJZdcEitWrIjt27dHc13z97KXvSyuueaa/V4/3S06SKcr5ToCBAgQIEBgj4D+od2zoM9r5+VqAgQIECBAIJ2APi+dfe0rG7jXnrD97RVwkHoYCBAgQIAAgbYC+oe2Ymmul1Mad6sSIECAAIGSBfQPJaeXd+0G7nnno7oOBRykHWK6FQECBAgQ6ImA/qGMoOVURk6qJECAAAECOQnoH3JKo65aDNzrytNuhgg4SD0eBMYjsGvX7tiwYWts27Yz5s+fG0uXHh5z5swez2LuSoAAgQkL6B8mDD7icrXltDt2x5a4K3bGzpgbc2NBLIzZ4f+2jvh4eBkBAgQIENinQG39g5jzETBwzycLlYxZwEE6ZmC376VAM2xfuXL9YNi+568Zuq9de6qhey+fCJsmUJ+A/qGMTGvKqRm2XxVXxo7YsRd/XsyL02KZoXsZj6MqCRAgQKAQgZr6h0LIe1OmgXtvorZRB6lngED3AldccXusWXPTlBuvWnVcLF++qPsF3ZEAAQITFtA/TBh8xOVqyum22Bw3xg1TJI6PE+LIOGpEIS8jQIAAAQIEnihQU/8g3bwEDNzzykM1YxRwkI4R1617K3DRRdfHunVbpux/2bKFsXr18b11sXECBOoR0D+UkWVNOW2K6+LOuGMK/BGxKE6Mk8oIRJUECBAgQKAAgZr6hwK4e1WigXuv4u73Zh2k/c7f7scj4BPu43F1VwIE8hHQP+STxbBKasrJJ9zLeOZUSYAAAQLlC9TUP5SfRl07MHCvK0+7GSLgIPV4EOhewHe4d2/qjgQI5CWgf8grj/1VU1NOvsO9jGdOlQQIECBQvkBN/UP5adS1AwP3uvK0GwN3zwCBiQs0Q/eNG7cNfji1+cHUJUvm+8HUiadgQQIExiXgjdi4ZLu9b205NUP3rbFl8MOpzQ+mHh4L/GBqt4+MuxEgQIAAgaitfxBpPgIG7vlkoZIxCzhIxwzs9gQIECBAoEIB/UMZocqpjJxUSYAAAQIEchLQP+SURl21GLjXlafdDBFwkHo8CBAgQIAAgbYC+oe2Ymmul1Mad6sSIECAAIGSBfQPJaeXd+0G7nnno7oOBRykHWK6FQECBAgQ6ImA/qGMoOVURk6qJECAAAECOQnoH3JKo65aDNzrytNuhgg4SD0eBAgQIECAQFsB/UNbsTTXyymNu1UJECBAgEDJAvqHktPLu3YD97zzUV2HAg7SDjHdigABAgQI9ERA/1BG0HIqIydVEiBAgACBnAT0DzmlUVctBu515Wk3QwQcpB4PAgQIECBAoK2A/qGtWJrr5ZTG3aoECBAgQKBkAf1DyenlXbuBe975qK5DAQdph5huRYAAAQIEeiKgfygjaDmVkZMqCRAgQIBATgL6h5zSqKsWA/e68rSbIQIOUo8HAQIECBAg0FZA/9BWLM31ckrjblUCBAgQIFCygP6h5PTyrt3APe98VNehgIO0Q0y3IkCAAAECPRHQP5QRtJzKyEmVBAgQIEAgJwH9Q05p1FWLgXtdedrNEAEHqceDAAECBAgQaCugf2grluZ6OaVxtyoBAgQIEChZQP9Qcnp5127gnnc+qutQwEHaIaZbESBAgACBngjoH8oIWk5l5KRKAgQIECCQk4D+Iac06qrFwL2uPO1miICD1ONBgAABAgQItBXQP7QVS3O9nNK4W5UAAQIECJQsoH8oOb28azdwzzsf1XUo4CDtENOtCBAgQIBATwT0D2UELacyclIlAQIECBDISUD/kFMaddVi4F5XnnYzRMBB6vEgQIAAAQIE2groH9qKpbleTmncrUqAAAECBEoW0D+UnF7etRu4552P6joUcJB2iOlWBAgQIECgJwL6hzKCllMZOamSAAECBAjkJKB/yCmNumoxcK8rT7sZIuAg9XgQIECAAAECbQX0D23F0lwvpzTuViVAgAABAiUL6B9KTi/v2g3c885HdR0KOEg7xHQrAgQIECDQEwH9QxlBy6mMnFRJgAABAgRyEtA/5JRGXbUYuNeVp90MEXCQejwIECBAgACBtgL6h7Ziaa6XUxp3qxIgQIAAgZIF9A8lp5d37Qbueeejug4FHKQdYroVAQIECBDoiYD+oYyg5VRGTqokQIAAAQI5CegfckqjrloM3OvK026GCDhIPR4ECBAgQIBAWwH9Q1uxNNfLKY27VQkQIECAQMkC+oeS08u7dgP3vPNRXYcCDtIOMd2KAAECBAj0RED/UEbQciojJ1USIECAAIGcBPQPOaVRVy0G7nXlaTdDBBykHg8CBAgQIECgrYD+oa1YmuvllMbdqgQIECBAoGQB/UPJ6eVdu4F73vmorkMBB2mHmG5FgAABAgR6IqB/KCNoOZWRkyoJECBAgEBOAvqHnNKoqxYD97rytJshAg5SjwcBAgQIECDQVkD/0FYszfVySuNuVQIECBAgULKA/qHk9PKu3cA973xU16GAjpuhaAAAIABJREFUg7RDTLciQIAAAQI9EdA/lBG0nMrISZUECBAgQCAnAf1DTmnUVYuBe1152s0QAQepx4MAAQIECBBoK6B/aCuW5no5pXG3KgECBAgQKFlA/1ByennXbuCedz6q61DAQdohplsRIECAAIGeCOgfyghaTmXkpEoCBAgQIJCTgP4hpzTqqsXAva487WaIgIPU40GAAAECBAi0FdA/tBVLc72c0rhblQABAgQIlCygfyg5vbxrN3DPOx/VdSjgIO0Q060IECBAgEBPBPQPZQQtpzJyUiUBAgQIEMhJQP+QUxp11WLgXleedjNEwEHq8SAwHoFdu3bHhg1bY9u2nTF//txYuvTwmDNn9ngWc1cCBAhMWED/MGHwEZfra067Y3dsibtiZ+yMuTE3FsTCmB3+b/CIj5GXESBAgEDPBPraP/Qs5iTbNXBPwm7RFAIO0hTq1qxdoBm2r1y5fjBs3/PXDN3Xrj3V0L328O2PQE8E9A9lBN3HnJph+1VxZeyIHXtDmhfz4rRYZuhexmOrSgIECBBILNDH/iExeW+WN3DvTdQ26iD1DBDoXuCKK26PNWtumnLjVauOi+XLF3W/oDsSIEBgwgL6hwmDj7hcH3O6LTbHjXHDFLHj44Q4Mo4aUdLLCBAgQIBAfwT62D/0J920OzVwT+tv9QkKOEgniG2p3ghcdNH1sW7dlin7XbZsYaxefXxvHGyUAIF6BfQPZWTbx5w2xXVxZ9wxJaAjYlGcGCeVEZwqCRAgQIBAQoE+9g8JuXu1tIF7r+Lu92YdpP3O3+7HI+AT7uNxdVcCBPIR0D/kk8WwSvqYk0+4l/FsqpIAAQIE8hXoY/+Qbxp1VWbgXleedjNEwEHq8SDQvYDvcO/e1B0JEMhLQP+QVx77q6aPOfkO9zKeTVUSIECAQL4Cfewf8k2jrsoM3OvK024M3D0DBCYu0AzdN27cNvjh1OYHU5csme8HUyeeggUJEBiXgDdi45Lt9r59zakZum+NLYMfTm1+MPXwWOAHU7t9tNyNAAECBCoW6Gv/UHGk2WzNwD2bKBQybgEH6biF3Z8AAQIECNQnoH8oI1M5lZGTKgkQIECAQE4C+oec0qirFgP3uvK0myECDlKPBwECBAgQINBWQP/QVizN9XJK425VAgQIECBQsoD+oeT08q7dwD3vfFTXoYCDtENMtyJAgAABAj0R0D+UEbScyshJlQQIECBAICcB/UNOadRVi4F7XXnazRABB6nHgwABAgQIEGgroH9oK5bmejmlcbcqAQIECBAoWUD/UHJ6eddu4J53PqrrUMBB2iGmWxEgQIAAgZ4I6B/KCFpOZeSkSgIECBAgkJOA/iGnNOqqxcC9rjztZoiAg9TjQYAAAQIECLQV0D+0FUtzvZzSuFuVAAECBAiULKB/KDm9vGs3cM87H9V1KOAg7RDTrQgQIECAQE8E9A9lBC2nMnJSJQECBAgQyElA/5BTGnXVYuBeV552M0RgnAfprl27Y8OGrbFt286YP39uLF16eMyZM1seBAgQIECAQOEC4+wfCqfJqnw5dRvH7tgdW+Ku2Bk7Y27MjQWxMGaH3rZbZXcjQIAAgdQC+ofUCdS7voF7vdna2RMExnWQNsP2lSvXD4bte/6aofvatacaunsKCRAgQIBA4QLj6h8KZ8mufDl1F0kzbL8qrowdsWPvTefFvDgtlhm6d8fsTgQIECCQgYD+IYMQKi3BwL3SYG1rqsC4DtIrrrg91qy5acqCq1YdF8uXLxIFAQIECBAgULDAuPqHgkmyLF1O3cVyW2yOG+OGKTc8Pk6II+Oo7hZyJwIECBAgkFhA/5A4gIqXN3CvOFxb+0GBcR2kF110faxbt2UK97JlC2P16uPFQIAAAQIECBQsMK7+oWCSLEuXU3exbIrr4s64Y8oNj4hFcWKc1N1C7kSAAAECBBIL6B8SB1Dx8gbuFYdra5MZuPuEuyeNAAECBAjUK+CNWBnZyqm7nHzCvTtLdyJAgACBvAX0D3nnU3J1Bu4lp6f2VgLjOkh9h3urGFxMgAABAgSKEhhX/1AUQgHFyqm7kHyHe3eW7kSAAAECeQvoH/LOp+TqDNxLTk/trQTGeZA2Q/eNG7cNfji1+cHUJUvm+8HUVum4mAABAgQI5Ckwzv4hzx2XWZWcus2tGbpvjS2DH05tfjD18FjgB1O7JXY3AgQIEMhAQP+QQQiVlmDgXmmwtjVVwEHqqSBAgAABAgTaCugf2oqluV5OadytSoAAAQIEShbQP5ScXt61G7jnnY/qOhRwkHaI6VYECBAgQKAnAvqHMoKWUxk5qZIAAQIECOQkoH/IKY26ajFwrytPuxki4CD1eBAgQIAAAQJtBfQPbcXSXC+nNO5WJUCAAAECJQvoH0pOL+/aDdzzzkd1HQo4SDvEdCsCBAgQINATAf1DGUHLqYycVEmAAAECBHIS0D/klEZdtRi415Wn3QwRcJB6PAgQIECAAIG2AvqHtmJprpdTGnerEiBAgACBkgX0DyWnl3ftBu5556O6DgUcpB1iuhUBAgQIEOiJgP6hjKDlVEZOqiRAgAABAjkJ6B9ySqOuWgzc68rTboYIOEg9HgQIECBAgEBbAf1DW7E018spjbtVCRAgQIBAyQL6h5LTy7t2A/e881FdhwIO0g4x3YoAAQIECPREQP9QRtByKiMnVRIgQIAAgZwE9A85pVFXLQbudeVpN0MEHKQeDwIECBAgQKCtgP6hrVia6+WUxt2qBAgQIECgZAH9Q8np5V27gXve+aiuQwEHaYeYbkWAAAECBHoioH8oI2g5lZGTKgkQIECAQE4C+oec0qirFgP3uvK0myECDlKPBwECBAgQINBWQP/QVizN9XJK425VAgQIECBQsoD+oeT08q7dwD3vfFTXoYCDtENMtyJAgAABAj0R0D+UEbScyshJlQQIECBAICcB/UNOadRVi4F7XXnazRABB6nHgwABAgQIEGgroH9oK5bmejmlcbcqAQIECBAoWUD/UHJ6eddu4J53PqrrUMBB2iGmWxEgQIAAgZ4I6B/KCFpOZeSkSgIECBAgkJOA/iGnNOqqxcC9rjztZoiAg9TjQYAAAQIECLQV0D+0FUtzvZzSuFuVAAECBAiULKB/KDm9vGs3cM87H9V1KOAg7RDTrQgQIECAQE8E9A9lBC2nMnJSJQECBAgQyElA/5BTGnXVYuBeV55F7Obaa6+NCy+8MG655Za477774vLLL4/TTz99v7U31/z2b//24PotW7bEm9/85vijP/qj1nt1kLYm8wICBAgQINB7Af1Du0dAn9fOy9UECBAgQIBAOgF9Xjr72lc2cK894Qz3t379+ti0aVMcc8wxccYZZzzpwH379u3xoQ99KH7yJ39y8P+efPLJBu4Z5qokAgQIECBQo4A3Yu1S1ee183I1AQIECBAgkE5An5fOvvaVDdxrTzjz/c2aNetJB+7fv4WXvexl8ZKXvMTAPfNclUeAAAECBGoR8EZs9CT1eaPbeSUBAgQIECAwfgF93viN+7qCgXtfk89k3+N8I/bII49E8589f81Beuihh8YDDzwQBx10UCYCyiBAgAABAgRyFvBGbPR09Hmj23klAQIECBAgMH4Bfd74jfu6goF7X5PPZN/jfCN2/vnnxwUXXDBlpwbumYSvDAIECBAgUICAN2Kjh6TPG93OKwkQIECAAIHxC+jzxm/c1xUM3PuafCb7HucbMZ9wzyRkZRAgQIAAgYIFvBEbPTx93uh2XkmAAAECBAiMX0CfN37jvq5g4N7X5DPZ9zjfiD1xiw7STEJXBgECBAgQKEhA/zB6WPq80e28kgABAgQIEBi/gD5v/MZ9XcHAva/JZ7Jvb8QyCUIZBAgQIECAwD4FvBEb/cHQ541u55UECBAgQIDA+AX0eeM37usKBu59TT7hvh988MHYunXroILFixfHxRdfHKecckrMmzcvDjvssDjvvPPinnvuiUsvvXRvlbfeeuvgv37DG94QRxxxRLztbW+Lpz3taXHUUUdNeycO0mlTuZAAAQIECBD4roD+od2joM9r5+VqAgQIECBAIJ2APi+dfe0rG7jXnnCG+7v66qsHA/Yn/p155plxySWXxIoVK2L79u3RXLfnr/mE1BP/XvCCFwyum+6fg3S6Uq4jQIAAAQIE9gjoH9o9C/q8dl6uJkCAAAECBNIJ6PPS2de+soF77Qnb314BB6mHgQABAgQIEGgroH9oK5bmejmlcbcqAQIECBAoWUD/UHJ6eddu4J53PqrrUMBB2iGmWxEgQIAAgZ4I6B/KCFpOZeSkSgIECBAgkJOA/iGnNOqqxcC9rjztZoiAg9TjQYAAAQIECLQV0D+0FUtzvZzSuFuVAAECBAiULKB/KDm9vGs3cM87H9V1KOAg7RDTrQgQIECAQE8E9A9lBC2nMnJSJQECBAgQyElA/5BTGnXVYuBeV552M0TAQerxIECAAAECBNoK6B/aiqW5Xk5p3K1KgAABAgRKFtA/lJxe3rUbuOedj+o6FHCQdojpVgS+T2DXrt2xYcPW2LZtZ8yfPzeWLj085syZzYgAAQJVCOgfyohRTuPLaXfsji1xV+yMnTE35saCWBizw/+dH5+4OxMgQIDApAT0D5OS7t86Bu79y7y3O3aQ9jZ6Gx+jQDNsX7ly/WDYvuevGbqvXXuqofsY3d2aAIHJCegfJmc9k5XkNBO9/b+2GbZfFVfGjtix96J5MS9Oi2WG7uMhd1cCBAgQmKCA/mGC2D1bysC9Z4H3ebsO0j6nb+/jErjiittjzZqbptx+1arjYvnyReNa1n0JECAwMQH9w8SoZ7SQnGbEt98X3xab48a4Ycr//vg4IY6Mo8azqLsSIECAAIEJCegfJgTdw2UM3HsYel+37CDta/L2PU6Biy66Ptat2zJliWXLFsbq1cePc2n3JkCAwEQE9A8TYZ7xInKaMeE+b7Apros7444p/7sjYlGcGCeNZ1F3JUCAAAECExLQP0wIuofLGLj3MPS+btlBWm/yvkM8XbY+4Z7O3soECExGQP8wGeeZriKnmQru+/Xj+oS774UfT17uSoAAAQLtBPQP7bxcPX0BA/fpW7mycAEHaeEB7qf8Ztj+xjdeFTfffG88/PCjceCBB8Sxxx4SH/vYab5DfAKR+w73CSBbggCBpAL6h6T8015cTtOmanXhOL7Dfc89/yv+Kx6Oh6P575sfY/3V+LV4ejy9VX0uJkCAAAECMxHQP8xEz2uHCRi4ez56I+AgrTPqT396c6xatWEwbN/z1wzd16xZGr/yKy+qc9OZ7aoZum/cuG3ww6nND6YuWTLfP+zILCPlECAwuoD+YXS7Sb5STuPTbgbiW2PL4IdTmx9MPTwWzOgHU5tPzV8fm+L/xP+J78T3+rdD4kfj/43XzOje41NwZwIECBCoUUD/UGOqeezJwD2PHFQxAQEH6QSQEyyxYsXl8fnP3zVl5dNPPyI+9anTE1RkSQIECBCoSUD/UEaaciojp6bK5nvh/zVuiW/Ht36g6APj6bE0TvVjrOVEqVICBAgUL6B/KD7CbDdg4J5tNArrWsBB2rVoHvdbseKK+Pzn75xSzP/6X4viz/5seR5FqoIAAQIEihXQP5QRnZzKyKmpsvmE+4ZYH7vi4R8o+qB4ViyOY/wYazlRqpQAAQLFC+gfio8w2w0YuGcbjcK6FnCQdi2ax/0+/el/i7e8ZWPs2vWdvQXNmfPU+PCHm6+UOTqPIlVBgAABAsUK6B/KiE5OZeTUVNl8Rc3/F5fGfXHv3qIPiAPi2fFDcUKc6BPu5USpUgIECBQvoH8oPsJsN2Dgnm00CutawEHatWge92u+P/ycc9bFzTffFw8/vDsOPHB2HHvs8+LjH/9F3yOeR0SqIECAQNEC+ocy4pNTGTntqfKheCg+HX8VO2Pn4DvbD4wD49nx7DgtlvkO97KiVC0BAgSKFtA/FB1f1sUbuGcdj+K6FHCQdqmZ1738aGdeeaiGAAECNQnoH8pIU05l5PT9VXb9Y6zlCaiYAAECBFIL6B9SJ1Dv+gbu9WZrZ08QcJB6JAgQIECAAIG2AvqHtmJprpdTGnerEiBAgACBkgX0DyWnl3ftBu5556O6DgUcpB1iuhUBAgQIEOiJgP6hjKDlVEZOqiRAgAABAjkJ6B9ySqOuWgzc68rTboYIOEg9HgQIECBAgEBbAf1DW7E018spjbtVCRAgQIBAyQL6h5LTy7t2A/e881FdhwIO0g4x3YoAAQIECPREQP9QRtByKiMnVRIgQIAAgZwE9A85pVFXLQbudeVpN0MEHKQeDwIECBAgQKCtgP6hrVia6+WUxt2qBAgQIECgZAH9Q8np5V27gXve+aiuQwEHaYeYbkWAAAECBHoioH8oI2g5lZGTKgkQIECAQE4C+oec0qirFgP3uvK0myECDtJ6H49du3bHhg1bY9u2nTF//txYuvTwmDNndr0btjMCBAgQmJiA/mFi1DNaSE4z4ivqxbtjd2yJu2Jn7Iy5MTcWxMKYHfq+okJULAECBDIR0D9kEkSFZRi4VxiqLe1bwEFa55PRDNtXrlw/GLbv+WuG7mvXnmroXmfkdkWAAIGJCugfJso98mJyGpmuqBc2w/ar4srYETv21j0v5sVpsczQvagkFUuAAIE8BPQPeeRQYxUG7jWmak/7FHCQ1vlgXHHF7bFmzU1TNrdq1XGxfPmiOjdtVwQIECAwMQH9w8SoZ7SQnGbEV8yLb4vNcWPcMKXe4+OEODKOKmYfCiVAgACBPAT0D3nkUGMVBu41pmpPBu49egYuuuj6WLduy5QdL1u2MFavPr5HErZKgAABAuMQ8EZsHKrd31NO3ZvmeMdNcV3cGXdMKe2IWBQnxkk5lqwmAgQIEMhYQP+QcTiFl2bgXniAyp++gIN0+lYlXekT7iWlpVYCBAiUJ6B/KCMzOZWR00yr9An3mQp6PQECBAh8v4D+wfMwLgED93HJum92Ag7S7CLppCDf4d4Jo5sQIECAwH4E9A9lPBpyKiOnmVbpO9xnKuj1BAgQIGDg7hmYhICB+ySUrZGFgDdiWcQwliKaofvGjdsGP5za/GDqkiXz/WDqWKTdlAABAv0T0D+UkbmcysipiyqbofvW2DL44dTmB1MPjwV+MLULWPcgQIBADwX0Dz0MfUJbNnCfELRl0gs4SNNnoAICBAgQIFCagP6hjMTkVEZOqiRAgAABAjkJ6B9ySqOuWgzc68rTboYIOEg9HgQIECBAgEBbAf1DW7E018spjbtVCRAgQIBAyQL6h5LTy7t2A/e881FdhwIO0g4x3YoAAQIECPREQP9QRtByKiMnVRIgQIAAgZwE9A85pVFXLQbudeVpN0MEHKQeDwIECBAgQKCtgP6hrVia6+WUxt2qBAgQIECgZAH9Q8np5V27gXve+aiuQwEHaYeYbkWAAAECBHoioH8oI2g5lZGTKgkQIECAQE4C+oec0qirFgP3uvK0myECDlKPBwECBAgQINBWQP/QVizN9XJK425VAgQIECBQsoD+oeT08q7dwD3vfFTXoYCDtENMtyJAgAABAj0R0D+UEbScyshJlQQIECBAICcB/UNOadRVi4F7XXnazRABB6nHgwABAgQIEGgroH9oK5bmejmlcbcqAQIECBAoWUD/UHJ6eddu4J53PqrrUMBB2iGmWxEgQIAAgZ4I6B/KCFpOZeSkSgIECBAgkJOA/iGnNOqqxcC9rjztZoiAg9TjQYAAAQIECLQV0D+0FUtzvZzSuFuVAAECBAiULKB/KDm9vGs3cM87H9V1KOAg7RDTrQgQIECAQE8E9A9lBC2nMnJSJQECBAgQyElA/5BTGnXVYuBeV552M0TAQerxIECAAAECBNoK6B/aiqW5Xk5p3K1KgAABAgRKFtA/lJxe3rUbuOedj+o6FHCQdojpVgQIECBAoCcC+ocygpZTGTmpkgABAgQI5CSgf8gpjbpqMXCvK0+7GSLgIPV4ECBAgAABAm0F9A9txdJcL6c07lYlQIAAAQIlC+gfSk4v79oN3PPOR3UdCjhIO8R0KwIECBAg0BMB/UMZQcupjJxUSYAAAQIEchLQP+SURl21GLjXlafdDBFwkHo8CBAgQIAAgbYC+oe2Ymmul1Mad6sSIECAAIGSBfQPJaeXd+0G7nnno7oOBRykHWK6FQECBAgQ6ImA/qGMoOVURk6qJECAAAECOQnoH3JKo65aDNzrytNuhgg4SD0eBAgQIECAQFsB/UNbsTTXyymNu1UJECBAgEDJAvqHktPLu3YD97zzUV2HAg7SDjHdigABAgQI9ERA/1BG0HIqIydVEiBAgACBnAT0DzmlUVctBu515Wk3QwQcpB4PAgQIECBAoK2A/qGtWJrr5ZTG3aoECBAgQKBkAf1DyenlXbuBe975qK5DAQdph5huRYAAAQIEeiKgfygjaDmVkZMqCRAgQIBATgL6h5zSqKsWA/e68rSbIQIOUo8HAQIECBAg0FZA/9BWLM31ckrjblUCBAgQIFCygP6h5PTyrt3APe98VNehgIO0Q0y3IkCAAAECPRHQP5QRtJzKyEmVBAgQIEAgJwH9Q05p1FWLgXtdedrNEAEHqceDAAECBAgQaCugf2grluZ6OaVxtyoBAgQIEChZQP9Qcnp5127gnnc+qutQwEHaIaZbESBAgACBngjoH8oIWk5l5KRKAgQIECCQk4D+Iac06qrFwL2uPO1miICD1ONBgAABAl0L7Nq1OzZs2Brbtu2M+fPnxtKlh8ecObO7Xsb9EgroHxLit1haTi2wXEqAAAECBDIS2B27Y0vcFTtjZ8yNubEgFsbsmEw/rX/I6EGorBQD98oCtZ39CzhIPR0ECNQoYOCbLtXGfuXK9YNh+56/Zui+du2phu7pYul8Zf1D56RjuaGcxsLa25umHP70Ft3GCRDopUBz3l4VV8aO2LF3//NiXpwWyyYydNc/9PKxm8imDdwnwmyRHAQcpDmkMJ4aDBzH4+qu+QsY+KbN6Iorbo81a26aUsSqVcfF8uWL0hZn9c4E9A+dUY71RnIaK2+vbp56+NMrbJslQKBagen+g8vbYnPcGDdMcTg+Togj46ix++gfxk7c2wUM3Hsbff827iCtM3MDxzpztavpCRj4Ts9pXFdddNH1sW7dlim3X7ZsYaxeffy4lnXfCQvoHyYMPuJychoRzsumCKQe/oiEAAECpQu0+QeXm+K6uDPumLLlI2JRnBgnjZ1C/zB24t4uYODe2+j7t3EHaZ2ZGzimz9W/YZAuAwPfdPbNys6ftP6TWl3/MCnpma0jp5n5efX3BNoOf6b7KU7GBAgQ6ItAm39w2ebacfjpH8ah6p6NgIG756A3AuM8SA0c0z1GBo7p7JuV/RsGaf0NfNP6e/7T+k9q9XH2D5PaQx/WkVMfUp7MHtsMf9p8inMm1Rvqz0TPawkQmLRAm39wOalzdH8G+odJPx39Wc/AvT9Z936n4zpIDVzSPloGjvzTCqRd3fmT1n/PP3TauHHb4IdTmx9MXbJkvh9MTR9LpxWMq3/otEg3Czl5CLoSaDP8aTOcH7W+NvWMuobXESBAoEuBtmdjc85tjS2DH05tfjD18FgwkR9Mbfasf+gyeff6fgEDd89DbwTGdZAa+KZ9hAwc0/r7NwzS+hv4pvdXQf0C4+of6peb7A7lNFnv2leb7vCnzac4RzVrO7gadR2vI0CAQFcCJf2DQv1DV6m7zxMFDNw9E70RGNdBauCY/hFqhu4+YZomB//AKY27VQkQmJzAuPqHye2gHyvJqR8557bLSQzDJzHUz81VPQQIlC8w3X9wmXqn+ofUCdS7voF7vdna2RMExnWQGjh61Pos4N8w6HP69k6gHwLj6h/6oTe5XcppctZW+p7AJD7FOYmhvkwJECDQVwH9Q1+TH/++DdzHb2yFTATGdZAaOGYSsDKSCfg3DJLRW5gAgQkIjKt/mEDpvVpCTr2KO6vNjvtTnJMY6mcFqhgCBAhMUED/MEHsni1l4N6zwPu83XEepPff/1B84APXx+bN34gXveg5ce65J8TBBz+9z9z2ToAAAQIEqhAYZ/9QBVAmm5BTJkEoYywCXQ71m3ttibtiZ+yMuTE3FsTCif044Vhw3JQAAQIzENA/zADPS4cKGLh7QHojMK6D1Cfce/MI2SgBAgQI9FBgXP1DDynHumU5jZXXzSsR8Gn5SoK0DQIEOhPQP3RG6UZPEDBw90j0RmBcB6nvcO/NI2SjBAgQINBDgXH1Dz2kHOuW5TRWXjevRMD3wVcSpG0QINCZgP6hM0o3MnD3DPRVYFwH6UUXXR/r1m2Zwrps2cJYvfr4vnLbNwECBAgQqEJgXP1DFTgZbUJOGYWhlGwFNsV1cWfcMaW+I2JRnBgnZVu3wggQIDAuAf3DuGTd1yfcPQO9ERjXQeoT7r15hGyUAAECBHooMK7+oYeUY92ynMbK6+aVCPiEeyVB2gYBAp0J6B86o3SjJwgYuHskeiMwroPUd7j35hGyUQIECBDoocC4+oceUo51y3IaK6+bVyLgO9wrCdI2CBDoTED/0BmlGxm4ewb6KjDOg7QZum/cuC22bdsZ8+fPjSVL5secObP7Sj3xfTf+GzZs3eu/dOnh/CeeggUJECBQp8A4+4c6xdLsSk5p3K1ankAzdN8aW2JH7Ih5MS8OjwUxO7xvKS9JFU9SoPn/N1virtgZO2NuzI0FsdD/v5lkAGNcS/8wRtye39on3Hv+APRp+w7SOtP2bxjUmatdESBAYDoCk/gHrvqH6SSR/ho5pc9ABQQIEMhRYKbDcv9mSI6pdleT/qE7S3f6QQEDd0/ExAWuvfbauPDCC+OWW26J++67Ly6//PI4/fTTh9ZxzTXXxOrVq2Pz5s1xyCGHxNvf/vY455xzWtXuIG3FVczFvkO/mKgUSoAAgU7Vzp7CAAAe0ElEQVQFJvUPXPUP7WLT57XzcjUBAgQIjE+gi2G53z4YXz453Fmfl0MKddZg4F5nrlnvav369bFp06Y45phj4owzznjSgfvdd98dRx99dJx11llx9tlnD177pje9KS677LLB66f75yCdrlRZ11100fWxbt2WKUUvW7YwVq8+vqzNqJYAgeIEJvEJ6+JQJlTwpP6Bq/6hXaD6vHZeriZAgACBfQvM9JPpzV27GJZviuvizrhjSpFHxKI4MU4SX+EC+rzCA8y4fAP3jMPpQ2mzZs160oH7ueeeG1/4whfi9ttv30vSfLr9y1/+ctxwww3TZnKQTpuqqAsnNXApCkWxBAhMRGBSn7CeyGYKXGRS/8BV/zD6w6HPG93OKwkQINBngS4+md74dTEs72Jo3+csc9+7Pi/3hMqtz8C93OyqqHw6b8R+5md+JhYvXhxr1qzZu+fma2he+cpXxkMPPRSzZ+/7R34eeeSRaP6z5685SA899NB44IEH4qCDDqrCzyYiDLw8BQQIpBLwD/xSyf/PupPy90Zs9Jz1eaPbeSUBAgT6LNDVkLuL+3Q1/O9znjnvXZ+Xczpl12bgXnZ+xVc/nTdiCxcujBUrVsQ73/nOvfu9/vrr48QTT4x77703nve85+3T4fzzz48LLrhgyv/OwL34x2bKBpqh+8aN22Lbtp0xf/7cWLJkfsyZs+9/EFPf7u2IAIFUApP6hHWq/eW+7qT+gas3YqM/Cfq80e28kgABAn0W6OKT6Y1fV8Py5j5bY0vsiB0xL+bF4bEgZof3mzU8o/q8GlLMcw8G7nnm0puqpvtG7LWvfW2cd955e12a73E/6aSTBj+6+tznPnefXj7h3pvHyEYJ9FrAd4ini7/5hPWHPnRj7NjxcDz88KNx4IEHxLx5Bw5+P2L58kXpCuvRyvff/1B84APXx+bN34gXveg5ce65J8TBBz+9UwFvxEbn1OeNbueVBAhMFWj7nd5tr29rPu77t62npuu7+GT6Hg/D8pqejO73os/r3tQd/0fAwN2TkFRgOm/ERv1KmSduzEGaNGqLEyAwBoFJfcJ3DKVXcctm2HvssX8aX/vat+M733k8nvrUWfH85z8zbr75rM6HvlWAdbyJST3/+ofRg9PnjW7nlQQI/KBA208qt72+rXdz/yvj83Fv3BO749GYHQfEIfGjsSyWj/zJZwP876Uw7vza5u36egX0efVmm3pnBu6pE+j5+tN5I9b8aOqVV14Zt912216tN77xjXHrrbf60dSePz+2T6DvApP6Duu+O+9v/3/5l1+O17/+C/F//+9jey952tOeEn/2Z8vj13/9xdjGLDCp598bsdGD1OeNbueVBAj8oEDbTzy3vb6t97/FV2J9XBWPxqN7X3pAHBCnxmnxE9G+BzBgnpqAT6a3fSpdP4qAPm8UNa+ZjoCB+3SUXNOpwIMPPhhbt24d3LP5MdSLL744TjnllJg3b14cdthhg6+Oueeee+LSSy8dXHP33XfH0UcfHWeffXacddZZgyH7OeecE5dddlmcccYZ067NQTptKhcSIFCIgO8QTxvUz/7sJfFP//TVKUX8zM8cFv/wD2emLa4Hq0/q+dc/tHuY9HntvFxNgMD0BNp+p3fb66dXxfeu+lx8Nu6I26e8bFEcFa+I6b9H3XODcf8Dgrb7cz2Bvgjo8/qS9OT3aeA+efPer3j11VcPBuxP/DvzzDPjkksuGfxA6vbt26O5bs/fNddcE29961tj8+bNccghh0Tzqfdm6N7mz0HaRsu1BAiUIDCpT/iWYJGixp/6qU/Erbd+Y8rSixf/SNx002+mKKlXa07q+dc/tHus9HntvFxNgMD0BNoOpNteP70qvndV1wP3cf8Dgrb7cz2Bvgjo8/qS9OT3aeA+eXMrJhJwkCaCtywBAmMTmNR3WI9tA4Xf+Hd+Z2N8+MM3xeOPf28js2ZFvOUtx8WFF/584bvLv/xJPf/6h/yfhaZCOZWRkyoJjCrQ9itX2l7ftq6vxJdjffxtfOf7vlLmqXFA/MKIXykz7n9A0HZ/rifQFwH9Q1+Snvw+Ddwnb27FRAIO0kTwliVAYKwCzdBx48ZtsW3bzpg/f24sWTI/5syZPdY13fx/BPb8aOp99z2490dTn/e8/8ePpk7wAZnE869/mGCgM1hKTjPA81IChQi0/U7vtte3YWju/YW4Iu6Ne6P5r2fH7DgkDolfitNH+tHUcf8DgjZ7cy2BPgnoH/qU9mT3auA+WW+rJRRwkCbEtzQBAgQqFWiG7v/7f98Qmzd/I170oufE299+fBx88NMr3W0/t6V/KCN3OZWRkyoJ1CTw/7d357GSFPUDwGsxATkWkFsEEblvdUFBIkdU5FSu4AJyKhI0hFMNIJcoeBCBBORMuAS5FAmrckdOTQgsCQiagK7uct+iIPcv30rm/YZh3nvTb/q96e75dMIfZPtVV33qOzU1366uLjuhX3Z5TbLWFgKTJWD+MFmyypVwFwNDI2AgHZqu1lACBAgQIFCagPlDaZSTWpB+mlRehRMgQIAAgUYKmD80slsr0SgJ90p0g0pMhYCBdCqUXYMAAQIECDRLwPyhHv2pn+rRT2pJgAABAgSqJGD+UKXeaFZdJNyb1Z9aM4aAgVR4ECBAgAABAkUFzB+Kig3mfP00GHdXJUCAAAECdRYwf6hz71W77hLu1e4ftStRwEBaIqaiCBAgQIDAkAiYP9Sjo/VTPfpJLQkQIECAQJUEzB+q1BvNqouEe7P6U2vGEDCQCg8CBAgQIECgqID5Q1GxwZyvnwbj7qoECBAgQKDOAuYPde69atddwr3a/aN2JQoYSEvEVBQBAgQIEBgSAfOHenS0fqpHP6klAQIECBCokoD5Q5V6o1l1kXBvVn9qzRgCBlLhQYAAAQIECBQVMH8oKjaY8/XTYNxdlQABAgQI1FnA/KHOvVftuku4V7t/1K5EAQNpiZiKIkCAAAECQyJg/lCPjtZP9egntSRAgAABAlUSMH+oUm80qy4S7s3qT60ZQ8BAKjwIECBAgACBogLmD0XFBnO+fhqMu6sSIECAAIE6C5g/1Ln3ql13Cfdq94/alShgIC0RU1EECBAgQGBIBMwf6tHR+qke/aSWBAgQIECgSgLmD1XqjWbVRcK9Wf2pNWMIGEiFBwECBAgQIFBUwPyhqNhgztdPg3F3VQIECBAgUGcB84c691616y7hXu3+UbsSBQykJWIqigABAgQIDImA+UM9Olo/1aOf1JIAAQIECFRJwPyhSr3RrLpIuDerP7VmDAEDqfAgQIAAAQIEigqYPxQVG8z5+mkw7q5KgAABAgTqLGD+UOfeq3bdJdyr3T9qV6KAgbRETEURIECAAIEhETB/qEdH66d69JNaEiBAgACBKgmYP1SpN5pVFwn3ZvWn1owhYCAVHgQIECBAgEBRAfOHomKDOV8/DcbdVQkQIECAQJ0FzB/q3HvVrruEe7X7R+1KFDCQloipKAIECBAgMCQC5g/16Gj9VI9+UksCBAgQIFAlAfOHKvVGs+oi4d6s/tSaMQQMpMKDAAECBAgQKCpg/lBUbDDn66fBuLsqAQIECBCos4D5Q517r9p1l3Cvdv+oXYkCL7/8clp88cXT3Llz06KLLlpiyYoiQIAAAQIEmioQP8RWXHHF9NJLL6XFFlusqc2sfbvM82rfhRpAgAABAgSmXMA8b8rJh+aCEu5D09UaOm/evPyD2UGAAAECBAgQKCoQN+xXWGGFon/m/CkSMM+bImiXIUCAAAECDRQwz2tgpw64SRLuA+4Al586gXfeeSc98cQTafr06WnatGmlX7h1Z9QK+tJpeyqQf09Mk3YS/0mj7alg/j0xTdpJ/CeNtqeCJ9v/3XffTa+88kpafvnl03zzzddTnZw09QL9zPMmO4amXqOcK3IZ3ZFNdxsuXIqMPuJFvFQhXszzivSCc4sISLgX0XIugTEE7P012PDgz3+wAoO9uvjnP1iBwV5d/A/WvwlXF0OjJ31iG6XYrsd2jO81EjNipsjYJ17Ei3gpIiBe+tdSQhUEJNyr0Avq0AgBE6nBdiN//oMVGOzVxT//wQoM9urif7D+Tbi6GJLcKBrHYkbMFIkZ8SJexEsRAfHSv5YSqiAg4V6FXlCHRgiYSA22G/nzH6zAYK8u/vkPVmCwVxf/g/VvwtXFkORG0TgWM2KmSMyIF/EiXooIiJf+tZRQBQEJ9yr0gjo0QuD1119Pp5xySjrqqKPSAgss0Ig21akR/AfbW/z5D1ZgsFcX//wHK+Dq/Qr4DHcX5DJ6ZLERM0XGHfEiXsRLEQHx0r+WEqogIOFehV5QBwIECBAgQIAAAQIECBAgQIAAAQIECBCovYCEe+27UAMIECBAgAABAgQIECBAgAABAgQIECBAoAoCEu5V6AV1IECAAAECBAgQIECAAAECBAgQIECAAIHaC0i4174LNYAAAQIECBAgQIAAAQIECBAgQIAAAQIEqiAg4V6FXlAHAgQIECBAgAABAgQIECBAgAABAgQIEKi9gIR77btQA6ZS4I477kg/+9nP0n333ZeefPLJdO2116Ydd9xxpAr77rtvuvjii99Tpc985jPpz3/+81RWs7HXOuWUU9JvfvOb9Ne//jUtuOCC6bOf/Wz6yU9+ktZYY42RNr/++uvpyCOPTL/61a/Sa6+9lj7/+c+nX/ziF2mFFVZorMtUNawX/y222CLdfvvt76nSV7/61XTFFVdMVTUbe52zzz47xX9z5szJbVxnnXXScccdl7bZZpv8/2J/crt+PH+xP7n+naXHeHT00UenQw45JJ1++uk+A1PLX4urxXd/zNlivhbjZcTJ5z73uVHr/utf/zode+yx6bHHHkurrLJK+tGPfpR22mmnkfPffffddOKJJ6bzzjsvvfjiiynmd2eddVYuu05H2S5NmfsWcfnLX/6Sv3/j98A///nPdNppp6VDDz30fWFQpMyqxlCRNvTicsIJJ+TPUfux7LLLpqeeeqqqBF3rVcTl/PPPT5dcckl66KGHclkzZsxIJ598cvr0pz891ONLLy5NGV+io4vETPzejRh59NFH05tvvplWW221dMQRR6S99tprqGOmF5cmxUytBkWV7Sog4S4wCBQQ+MMf/pDuvvvu9KlPfSrtsssuXRPuTz/9dLrwwgtHSp1//vnTEkssUeAqTh1NYOutt04zZ85MG220UXrrrbfSMccckx588MH08MMPp4UXXjj/2UEHHZSuv/76dNFFF6Ull1wyT05eeOGF/KPoAx/4ANw+BHrxj6Tj6quvnn7wgx+MXClujiy22GJ9XNmfhkDEdcTwqquumkHi5l4kk2bPnp0TPmJ/cuNkPH+xP7n+7aXfe++9abfddkuLLrpo2nLLLUcS7j4DU9cHVb/SlVdemRMTkeDYdNNN07nnnpsuuOCCPF/46Ec/+r7q/+lPf8rJ+JNOOikn2WNBRSRU77rrrpxYjyNu8EcSPuYX8T33wx/+MMVCjL/97W9p+vTpVSfJ9ZsMl0hu1H3uW9QlxqCrrroqJ04PO+yw9L3vfe99CfeiZVYxgIq2oReXSLhfc8016ZZbbhlpcsxtll566SoSdK1TUZc999wzj0OxUOiDH/xg+ulPf5oXEMUNio985CNDO7704tKE8WUiY+8f//jHfGN3zTXXTJFLmDVrVv5N+7vf/S596UtfGtqY6cWlKTFTmwFRRccUkHAXIAQmKDBt2rSuCfeXXnop/fa3v51gqf6siMCzzz6blllmmbyierPNNksvv/xynrBfeumlKVZVx/HEE0+kFVdcMf3+978fmaAUuYZzRxfo9I8zI+n4iU98YiQBxm9yBeJmXiTdd911V7E/udRdS2/5f/3rXxf7U+T/n//8J9/0jkRqJDxb443xf4o6oCaXiSR5xEk8mdI61lprrfxUYjwd0XnEnOHf//53ioUVrSNuMn/oQx/KT8zF6vbll18+J1UjuRpHPFUUK3MjEX/ggQfWQqZsl2h0JDfqPvct6tLe2R/72MdyXHSucO+nzKoEUz9tGM0lEu7xO+mBBx6oSjML16Mfl7jY22+/nceWM888M+29995DO750wne6NGV8iXb0GzNRRnynbbfddvnG8LB+J3X7sLa7NClmCg9M/qCSAhLulewWlaqDwGgJ95hExp3oxRdfPG2++eZ5NVQkhR3lC8RjdvGIXaxyX3fdddNtt92Wt5CJFe0xkW0dG2ywQf6R3fkIa/k1Gq4SO/2j9ZFwjxU7MRGMRERsd3L88cfXZvVfXXowfpRcffXVaZ999skr3ONRbLE/db3X6b/22muL/Snij5iPGx2xhUP7DT7j/xR1QA0u88Ybb6SFFlooj5HtW8LE9kOR5Ovc9iyaFKveY6Vy/Nc6IsZiG5rYMuTvf/973mbm/vvvT5/85CdHzvnKV76S53ud2wlWkWkyXFrJjTrPfSfi0t6/3RLL/ZZZhfjptw1jJdxjoUA8+bjAAgvkRGRsnfHxj3+8Cs0etw79usQFXnnllfzbMMao7bfffmjHl07sTpcmjC/Rhn5jJn5TxRzny1/+cr5Z9cUvflHMpJR/a3a6NCVmxh2InFAbAQn32nSVilZNoFvCPR4xXGSRRdJKK62U/vGPf+S9QGPrk9jOJCaVjvIE4ks2fujG43Z33nlnLvjyyy9P++23X1511n5stdVWaeWVV86PlDvKEejmHyXHfoxhvdxyy+W9Ko866qi8BcrNN99czoWHvJS4ubTJJpuk//3vf3msiZjfdtttxf4UxcVo/mJ/ajog3gURN7Fj24J4LL894W78n5o+qMNV4sm22KYhtgCMLRxaRyT1IjEeW8B0HrFQIraK2WOPPUb+qT2m7rnnnrwlxOOPP55XureOb37zmzkhf+ONN1aeZjJcotF1n/tOxKW9s7sllvstswrB1G8bRku4x1Mkr776at6WKbYiiieV4t1MsVgjtoKs+tGvS7Tv29/+dh4zYp4c32XDOr509nWnSxPGl2jDRGMmntyL77L4XRvbLsWTffvvv39mG+aYGculKTFT9XFQ/XoXkHDv3cqZBN4j0C3h3kkUL+qK5HskCXbeeWeCJQrEpCz2sYv9VVsvRB0t4RIrAWJl2jnnnFNiDYa7qG7+3UTiZtOGG26YbzrFI3+O/gRilcy//vWv/Ph+vOAv9iSO1ZqxarPbzSax359351+P5h8r3DsPsV+u/dy5c/NYctNNN6V4aimOXhLuPgPl9kMdSmslNyIhETcoW0fcrIkt5yK513lEwj2S8bvvvvvIP1122WUptouKG5yt5EaU/eEPf3jknAMOOCBFbN5www2Vp5kMl26NrtvcdyIu7e0eK+FeJAarFkCT4dKtjf/973/zHP273/1uOvzww6vG8L769OsS+7f/+Mc/TrEX9frrr5/LH9bxpR23m0sTxpdow0Rj5p133skr2WMrvVtvvTVvJRMr3GPuM8wxM5ZLU2Km8gOhCvYsIOHeM5UTCbxXoJeEe/xFbHnyjW98Y2TPT479Cxx88MF5whEvK4vV1K3DlgL92/ZSwmj+3f42VsLH0x3t++r3cg3n9CbwhS98If9Qjf2HbSnTm1mZZ7X8uz09I/bLlE55zI/tQdpffh1b+8R38XzzzZdXC0Z/2FKsXPc6ljaRx/dtKTOxrXZGi486zX0nEi/t7balTPcoGG2Fe7ez48ZoPA3Z/s6Fqo49/cTLqaeemlf0xwtj4wZy6xjWLata7R/NpQnjS7Shn5hpN4h8QtzgjfnOsMfMaC5NiZmqjn/qVVxAwr24mb8gkAV6Sbg///zz+VGw8847L78Ux9GfQCSwItl77bXX5pUh8YOu/Wi9NO+Xv/xl2m233fI/xUqrWAHvpan92cdfj+ff7QrxuOx666038mLb/muhhHaBSLLHS4HPOOOM/NJUsT+18dHyj60oOg+xX25fxN6usXVH+xFPday55pr5hnZ8DnwGyjWvc2mxL/SMGTPyI/itI55Eia3oRntpasRYzBVaR7yDJPZnb39pauzxHitx44gkSuzDXLeXppbp0i1G6jj3LRov7e0e66WpRayr+HmbDJfOdsZ2GbFwILZnOu6446rI8L46TcQl9q2PZHskSzfeeOP3lNl6AeawjS+BMJZLU8aXaMdEYqaz/fHE1WOPPZZ/Aw9zzIzl0qSYqcVgqJLjCki4j0vkBAL/LxCPdMWLIuOIl2b9/Oc/T1tuuWV+gVv8d8IJJ6RddtklP248Z86cdPTRR+ftHx555BEvjSwhkL71rW/lvaqvu+66tMYaa4yUGC9eWnDBBfP/H3TQQWnWrFl5L9bokyOPPDLFj7/Y3qF9ZWQJ1Rm6Isbzj0lgPIIfe4ovtdRS6eGHH05HHHFE7pvYc5l/fyET40kkgCKxGImh2KoqHkuOrQxidZjY7893vL8eyz9e9ib2xxMs/9/bt5Qx/pfvW+cSY1/xvfbaK28lF9vKxMKHeMdI7BMdW/3FIohYENFKvsfj+Ztttll+R0Ak5WOe8f3vfz9vWxeJkjgisR7nX3jhhfmGf+wJH4mP2BN++vTpteAq2yXmxU2Y+xZ1iZstMceJI+Y8e+65Z/4v3q0SK7XjGK/MOgTMeG3o/Bz14hLz8h122CG/qPiZZ57JSejYGi/ekRKfzTocRV1iu5R4r1f8hol3QbSOiJf4b1jHl/FcmjK+9DIedH6W4rsmnoKIm1HxuYqbwbG4IJ4CiZXuwxoz47k0KWbqMBaq4/gCEu7jGzmDwIhA/LCKBHvnsc8+++QvwB133DHNnj07768cSfc4N/ZbiwSZo3+BeKqg2xE/fvfdd9/8T7HX6ne+8508qX3ttdfyNhuxwk0fTL5/POb4ta99Lb8EKiY8Yb7ddtul448/Pt/8cPQnECtbYg/HeGojbjLF3p8x+Y5ku9jvz7aXvx7LX+z3Ilj+OZ0Jd+N/+cZ1LjG++yOhE2Pmuuuum0477bScVI8jYidWJrc/nXLNNdfkJHvrUf1Ivre/fydWFJ544on5BezxwvZIxJ911lm57DodZbrEPKspc98iLrGopn1Lw1b/b7755vkmTOsYq8y6xEzZLjNnzsxbQj733HP5qaRY7R2/lbq9C6XKRkVcYqzpfEIr2hbz47hhFccwji/juTRpfIk+LhIz8V0UN3bmzZuXFy7F03yHHHJI3kKydQxjzIzn0rSYqfIYqG69CUi49+bkLAIECBAgQIAAAQIECBAgQIAAAQIECBAgMKaAhLsAIUCAAAECBAgQIECAAAECBAgQIECAAAECJQhIuJeAqAgCBAgQIECAAAECBAgQIECAAAECBAgQICDhLgYIECBAgAABAgQIECBAgAABAgQIECBAgEAJAhLuJSAqggABAgQIECBAgAABAgQIECBAgAABAgQISLiLAQIECBAgQIAAAQIECBAgQIAAAQIECBAgUIKAhHsJiIogQIAAAQIECBAgQIAAAQIECBAgQIAAAQIS7mKAAAECBAgQIECAAAECBAgQIECAAAECBAiUICDhXgKiIggQIECAAAECBAgQIECAAAECBAgQIECAgIS7GCBAgAABAgQIECBAgAABAgQIECBAgAABAiUISLiXgKgIAgQIECBAgAABAgQIECBAgAABAgQIECAg4S4GCBAgQIAAAQIECBAgQIAAAQIECBAgQIBACQIS7iUgKoIAAQIECBAgQIAAAQIECBAgQIAAAQIECEi4iwECBAgQIECAAAECBAgQIECAAAECBAgQIFCCgIR7CYiKIECAAAECBAgQIECAAAECBAgQIECAAAECEu5igAABAgQIECBAgAABAgQIECBAgAABAgQIlCAg4V4CoiIIECBAgAABAgQIECBAgAABAgQIECBAgICEuxggQIAAAQIECBAgQIAAAQIECBAgQIAAAQIlCEi4l4CoCAIECBAgQIAAAQIECBAgQIAAAQIECBAgIOEuBggQIECAAAECBAgQIECAAAECBAgQIECAQAkCEu4lICqCAAECBAgQIECAAAECBAgQIECAAAECBAhIuIsBAgQIECBAgAABAgQIECBAgAABAgQIECBQgoCEewmIiiBAgAABAgQIECBAgAABAgQIECBAgAABAhLuYoAAAQIECBAgQIAAAQIECBAgQIAAAQIECJQgIOFeAqIiCBAgQIAAAQIECBAgQIAAAQIECBAgQICAhLsYIECAAAECBAgQIECAAAECBAgQIECAAAECJQhIuJeAqAgCBAgQIECAAAECBAgQIECAAAECBAgQICDhLgYIECBAgAABAgQIECBAgAABAgQIECBAgEAJAhLuJSAqggABAgQIECBAgAABAgQIECBAgAABAgQISLiLAQIECBAgQIAAAQIECBAgQIAAAQIECBAgUIKAhHsJiIogQIAAAQIECBAgQIAAAQIECBAgQIAAAQIS7mKAAAECBAgQIECAAAECBAgQIECAAAECBAiUICDhXgKiIggQIECAAAECBAgQIECAAAECBAgQIECAgIS7GCBAgAABAgQIECBAgAABAgQIECBAgAABAiUI/B8uIcd76T9AFgAAAABJRU5ErkJggg==\" width=\"1500\">"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {}
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "It can be seen from the above plot that\n",
    "* The actual impact of rank is in line with the intuition - the smaller the value the better the result.\n",
    "* It is interesting to see that the optimal value of reg is around 0.1 to 0.15. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Get the best model."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "source": [
    "als = ALS(\n",
    "    rank=best[\"rank\"],\n",
    "    regParam=best[\"reg\"],\n",
    "    maxIter=15,\n",
    "    implicitPrefs=False,\n",
    "    alpha=0.1,\n",
    "    coldStartStrategy='drop',\n",
    "    nonnegative=False,\n",
    "    seed=42,\n",
    "    **HEADER_ALS\n",
    ")\n",
    "    \n",
    "model_best_hyperopt = als.fit(train)"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Tuning prameters against other metrics can be simply done by modifying the `objective` function. The following shows an objective function of how to tune \"precision@k\". Since `fmin` in `hyperopt` only supports minimization while the actual objective of the loss is to maximize \"precision@k\", `-precision` instead of `precision` is used in the returned value of the `objective` function."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "source": [
    "# Customize an objective function\n",
    "def objective_precision(params):\n",
    "    with Timer() as time_run_start:\n",
    "\n",
    "        rank = params['rank']\n",
    "        reg = params['reg']\n",
    "        train = params['train'] \n",
    "        valid = params['valid'] \n",
    "        col_user = params['col_user'] \n",
    "        col_item = params['col_item']\n",
    "        col_rating = params['col_rating'] \n",
    "        col_prediction = params['col_prediction'] \n",
    "        k = params['k']\n",
    "        relevancy_method = params['relevancy_method']\n",
    "\n",
    "        header = {\n",
    "            \"userCol\": col_user,\n",
    "            \"itemCol\": col_item,\n",
    "            \"ratingCol\": col_rating,\n",
    "        }\n",
    "\n",
    "        als = ALS(\n",
    "            rank=rank,\n",
    "            maxIter=15,\n",
    "            implicitPrefs=False,\n",
    "            alpha=0.1,\n",
    "            regParam=reg,\n",
    "            coldStartStrategy='drop',\n",
    "            nonnegative=False,\n",
    "            seed=42,\n",
    "            **header\n",
    "        )\n",
    "\n",
    "        model = als.fit(train)\n",
    "    \n",
    "        users = train.select(col_user).distinct()\n",
    "        items = train.select(col_item).distinct()\n",
    "        user_item = users.crossJoin(items)\n",
    "        dfs_pred = model.transform(user_item)\n",
    "\n",
    "        # Remove seen items.\n",
    "        dfs_pred_exclude_train = dfs_pred.alias(\"pred\").join(\n",
    "            train.alias(\"train\"),\n",
    "            (dfs_pred[col_user] == train[col_user]) & (dfs_pred[col_item] == train[col_item]),\n",
    "            how='outer'\n",
    "        )\n",
    "\n",
    "        top_all = dfs_pred_exclude_train.filter(dfs_pred_exclude_train[\"train.Rating\"].isNull()) \\\n",
    "            .select('pred.' + col_user, 'pred.' + col_item, 'pred.' + \"prediction\")\n",
    "\n",
    "        top_all.cache().count()\n",
    "\n",
    "        rank_eval = SparkRankingEvaluation(\n",
    "            valid, \n",
    "            top_all, \n",
    "            k=k, \n",
    "            col_user=col_user, \n",
    "            col_item=col_item, \n",
    "            col_rating=\"rating\", \n",
    "            col_prediction=\"prediction\", \n",
    "            relevancy_method=relevancy_method\n",
    "        )\n",
    "\n",
    "        precision = rank_eval.precision_at_k()\n",
    "\n",
    "    # Return the objective function result.\n",
    "    return {\n",
    "        'loss': -precision,\n",
    "        'status': STATUS_OK,\n",
    "        'eval_time': time_run_start.interval\n",
    "    }"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "### 4.2 Hyperparameter tuning with `hyperopt` sampling methods"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Though `hyperopt` works well in a single node machine, its features (e.g., `Trials` module) do not support Spark environment, which makes it hard to perform the tuning tasks in a distributed/parallel manner. It is useful to use `hyperopt` for sampling parameter values from the pre-defined sampling space, and then parallelize the model training onto Spark cluster with the sampled parameter combinations.\n",
    "\n",
    "The downside of this method is that the intelligent searching algorithm (i.e., TPE) of `hyperopt` cannot be used. The approach introduced here is therefore equivalent to random search."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "source": [
    "with Timer() as time_sample:\n",
    "    # Sample the parameters used for model building from the pre-defined space. \n",
    "    sample_params = [sample(space) for x in range(NUMBER_ITERATIONS)]\n",
    "    \n",
    "    # The following runs model building on the sampled parameter values with the pre-defined objective function.\n",
    "    results_map = list(map(lambda x: objective(x), sample_params))\n"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "source": [
    "results_map"
   ],
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": [
       "[{'eval_time': 9.468051671981812, 'loss': 1.027085217204854, 'status': 'ok'},\n",
       " {'eval_time': 8.947720766067505, 'loss': 1.017764730532703, 'status': 'ok'},\n",
       " {'eval_time': 9.599841117858887, 'loss': 1.2995721337596726, 'status': 'ok'},\n",
       " {'eval_time': 8.645057439804077, 'loss': 0.998792404289471, 'status': 'ok'},\n",
       " {'eval_time': 9.246882200241089, 'loss': 1.159882028048988, 'status': 'ok'},\n",
       " {'eval_time': 9.159096479415894, 'loss': 1.3707996773718212, 'status': 'ok'},\n",
       " {'eval_time': 9.555922508239746, 'loss': 1.4255606225971154, 'status': 'ok'},\n",
       " {'eval_time': 9.555083751678467, 'loss': 0.9974149852593205, 'status': 'ok'},\n",
       " {'eval_time': 9.14759874343872, 'loss': 1.0233910316377184, 'status': 'ok'},\n",
       " {'eval_time': 9.288854122161865, 'loss': 1.1079683856151636, 'status': 'ok'},\n",
       " {'eval_time': 9.035703420639038, 'loss': 1.4973257401273627, 'status': 'ok'},\n",
       " {'eval_time': 9.43152904510498, 'loss': 1.3660611992616116, 'status': 'ok'},\n",
       " {'eval_time': 9.249063491821289, 'loss': 1.3212144812805433, 'status': 'ok'},\n",
       " {'eval_time': 9.086166143417358, 'loss': 1.1874756727128037, 'status': 'ok'},\n",
       " {'eval_time': 8.880879878997803, 'loss': 1.017539254275622, 'status': 'ok'},\n",
       " {'eval_time': 9.382610559463501, 'loss': 1.0726440276761462, 'status': 'ok'},\n",
       " {'eval_time': 9.176624774932861, 'loss': 1.4048578426830673, 'status': 'ok'},\n",
       " {'eval_time': 9.292484045028687, 'loss': 1.4308198992737957, 'status': 'ok'},\n",
       " {'eval_time': 8.882294178009033, 'loss': 1.2712116101690774, 'status': 'ok'},\n",
       " {'eval_time': 9.89118218421936, 'loss': 1.0887572322216503, 'status': 'ok'},\n",
       " {'eval_time': 9.102333545684814, 'loss': 1.449849882363006, 'status': 'ok'},\n",
       " {'eval_time': 9.27437686920166, 'loss': 1.0465564408902348, 'status': 'ok'},\n",
       " {'eval_time': 9.215168714523315, 'loss': 1.248519446580608, 'status': 'ok'},\n",
       " {'eval_time': 9.225409269332886, 'loss': 1.0265498645574211, 'status': 'ok'},\n",
       " {'eval_time': 9.08506464958191, 'loss': 1.254533287299843, 'status': 'ok'}]"
      ]
     },
     "metadata": {},
     "execution_count": 30
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Get the best model."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "source": [
    "loss_metrics = np.array([x['loss'] for x in results_map])\n",
    "best_loss = np.where(loss_metrics == min(loss_metrics))"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "source": [
    "best_param = sample_params[best_loss[0].item()]"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "source": [
    "als = ALS(\n",
    "    rank=best_param[\"rank\"],\n",
    "    regParam=best_param[\"reg\"],\n",
    "    maxIter=15,\n",
    "    implicitPrefs=False,\n",
    "    alpha=0.1,\n",
    "    coldStartStrategy='drop',\n",
    "    nonnegative=False,\n",
    "    seed=42,\n",
    "    **HEADER_ALS\n",
    ")\n",
    "    \n",
    "model_best_sample = als.fit(train)"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 5 Evaluation on testing data"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "The optimal parameters can then be used for building a recommender, which is then evaluated on the testing data.\n",
    "\n",
    "The following codes generate the evaluation results by using the testing dataset with the optimal model selected against the pre-defined loss. Without loss of generity, in this case, the optimal model that performs the best w.r.t regression loss (i.e., the RMSE metric) is used. One can simply use other metrics like precision@k, as illustrated in the above sections, to evaluate the optimal model on the testing dataset."
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "source": [
    "# Get prediction results with the optimal modesl from different approaches.\n",
    "prediction_spark = model_best_spark.transform(test)\n",
    "prediction_hyperopt = model_best_hyperopt.transform(test)\n",
    "prediction_sample = model_best_sample.transform(test)\n",
    "\n",
    "predictions = [prediction_spark, prediction_hyperopt, prediction_sample]\n",
    "elapsed = [time_spark.interval, time_hyperopt.interval, time_sample.interval]\n",
    "\n",
    "approaches = ['spark', 'hyperopt', 'sample']\n",
    "test_evaluations = pd.DataFrame()\n",
    "for ind, approach in enumerate(approaches):    \n",
    "    rating_eval = SparkRatingEvaluation(\n",
    "        test, \n",
    "        predictions[ind],\n",
    "        **HEADER\n",
    "    )\n",
    "    \n",
    "    result = pd.DataFrame({\n",
    "        'Approach': approach,\n",
    "        'RMSE': rating_eval.rmse(),\n",
    "        'MAE': rating_eval.mae(),\n",
    "        'Explained variance': rating_eval.exp_var(),\n",
    "        'R squared': rating_eval.rsquared(),\n",
    "        'Elapsed': elapsed[ind]\n",
    "    }, index=[0])\n",
    "    \n",
    "    test_evaluations = test_evaluations.append(result)"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "source": [
    "test_evaluations"
   ],
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Approach</th>\n",
       "      <th>Elapsed</th>\n",
       "      <th>Explained variance</th>\n",
       "      <th>MAE</th>\n",
       "      <th>R squared</th>\n",
       "      <th>RMSE</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>spark</td>\n",
       "      <td>133.150098</td>\n",
       "      <td>0.289853</td>\n",
       "      <td>0.776693</td>\n",
       "      <td>0.252805</td>\n",
       "      <td>0.976022</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>hyperopt</td>\n",
       "      <td>235.779974</td>\n",
       "      <td>0.299736</td>\n",
       "      <td>0.790172</td>\n",
       "      <td>0.240981</td>\n",
       "      <td>0.983563</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>sample</td>\n",
       "      <td>230.902271</td>\n",
       "      <td>0.287638</td>\n",
       "      <td>0.791199</td>\n",
       "      <td>0.232688</td>\n",
       "      <td>0.988922</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Approach     Elapsed  Explained variance       MAE  R squared      RMSE\n",
       "0     spark  133.150098            0.289853  0.776693   0.252805  0.976022\n",
       "0  hyperopt  235.779974            0.299736  0.790172   0.240981  0.983563\n",
       "0    sample  230.902271            0.287638  0.791199   0.232688  0.988922"
      ]
     },
     "metadata": {},
     "execution_count": 35
    }
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "From the results, it can be seen that, *with the same number of iterations*, Spark native construct based approach takes the least amount of time, even if there is no parallel computing. This is simply because Spark native constructs leverage the underlying Java codes for running the actual analytics with high performance efficiency. Interestingly, the run time for `hyperopt` with TPE algorithm and random search methods are almost the same. Possible reasons for this are that, the TPE algorithm searches optimal parameters intelligently but runs the tuning iterations sequentially. Also, the advantage of TPE may become obvious when there is a higher dimensionality of hyperparameters. \n",
    "\n",
    "The three approaches use the same RMSE loss. In this measure, the native Spark construct performs the best. The `hyperopt` based approach performs the second best, but the advantage is very subtle. It should be noted that these differences may be owing to many factors like characteristics of datasets, dimensionality of hyperparameter space, sampling size in the searching, etc. Note the differences in the RMSE metrics may also come from the randomness of the intermediate steps in parameter tuning process. In practice, multiple runs are required for generating statistically robust comparison results. We have tried 5 times for running the same comparison codes above. The results aligned well with each other in terms of objective metric values and elapsed time. "
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "# Conclusions"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "In summary, there are mainly three different approaches for running hyperparameter tuning for Spark based recommendation algorithm. The three different approaches are compared as follows."
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "|Approach|Distributed (on Spark)|Param sampling|Advanced hyperparam searching algo|Custom evaluation metrics|Custom data split|\n",
    "|---------|-------------|--------------|--------------------------|--------------|------------|\n",
    "|AzureML Services|Parallelizing Spark sessions on multi-node cluster or single Spark session on one VM node.)|Random, Grid, Bayesian sampling for discrete and continuous variables.|Bandit policy, Median stopping policy, and truncation selection policy.|Yes|Yes|\n",
    "|Spark native construct|Distributed in single-node standalone Spark environment or multi-node Spark cluster.|No|No|Need to re-engineer Spark modules|Need to re-engineer Spark modules.|\n",
    "|`hyperopt`|No (only support parallelization on MongoDB)|Random sampling for discrete and continuous variables.|Tree Parzen Estimator|Yes|Yes|"
   ],
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "source": [
    "# cleanup spark instance\n",
    "spark.stop()"
   ],
   "outputs": [],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "# References\n",
    "\n",
    "* Azure Machine Learning Services, url: https://azure.microsoft.com/en-us/services/machine-learning-service/\n",
    "* Lisa Li, *et al*, Hyperband: A Novel Bandit-Based Approach to Hyperparameter Optimization, The Journal of Machine Learning Research, Volume 18 Issue 1, pp 6765-6816, January 2017.\n",
    "* James Bergstrat *et al*, Algorithms for Hyper-Parameter Optimization, Procs 25th NIPS 2011. \n",
    "* `hyperopt`, url: http://hyperopt.github.io/hyperopt/.\n",
    "* Bergstra, J., Yamins, D., Cox, D. D. (2013) Making a Science of Model Search: Hyperparameter Optimization in Hundreds of Dimensions for Vision Architectures. Proc. of the 30th International Conference on Machine Learning (ICML 2013).\n",
    "* Kris Wright, \"Hyper parameter tuning with hyperopt\", url:https://districtdatalabs.silvrback.com/parameter-tuning-with-hyperopt"
   ],
   "metadata": {}
  }
 ],
 "metadata": {
  "celltoolbar": "Tags",
  "kernelspec": {
   "display_name": "Python (reco_pyspark)",
   "language": "python",
   "name": "reco_pyspark"
  },
  "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.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}