{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "This notebook mainly implements the **Structured Low-Rank Matrix Completion (SLRMC)** of the following paper:\n", "\n", "Jonathan Gillard, Konstantin Usevich (2018). Structured low-rank matrix completion for forecasting in time series analysis. arXiv:1802.08242.\n", "\n", "Technical highlights:\n", "\n", "- Embed a time series into a Hankel matrix and the missing data (to be forecasted) are stored in the bottom right-hand corner of this matrix.\n", "- Consider a matrix completion problem for Hankel matrices.\n", "- Choose a proper weighting scheme for the known observations." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from numpy.linalg import inv as inv" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The solution of problem\n", "\\begin{equation}\n", "\\begin{aligned}\n", "&\\min_{\\boldsymbol{x}\\in\\mathbb{R}^{(m+n)}}~\\|\\mathcal{S}(\\boldsymbol{x})\\|_{*} \\\\\n", "&\\text{s.t.}~\\boldsymbol{x}_{[1:n]}=\\boldsymbol{y}, \\\\\n", "\\end{aligned}\n", "\\end{equation}\n", "is given by\n", "\n", "1) Construct the scaled initial data vector $\\tilde{\\boldsymbol{y}}$ as\n", "$$\\tilde{\\boldsymbol{y}}=\\left(y_1\\cdot\\operatorname{exp}\\left(\\frac{\\alpha}{2}\\right),y_2\\cdot\\operatorname{exp}\\left(\\frac{2\\alpha}{2}\\right),\\cdots,y_n\\cdot\\operatorname{exp}\\left(\\frac{n\\alpha}{2}\\right)\\right),$$\n", "for a given $\\alpha$.\n", "\n", "2) Solve the problem\n", "$$\\begin{aligned}\n", "&\\min_{\\boldsymbol{x}\\in\\mathbb{R}^{(m+n)}}~\\|\\mathcal{S}(\\boldsymbol{x})\\|_{*} \\\\\n", "&\\text{s.t.}~\\boldsymbol{x}_{[1:n]}=\\tilde{\\boldsymbol{y}}. \\\\\n", "\\end{aligned}\n", "$$\n", "\n", "3) Scale back the weighted values\n", "$$\\hat{\\boldsymbol{y}}=\\left(x_1\\cdot\\operatorname{exp}\\left(-\\frac{\\alpha}{2}\\right),x_2\\cdot\\operatorname{exp}\\left(-\\frac{2\\alpha}{2}\\right),\\cdots,x_n\\cdot\\operatorname{exp}\\left(-\\frac{n\\alpha}{2}\\right),\\cdots,x_{m+n}\\cdot\\operatorname{exp}\\left(-\\frac{(m+n)\\alpha}{2}\\right)\\right).$$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "def scale_to_data(vec, alpha):\n", " dim = vec.shape[0]\n", " new_vec = np.zeros(dim)\n", " for i in range(dim):\n", " new_vec[i] = vec[i] * np.exp((i + 1) * alpha / 2)\n", " return new_vec\n", "\n", "def scale_back_data(vec, alpha):\n", " dim = vec.shape[0]\n", " new_vec = np.zeros(dim)\n", " for i in range(dim):\n", " new_vec[i] = vec[i] * np.exp(- (i + 1) * alpha / 2)\n", " return new_vec" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "What is structured low-rank matrix completion?" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "def hankel(vec, window_length):\n", " column_num = vec.shape[0] - window_length + 1\n", " hankel_mat = np.zeros((window_length, column_num))\n", " for i in range(window_length):\n", " hankel_mat[i, :] = vec[i : column_num + i]\n", " return hankel_mat\n", "\n", "def hankel2vec(mat):\n", " dim1, dim2 = mat.shape\n", " new_mat = np.zeros((dim1, dim1 + dim2 - 1))\n", " for i in range(dim1):\n", " new_mat[i, i : dim2 + i] = mat[i, :]\n", " return np.true_divide(new_mat.sum(0), (new_mat != 0).sum(0))\n", "\n", "def svt(mat, tau):\n", " [m,n] = mat.shape\n", " if 2 * m < n:\n", " u, s, v = np.linalg.svd(mat @ mat.T, full_matrices = 0)\n", " s = np.sqrt(s)\n", " tol = n * np.finfo(float).eps * np.max(s)\n", " idx = np.sum(s > max(tau, tol))\n", " mid = (s[:idx] - tau) / s[:idx]\n", " return (u[:,:idx] @ np.diag(mid)) @ (u[:,:idx].T @ mat)\n", " elif m > 2 * n:\n", " return svt(mat.T, tau).T\n", " u, s, v = np.linalg.svd(mat, full_matrices = 0)\n", " idx = np.sum(s > tau)\n", " return u[:,:idx] @ np.diag(s[:idx]-tau) @ v[:idx,:]\n", " \n", "def SLRMC(vec, window_length, forecast_horizon, rho, maxiter):\n", " mat = hankel(np.append(vec, np.zeros(forecast_horizon)), window_length)\n", " pos_missing = np.where(mat == 0)\n", " dim1, dim2 = mat.shape\n", " T = np.zeros((dim1, dim2))\n", " Z = mat.copy()\n", " for it in range(maxiter):\n", " X = svt(Z - T / rho, 1 / rho)\n", " Z[pos_missing] = (X + T / rho)[pos_missing]\n", " T = T + rho * (X - Z)\n", " vec_hat = hankel2vec(X)\n", " return vec_hat" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Toy example: forecast $\\boldsymbol{y}=\\left(1,2,3,4,5,6,?,?\\right)$ where we assume that the last two values are 7 and 8 (i.e., ground truth).\n", "\n", "Ready? Let us do it!" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Original data sequence:\n", "[1 2 3 4 5 6]\n", "\n", "Forecasted data sequence:\n", "[1.00016074 1.9998642 2.99993493 4.00004691 5.00041026 5.99939527\n", " 6.99627719 7.9884971 ]\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "/Users/xinyuchen/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:13: RuntimeWarning: invalid value encountered in true_divide\n", " del sys.path[0]\n" ] } ], "source": [ "y = np.array([1, 2, 3, 4, 5, 6])\n", "alpha = - 1 / 2\n", "window_length = 5\n", "forecast_horizon = 2\n", "rho = 0.1\n", "maxiter = 200\n", "y_tilde = scale_to_data(y, alpha)\n", "x = SLRMC(y_tilde, window_length, forecast_horizon, rho, maxiter)\n", "y_hat = scale_back_data(x, alpha)\n", "print('Original data sequence:')\n", "print(y)\n", "print()\n", "print('Forecasted data sequence:')\n", "print(y_hat)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This rather obvious sequence can be easily forecasted. But how about more complicated sequences?" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [], "source": [ "import scipy.io\n", "\n", "tensor = scipy.io.loadmat('../datasets/Guangzhou-data-set/tensor.mat')\n", "tensor = tensor['tensor']\n", "data = tensor.reshape([tensor.shape[0], tensor.shape[1] * tensor.shape[2]])[0, - 9 * 144 :]" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "code", "execution_count": 25, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAzIAAACPCAYAAAAhie4bAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztnXWYXNX5x79n49nduLsjgSBJA4QWiyJFiru18IMUKxQpUJwiLW5NgRYJEDRAISQhCSVFAhHiECPustlNsvHz++Odt+fcO3dm7vjcue/nefa5M3fs7JVzXn+V1hqCIAiCIAiCIAhBoiTfAxAEQRAEQRAEQUgWUWQEQRAEQRAEQQgcosgIgiAIgiAIghA4RJERBEEQBEEQBCFwiCIjCIIgCIIgCELgEEVGEARBEARBEITAIYqMIAiCIAiCIAiBQxQZQRAEQRAEQRAChygygiAIgiAIgiAEjpp+3qSUWgygCsAeALu11n2UUk0AjADQCcBiAGdprTfF+55mzZrpTp06pTFcQRAEQRAEQRCKmSlTpqzXWjdP9D5fikyEY7XW663ntwIYp7V+SCl1a+T5LfG+oFOnTpg8eXISPykIgiAIgiAIQphQSi3x8750QstOAfBK5PErAE5N47sEQRAEQRAEQRB841eR0QDGKKWmKKWuiOxrqbVeBQCRbQuvDyqlrlBKTVZKTV63bl36IxYEQRAEQRAEIfT4VWSO1FofCuB4AEOVUkf5/QGt9TCtdR+tdZ/mzROGugmCIAiCIAg+2LMHuOYa4Oef8z0SQcgPvhQZrfXKyHYtgA8A9AWwRinVGgAi27XZGqQgCIIgCILg5LvvgGeeAS66KN8jEYT8kFCRUUqVKqXK+TGAQQBmAfgIwMWRt10M4MNsDVIIFloDw4YBFRX5HokQBtasAaqrySKpFPCf/0S/54YbgG7dcj82QRCEbKI1bXftyu84BCFf+PHItATwX6XUdADfAfhEa/0ZgIcADFRKzQcwMPI8sCxYADz9ND3WGhg/nly2jP1YcMITKQBMmAA8+ihw5ZXA0KFm/4YNuR+XUPxUVQGtWgFXXEHXHkBKtJsnngAWLszt2IKE1sCXXwKzZgFvvZXv0QiCkIhhw4AxY4Ddu/M9kmCycSOwYkW+RyFkgoSKjNZ6kdb6oMhfT631A5H9G7TW/bXW3SPbjdkfbvY48UTg2muBtWuB994D+vcHnn0WGDmSrLw1a0oMqheXXgqURK4irYHjjgNuiRThXrSItqNHA82aAcOHA7Nn52ecQUNr4OOPgeXLzb7Ro4E//jF/YypE7r+ftu+9B+zcSY/feAN46SXv99tWy+XLgblzszu+oPD228DRRwMHHgice67TOCEIQuFx5ZXA4MHA1q35HkmwqK4Gpk4FTjoJaNcOkI4gwSed8stFxapVtP34Y+DMM+nxpEnOuNOpU3M/rkLnX/+i7VFHAffc43yNvVizZtH2gguAww/P2dACzUsvASef7DymF1wA/PWvwPTp+RtXobB3Lyl1jzxCz+vXB5YtM6/feivw7rvADz84P7fJatl71VXAGWdkf6xB4IsvnM83xW1tLAiZQRTm1LCP2wknmMdDhwJ/+INEkMRi925gwACgd2/gm29o3y9+Abz2Wn7HJaSHKDIR6tSh7d//bvZ98QWFrjAcHjV8OHDEETkbWiCYODFakfnhB2DQIOCjj8y+LVtyO66gsGoV5XG88QbwzjvA735n9jP77kvbTz/N/fgKjffeI6UOIG/qhg1Oy9rmzWSQOOQQYN48s39jxG+sNS1kixaFT5hizxWjNSUM26xaRWEXJ54IrF6du7EFnaeeAhYvzvcogkOnTsDll1O0g9bA9u3RxgchGi8vzKRJwHPPAY8/TqHyQjQvvwx8/XX0flvuE4KHKDIRWJj5/nugrIzCzFaudL7nu++ADh3IMv7tt8COHbkfZ6GwebO3UnLaaSSMX3ghhfGMHUux9wCw//60laREw9q19Hf//ZTHcf75wFlnmdeXLjWPlaJt2K3ln3/uPEbnn0/bMWPMPvsa22cf87hvX/Lm/PwzKT/bt4crf+v778lo88UXZs77+utob/PKlcCHH5LS3Lo13e9CfNatA667jpRnITGLFtH89vLLQJcuJExefz0dP8ldiE+ilnyVlbkZR9D4z39Ihuvb17l/2zba7t1Lsp0QLEKryLAw064dWXdtYWbffYEWHu09X3rJGb4S5smiUSMSEGvVouennQbMnAm8/jrF2L/6qvP9TZoAV19Njze6sqkqKsKZf7R5M9CyJdC1qzmOTOPGFNa4eLEROPkaDfN1BwA//uh83r+/eTx9OgngsaiqoryY7783++w8pGJn9Gjannwy0KABhYa+8AJQt64xNAAUdlddbZ6/+WZOhxk4NmwADjuMHldUhNvI5RcuzsFMnEheBSD6Hg87b79NRoVp0+heXL8+/vvDvkbEYvZsoGdP8jTbbN1Kxuszz6Rom1Gj8jM+ITVCq8jUqEFWoBUrTJw8Ky+9e5vHnTvH/o6wTxYrV5Ll+5FHgPffBw44gHIVGPYgAEB5OSkzQLQic+CBdC7CBoegbNkCPPmk87Xf/AY4+GASvI8+GhgyxISZhd06XrOmebxkidPo0KsXXU/x+OorpyJjGyeKnZ9+om1VFV13l15KxofTTiMlh0PMFi4EbrqJjnV5OXlVpeqbN3v2AOec4zTGjB1LOVyi0Dj58ktTFc/tdSkpMWuEFIVxcvbZJHwfeihw3nmxFRk2iIVdNvFizx5SkHv2BO64wyjNAIUfP/00yTEAMGNGfsYopEZoFRmbpk1JAH/5ZSrleuutQPPm9FrNmhT64xWT+tRTJmzq228pPjUM/OY3zuctW3q/zw7pWbKEjjNAyYl2GBBbxMOUq3Dqqd55VpzIf/zxtHiVl5OlcvRoE1IW9kXKzltr1oy2//iHycXq0AFo04aEy8svp309e5rPfPst5SHtsw8JT7ZSU+xMmRK9r0cPOn7t2lHi6403mtd27wZ++UuyAnfrFq571Oaf//Q+dgDwwQcU7mjz619TDtc772R/bEHi6KPJYw9Eh24rZZLUbSEToP1hvfa8Qrjt0FobXnNHjiTvjWBYt44MC50707zft6+R89y4ja1CYRNKRcaOn+/dm274L74gi8eqVZSAyFbemjXpYq9fHxg3zim0P/UUTcwACaV235RiRWtauG26d/d+76hRwGOPmedsbVu8mMpGau0MX7EfFzsffuj9//bqRQvX6aeTMN67d/R7wuyR0dopULIH8Le/JeERIG/rkiUkfLdtS/t69iTluXFjEkqXLgX+9Ce6f997L7f/Q76YOJHC6mrXpuennkr9siZOBEpLzfsuucT5ubvvNo/DqERrDVx2GdCnj3P/t9+SoeHMM50Cke0xDPO9Go8dO6IVmZIS43V+/31qdMvUrGmMEmHDq0R8rKI5XbvS9tVXycsqGHjuatDA7FuwgNYCXicYbh0hBINQKjJ2PkyvXkD79mSJtOGFnauZAdQjZezY+N9d7FYjL5f2fvt5v7dTJ6rE9be/USlc9sgw/fo5Q9HCknQdK5SJvVS2UOllMaqspOtMa7IchalK0ltvASNGJH4fC5OsPJeXAwMHkrcGIKvchReSVW7ePBJS//a37Iw53+zcSWGKRx1Fz085hbZ79wLHHhudD2iH0z71FB2j11+n52GsYGYLjSNGkJdq504yXvFrxx1n3mMbdsKY++eHu+6iEE8brUmR6d6dinC0akVGh7176fV//jP34ywEJk703v+rX9G2Tx+jwNj37pIl2R1XUPjpJ4r6YE9+ebl5rUEDkv+WLydZhbErXQqFTygVGVsY79bN+z09e1L8+PDhzv0HHhjdQMlWXoq9ORXHNV9wgdnXqFH8z/zhD8bDYOOuDhIGd67W1IjLi4EDo/c9+ihZ1tavp/yF888nReaUU8iC2bSpEVDDgB0uwY1X48Fdr1k55O0ll1AoS9Om9J4pUygnpBgZOdLZe+j442nLAqKb0lK6tq68ErjmGtrHnug1a0gguOuu4p/rGFt5O+ccqjDoLhDRsyd5uACn8WH+/OyPLyjYvToefpjm+/btzb5XX6Vri4smAOQJs0vjjh+f/XEWGm4ZhOFy/JdfTp6FTz8lLzNjG2HDitZ0nLp29VZkbOxrcf782POjUHiEXpGxL16bmjUpZ8au5BPrM3YSrB2/X4zwAv7731NcPTck9IM9sZ5zTnRuQhgUmfHjKZGQXdnsOfC6zgCgY0cKs2jalLyGTZqQpe3jj817wpKsvnChCT1p2RJ46KHEn6lRg7bsdeDFiXNreFvMsOW7XTuan448kp6zZ8aLPXuomhnTqhVtb7iBLMD33ktl1sOAHeIE0HFxW2w7d6YQxZ07nYrMqlV0r4YpbNaLuXOdzaWZ+++P3sfXJ8PVLgGqUBimAgo//kil0TnqoVMn2l52GVCvHj2uW5e2xx/vvPZ27KB7ftky8jYXe7SIF+wR3bnThHnaoWU2Z59NhuqOHel+lRLgwaFm4rcUH6zInHYaXbzJ4r4R7FCCqirqu1CssCLTrp3TcpYsr75KFVZ+/JE8DRddVPyhZTt3mqooDz9MXq2DDwYeeAA46CB/3xGrh8zOnRSC0KFD7JylIDNlislR6NEjOhk4FldeSdcVhw1wMjErMO5wx2JjyRIKWezc2cR99+hBx6RxY//fwx4Zu99MWCpLuRWZNWuMV4s58EDyYpWUOJXjadOo1PXJJ8cvC16MPPss9e1Qyjsk+de/pvDOiy927j/00Pjf+913Jqyq2Ni4keaqp58mOYMLmFxyCXmgDzyQGvm2aEGe+dq1TfEEL1asMOG0p55qQtDCgh2+yBUZY3lk2rQhI+OECRQqOnUqpRJ060YFT5RyVmIVCodQemRYYH7mGZP4mgxsAfGi2D0y33xDXgG20CbL1KkUd89lIvfZx8SXF7tH5oYbqLJd3bpkWWzQgIohDBoUu/Kbmwsv9N6/aRMwYAAJqcVGRYUz0fqooxKHMzL16gH33WdysRIpMsVmtbz+ejIWuP/PJk2SW5S9FD53iG2xUV1Nhi4OZ1q0KNrYcsUVZPW1DRF8bfXoYUIbWSANE7//PVVte/tt75CwW26ha/DZZ5372esQi2LOX3j8cTLy8TGZOpWOBysgxx5La29JCc2Bjz4aHULG4chuwthI2U7a5+smliLDHHEEGWpvuYXC9o4+muTEk0/O3jiF9AilItO7N8V4pxNWMmaMd5J7MSsyWtP/PXCgCdlJlkMOMXH3TPPmtKBx2FCxwuFgnMi6eXPylsVBg7yvOzuOf+nS1MdYiNihcw88ADzxROrfxYpMw4a0dc8BDzxAC1gxVJvSmgwPgH/FLxYlJc7CHEDxhjRqDdx+O/2/b79NxofSUrJsc/EIgHL//v73aMGbwxjd5dWLTUmOh591kO+9q6+mqqFM8+akIMbKgSvmvCM2rN5+O1UL/OQTWjN/8xtaP667LvF33HRTdOVBIJyFOuy1kKvkJVJk6tYF/vxn03MLoHXj3//O/PiEzBBKRaZPH5okUvHGMAMHRneHBagEc7EuWIsW0WR47LGZ/d7atUmwL/YqK5kKOWRvls2cOeZxMVUxq652WmBvvtlZ1S1Z+valLXsYbE9DWRlw552U9xX0MCCtSQjnsChW4NLhv/91Pl+zpjjnurFjgQcfdO7r2zfaePPXv3p//txzSflxKzLFoBz7xU9Jc9uIcPTR1BahqoqMWn//u3cOXFlZcXtkOO9Fa+Cee6gq3qGH0jE56SQyKPjhkENoe+mlZl8YFRm7YTLnvPhZP9zhjkJhE0pFJlP06kXbr74ytd7nzIld4z3ocEfmww/P/Hd37Eg5Htu2Zf67CwWu8pRumElNj8y2f/zDPHb3ZwgyPXsCZ5xBj3/6yft/T4annqJqeWxFb9KErufXXnNWkwu61XfECCqowWQiQfqQQ0zISocO9J3FKJyPGkVWWdvz16+feTxnDjXAjBWa17w5cNVV0aXTi81TGo8PPgC6dIn/HreXsGZNUlRsbGXnq6/oHv3gg9RyW4MAhyLasFKSDO3bkzJ0/vlm3+9+R4YbuzhRsbJhA/3/S5YABxxA+1asoOvLjzIYy8gdpubJQcK3IqOUqqGUmqaU+nfkeWel1CSl1Hyl1AilVBr+jWBywQVk/e7Xz+muLLamcbt2UbzzHXfQ4m13Sc8UdetSCUm7lnuxsWwZHUdu3JgqLMz37GlKL0+YYF4vJkXG7sORKHbeD3XrOotUlJRQ+NUFFzirIwXd6utOrs50pScOb/zTn4rPkzphAlXOuu46k79ne1f2249y3BLhDltcupSMGRzuV8ysXOmdr7d5M3k9W7TwF548c6Z53LWrme/efjszXsZCw0t2SEWRYdw9ojZujF3OuRhYtoxy95o1o3t4wQIqkACQkpgorMzmhReii3qwR18oLJLxyFwHwO4x+zCAx7XW3QFsAhC6vrtKkScBcN4gxWalHDnSJB8+/XT6VnEvOOzKjktNl9dfp5yHQuDJJ2mRshuWpQoLkcOHU2z5zTc7Xy8mRYapXz+9UFA//OpXZMU7/vjgKzLuHAW/VfESwaFkbG1//vni6ri+YwdVY2Nlt1072qbihXYrMkuWUAhVv37+q+4FlTVroguYdO1KBU7uvTe6Elws7KIy5eVOT8zIkemPs9CwFZmPPgLOOy+9kGSvIjLF3B+lQwfT3JwNBt27m7UjGWPYlVdSb56vvnJed8UYTht0fCkySql2AE4E8GLkuQJwHIB3I295BcCp2RhgUGjQwFiLMu2R2bkzfxW9tm51CkXZKnv59NO0zWRfjwsvJC+S30UzW8yZQ9WjADPJpsOzz1IM+kEHkTL98MNUhGHwYJqo33+fQl+Czvbt5rE7TCebdOkC/PCDs4mkX8aONWU+88nq1XRtVFZSyCbfX+nCizj3QQK8c7ZSYfjw/Hsr5s4lyy2HDZ9/PiWdp1Km257LSkupGtWUKfT8xx/THyszcSJwzDG0ThQCjz9OlvEWLYATTqB9y5c7S3cnw7BhJMzXq0ehoFxS94wzvEOxgsi0acC111KFRgD4+mvy3A8fnl7JX75uDzrIeMDsSl7FhJ0natOxoylUksr626+fs8/buHHJf0cQ2L69cOaQZPHrkXkCwM0AWJdvCqBCa83TyHIAbb0+qJS6Qik1WSk1ed26dWkNttDhBNGKiujE2FRZvpzqmOer30VZmdPiyjXpM03TpuQKzoTCVlHhrHiT79KndnhUoh4JfigtpSo2NgMHAp99Rtann3+m5+la3nbvpjyLfFmPbc9SpjwKfuCGtwcfHB1Pvnu3ETa8GDSIrPn5rsC3Zg3NG+Xl1APBXXEsVdgzMXiw2deuHSnPqSh+zOzZFN53+unpjS9dZsygLSsygwb5a7zqhV3h7KmnnApuJkOjLruM+rVk0pudKnv3UjU3gLwB77xD3s22bWM3IkzE735HcwEL9AcfbF7LlEL45ptkhMgXv/89GRtGj6a5zl0oIlVq1CBZ5PPPKW+kd2/KB3zhheIrx2yHIXJVSoAUGZ6zUzUkchEGgNbWYkNrUtaGDMn3SFIjoSKjlDoJwFqt9RR7t8dbPR1uWuthWus+Wus+zXNpVs0DfPM88gh5Luzu66mwdy8JVVzm1LZQ5wIvQTjdMq7xaNw4M5Prs8/SOWDSEbAyAQvDf/5zehW3/GDnL6Xb52P6dEoa/+1v0/seN1qTohnre6uqqO8EN189/HBnMYNsYyvrd97pfO2qq+g69RJC7fuzTRuylOeL1av99yZKhptvpuuiTx8qmgAAL75ICkg6CdjsQUxV2E0XremPK/4lSlT3Q82apMzUqEHVo+zl79NPyUOTCTiseelS8lbkK/Rl0SJn7lC9eqRAZ7pBb/36xvr+6aeZ+c7zziOlNV/wtbF+febvgSOPJO9gw4ZUCQ2geez++zP7O/nG7vHUqRM1AAVM+D+QuiLjzjEsFk8g88UXZACdMCH/RrhU8OORORLAyUqpxQDeAoWUPQGgkVKKsyXaASjCyPzkYEWGwyPSrXzkFsDdzdiyjfuCjtdBOBM0bkwhPWvX0vOKClMNLhlsi/lhhzktNflg0SLybN19d/Z/y87V+uST1L5j2TKKD+YmlJn2Bi5eTIrmSy95C10XXEACEYeivPxydNJqNrEVmTffpPHcey89f+UV2nIy/dKlpNScfLKzshVgwojywZIlpExlmpIS46047DBnIvL8+amH1XJe4cqVRoHNJXfeSZ6ljz6ia83dZDBVmjWje18pp5UYSG1u84IF3/vuI+8bV5fMNe++6+wHk80ogh49qJnyPfcUR6VLu9JpMgnpyXLiiRT6BwDFFiBjR3M0b05heVOmOBXDVJVqNlJxeB7LKEFn+3ZTIAGgMOFMzUu5JKEio7W+TWvdTmvdCcA5AMZrrc8HMAFApDAqLgYQ8M4L6cMLFV/06XpQ3BNNrhSZYcPIQmWH1fToAbzxRnZ/t3Fj2nLoQNeu5O7UmpL2/VZb+fFHsgYuWkShXNOm5Tf2c+FCsvCmE+vsl6uvJvdw166pNfDavp0E+WHDzL5ML6x2RS2vMC2OgWeDACdc5wr+Pb6fhw+nBrqACdFauZJCUTp2pH4iH39M1xlgrMQrV+YnHn3NGlIWM5GPlQgW+Hv3Jg9uquFNrMhUVZnQvlyxcyfw2GN0vqZMyezvN2tmvLB87QwbRgJDJvNkABMCaoey5pKNGympurKSPKpnnZW936pRg5pGbtuWfllrLosPUIXOfMA9TgBjKMgW119PXu5Ui8IsX1545ek3baJ7mKlfn/44lPuaa4BTTvHfh8cN5wZzQ1I2uGzZQoJ/IeRFpsIdd9B6f9995MHfts1UagwS6fSRuQXAH5RSC0A5My9lZkjBxV0D354orruOrCHJ9JhhC8Nzz9E2F4rM0qVUAvnNN81Ne/LJuckzYQF31SoSKPj/nzmTbrgLLvD3PXPmUL+Bzp0p2bSqio59LkMuNmygULKvviJhNhOhKn5o1Yr6YFx+OXk0Vq0CBgygfAk/YXt2OCSHWmTa+mRbztzfPXKkuc7ffpusadm0UHrRsSMpI+6eAXv2GKF01SqT9PnOO7R9+mlaLI85hiyCjz5KCmWuQ31YAXR7iLIBe2C4V4MtFCaDu9LjYYc5BZNs8tln1HiVy7RmsqpT8+bm+uU4+333pX48s2dn5jfcBq98hZZt3EjGqPJyapqcbcMNV4BMt/y3XQyGw7hzzcqVxuuci2aMnTolp/D+7W8ULQGQou9VWjsX7NkTfb3v2UPRIvb65q5w+dRT6VW5GzyYvp+jUu68k45DmzZkbLVL+gcJ9sRUVdF9kI2KtLkgKUVGa/2F1vqkyONFWuu+WutuWusztdYZ7lQQPNza/ogRpLjs3Us30qefkpDpFxboeNJw94bINNOmUaUUjgfliWvkSHLjZxvbBWwn6//tb/4+v3s3JZYuWkQWYsAk5n3+eW6ql2lNyZVPPmlCPebMIYE2l3BzxwMPJIF74UJnr5lY2Ja2U08FLrqIrE2ZTAy1v8telKqqTMNFJl89mY4/npQ/mzffNAvkypXmeuLtJZfQvVKvnrNsbDYa5H76KRV80JqU/mnTyBty+eXkIapVKzOFJRLB54dLgnspMo8+agqh2KxdayzRbkXmu++Am27Kzj27cqXJcZo+nZRPwHicM5mj8Oc/01wAkCdm0CAK2ezRg7xmmfAUuwW7fJVf37TJWeAg23AI6JAh6UU/2NfYJ5/QOpvLHjXLltF9dNttNFfsu2/2f7NzZxMWm4iqKroXe/d23qcLF1II8IgRZHTKhQJ9992k8Nmy0LXXUpEEm0yFhto0amTCdT/7jLZ2Rdfbbou+DrdsyZ+H1A87dlB4cIsWJuogkGitc/bXu3dvXex89x2njdLfq69qvWyZeX7xxf6/69576TNLltD2+ee1rq7WeseOzI971y6tO3bUunVrrT/+WOurrqLfPP30zP9WLDZv1vr2282xOuEE57EEtN661fuze/ZofeSR5n3jx5vXRo+mfV98kd3xDxoUPV7+e/TR7P62m717tX7kEecYbr018ed++1utGzXS+umn6Tq75Rb67B//mLmxPfecGdP775v9Z50VfdyaNMnc76bC2Wd7n89779W6Xz/z/KCDnJ+zr4WnntJ69uzMjou/+x//8B7f4Ydn9vdiUV5OvzdyJG1HjIg91nj7Bw50jn/IEPP/ZZJNm+h7e/fWetQorW+8kZ63b0+vjxlD83W2+ec/6XfnzaPnZ52l9TXXJP89ixdrrZTz2DVoQOtErti719zT/frl7nd37TL/87//ndxnKyu1btxY63ff1fqDD+g7atbUulMnrUtKtN5/f613787OuN28+CL9/syZufk9rbUeNszIFm727nU+nzrVOZfx4xtuoG39+uYcbNqU3XEfcID5/XHjtJ440Xv+u/TS7Pz+nj3mNxo10rprV+fvPv208/033aR1y5b0uLpa6+3bszOuVGndmo6V+5wXCgAmax+6hSgyWWDKFFqg6tYlAecXv6AjXaOG1r16+f+e668nQWHHDvr8nXdqXbs2PX7wwcyOedSoaMFy+fLYikO2WL7cTAqffUaLuz1RTJsW/ZmVK53vadtW66oq8/rPP9P+YcMyP97HHtN61ix6HEuJ4Uk3HzRsaMaw//5af/651t9843zPxo1aH3ig1hMmaD1ggNZ9+5rXVq+mzx5/PE3iNtOmaX300c5j7Yf77zdjeuEF2ldZSfdLv370+O67tb7vPq3nz0/2P84sO3bQMRs6NP75ffJJ5+fc7z/rrMyOK95YAK3//OfM/l4sjjiCfm/uXNq+9JLWjz9uBKTqan+KTN++5vmFF9KC36JFcoYfP8ya5X28tm3L7O8kggWwUaPoeaxjlIjbb6d15fTTnf9PLg1QPL8CWh97bO5+V2utzz3X/DbPwz/+mHhOmjyZPlNSQsInoPUllziPYaaV6FgMGaJ1hw65FSbHjDFKyP7703HUWusrrtC6SxdSEvfu1fqNN+h+9rpnjjvOe3824fkm0d8jj2RvDPwbEydqvWYNrZ28z22MGDBQ/B/sAAAgAElEQVSA9ldXa92/Pz1euDB7Y/PDwoUkM23fTuO5++78jiceosgUAH36OG+uU07Ruk4dmiT8cOGF5CXRWutf/pImu2xNGPfdR99ZWZnZ702W3bvN/7dhg9arVkVPUm+84fzMv/9tXnv+ea23bHG+vmcPHfebbsrsWFlIa9yYnseaVL/8MrO/mww8kdqWLPe1wwv5iSeSZdotdJ90kv6fMmPzq185hTG/3HgjWT9ZKf/2WyNY2Ip0oWF7YPh48eOPPnK+d+dOp+fpgAPif/fo0Vr/97/+x+J1nTVvTp6Yv//d/xyTLuvWaf2f/9CCDmh9+eVmPHv3krLLz93CJe/fs0frffcl4dtWsk86Sev99svseMePz73w5QUbX559lo5TquM45BCtjzpK6xUrSHn+4Qc6lnXqRBsesgULxYDWPXrk5jdt+LdffpkMiABFFNhUV2u9aJF5/sknzvNfWqr122+b5506kRFx7tzMj3fqVK0feojOz8aNpEzdfnvmfyce8+dH3wPr15vHa9eSF9l+/eijnc/r1cv9vTR4cOx1tmlTUqo/+yy73jT+veXL6bl9nE491fnedu1o/9Kl5j29e+d3neNxzJhB2+HD8zeWRPhVZNJJ9hcSYHeDBSjufscO/zGTGzeamONLL42uzqJ1+mNkZs6khPRcJ1a74fKGAP3vrVpR/C33rACo/8iaNSa+3E4qPfTQ6F4tJSWUdJvpYgmcqO7O47DzPPbuNRVP8sGbb1IPlqOOiv0erhL2yScUr33yyc7X+ZyMGkUVTThZn+OQk02Q3bSJKqTceCM9f/ZZE1uc7T476eCuAnb99eaxOzegVi2635l582L3HvjuO0omPfpof+Pwyqt46im6Hr/5BrjiitwlbTZrRtcWnze7d86hhzobXMbqT7BxI8XeN25sGm4ClEsydy5dGxs30jFasCC98XoVrsjHnNeyJd1XK1Z4558995ypgheLNWvoPYMHU+z+M89QM8WhQ2mdyUV53YoKZ/+VdCuIpcKsWbRduBB4+GF67K6ed845tL7xveO+Fvfscc7TX3xB5+eww5wVxdJh+3Zaqw49FLj1VsrPWr+e1gjOMcsVXpX57L4ymzc7891q1wb+9S/n+6urvb87m5XfvGSeHj0oL3bxYipiMHiwU47INJx/2ro1bfffnwrDDBxI192YMbS/qsqUkl+92oxpyhST38hkssCIX7g6nt1gNqiIIpNF3A3pWLHxq8hs2WIWWa/E3UwJ5j/9BHz4YW67p8fjl790CnVnnumsCrJtGyk499xDJUdPPNG8FqsJYIMGmU8cZ6GoVi3aclGG116jSe2NN3JTcjkezZqR4te2rXO/neRpVxHr1y+6X5DdBHXCBJO4zd/Bzen8snkznY/77yfhddcuo8jYHZQLjQEDaNuhAxWjOOYY85pXkrPdAHHnTu9SzPPmmWt7zx5/yd92fyn+XbfRJNfweZs3z+z74Qfn//zyy6RkuYWR//6XhG53v2QWFNavJyFqzBgqHJAObkXmscdMv6JcUlJC/9/kyc714JRTSDAfOpTm/Ndfj+4LM3kyzStcUGLwYOfrnASfC6XitdfM4y5dqGN8runZkwpzPPAA8M9/0j63MP1hpDkEC5asyPDctX27MZo99BBVLhwxgtYMFkzTYc2a6Llt8mTTA4dLc+eKOnXoeNkKzRNPmMeVlc4mkGecQUrC2rV0jBivIjbZKko0aZJ3IYuyMipe4K4amy3efpuMd3Zxpz596JqZO5fux+nTnXPhTz/R/G4rDTNnktGwe3dScrJZwvm775zn1yZfFegyiSgyWcRemF94wfRJ+fBDKtWaiJ07TZUkL0HFTxUqPzzzDFmLuetvvpk40dlYLRYzZ5qO4EwsRaZhw+jKSOnCFs/t2+l8VldT5arSUlr8st1ANBncVmfbmmb3cvnww+jqe4895rxe2YLMXdDHjEnOO7hzJy2kJSUk/FRWGkWmbl3/35NrBgwgj9Rrr5GwY3s9vBSZ0lLg/PNNNRivRmNffklbXoT9VIf7+mvaLl1KAtK4cc6O6vnAvmYuuST6daXomP3jH1Txx/ZO/eUv9PrQoc7PcEPFBQuMUJBqHwjGXQVt0KDo6nS5okkTuneGDDH7PvqIDDfMhRdGzyPuRrd2U1LAKDKjRgG//312K3CxknX55eQRyUXpYC/Yk7d3L3mnJk4kZWXkSOf8xnPWypV0/Pfbj66BF1+k/WeeaSpmDhxIhiA/a1Ei2Gtk84c/GG9PPjzRf/pTbC/w5s1Gkbn3XrpvAZJpuOw1AFx2GW3tazgblQZ//pnO8axZNKfa63yu14x69bx7m9ktFr76ytknis//0KGk8AB0X3/8sfEyjx+fnfECwJFHUluN//zHuf8vfwluyWUbUWSyCNeFP+II6pTOlu3nnqOSgYmsr7YiU7s2CQLnnWdez1TDsU8/pX4r3EehUJk+3fS2AUjQcS8AsSa1Bg1octaaPBTuMr+pYFt3r72WlINC9SiwAMjeI1uR2bSJrL/Dh9PC7aZJE6eQyQvV2rWknM+ZQ6F/Dz3k9N7Ewr6u2VMWBEWmbl1SGrzC9NhI4eb1100InZfnato0UjK5AakfL+uMGbSQt29Pi1ChNTC76abofb/5jXl8wgnknWG++44Wd7fXkBWZY48F/v53erxgAR1Hu+xpMqxdS5Z3Lrmc6+abNuwVcFuwvTx3I0eaUEbb2zBnTrRyx4LSXXdR2Ga2mrIuXkxK9YMPGkUgX7z6qnnMQnWbNjTP2/cm95xZtcqU0h09mhQxN0qRtzRRiJ8fbOH+pZdIwdyyxfTxyLVHhrFDsH7xCzMPbd5s5JOBA53js402++9Pa8knn5iw0mx0vbe9YrVr07XH+3LdMDkWdojlV18BV11lnnPIY7dupn8U97tjPvww+w1ZOYqgpIQ8NLfemt3fyxWiyGQRVmTY+ugW8hI1RNu1ywieAFmKXn89MxYiRmuydhRKWFk8evVy5h3Yky0APP987M82bEgC87RptJCMHJl+jpF7wt62rXAVmYsuIkv5fffRc7ci07evU0l2Y4fITZxIn9m6lb6zdm3g3Xepjr5tyYuFfV0HSZGJh7sBm015OS22Xh6ZOXPIgMAKpB9FZt06Z5+aQqNTJ+fz9u1JcXnoISN0XHml8z0sfNuwImMzezaFErEikixr19K8/N57JAxlsmdMpvAycJ12GvWjWbbMeBVeeME7t8KtVK9enfEhAjDr17HHZuf7k0EpWscmTYo/j/GxW7XKhC7GY//9KSwoVn6bXzjM7+CDyYvBzZ05VCpfiox9H+2/v/Hs2vO4uyeLrciUl9OcXVJi5J1seGRYSerfn8Zcty49fvBB0zA83xx6qEkBGDmSjC3uSIhWrZzG19/8hgxTABkD7bzLTMJKO3P99U6jcNARRSaL8A3Prv3SUqcFxN3EyY1tuWaUInfwzTdnpulTdTUJ9PlO8veL7VLevNmZ3/F//xf7cw0b0oLEjTKB9Bd4d65TIXtkysspfpwt0BybrTUpJbE8CjZjx5JlqbKSFGqArOhHHmkKBgCJY6Tt67q8vDgUmUR07uxdFKGyko49C+1+FJn16709Z4WCPc99/jmFgDZoQIYYO6fCxkuRcf+PrVoZL0aqYbWsyNSo4f2buWTUKJNw65fnn6d559hjo5VBm2efNY9jFVlIFzaGFEqBjk6dyCDTv78JhXV7+ezQMj+KTM+eNF/ZjYJTYcYMMmSyd6dhQ9ryucmXInPaaRTedu65FGbE44qnyNjKv/24QwfyECebM+mH+fPpvH7+uVG+SkrIeOZl8MgHSlEi/6GHmvV18mQy4pSWAo88Et1YvF8/MmRx+Ogrr2THK+NWZLLRMDSfiCKTRVg47N6dtkqZiQKgZLt47NoV29Jbvz5NNOnGP3PX8VwlyqULx38DNNn6LXjgZXlNJ+RC62hFVOv8LUh+4fGxELJ1K11DfhSZAQMoTKp2bZOYWF5OAqatUCaqllRMHplnn/VnRWvWzPu4bN9Oyi8bPezjGIsNGwpbkQHM+T3gAOecd8wxTgG8Xz/aegmVbgHlnHOcz2+/PXlLOSsyhUDv3s6QLK1NGGIsHnqIhCO318vN1Veb642F5auuymzxEZ5DCtF4c8klZNjiwiQAhYktWULHefXqaOHOi1/+kq7l225LfSxz5lD1SPsa53si3x4ZgI7DG2/Q+HidtJP93UKvHcpoG0Dr1iXF/PvvMz/GhQvzl8uWLHy8OnWiRPpLLyU5649/NPffyy9TFAwn2r/+Op2DrVtN3mQm4XPG51cUGcE3PXqQi5HjTgFzQTVsSBd3vIV4505naJmNWyCNxe7dpDDFqtgVNEWmSRM6LuecYzwyXbqY/yMWtjDFLFyY+jjWrfNWhApxUbfh8VVXU8WUm2+m535yWwCyZHfrRpYngK6bxo2dAniiGGl3jkxVlbmOg6bIXH018Pjjid/XvLm3p6q6mv7nZD0yhWKFtPnoI5OIPmIECY5eCpftCXnvPSomccUV0e+zjTjLllFum82DD1I4TDJVkgpJkQGiLbTxPOPnnENC+M6dzoTrWDRpQusHe54zXVGskBWZsjLyXl10EQmNl19OBsXFi+ke27XLn0ema1cSQD/+OPWQKfYC3Xmn2cfzbb49Mm5q1aLzaYdtJwqbtTniCMojzGTi+qZNdM68KqQVIqwk2Mn/bi69lKo6smJTUkJh83Xreld2TJctW4BTTzUhyaLICElxyilObwAvumxdYEu0F16hZYxfRWb0aOCOO2JbjTlpNiiKDECTLVch27iRwgcShTd4eWRefDH1sAtWDP/yF6fwXYiLug2Pb/FiWpw5r8iPR4bp0MHke5SXk8BkV9tK1iOza5epKBc0RcYvzZqZnhE2HI7I5yXefACQYWLTpsL0yPz615TID5CS/O233v0cbIGkVSuqpmN7Wm3eeYdyMdq1M55tm/nz4+fG2VRX03xXSIoMz0sc8sqCodd8zMcW8KfIlJRQKK57jsuUkMThM4U+5116Kc31HTuSR2bSJNrvR5EBqBfS3r2phzPyWmHnodatS3NgoSkygPGSx/LIuN9rw5VPP/3U7OPrJFU4hDueYlBI8PFK9pw2akRFdRYtSpxn6peqKpojq6poTmFDejzlNIiIIpMnOM4z3gUbL7SMF49EkwRr/D/84P160DwyDCsyq1ZF95/wwu1xqFmTEgjdPRj8ws3CevRwVv4o9EWdJ1d3CddYZau9sJscsiJjC+iJFBm3RwYgS3lJSXGUgvSieXMK4bNLwQJGkWGBP1Go6KZNJIgWoiLjF264euSRid97xhmm9HysudBv40d+XzLXei5YtcoUcOF52Mszc8AB5nGi0DKmWTOyZv/tb2ZfpoSkQvbIeMFlqtmz51ehZYUn1fL9bCy0hX4OM+d5s5COYVkZrW9+FBm3AbFpU1Ky2Qs4Ywa9h/v4pAKXqM5nhcFk4OOVilGOi6EkMmj55ZxzSEb5+Wc6r7zOiEdGSItPPqGEV55E4y0qfkLLEikybCmP1YQzqIpMgwZ0fH76yV9DJ1vZWbfOhPQlqhwXC1ZkSkudYWuFZFnzghdMd7dqv9ZJwBniU1YW3UMlUXiU7ZFhgXzxYpr4891ANFvw/+kOg2JFRilS5BLlfHAIn1ffmqBQty5dI7bV1i9eJeLZSKMU9caIBc+B7gTwfNOqVWwFxp5bWrWiZODyckpC90PjxrQG2CWxE3nx/RK0cNDTTyerPgvZfsMz+dwkCl+OBXtk3N4LNq7Vr19Y815ZGf2vfkLLvPo6tW5tjjHf42PHpj4ebmJaaPdtLFhJSEU55XspU/eo3UA01IqMUqquUuo7pdR0pdRspdQ9kf2dlVKTlFLzlVIjlFJF5qzKDiecQKUX+UJKpMgkCi3zq8jEujGCqsjYbmZ3nLkXtiJjW7NTvaFtRcb29hSSZc0LHp9b2UhGkbFhj4xNoh4f9nXNIZazZwdHIEoFvv5s78GePXQs+JzUqJHYI5OvTuCZpkmT1Eoff/kl9ZPiansXXUSVoPh6/stfYn+WK0YVcql5ex5esMBZLat5c/p/N2/2HwrKioxNJhUZVsKDgu3J8qvI8L2WriLjXmNZSS20e7m0lP7XeB6Za66JfQ3a1QU5jzSdynYrVpgwySDA61g6ikymPDJ2BbTSUqN4hk6RAbADwHFa64MAHAxgiFLqcAAPA3hca90dwCYAHi2lhFgkumC1ju4jY+NXkWELrlesOmAm56CUX2YOO8w8TtYjY5NqrKityNiCfKErMnzduBWZVBea8vLoBS3Rgm9f16zIrFhR3IoMH1/7fuV7n6+ZmjUTKzJBC+fJNI0amZLFjRpRFbQtW4DPPkv82WnTSGEv5B489jzctatz3uKwy2QUB3f+GmCuu0WLSJi2G/clQyGXm4+FrTz79WrWqEH/J8/5ycL9RNzeCx5LoSkydmhZzZreXpennopdYbFVK+ORYUU8nRLgK1bQfRuUsON0QssyrcjY683y5SFWZDTBokmtyJ8GcByAdyP7XwFwalZGWKQk8sjs2UPKTLrJ/ryIbd/uneQZVI9M587kWSkrc/aGiYVbkeG4+1QTX3lRq1/fadkrxGpSNrE8MsnCE2JZWXS+RiJFxvbINGxozk0xKzJe97tbkalRI3FoWdgVGTcHH0zbkSMTv3fZssJPGPaah4cONWWqk8VdURAw19Dw4eQtSLWaWRAVGfaClJYmJ8xxuFUqVFZ6Gwp5X6EpMuyRiRcREo/Wrema277d5AS6Q5mTwW/z0kIhndAyv0Vf/LJtG/CLX9Dj1q1DnuyvlKqhlPoBwFoAYwEsBFChteZldzmAgEQwFgaJFBl2Caab7M+L2N693u8NqiKjFNXnX7/e34Lk9jhMmEBVliorU1ugbI+MLcgXUkUkL3jR3LmTLILdugEnnZT894wbR4mE5eXORaZ27cShZW5PIyfhFrMiw/ex3bXdnWPgJ7RMFBknfL9NnWr2VVTQcfzTn0x8PUBzbaFfY14C7zPPOBvOJkPjxtFrDF9D3K8ilRA//p6gXYf8vyZrcGLhPhUqK72PMa+5qR7/bMFK244dqVnu2ViwcGFmFJlkQikLAZ7r0/HIZCr8c9s2Kqry/fdUvTa0OTIAoLXeo7U+GEA7AH0B7Of1Nq/PKqWuUEpNVkpNXue3vEwISORCZIEnU8n+gLeAWVFBLttCX+C9aN489RuyRQvqAQKkJiTEUmT8VFDLJ7VqmYm2fn3gxx+p/0eyHHMMNXlTyrkQd+qUnEcGIIUSyE5H6ELBy3DhVkqSCS0rNCtuvmBh0C5msmED8N//Ur7M0KFmf6qCWS5JJ5fAC6/wKb6GeG3YsiU5z/SqVdSM9Ouvg6fIsEcm2WIZ6XhkqqrCpchwzuq8eabSWzoRALEUwUKF76V858hoTfJh/fpAnz50LkMbWmajta4A8AWAwwE0Ukpx1GI7ACtjfGaY1rqP1rpP80KX8nJIIo9MooohfhUZu+mjlyKzeDHV1w9SwmaqPPecswzkUUfR8R03LvnvipXsH0vxLCTY6stlf9M99/bnO3aMv+Bz7pd9XZ93Hm3tvKdiI55HRkLLUoeFQVsQr6hweqKZIHlk+vfPzPfZluxrrqEtX0NsLd+7NznB6fnnqRnpkiXBuw5ZIE42J5TzRubOBS65JPF9ahMrtKxQFZnSUpMjk4rAy/2efvzRKDIVFamHcQdNkeE5J5XwrUwqMrt2kWHMNnqxIhOUfCO/+Kla1lwp1SjyuB6AAQDmApgA4IzI2y4GkEal8PCRbmiZH0WmspIm3kMPpedeiszChYUfN54prrrK9LAA6Bh26ECx88mydStNBrVrxy6kUKjYikymado0viLDAoCt8HFCcirleIOCH4+MhJYlT40a0d6pigrqSwQ4BcggeGQaNiTP5D/+kZnvswVAbqLJQpLdFyUZb8NPP5nHQbsO/XRd94K9FKedRiWw7Wpyidi2zdvTxopMoYV1l5XRPFNdnZow3qABRSbMmEFCfcuWNO8nKpawdy/lbbnDqmJ5tAoVVthSkQsymSPj5b1nRcbdmDno+PHItAYwQSk1A8D3AMZqrf8N4BYAf1BKLQDQFMBL2Rtm8cETaqqhZX409++/p5tqwAB67qXILFoUHkXGiyZNYldficfWrZkPA8kVvHBmIzyprCx+jkwsBb1Ro2D3RkmEH4+MVC1LDVZWeE6sqACWLqXH9nEKgiIDAPvtlznPkT1HcQhsdTWtC5s3mwpuySgythAftBBH7qvB1RL9wjky7j5QfojlCeTrttAMYXzNbNqU+v3StKkJ9+SS1+5mwG6++opaU1x8sdmndWyPVqHCSkIqkQ6ZzJHxKtXP11qxKTIJHUxa6xkADvHYvwiULyOkAF+wqYaWJVKEAGD0aFKEhgwBHnkkWsDcsoViV9lSF0aaNEktftdtZXvsseAkJGbDI/Pzz7TIDx8eXyhKpKAXK349MvFCVrZsAW64wfkZga7nNWvIuzpvHglM7GW157ygKDKZxLb224pMdTUZFdq2pVK5fksLa+1UZArNm5CIq6+mBqrJlpxmjwyvFcmUYt6+3fu64/Wj0MK6+Zxu2JD6/dK4sQlr79QJmDSJFCPuXO/miSeAb7+lx++8Y/Zv3UrXXBA9Ml5lqxORydAyL0UmzB4ZIQv4DS2LJfDVrEl/8S74sWOBX/7SVJVyKzJsISn0ksHZJFMemRtuoNjpIJANj0ynTpRQWF4eP3k4UchkseLlkeHrzrbMxvPIDBtmHodNII8HH78OHWhbUWHmNrvYSSyBspix5yie56urTVgZC5Z+PTJbttAfH3M7VDcIdOhAPYeSXfMaNHCG4iWryHh5ZFioTEXgzSZ8btetS0+R4fBOPx6ZG24ARoygx7bSws1Eg6TIZMIjky1F5q9/pVSDww9P//sLiQK7hcJDusn+AF30sS74PXtMfgxPTDwpMEGcJDKNV58FPwQ5tCybOTJlZaTExFrow+qR4fu4spJCh0aPptDPevVMlZ9EoWW2UFFoVtx8wtdz27YkFFZUGMHcFp7C6JGJFVp23330uG2kaYJfRYYVw3vuAV57zRTqKHbatHF67jOhyLDAW2iKDEcWLF+e+hpnRydweX23IrN7N3D88cAXXzj32+sSG1+DJKOkc14zmSPjpcgccggwZUrwPKmJKLLaBcEhkebtx3IdT5FZtowW7n32MZOKeyIRRYY8MhUVNPkkM/EEWZHJZo4MW3iXLjVNR23C6pEpKSHlbdo0quZz0010Hvr0MRVkEoWWBSlOPJfw9dyoESXLV1QYQZPnPK3DqcjECi17/nl6nGyODCsyHToAp5+emTEGgfbtnc8TVQu1iXXdsde60IwSdq5iuopMs2am2JDtHQWoYupnn1GpdBvbuBtEGSWd85rtHJlipcBsAeHBr0cmnuU6niLDlWX22Ye0/Dp1onNBgjhJZJomTUziq5uFC2my9SLIikw2PTI9etB23jzv18PqkQFIeeNeOR06UK5Bz57m9UShZcUW15wp+D5s2JCUGVuRYeFp9266z8OmyNhzVHk5Kc1Llph9nB+ZrCITlHzATOHO7UgmpyiWR4aL7PTqld7YMo19btNVZH7xi9iGVL4O3dde0BWZY4+l7cEHJ//ZmjXJ6JUJjwxfo2HIpxRFJk+wRTrd0LJYn3/3XfrsAQeQZaBp0+gQqiBOEpmGrU9r1kS/1q1b7EIIQVZk2FKTjSph3EOAFZnXXwduvtm8HlaPDEBCNCehN2tGhgUO7QESh5ZlqttzsdGiBW3btIlWZDZupGuO58mwKTJ2vwilyEPz/ff0/MknTUVLv4oMC6OiyPj73O7dZIDwUmSGDKEE92QLD2SbTCgyvMb07m36rLkVmUWLnM/ffBO48UYS4tmrwTJKkLzR555L+UGp9EVTiq6VTMz1XDWOcweLGVFk8kRJCQlzsVzUiZL9gdgeGa7HftFFRlj1qs4ligxZjADgP/9J7nNBVmS4L+2ll2b+uxs2JMFy3jwSIi+8EHj0UfN62D0yDCs0bdqYfYlCy0SR8eaBB6jL/G9/61RklCKBaN268CoybsrKgOnT6fFJJ/lvrMyIR4bwe7x4fY513R12WOGFlmVCkWnZkraDB5MyXVYWHVrGgjbTvj3dv1qbeTCIOTKAWWNTgQvmpMvcuXQ8+VwUM6LI5JEOHSh0aebM6Nf4Qo43kcRSZJYvJ6GHhXTAuzqXKDLAvvtSN/rPP0/uc0FWZO66C5g9m5LOs0GPHqTIcPlNwFjYwqzI2MIMe6xsj0yi0DJWZJ54IvNjCzINGwJHHEHHr1EjChPdsgXo2pVeX71aFBmmtNTci40bm7ATUWTi457r/Xpk+LrLVF+gXFCrlvl/U13jrr8emDyZqqYCdL24PTLc64lp1iw65D6MMkqDBtGFmVJh7lySbwpNUc4GosjkkR49qGZ6r16mhjrDORsNG8b+fCxFhoUkzlcA4oeWBcltm2mUojwi96SaiCArMvXreyfiZ4oePSj/w/YA8sLEAlNQj1062B6ZFStoa3tkEoWWbdtGwvq112ZnfMVAw4Y0z1VXmxwEW5EJkkCZDezk//JyElpr1vTv7du0iebMMK8ZgOlv8tJL5PGLBa/PQbvuWMFNtbpVrVoUVsY0ahTtkeHmpEybNrEVmTBdb+Xl8ZtK+2Xp0vA0OxdFJo9w2VXACDYM38DJKjJaA1deSY85XwGIHVpWWlp4nYVzTevW0ZNqvGQ7rUmoDEM1kFTo0YOER7tQAl/PPEGHaWFivLwBbo9MotCyevXCYWFLlUaNoru3r1kjHhmGBdN69Uz+TP36/j0yFRW0JhVayeBcwiHhCxZQOOOBB5KQ7pWvmii0rFBhRSZTBicO+bRZudJ4TQFaE7wUmTp1gnf80iFTHpnKyvB4skI8HeUf2yrudlWzR7VesaoAABhfSURBVCbeheilyFRUUBJd48ZOa2+s0LKwXOjxaNOGBG+7KlS85l3V1aTMhNGr4Ae2Ak2ZYvbxxMwhk8VWx94P7gIHdeo4Q3T8hJaFoQJNOnCMPeDtkQmTQMScey7lqgHmvrMNCfXq+VdktmwJ75rx1ltAv36Uc7B1q6m6tWYNra8nnBD9maB6ZDi3NlNrnFdo2YoVVPAAoFwaILotRWVl+IxemVJk7Ma1xY70kckjRx1lHrurZm3eTItOPG9JnTrRigy7uZ9+2mm5bdiQFnO7pv3mzeFdlGxatyZL+IYNJkkvniIT5vAoP3C/CrsEs3hknBWkAPLG2PdozZrxPYGiyCSGKyQBdB3WqkXXXJgVmTfeMI9ZkbHn/fr1/YeWVVWF0wgBAGefTX/77edUZJjx46M/E9SQRp7DM+mR+eEH87yqigTtjh0pBKppU9rv9shUVYVPRslEaNnu3bSWhGWdFY9MHrFDv9audb62eXP8sDLA2yOzfj1teSJi+LvsfinikSFat6btiBFmXzxFhr1nosh4w9ce9zICzMQcZo+MW1i0PaaAv9AyCWeMj63IlJbSdRZ2RcbGS5FJxiNTVRUe4SgWpaV0vPzkVQY1tIwVC66emi7uHBkO/2zThqqV8bzGx+mMM4Bx48Ipo2TCIxM2g6EoMnlEKWDVKnK7enlkEt3AoshkhoMOou2TT1K435gx0YmJNqLIxIevPfuadntkwnjsuDohe/3s/BggcWjZtm3ikUkE95QB6DrkUqZBFSgzjVdoWTI5MqLI0PFij4ytOAMmrJEJamgZz01+q7MlolEjunY4fNtWZGz4/pwzh3ochVFGEUUmeUSRyTOtWpFnZtUq537xyOSO7t2B3/2OkjfnzqWmXHaiuttKzt6aROcnrLA1DzA17O0cmfr1w1lg4oorSJA891x6zsnojJ+GmKLIxMcO1z3iCOOR4fzAsJUNdhNLkUkmtCwswlEsSkvpOHz/PdC3r/M1dwWzoIaW3XUXcN11wAUXZOb7+Jphj3wsRcZ9nDZtCt/11qABtSmI1ezcD6LICDmnSxdnzw0g84qMV3ddUWQMdtOoPXuAL74wz93HmD0NYWg0lQq1a5sJ9PDDaWt7ZMIyubp54QVamH//e+DBB4E//tH5ut+qZUJsSkuBe+6h41unDgnuW7ZQwj9AhqMww4qc7XFONtk/jGGhNvXrUyGTOXPM/MbYBjAguJ7Ahg2pX1WmQllZzmAB26v8PBB9nKZPB/r0ycwYggKvj+l4ZUSRcaGUaq+UmqCUmquUmq2Uui6yv4lSaqxSan5kG3JbV+r06EEToK2B+1VkuIIWs24d7XdPQOKRiY8dkjJ3LvD22+a5KDLJw7HVRxxBW55Yw5wsrBR5Xbp3B267Lfr+ThRaFuTeRbnkz3+m4wuYxNlVqyjxn6sxhZUzz4zeJ6FlyWHfgyecAHz4oTEcurvV83ob9mPGcgYL5ytX0jFxHxcvhe/ss7M7tkKDc6ft4gjJErZcVD8emd0AbtRa7wfgcABDlVL7A7gVwDitdXcA4yLPhRTo0YOUEdsr40eRKS+nz9lhAatXk4Dt7jXhVmS0FkXGxksp6d+ftm5FhgszcK6DEA1b0Y45hgR0O7Qs7It6LBKFlvmZEwQntkemVSvpwdOhA1XXsiuZ1asnoWXJwIrMBRcAhx0GnHwytTwAjEfm55+Bd98VTyDjpci4vTGAtyJjF0UKA7/6FRldxo5N/TvC5pFJWH5Za70KwKrI4yql1FwAbQGcAuCYyNteAfAFgFuyMsoiZ7/9aDtzpukt40dosScH9sAsW0ZVQNy4FZlt2yjxThQZwh2K99prpOyNG+ftkWnShCYbwZvPP6fwgS5daDKV0LLEJAotE0Umedgjw4qMABx7rPO5X4/Mzp30F/b7l9daOxewvJyesyJz4IHkQb3mGlpjw15tkK8ZFrCTUWTcZeuLndJSUt5YOU6FsCkySeXIKKU6ATgEwCQALSNKDis7LWJ/UojHAQdQONikSfR8504Snv14ZABnzfHly4F27aLf26ABWSM5NpoFS1FkCPek2qJFdHMuZvFiCStLRJ06piGhXYVl3TqnACAY7NCyvn2Bnj3Na+xBFUUmOdgjs3KlKDKx8KvIXHUVbcNYqMOG1wX3vdi+vcn94Gpfixeb8v5hxq9HJiyCdyK4xHeqhE2+863IKKXKALwH4Hqtte80JKXUFUqpyUqpyevcJT0EAJQc3acP8O239Jy9Jsl4ZAASdpYv9/bI1KhBwjlXR2OFRgQjYt99gfnzTff1li2disymTeSdmTkT+PRTp5ApxKdBA6Nsx7o+BWdo2fffUzIxs3UrvSb3a3KwN3DhQqBr13yPpjBhZY9L48bizTdpG/a5j+9Rd4Wtxo2jy/bPmycKNOCUVbSOrci0bg189FH8Pm5hgEt8pwofP3d58GLFlyKjlKoFUmKGa63fj+xeo5RqHXm9NYC1Xp/VWg/TWvfRWvdpLkkFMenWjcLCAP+KjNsjs24dFQyIJSi2bWssRnPn0rZHj9THXGx060beMMDpkdmxg47TgAHAd9/RPk4mFhLDHpmqKrq2vTyGQvzQMr9zguCkrIzu323bZK6LRfPmJJx7CY8VFaZwR4cOwEknAb/+dW7HV2jwGuEOLWZFxk74nz9fFBnAWbVs40a6J70UGYCur4YNqYgCyylhI12PzKZN9B1hCX/3U7VMAXgJwFyt9WPWSx8BuDjy+GIAH2Z+eOGheXNKItc6eUWGPTJ807v7UzC2IjN9OlBSYnJyBKJ3b9o2b24UmfXrTVnrefNoK14F/7Aiw9eeKDLeeFUt+/JLWtBEkUmN884zj8OWNOwXrti41mWK3L2bhPMrrqDnq1YBnTvndmyFCCt27L1nGjUCZs+m5H9m715ZKwCnrBKrh4ybk0+mSIkwkgmPTFi8MYA/j8yRAC4EcJxS6ofI3wkAHgIwUCk1H8DAyHMhRVq0IEsPW60B/6Fl7JGZPp22Bx/s/X5bkZkxA9hnH+lL4eazz0h4rFnTlC6cPdu8Pm8e5RqFvbFeMnB4z/Ll9FwUGW84tMy2xB19NOUmiCKTGt26AU8+SY8PPDC/YylUWJHhsvIMryf/+hddk5WVku8BxPfIAMCsWc79ovyR0le3LgnYfJ2Jpyo26XpkRJFxobX+r9Zaaa17aa0Pjvx9qrXeoLXur7XuHtluzMWAixWOulu7NvXQshkz6HtiTRBt2wIbNlCpzRkzgF690h93sdGsGZU/BEzPCV7QAQoVaNQofJVU0oE9MqLIxIdDy9j7x0ydKopMOlx7LeW5SYEOb/i4uD0yEyfStkMHk1spiozpdj9woHO/bdyyleZOnbI+pEDQuDHw+uvmuEnRl9hkwiMTJmNrUlXLhOzBisy6dakn+y9cSOETsXolsAD53nsUxyuKTHy8FJl586JLNQvx4WR/VmQShRSEFQ4tcysyJSXA0qX0OEyLUyYJWmf1XBLLI8Oe6D17gOHD6fE+++RuXIXKr35FIeDuEG67UaYdXiYeGaJJE9NXh58L3mQiR0Y8MkLO4cVk3Tr6AxILzPXrk5DjtyJU27a0vfBC2h50UOrjDQNlZeR5+ekns2/XLlFkkqVhQ7pGlyxx5h4JTji0zK3IKAW8+CJVixJBUsg0zZrROvLBB6YPCkDeZ4DCke+6Czj3XOCII/IyxEDA3dQvuQT43e/MfvHIEG7FRTwyseGS6IkqCcZCQsuEvMA39fr1ZBmrWzdxTXWlTP5BvNLLDCsyAOXRHH98+uMuZuxcGDscSibg5OAFbOZMCSuLR6zQMqWoFPOgQdKZXsg8NWoAHTsC48cDhxxi9i9Y4OwZM2xY7scWJDiku39/59oteaiE25vsLpYgGNi75+5h54ddu8gYLoqMkHN44tu6lWKVW7b0J7Rw5+oNG+iijyco2orME0+QFU6ID4ftDR1qJt4OHfI3niDCiszUqaLIxINDy9zttlasoHlBcjyEbMFeUi7BvGEDXXd/+IN5Dxc/Eby57DJg1Cjg/POlsaMXIm/4p3592qaSJ/P+++QdHDw4s2MqZOTSKhBYA9+yhTwyfoUWTqTmHjTxBEW7y+vhh6c2zrCxYwdtTzzRTC7SjyI52IO1axdw1FH5HUshwwUkOOn62GNpy4qNKDJCtmDLL1+Db79N2/POA0aOpAatQnxKSoAhQ8gAWasWebeefz7foyocNko5KN+wPLh2LRWO+Ppr/5+dNYuuwRNOyM7YChFRZAqEunVpIty6NTlFhj0yfhpcKkWemDFjJPk1WXr2NJ4FUWSSww7Fu/rq/I2j0OEwntWrKZdo/Hi6XxnOoxOETHPHHbTl0ORvviGj2EEHAaecAvTpk7+xBZWpU4H/+798j6JwEEXGP2w0nTCBFJNbb/X3uccfB+6/n+S7MHnAQvSvFjZKkRbOioxfoYUrQs2aRda0RMnA110XXTZSiM3nnwOvvkqTwrnn0j5JuE4OO8mTJ2ghGo6lX7rUFJSwwxjFIyNki8suAy6+mKpZXnkl8NprFIosOVlCprj3XgpP/PxzYMqUfI+msOHomYULnc8TwaGgYTNUSzeMAqK0lJSSdeuS88isXEmJ1PvuKwl0maZ/f/P43nup4luXLvkbTxCR4gj+4ByExYtNv46OHc3r4pERsgknB3NSvzQsFDLJaaeZCqtCfLhk9/jxtE025ypsiox4ZAqIsjLKddmzJ7kcmaoqalgmSejZpaREvDGpwP2QwuTqTgVbkfHyyEgzQiGbuPuWiQFCEPIDKzIzZtB2587kPi+KjJA3SkuNKzEZj0xlJYWkSVUZoRBRCnjlFQp/FGLDVrfdu02p0qZNKYH49ddNIrYgZANuwMd9T7jQiSAIucWtiLhL8icibEZDWRoLiLIyo4En65EpK3N2FhaEQuKii/I9gsLHNkSwIqMUlXQVhGzDpeZ79iSvIDd4FAQh9xx9NKUZ7Luvsym3m717gcmTgb59zT6tsz++QiJkelthU1pqLkC/8fDl5aYbuCgyghBcbEUmTM3MhMLgzjsp6f+ll0gouu++fI9IEMLLhAlk2G7WLLq3mM3ddwOHHQZMm2b2hU2REY9MAWErIn49MizwbN8uoWWCEGREkRHySbt2pMQAwKRJ+R2LIIQdpagkf8eO1E+mstK7ehn3fFq0KLfjKyTEI1NAsCBTu7YJLUmEnQwsHhlBCC5eoWWCIAhCeOnVi7YzZ3q/znnVc+aYfWHzyIgiU0CwINO2rf9kLU7MBESREYQgY5fYFI+MIAiCcNBBtJ0+Pfq1XbuoOAxA/XkYUWRcKKVeVkqtVUrNsvY1UUqNVUrNj2zFfpgBuGdEMiFitkdGQssEIbjYhghRZARBEIR27WhtmD8/+rXFi83jL780j0WRieZfAIa49t0KYJzWujuAcZHnQprstx9tk2kaZQs/4pERhOBSq5Z5LKFlgiAIglJksF6yJPo1L+UmjCRUZLTWXwLY6Np9CoBXIo9fAXBqhscVSlJRZABTc1wUGUEINvvvT94YbogmCIIghJuOHakE85gxTm/LggW0ZdnvzTdzP7ZCINUcmZZa61UAENnGLBaslLpCKTVZKTV5XbwacgK6dAEGDQLeeiu5z7VvT1sJLROEYDN1KpXarFcv3yMRBEEQCoGOHSmZf/Bg4OWXad+kSRRaVrMm8OOPwOrVQI8e9FrYQsuyXn5Zaz0MwDAA6NOnT8gOb3LUqAGMHp3859q1I818z57Mj0kQhNzh7ugsCIIghBu7qNMbb1BI2cMP0/PWrUkGBIANG2gbNkUmVY/MGqVUawCIbNdmbkhCslx7LW27ds3vOARBEARBEITMwZ4WABg/3igxgLN5et26uRtTIZGqIvMRgIsjjy8G8GFmhiOkwmmnkTdGFBlBEARBEITiYd99Y7/mpciIR8aFUupNAN8A2EcptVwpdTmAhwAMVErNBzAw8lzII377zgiCIAiCIAjBoFs35/MDDgAef5web9li9teMJIuETZFJmCOjtT43xkv9MzwWQRAEQRAEQRAi1K4N3HEH5cq89x7w7LO0/4YbgG3bzPu44NP55+d8iHkl68n+giAIgiAIgiCkxn330fbyy82+554DjjnGPK9fH9i0CSgvz+nQ8o4oMoIgCIIgCIIQIK66Knpfo0a5H0e+kcwKQRAEQRAEQRAChygygiAIgiAIgiAEDlFkBEEQBEEQBEEIHErnsE6bUmodgCU5+8H4NAOwPt+DEJJGzlswkfMWPOScBRM5b8FEzlvwkHOWXTpqrZsnelNOFZlCQik1WWvdJ9/jEJJDzlswkfMWPOScBRM5b8FEzlvwkHNWGEhomSAIgiAIgiAIgUMUGUEQBEEQBEEQAkeYFZlh+R6AkBJy3oKJnLfgIecsmMh5CyZy3oKHnLMCILQ5MoIgCIIgCIIgBJcwe2QEQRAEQRAEQQgooVRklFJDlFI/KaUWKKVuzfd4BEIp1V4pNUEpNVcpNVspdV1kfxOl1Fil1PzItnFkv1JKPRU5jzOUUofm9z8IN0qpGkqpaUqpf0eed1ZKTYqctxFKqdqR/XUizxdEXu+Uz3GHGaVUI6XUu0qpHyP33RFyvxU2SqkbIvPjLKXUm0qpunKvFR5KqZeVUmuVUrOsfUnfW0qpiyPvn6+Uujgf/0uYiHHeHo3MkTOUUh8opRpZr90WOW8/KaUGW/tFzswRoVNklFI1ADwL4HgA+wM4Vym1f35HJUTYDeBGrfV+AA4HMDRybm4FME5r3R3AuMhzgM5h98jfFQCez/2QBYvrAMy1nj8M4PHIedsE4PLI/ssBbNJadwPweOR9Qn54EsBnWut9ARwEOn9yvxUoSqm2AK4F0EdrfQCAGgDOgdxrhci/AAxx7Uvq3lJKNQFwF4DDAPQFcBcrP0LW+Beiz9tYAAdorXsBmAfgNgCIyCfnAOgZ+cxzEYOeyJk5JHSKDGgyWKC1XqS13gngLQCn5HlMAgCt9Sqt9dTI4yqQUNUWdH5eibztFQCnRh6fAuBVTXwLoJFSqnWOhy0AUEq1A3AigBcjzxWA4wC8G3mL+7zx+XwXQP/I+4UcopRqAOAoAC8BgNZ6p9a6AnK/FTo1AdRTStUEUB/AKsi9VnBorb8EsNG1O9l7azCAsVrrjVrrTSCB2i1kCxnE67xprcdorXdHnn4LoF3k8SkA3tJa79Ba/wxgAUjGFDkzh4RRkWkLYJn1fHlkn1BAREIgDgEwCUBLrfUqgJQdAC0ib5NzWTg8AeBmAHsjz5sCqLAmf/vc/O+8RV7fHHm/kFu6AFgH4J+RkMAXlVKlkPutYNFarwDwVwBLQQrMZgBTIPdaUEj23pJ7rvC4DMCoyGM5bwVAGBUZL2uUlG4rIJRSZQDeA3C91roy3ls99sm5zDFKqZMArNVaT7F3e7xV+3hNyB01ARwK4Hmt9SEAtsKEungh5y3PRMKKTgHQGUAbAKWg8BU3cq8Fi1jnSc5fAaGUuh0UAj+cd3m8Tc5bjgmjIrMcQHvreTsAK/M0FsGFUqoWSIkZrrV+P7J7DYewRLZrI/vlXBYGRwI4WSm1GORCPw7koWkUCX8BnOfmf+ct8npDRIdgCNlnOYDlWutJkefvghQbud8KlwEAftZar9Na7wLwPoB+kHstKCR7b8k9VyBECi2cBOB8bfqWyHkrAMKoyHwPoHukykttUKLWR3kek4D/5VW8BGCu1vox66WPAHC1losBfGjtvyhS8eVwAJvZbS/kDq31bVrrdlrrTqD7abzW+nwAEwCcEXmb+7zx+Twj8n6xVuUYrfVqAMuUUvtEdvUHMAdyvxUySwEcrpSqH5kv+ZzJvRYMkr23RgMYpJRqHPHGDYrsE3KIUmoIgFsAnKy13ma99BGAcyLVATuDijV8B5Ezc0ooG2IqpU4AWYxrAHhZa/1AnockAFBK/RLARAAzYXIt/gTKk3kbQAfQQn6m1npjZCF/BpT8uA3ApVrryTkfuPA/lFLHALhJa32SUqoLyEPTBMA0ABdorXcopeoCeA2UA7URwDla60X5GnOYUUodDCrQUBvAIgCXggxccr8VKEqpewCcDQpxmQbgt6D4e7nXCgil1JsAjgHQDMAaUPWxkUjy3lJKXQZaBwHgAa31P3P5f4SNGOftNgB1AGyIvO1brfX/Rd5/OyhvZjcoHH5UZL/ImTkilIqMIAiCIAiCIAjBJoyhZYIgCIIgCIIgBBxRZARBEARBEARBCByiyAiCIAiCIAiCEDhEkREEQRAEQRAEIXCIIiMIgiAIgiAIQuAQRUYQBEEQBEEQhMAhiowgCIIgCIIgCIFDFBlBEARBEARBEALH/wOi7AHQN2wZYQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig = plt.figure(figsize = (14, 2))\n", "plt.plot(data, color = \"b\")\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "def compute_mape(var, var_hat):\n", " return np.sum(np.abs(var - var_hat) / var) / var.shape[0]\n", "\n", "def compute_rmse(var, var_hat):\n", " return np.sqrt(np.sum((var - var_hat) ** 2) / var.shape[0])" ] }, { "cell_type": "code", "execution_count": 27, "metadata": { "scrolled": true }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAzIAAACPCAYAAAAhie4bAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsnXd4FNXXx7+TRkISWhJ6Cb2EDtK7ICCIoKIoiICIgg0bIOJrARVFxR8WEARFQEVFaYIU6Sgt9A4JNQQSUkghfc/7x9m7M7PZZHtL5vM8eXZmdmb2zmTm3nu6RETQ0NDQ0NDQ0NDQ0NDwJnzc3QANDQ0NDQ0NDQ0NDQ1r0QQZDQ0NDQ0NDQ0NDQ2vQxNkNDQ0NDQ0NDQ0NDS8Dk2Q0dDQ0NDQ0NDQ0NDwOjRBRkNDQ0NDQ0NDQ0PD69AEGQ0NDQ0NDQ0NDQ0Nr0MTZDQ0NDQ0NDQ0NDQ0vA5NkNHQ0NDQ0NDQ0NDQ8Do0QUZDQ0NDQ0NDQ0NDw+vws2QnSZIuA0gHUAAgn4jaS5JUCcBKAJEALgN4lIhSijtPeHg4RUZG2tFcDQ0NDQ0NDQ0NDY2STHR09G0iijC3n0WCjJ7eRHRbsT4NwD9ENFuSpGn69anFnSAyMhKHDh2y4ic1NDQ0NDQ0NDQ0NEoTkiRdsWQ/e1zLHgSwVL+8FMBQO86loaGhoaGhoaGhoaFhMZYKMgRgsyRJ0ZIkTdBvq0JE8QCg/6xs6kBJkiZIknRIkqRDiYmJ9rdYQ0NDQ0NDQ0NDQ6PUY6kg05WI2gIYCOB5SZJ6WPoDRLSQiNoTUfuICLOubhoaGhoaGhoaGhZQUAC8+CJw6ZK7W6Kh4R4sEmSI6Ib+MwHAnwA6ALglSVI1ANB/JjirkRoaGhoaGhoaGmoOHAC++goYPdrdLdHQcA9mBRlJkoIlSQoVywDuA3ASwFoAT+l3ewrAGmc1UsO7IAIWLgRSU93dEo3SwK1bQFYWayQlCdi5s/A+r7wCNGjg+rZpaGhoOBMi/szLc287NDTchSUWmSoA9kiSdAzAAQB/EdHfAGYD6CdJ0gUA/fTrXsvFi8CXX/IyEbBtG5tsBcplDTWiIwWA7duBOXOAZ58Fnn9e3p6U5Pp2aZR80tOBqlWBCRP42QNYiDbmiy+AmBjXts2bIAJ27QJOngR++cXdrdHQ0DDHwoXA5s1Afr67W+KdJCcDcXHuboWGIzAryBBRLBG10v9FEdEH+u1JRHQvETXUfyY7v7nOY9Ag4KWXgIQEYNUq4N57ga+/BlavZi2vn5/mg2qKsWMBH/1TRAT06QNM1Sfhjo3lz02bgPBwYMUK4NQp97TT2yAC1q0Drl+Xt23aBLzxhvva5InMmsWfq1YBubm8/NNPwOLFpvdXai2vXwfOnHFu+7yFX38FevYEWrQAHn9crZzQ0NDwPJ59FujfH8jMdHdLvIusLODwYWDwYKBmTUCrCOL92JN+uUQRH8+f69YBw4fz8v79ar/Tw4dd3y5P54cf+LNHD+C999TfCSvWyZP8OWoU0KmTy5rm1SxeDAwZor6no0YBn34KHDvmvnZ5CjodC3WffMLrZcsC167J30+bBvz+O3D0qPq4FEXJ3okTgUcecX5bvYEdO9TrKcWWNtbQcAyawGwbyvt2//3y8vPPA6++qnmQFEV+PtC3L9CuHfDff7ztnnuAZcvc2y4N+9AEGT1lyvDnt9/K23bsYNcVgXCPWrEC6NzZZU3zCnbvLizIHD0K3HcfsHatvC0jw7Xt8hbi4zmO46efgN9+A555Rt4uaNKEPzdscH37PI1Vq1ioA9iampSk1qzducMKiTZtgPPn5e3JersxEQ9ksbGlbzIlLFcCIg4YVhIfz24XgwYBN2+6rm3ezrx5wOXL7m6F9xAZCTz9NHs7EAHZ2YWVDxqFMWWF2b8f+OYbYO5cdpXXKMySJcC//xberpz3aXgfmiCjR0xmDh4EQkLYzezGDfU+Bw4AtWuzZnzfPiAnx/Xt9BTu3DEtlAwbxpPxJ59kN54tW9j3HgCaNeNPLShRJiGB/2bN4jiOkSOBRx+Vv796VV6WJP4s7dryrVvV92jkSP7cvFnepnzGGjeWlzt0YGvOpUss/GRnl674rYMHWWmzY4fc5/37b2Fr840bwJo1LDRXq8bvu0bxJCYCL7/MwrOGeWJjuX9bsgSoV48nk5Mn8/3TYheKx1xJvrQ017TD29i5k+dwHTqot9+9y586Hc/tNLyLUivIiMlMzZqs3VVOZpo0ASqbKO+5eLHafaU0dxYVKvAE0d+f14cNA06cAJYvZx/7H39U71+pEjBpEi8nG0VTpaaWzvijO3eAKlWA+vXl+yioWJHdGi9fliec4hktzc8dAJw9q16/9155+dgxnoAXRXo6x8UcPChvU8YhlXQ2beLP/xsUjeeDluCHH4AFC4DAQFnRALDbXVaWvP7zzy5tpteRlAR07MjLqamlW8llKSI5h2D3brYqAIXf8dLOr7+yUuHIEX4Xb98ufv/SPkYUxalTQFQUW5qVZGay8nr4cPa22bjRPe3TsI1SK8j4+rIWKC5O9pMXwku7dvJy3bpFn6O0dxY3brDm+5NPgD/+AJo351gFgbAgAEBoKAszQGFBpkUL/l+UNoQLSkYG8L//qb976CGgdWueePfsCQwYILuZlXbtuJ+fvHzlilrp0LIlP0/FsXevWpBRKidKOufO8eecu5PwTc7T+G3sX1i+nBURmzbJLmYxMcDrr/O9Dg1lq6qW9c00BQXAiBFqZcyWLRzDpQk0anbtkrPiGVtdfHzkMUJLCqPmscd48t22LfDEE0ULMkIhVtrnJqYoKGABOSoKmDFDFpoBdj/+8kuexwDA8ePuaaOGbZRaQUZJWBhPwJcs4VSu06YBERH8nZ8fu/6Y8kmdN092m9q3j/1TSwMPPaRer1LF9H5Kl54rV/g+AxycqHQDEhrx0hSrMHSo6TgrEcg/cCAPXqGhrKnctEl2KSvtg5Qybi08nD8XLZJjsWrXBqpX58nl00/ztqgo+Zh9+zgOqXFjnjwphZqSTnQ0f/qDfe864AAaNeL7V7MmB76+9hrvUw8xqJCfiG5dCQd/voAmDfJK1Tuq5Pvv5XtnzJ9/srujkgce4Biu335zftu8iZ492WIPFHbdliQ5SF05yQR4e2l99ky5cCtda5WIMXf1arbeaMgkJrJioW5d7vc7dJDnecYYK1s1PJtSKcgo/efbteMXfscO1njEx3MAotDy+vnxw162LPDPP+pJ+7x53DEDPClV1k0pqRDxwK2kYUPT+27cCHz+ubwutG2XL3PaSCK1+4pyuaSzZo3p623Zkgeuhx/myXi7doX3Kc0WGSL1hFJYAMeP58kjwNbWK1fYBaNGDd4WFcXCc8WKPCm9ehWYPp3f31WrXHsN7mL3bnarCwgAaoHNUOMfuIXdu4HgYHm/MWOAyriFGDRAIipjw98+uIBGmIwvkHH0YqnL2EEEjBsHjG1/nAeLlSuB3btxdeJHCA1ldxTlhEhpMSzN72px5OQUFmR8fGSr8x9/cKFbgZ+frJQobZhKEV/UK1i/Pn/++CNbWTVkhAKwXDl528WLPBaIcUIgSkdoeAelUpBRxsO0bAnUqsWaSCViYBfZzACukbJlS/HnLulaI1Mm7aZNTe8bGcmZuD77jFPhCouMoEsXtStaaQm6LsqVSViplJNKUxqjtDR+zohYc1SasiT98gvPI80hJpNCeA4NBfr1Y2sNwFq5J59krdz580D79vyclkRyc9lNsUcPXn/s/nREgF/kGn4J6njA3Fw0ffk+rMDIQueZgykIbduQJZ1SREYGMAHf4jhasW/PiBFAjx6ovWA6pAyeHfXpI++vVOyUxtg/S3jnHXbxVELEgkzDhpyEo2pVVjrodPz999+7vp2ewO7d6nU/5GEVHsLHDRYB4L5LCDBKV/grV1zUQA/n3Dn2+hCW/NBQ+bty5Xj+dz0mB++Ol4MllZkuNTyfUinIKCfjDRqY3icqios9rlih3t6iReECSkrhpaQXpxJ+zaNGydsqVCj+mFdflS0MSoyzg5QGcy4RF+IyRb9+hbfNmcOatdu3OX5h5EgWZB58kDWYYWHyBLU0oHSXEIVXi0NUvRbCofgcM4ZdWcLCeJ/oaI4JKYmsXq2uPTQ86rS8kpCg3nnDBvhu24K++AeJoXXZDH3kCJJb95b3WbUKn756o8T3dYKbV3PxLZ4DAOxHB9yqLAdhNQFHpUc1I8zs9BfKIBu1K6ZjENYDAC5ccH17PRVRq6MMsrHs4zgkJ/MkUvDjjzzZFEkTALaEKVPjbtvmmrZ6EsZzkPfwDh7Cn5hycQIikICJI1JwsdkQxA57DXM2NEMvcBYFpRK2tELEyZvq1zctyBh49lm8810tlEE2AH5vhQCt4fmUekFG2ZEq8fPjmBllJp+ijlEGwSr990siIp7lhRfYr14UJLQEZcc6YkTh2ITSIMhs28aBhMKULSwHpp4zAKhTh90swsLYalipEmva1q2T9yktweoxMbLrSZUqwOzZ5o/x9eVPYXUQg5OIrRGfJRmh+a5ZUz9RLOAN8Q2785c//MAPFJEq3WDEkC5samjdGnde+j8swni8W2cJAODXudfx00+uvhL3kHGQfXs+wHR0wn4MK/gdCR047VEzsFA44MI8zNg3GJnTZuHdy2OwHg+gI/YhPp5vbWlymzXFmTOchbERziEbQYhDTQCEWbMK79u1q3pdZLsEOENhaUqgcPYsp0YXXg9RtdMxDbORUpa1gp2wD+12fg6sW4e6f34Ov/NnsBpD0Q+bkZPD7/y1a2xtLuneIqYQFtHcXNnNU+laZuDXXwEADzQ4izp1+H3VUoB7EUTksr927dqRJ/Dbb+yYM2wYUU6O9cdnZQnHnsJ/5845vr2exPz5fJ3Xr9t2vLhPubm8fvYs0Y8/8rbffnNcOz2RnByiSZP4Wpcv58/27Yk2bSK6edOyc4waZfq5y8kh2rqV6Px5516Duzh0SL7WRo2IUlIsO+7uXaIZM4gyM3n9nnv4HCtX8vqaNer7WNK4fJmoSROiunX1G27fJmrXjvIj65GubVv1xf/6K5G/P/8BREeOGM6TkMCbumAPEUD9sIleftk91+RSLlyg+DYDiQBqilOGW+WLPEpHMC3BGAJ0lF2rQaGXMgHh9Jn0GtXHBRoyxN0X4nq++opo+HCiRx8l6tOHb8tRtDTcn8f7J5FOV7gv27+/6PEVINq1y91X5jySkohGjya6c4fXP/6Yr1l8ju8TQwRQwVffkM7Pj7Z2nE4FzaJUNygZFegQ2ha6bxcvuvfa3IGYWwBE06fz59mzJnasUoW/7NyZrg19gcKRQKtXEy1eTLRzJ1FBAZFO5/Lml3oAHCILZItSKcgsWMBXHhdn+zmK6mQPHXJcOz2R0aOJKlUiys+37fjDh4nmzVNvu36d792339rfPk9GCDGBgUTx8UTlylk/KG/aZPq5u3mz5E7GU1LU1zp+vO3nEnP3f/7h9T171OcuaYPV0KFEAcim+1rd5NG4XTtZkl63Tn3x/frx55YtRKmpqvMUFPBXzXCSCKD3MYO6dnXTRbmIu7fS6Fq5pkQAHUVLij2fR0lJ8u2aj2cp17cMxf/wt+o+XqrVjR7HCjpVprVhm4QCoitXiF591fbO08uoh4vUB1upPQ7QBgygKJxQ3acj3/Fg+fXX6sfw1k1dsYLMd9+5+cKcyIwZfI0ffsjrjz1GFBlJ9PvvvH3FKwd5Ye1aog4d5JsyaxbRoEFEGzbQnn7vUgEkqoyb1B07qRHOEkB08KB7r80dvPuufIseeYQ/016YRjRuHNHbbxOdOME71qypesj2lOlNjRvLm3x9iQYPdu+1lEY0QaYYDh4keucd26wxgs2biZo2LdzJbt/uqFZ6HjodUdWq3Lk6kpwcIkniTqckU6uWY4QNU8/d0aPy8pUrjmmvp3D8uHxtH3xAlJFh+7latVIrHM6eVd/HmTOJpkwpNI/3SnQ6VjR+hlf44oYPly9SUFBA9MYb6ptQhLmrbFmiaogz7PdGpZI5o9TpWHv7Mx4jAmgOXqMKZXMM8oe4TR+PPqm+b+fPE506Rf/7H68u6LbM8F0vbCNdS7014tgx916gC0hLI7qDUCKAtqMnEUB/4kEigEaC78v1L2QT/I4dfGsexJ9EAL325C2aOtW0IDN1qvuuy9m8/z5fYxXcpK8nnqCQEPYc0RXoaN1aHRVs1Guy9uxhoVjclJgY+STR0aSTJLoCecBZj/vp0JSV7rswNzFunHyLunQhaoDz6oepZk1+H5XbJk4kAmgSvqLySFF9peFaNEHGBbz+euFOtlmzkqfVFVy8qB+gFzj+3NWqEY0d6/jzehJKBZo9tGxZ+Ln76Sd5eedOx7TXE7h7V9ZGAkR5efad75ln+DyXLvF6YqJ87pAQeXnpUrub7lZ0OqJffuFr6Yy9lOIfLl/crVuFD/j2W/6uefMiz3n4MFEQMg3n+Vfqou7r7JEwPQi2euroJirTbnQlQEe9e8vfi9uo0xFrdQE2VetJSCD65huib7/OoymYTWkIobsIlA/cvdv1F+Vivv+eDNd7EfVUnVUl3CYCKOe+wdyZ6V/GvDyivAeGkdIkYdzPiYl9SWXOHL5OIfR9jsm0YfDXPLEYOVLu6E+fJrpwgejee036eac99bxpKbCU0a9PPlWuzJc+IeIPikYbXjlyhKhNG/m+BAcTHThAlJzMWiz99hSUpyqIL623z+1ogowLEP6Xe/cSnTkjvxNpae5umXOYNYsM2n9H06kTUYMGcixDSSRK78q8dq195zEObQCIeveWl3/+2THt9QTq1pWvyxHxZ1lZRPv2yesFBfzsLVtGNGKE/FszZtj/W+7k55/Vz8dzzXebFVTozh2W7Iph2FB1UMOdfaf5i2V664NJB3TvYvJkoqgA1tw+g28JIHrrLfn706c5Hs2ATmfSXWzVKr4l9XBRFRtCa9Y4/yLczJAhpHpO8uBLBNAV1CKA6AQUcR1KyeTpp3nbffcRLVtGVcLyDLvt3Su/o48+6rZLcyoffcTXF4tI04KI8I8ypYxQsn49EUA3ewyn/9DRcHyVijmlIlbm9m0i3QsvEgG0M+Jhg2B4QWrAJncivodDhxJ98QXRjRvqEzRvbrhn+3EPVcd1AljW0XAdDhdkAPgCOAJgvX69LoD9AC4AWAkgwNw5Spogo9NxMC2RHOcB2B4I76nk5hI9r1fwSJL9WnFT9Oql15pMcPy5PYVy5YheeMH+8wjLTlQUUY8ehce6zz6z/zc8BeV12eMKagm7dsm/5e0TpS+/VN+7du2ILQHXrtl13mHDyHDSJFSkU7Xuo8sx+fxwA0QLFzqk/e6kVSuiT5ouIQJobEcO8F+/3vrz7Nyp/h9sXRTLk6npSxzfaA+jfXuiNP+Khou/F1uoNQ5T+oHT9PbbRJ+XfUu+MWXKyBPzQYNUNy117hLDasLROPruU9nVpySGGr35JpGEAspEEM3Fy/Q8vqSsvoOJFi0iCgqS743IllMUeXlEc+fSyT0pVB3XaTP6EgE0ABvovfdccy3u4OpVDh0QSUmM/6IqJ1h2otu3afkHl+lkxa5EAP2N+zSrjBuwVJCxJv3yywCUNWY/BjCXiBoCSAFQ6uruShKnxwXUuclLWjXn1auBr7/m5S+/VFeudhTVqvHnuXOOO+fy5cAHHzjufPbwv/9x/RdlwTJbEak4V6wAduwApkxRf29cMbskULYsV6R3Jt2781A1cKD3F0QzTgPfqhWAbt04H6sdsA6L+R5jUf/aDvw2eKlcNltZ6McLyckBTp0CegfsBSpVgq5REwBAp07Wn8s4tXdMKlcEnv9hMvbvt7elns2tm4TAgkwsxWi0wWH8g75Ir98GIfc0xfvvA68kvMn507/6im/6U08BCxZwsSwF5Y/L1SAjWtfA2I/kaqOrV7vsclxGWhpQC9dQFlnoM6kpUp54AWU2rwPGjwf++kve0d+/+BP5+QGTJyOiYQXcQA0Mwl9IQASexuISXR+ldm0uUzAe3yETcrXt7VJvPIUfUK6+iQrTpggLw8jpdRD1+/sAgPZBpwxfKftADc/AIkFGkqSaAAYB+E6/LgHoA+B3/S5LAQx1RgO9hXLl5MKEYkx3FLm57quxkpmpnhR17+6c3/nyS/50ZF2PJ58EZswAbt1y3Dlt4fRpYPJkXr7nHvvP9/XXwKpVPDmVJODjj4HNm4H+/YHISK47s3Wr/b/jbrKz5eUIC8cfR1CvHnD0qLqIpKVs2VJoLuYWbt7kZyMtjSuDi/fLXoiA1XgQx9qOxUU0QBnkYuTlWSxdd+vGUoClLFoEfPSRYXXFCuC//xzTTls5c4YLpDZK3AN06YInRvlg6lSu42Qtyr4sOBhY/Gso8uCHMCTh7FnHtXn3bqBXLx4nPIG5c4G063fgr8tFbpNWOIo2uH6d66EYCA7mirbPP8+Vp//+G5g4EUhMlPfx9wd27sTChbKiyyfpNvbu4ZnkI4/IBW+9nSNHgJdeAlJTgZY4DgBo+XgUVqzg9xgA0Ls3FyV6912Lzyue22atArBKegQDsRHXzpfMokanTwMAoTe2oRv2YAv6oQP2Y8sPcXio/Db8iKesH3/79AHeew9hWXEYieUAuEZwSSQ723P6EKuxxGwDFljaAegFYD2AcAAXFd/XAnCyiGMnADgE4FDt2rVdZZFyCyKV68aNjovnvHZNznblDoyts5bW77CFrl1JFVRrKykp7AYr2uxubxe9uzIBzo+HFhl0AY7/sIe8PM5Qp4wpcSUxMfK1uLIOx+zZ8u8a+5Pn5RX/DojjjF2uXc1jjxE1bOj48wof/oMHiQZgg+GCTz32LqXc+xAHJRtz7hzvt22beru4WcnJdFKfAKxaNce32RpWf3CS/sZ93JjZs+06V16efImLF/PnTVSmBZhA99/vONeoBvoSNsePO+Z89iDSdD8LLjiW/dMq87WtlixRDzL/938ca/Xhh4bngzIyDN/nPfyoYVeRPddefvqJM5G6iy5d+HrCw4m+rvIe+3A7KNh2zx4OfctYxWnCB2EdzZ/Pt7Uk8csvRA/jN8NzMgWzDe+FeF6WLbPhxPrO6SSa0WR8Tj9ilMPb7m50Oo5HdcT8y5HAUa5lkiQNBpBARNHKzaZkoiIEpYVE1J6I2ke4Uq3qBsqX589PPmHLhbL6ui3odECtWnLldqWG2hWYMkFXqOC836tYEUhJsf88X3/N/wOBLZp1RxITw5//93+siHQmUVHy8qFD9p3r2DFg5Ur2anAkRKyMLeq86enAtm3A9eu83qkTK+9dRe3a8vLbb6u/mziRn9OCgsLHKd/P6tVZU+4ubt4EqlRx/HmnTOHnon17YPYvsp/kiyu7YfOhMFaL9uzJLkOCDRv486ef5G3KzqVSJYMF0WTVbRdABFCBDm0XPIP+2Mwbe/Wy65x+fkClSoCvLxsdIiKA2whHFdzChg3ATwsz7G84ZLfmq1eBvXvd5/oSGwvcey8vj8d3SKjdHmUefRANGxZ/HLp25c8ff+SLmD4daNwYaNeOtx8+DCQkGHb3W/Urzm+5AkB+tOzliSeA++5zzLlsQUyNWtzehieSvgQaNFD7q9tB165sHQwe1At5QaF4EGswcSIwa5ZDTu8xJCUBXbHXsB7UpQ0A2f0fsNEjIioKW/vPQRROYy5exZNYjvxzMXa21rPYsQO4dAnYvh2Ij3d3a6zHEteyrgCGSJJ0GcAvYJeyLwBUkCRJREvUBFACPfOtQwgywj3iwgX7zmc8AU9Ksu981mL8QD/+uHN/r2JFdukRY1ZqKrt6WEtqqrzcsSNw4oRj2mcrsbFASIhVHgE2oxz7lC7V1nDtGvDsszxZBWxzqymOy5dZ0Fy82PSka9QonhAJV5QlS4DKlR3bhuJQCjI//8zteZ9dpbF0KX/evs2fV6+yUDNkCNCli/o80dFwG1eusDDlaHx8gJYtebnVsHqG7QfQAbF39A/Krl3Aiy+yT8w33wBxcbw9KIg7xR9+AD77THVeEVd444YswLqSt98GJkX8hlrX/sProd+yf2DHjnafNzyc331J4vHhNJqhJY6jDQ5jxEsRwMyZdv+GEP5mzmTvvl9+sfuUNvH77zwhCkAOWuI4klr1YSnOHI0asQ/kk0+y5q5MGd6uFGSEf/D8+YCPD+pvWYDGjYH33gPu3nXK5biUjAygDLKxAfejQv5t4NVXHf8jZcrAf2A/jPNdio7Yp/LisxpXT0YsIPfyDTwBWVky9cfmiI5WK0fMCtVFcKTRY9iF7jgG7vyyvvnenqZ6DNnZPN6z0pPg70c2zbncjVlBhojeJKKaRBQJYASAbUQ0EsB2AI/od3sKwBqntdJLEIKM0Mzaa0Ex7mhc1XcsXMgaqhiF0qFRI7VC1RlUrMifrVvzZ/36QLNmPNn94AP2obeEs2d5zhQbC7Rty/7H7vT9jInhuAvJlB3TwUyaBAwYwPdu/Xrrj8/O5on8woXyNgcpBg0IIQBQC52CvXqlmlAI2BmfbjXi98T7vGIF8M47vFxWHz964wbHw9SpA3z6KVtfRZy70BLfuMHPoKu5dYuFRUfEYxVLQAAGtEtEKxxF43ahuA0TEu/zz7MwAwBr1nBHMnZsoQwV6cl5/JnOc1lXkpsLfP450DnlL9xEFexsON5hNy88XLbCli0LRKMd6uES5kqvwj8/m820WTbGLHz5JTBunGFVJBC4dMnORttIcjIn5EjecQIByEOTke0sP9hUJxMWBlStys/K7/pw3HvuAQYOhM9vK/HWWyzEXL1qQ2PXrjU8l5mZ8ua8PBvO5QDi4oCnsBSByMEfAxYCzz3nnB/q1Am+BXnYh85odGSlTaeI3xvLD7bIAORukpORkqRD2pdLUQWy5S6oXjW0bcvLL74IPPggK2JsoeWgWuiJXfjx1WP4A8MQtHQ+rl8jZGSwstUT4iJtYcYMHu9nzgR6VjqJnLBq6CNtd3ezrMbGfysAYCqAVyVJugggDMBixzTJewkJUa8rs0e9/DIwaBBrXixFBPiLeYArBJmrV4FXXmFN9Msv87YhQ7jfdzZighsfz9pscf0nTvALN2qUZec5fRp3eJ3fAAAgAElEQVQYPJgzhN1/P0+OBg1yrctFUhLPUfbu5clsvXrmj3EEVasCGzcCTz/Nisz4eKBvX/ZUsMRtT+kOKVwtFF4dDkGZuML43KtXy8/5r7+yNs3RgpQ56tRhYeTgQfX2ggJ5UhofLwd9/vYbf375JQ+WvXqxq8icOSxQutrVRwiAxhYiZ3AtKxzH0QrNmwPJqMQbH3lErXUQKnMx45w0CXjzTZ4IzZsHAMi7rU712LEjCxeu4O+/WZboUiYah9AeBWTPsKgmIkJ+foOCgP/QGQDQk3Yi3rcGf2H8oFnKd98Bv/yCxAT1A+Yu17LkZFZGBR/aCQCQunW1/6QiY4Kw4NWowQ/HpUuoX5UlkCtXFPvv329eg7hwIb+ozz8PQJ0MRrhxu5qEuDx8LE3DTvRAs0/HmT/AVoSWEMCMkyP45iUmsgm8GK3LZ5+xtwQATOzGyQgMgsycOcCffzqrxSoKCowUvHFxQFgY1necib7Z63AQ7fEA1uK71l+pNIfz5tmX5a5/fx4/H38c+Af3wu9OMjrUjkf16qxsVRlvU1JY8+UuqdgKhPt5ejrQIHk/pFu3+B3zNiwJpHHUX0mrI2MKZcxieDhRerocAAkQ/fqr5ef66is+ZutW64+1hcOHudCyr6/6OlTVu53Ic8/Jv3nvvfLy6NHycnHk5clxxSJONztbPjY+3vnXoNNxogdR8Fv8vfaa839biQhwDAuT27BqlfnjPvhA3v+bb+R778jAUGWxRmVSjLS0wskl3Jm3X6eu/UjLlhFFRvLyokVEY8bwcs2a/JmeLh/booV8nDMK5P71F9d00emIDh3id/fsWaJx4zhphr8/F/90NuLaZ88mehqL5BeWyPBCn241gjYNnlfoRbh1iyjpC64q/MKAC4b7FYi7VA6pJElEN286vs1xcXKg/dGj/JuRiCWdjw+9i/+jnj0d91vR0UR//83Lx44R9e9XQAXNuBjkI9Lv/OMffMBFMq2JNk9K4oBwgKLC4lXP6aRJjmu/NTzyCFHTpkR0//1EjRs75qSHD3OF0mnTOIsOEdHvfN9urD9EgI4Aoqz4FLkYmXj+lMTFEQ0ezOcx6lz+/VdenTePA+NdWaPm6lWi7uCiQ1k//eHcH8vOJpo0iTb14OrWBR/OVhedMpGlRPTLPj5c9H4qPjLsv6yOXA/IFfOEGTP45xITiatelimj+n/OwnQCiEaOdM7vx8UR9cdGIoC6YZfip3U0bZq+zxW1kLZto/R0othY57TFEXTqqKMnG+2jKz6Rrp/wWQAcXRDTEX+lQZA5cEDdT/74I2ceE+tPPWX5ud5/n4+5coU/58/nF8UZhQHz8ojq1OGsQevWEU2cyL/58MOO/62iuHOHK2iLe3X//YUntZmZpo8tKOAJnNhPmSBp0ybetmOHc9t/332mJ+IA0Zw5zv1tY3Q6ok8+Ubdh2jTzx40fT1ShAo9tOTlEU6fysW+84bi2ffON3KY/FOP2o48Wvm+VKjnud23hscdM/z/ff1/ONARwEUUlymdh3jyiU6cc2y5x7kWLTLevUyfH/l5RhIby761eTTQQf8kXTMSag5kzyRd5BOj4pVRIVwDRA1hDBNC50HbUFof43QVPSCUU0KJFjm1vSgr/brt2PC9+7TVe3xA0jCgkhHYtv2JvzVDz3LhBf3xwmgCinNr15WrtQOEZdFYWVx437vjWrjUc0x9/UwRu0Xb0pNmYQuXKuUaIFeh08jvdvVMuUUgIDyDO4swZw7UnIJza4wBduW+8+gXIzlYfM3SoyRfl95X59OefvOrnx4oKHx9OvucqYea774j+D++STpJYUnABCxcS/YtOVBBUVn1P7rmHdNGH+QHST2gPH1b3ZWvwABFAsdW6qI5dv965WU2JiJo3l3/y9Buc7W4j+ivehY0EEI0d65zfLyggaoDzRABNDPqe6tcnqo7rdBuV6FO8Sl9+SapJyOuvE1WpwsdmZRHlzf6UUz96CO+HGk0SHnjA3U1SoQkybiQ6muj8eaLAQJ7g3HMP32lfX7Z4WMrkyTxRyMnh499+mygggJc//NCxbd64sfDE8vr1ogUHZ3H9uvxO/f030Ysvqt+zI0cKH3PjhnqfGjXU2vFLl3i7M9Iwf/45Z2ckKlqIAYj++cfxv20J5cvLbWjWjK17//2n3ic5mS0I27cT9e1L1KGD/N3Nm3zswIGF0zkfOULUs6f6XlvCrFlymxYs4G1pafy+dOnCy+++SzRzJtGFC9ZesWPJyeF79vzzxf9///c/9XHG+z/6qGPbVVxbAM5g6wo6d+bf47mljtZP2UlzPyugK1f4+6wsuU2mrkFoosXfi/ifYXlYhW1WKX4sQaR5Vv5Vh77TmTHDsT9WDLt380/Gd32IJ/6iMevXyzvdvi3n3l++XN6elFToIv5s9IZhuSwyXKqAEv1raxyma0H6XNC//+68H1S6OAB0EfWIALrT50Gi118nAihz4w7Wzq1dy/njAdbQAERbthiOrSQlGwwSY8YQheKO4dSOFqKLYsAAoo1Bw0jXqJFrfpDY+PcSvjDchzzJT21GBkjXty/9vDyf5s6VN3fDLiKAZuIt6ts7nzajr+FLVlg4t92ivwnEXTqFpnQFtQjQ0UD8RfvQgYKRTgAr8ZyFH3IpGwGU1rQDZY5/iT6p+pnhHrz2bLrqXe6rvz1ZWUQvtFeY/kQH6QZiYnjOlJ1NdA4NuT1du7JJ39UTPjNogowH0L69erx58EG2hOblWXb8k0+ylYSIqFs3otq11edzJDNn8jmd4QZjDfn58vUlJbFS13ji8dNP6mOUdVrmzy9cq6WggO/76687tq1iklaxIq8XNanctcuxv2sNoiNVarKMnx0xkA8axPMm40n34MFkEGaUdO/O24XHh6W89hprP4VQvm8fu0cBakHa0+iiVkAaPAgAni8pyc1VW56aNy/+3Js2cb0HSzH1nEVEsCXm228t72PsJTGRaOdOdhMDiJ5+Wm6PTsfCrlg3FngBopY4qrqIvBBZ8v6x8Sx2VXIg27YVvm8d8R8v/PWXY3+sGITyZf+g99SNEX5tmZlEw4erOzaB6Kx796YTZe8pdEFvVl1Cz/stoIJbiS65ls2b+adX4HG5Hc62LHTsyL/z8MOG3/x5bjxdjE6lHPir70lwMOnKl6fLR1NImNuOTv6eCKC6iKFgpFNwMNHfn54gAuhR/EKRkaxEPHPG8U0/fJhdMQsKWInUULrA7Rw+3PE/VgQXLhA1xDnDPfJBPt2+nkVb0Ud17zpgn+rR/BSvUjYCKBjpFBREBOjoOXxjuJfOFmT6640vPbGdCKDxWGhoX1gYC9V//+1caxpANBNvFXrvCKB5bRQ1kVaupK5VL9IQrKarV4nWYrDhu5S6rYnu3nVeI820HyA6fkxHCQinuFYDC1swPQRLBRnHRTVqFKJZM/X6wIFATo7lWWWSk7kOAcCJfoyzsxDZ30bBiRMckO7qwGpjlNk6K1Xi4PVffwX27ZO3jx/PAZoiE5kyG1jbtoVrtfj4cNCto5MliED1tDT19mHD5GWdjmsKuYuff+YaLD16FL2PyBL2118c7DpkiPp78T/ZuJELHYtgfZEl1doA2ZQUrnHy2mu8/vXXcnyus+vs2INxIqvJk+Vl8Z4K/P35fRecP190FfIDBziYtGdPy9phKgPfvHn8PP73HzBhAtcwcQXh4fxsif+bsnZO27bAww/L66bqE6RCLky1pfkr8MuQg/7b+RxG2TOHkJ3Nz1z//sDFi/a111TiilpB+o7B0XnGi6FKFX6v9tQcYdj2Dt4Fdu7kTGTBwZxFYto0/jIhgR+iqVM5uLp7d9xasRUt7h7Atv4fc7GQtWsBX198eHMcvsp/DtlvzDBdDMyBpKbKSUEaQl9v4PBhOeWfs9iwATh+nIP29Zy8XRUffVMej+Nn5PqUkffNzMSBoF6IbF0BuZU5JWFCHqfI/B5jkYFQlMnPRJcUzlf/DBZhxw7+/3TsKGcPt5fsbB6r2rblf+uxY5zgZg7pO0JXZOfQU6sWcAGci/gMmkAHX8z6NBB98Q8qIAVpD40BAISAsxMFBHDW9EH4CzvQC5kI0Sfbk7ALPLhMwSeoiGSnxriLOU8wONHDUbRGo0acp+DyZSAykvsJS7J+28rgwcBfGKTaltmE06O9eERO1JCVfBff3RyENRiKgNdewANYjxmYiUFYjwqXjoJekVNsO/k1NQHhg1YrEYHboOHD5cHcW7FE2nHUX2mzyLwhW/sNmnmAta+W0LMnUY8evKzUbIq/RAcp3M6eZYvFsGGOOZ+9dOtGJgNuja9/+nTZa0D8Xb5s+pzNmjk+3ufgQf7NwEBeb9SI1zMyiPr1K2w5cifKIH5ArbFSxnN06VLYheypp9THinjt3r15ffJk69ry8MMcFFxQwNasESPkhBbutF6ZY906bmPt2hw7pKzcfvp04f0VxcgJ4EQUxpw7p97Hkvg3ZRxepUr8uXWr/ddnD0bePib/pk0jeuYZOZYUYFceAmiV9BD9POQneWfFQ3n9XAZ9pvfemDDBvnbOm6du0+efE938ZCmvuNiPsWZNvsyYD3+hdjhIdRHDLj76xh1CW/r1m0TKDq7IEfyKILLpmGW4hkOHFCetUkV9gY72zTNC3E8/5FIu/Oj4/VOd+nuFyMkh8venDcGPEMDxLQDRfZ3TiJYs4Y4GoDfwMQHsVkNEtPRptUvjOHxnSJ6QEVKZqKDAEFu5ZIn9zRQuusq/hQuJTuxMomwE0JUeI12bYYB4TGhfPY7KI6VQ2878wAPrQPxFANETTxBHrQN0eMwXhv3q1+fP5XiCCKALqE/JC1ZyUKiDg8b37ZM9C4ZjJRFATXGK2rZ16M+Y5e5doutnjDLTzJ9P1ytyEg9hebky9WvKh49qvy4t0ggg+gqTSCdJ9OjANGqg98bcv995bd6/n2juXF5uB564pCOYzqAx5WW5yHxvA9AsMu5HVOsFgAUL5Dopa9ZwqlZz5OayJgQobN0BuAqrI/jqK9YWv/eeY85nL7t3c2E1c5w4AUNFcEFR1czLl5eL7jkKkQYyO5v/n1lZwJgxrEzdvNn5BUStwdjSpixdoazlsmZN4Vz7n3+ufl5FGufLl/lz82brrIO5uawA8vFhK2BammyRCQy0/Dyupm9ftkgtWwbMnq22ehhbZAB+DkaOlGvQmCo0tmsXf4rU7ZakyP73X/68epUtk//8I1dUdxfKZ2bMmMLfSxLfs0WLON2xsE6loxxGRR3BGN/l6D73IWDiRM71qjB/XT8Yj/PnC/+OLShT7QJsTajiq8/77kKLDMDPzObNQMfPH0M02uMS6qEOXcbSipPREsfQHtF4dFI40jN9OAf/r78ajj2ADoblNm0UJ50/Hzde+ABdsQcFvv7A0qXQfTzHdMEmByCKb0555BL8kY8Wj5kYqJxJQACQmoqVQzjdt07HhWA3/xeK+AFjsWv0d8h59iUsB+fuF33WtfQKqtMsxnjuxF59FcEZCcChQ+jXjy2OloxF5jh5svC2V18FfH/8HmWQi+sj3nCuGcEE06cDTfpUxx1UKPRdei5r6MsgB++/z+8tNm4EAPgNGmDYT5Qx+uW+7/E9xqABYlDxuceAN94o/LLZwaVLQKdOfB9HjgSql2OLzF2UdfmYERQE1GgSygPAo4/yfZkwAT++cRJlkG0oypl+6ip8ocOfGIoz9QfhEfyGsS+Fok4dYD0GQyLCrY3RBivztm3Oa3PXrlxWY+dOoCO44FQIMnG3Y2/4BbrIfO9MLJF2HPVX2iwyP/zAknbnzryuzF5mifa1XTvO3CWYPZs1I8pzOIJ69TwuWYVJjh0jevll+drHjCFVIGJx96N/f066oNOxH//Qofa3R/x/xV94uHOT9diDSOXtr3cfT0iQv2vYkKhtW6IVK4o+XpmKeNAg3hYczBYVgFOYfvQRJxcwR//+ckKB3r3ZAqfPqErHj9t+je7AkndZpC81laBj0iT2xf9Jb4ywJLvZuHFyJhxPQtwLUwH1ilAGAjiOR7neoIHRybKyKLkNm/y6Yrdhv759+R7ZGsv3zDNEVatyvCLAmRJp+nRW5RubIp1MRETh+wTIsWOGPk25sno1HWw9ngKQTYBpS2ByMu9aHinycYGBnGvagYgg/w8/JNlc+e+/Dv0NS1H2T+PGmb6vSuvKkwMTeUP16nIE+WOPcWBmQABn6yDu61q0sL99K1bIbVi8mOiFF3j5cp3udBDtVGnoXYnS2n7PPWwlAoi2zDtNBI4X+u8/4qC/5s2J6talmIs6wzF//smxogUFRLt35NMfUGSH27vXYe1csEA+7dixRLlzeUALR4LDE6lYTH6+yuok4jx9wIG+m30HEAF0P9Yb2r59O1vqw8DP3+Kmnxi+69SJb7Mz8JMNvTQKP8oryiQiHgg0i4z7qVyZP4X2sYKR4uPUqeKPz8tjX3vB1KnA8uWO0RAJiFjb0aqV487pLFq2VMcd3LmjjheYP7/oY8uXZ83/kSPA4sVcHIussCKYwtjf/u5d1tZ4IqNHs6Z85kxeV1pkUlKADh2AJ54o+nhlHNLu3XxMZiafMyCAC2+/+Sb/T3Jyim+L8rkuV857LDLFISynpggNBWrWNG2ROX0aaNFCrvtnSRxXYiLHjnkqkZHq9Vq1gCVL2CJTk0MU8Oyz6n3q1DE6SWAgEt/6HwCgGuTgmlOngKgoVWiEVSQkcL+8ahVr58uVA9/0sDD7zT0OwjgG6gDYOvU+3sa1tg/ii6hFyEUZLFgANG1a+Hhh+b+DChiIDYid9Cm/YJs2qas024kYv/q1uS1XTG7UyGHntwZJ4nFs//7i+zFhkTmfHI4R3eN4w7//cqXqH35gM9mIEbxcUIBmzYBz54qOb7MUEd/aujVbMUaNAvyQh2pxh7Ab3VG2rH3ntxXle9SsmWzZTcvmDq0Mcjh8Ys4cNoeMHYtKYfJgEBrKfbaPD1C5mi8ewp9Y9+EJ/rKYApvWImLv7r2X2+yfyxaZqe+WNRQMdzm+vqqBsW1b/tPBF1kIRJMCNsNdQl3DPlWrsqU+CeGIRV30Dj6A4/r6ovv2qeMuHUn16vJyKNJ5oWtXYPhw5/ygi/GMnruEItxNCgr4MzhYbT3etKn445WuZQJJ4qDgKVMcE5+VlcUTencH+VuK0nXszh11lfjnniv6uPLleUBq107edvOmfW0xTtqQleW5gkxoKPD99zypBORi60QslIjJT3Fs2QJ8/DELHsuX87YaNbg/FAkDAA5gLQ7lcx0aWjIEGXPUrWs6KUJaGt974dVkiSBz+7Ys+Hgiyn5u61Z2AS1XjhUxy5aZPqaQIAOgXONqAGRBpmpVOWGArW61QpDx9VX85u3bLncrA9gjpWVL8/vdi38QjkS8g/cxfz73O717FxYGlYii639jIA710AcVT53Ks1XxstlJVhYgQYfmL/Rk36PISLfcR0FkJCtk7r1XdoU1LlIuBJkbN4AydavLGpXgYLnz6dWLtTSxsYiK4v7qwgX72nb8OCsyjxzh9fLlge7YjYD8LPyHzm4TZIYN42QGjz8OfPSRnKMhNUt2LQukLBbs2rYFZsxg4V+Pcrl2bfa2OpDcgCcqDhRkLlzg/+vWrXrhK5MFmdffDnLnI6dCkoDoaL5Nd1EWtXAdADBtfiSCg4FPPgEaN5b3P4AOqHnjAFq0AEYPZ83i0qVwSrIEk4LMpk3Fa+C8CE2QcSJictiQk4NAktTJXD74oPjj8/KKfs7KlmXNtxCSbCWDk5IYfPQ9ndq15eU7dyzPRKbscAX29LNEhQVRIrhtQLIU0T5hkcnM5GfIEkGmb1/ONBYQwJm2ABZEqlZVC5QidqgoSpJF5uuvLdOihYebvi/Z2Sz8CqWH8j4WRVKSZwsygPz/bd5c3ef16qWegItETdWqFT5HpQaVkAt/1NRPCEaMUH//1lvWa8qFIKPi7FkO1nIx7doB330nrxPJmfyUZCAUSeB/+OzZwKFDha1exkyaJD9v8TcV5tQ7d8xr0CwkKwtogyMIjDnNGxQxPO5mzBhWbH34obytY0fgyhW+zzdvqid3Kpo358+TJ9GtGz/Lb75pe1tOn+bskcpnvHx54EV8iSS/KliHB9w6blSvDvz0E7dPjJNKQabq4g9YknjxRUCSVIZLpQI0MJAF8/+OBLLmRkhtDiAmBmjQQLFBuD94iBVVSZkyLMgAAKpUwejnyiIjg8OGhAFnyRLgevUO8L9xFXjjDSz9rSx2Tf4DlTNjDXGTjkTcpnLlWJDRST6eP1mxAs97CkoQjRqxC9PChfI28UCVL89CRHEDcW6u2rVMifGEtCjy81lgMk4RLPA2QaZSJb4vI0bIFpl69eTrKApT2UBjYmxvR2KiaUHIUy0yAtG+rCxOszxlCq8buz0Wha8vDyjR0bweEsJCkHICbirFrRKlRaZcOSA9XX6OvU2QmTQJmDvX/H4REaYtVVlZfM3WWmQ8RQupZO1aTuENACtX8sTRlMCltL6sWsXJJCZMKLxfQKAP9qAbnsEixO84h0m91b55H37IBgZzFkAlhQSZzEz2+VOaal2IUkMLFG8ZHzGCJ+G5uTxPNEelSjx+3LwJ4OOP8SeG4g7KGYK27SUrS5Fy+cSJwvnJ3UhICLsajx7Nk8ann2aF4uXL/I7l5ZkWngHImXWOHkX93DN44w1g3TrbY9eFFejtt+VtFSoALXACu/17IxtBHjOn9PfnMSI5UxZkyu7fxpH2JrJ4GD+vnTtz8pEbTfpw9Lq9Pnlgj4GkJKB+fcXGzEyPzdVfpgyQB/3ErYiU2mPHAq+v68Urn34KAOj+xcOIRX1MeIbsdns3JiMDGDqUlY4hyEBuQIjaX9zL0QQZJ/Pgg2prgBh0hXahOCu/KdcygaWCzKZNwIwZRWuN0/VWRm8RZADubEUWsuRkdh8w16eZssh8953p2haWIATDjz5ST769RZC5fJkHZxFXZIlFRlC7thzvERrKEyZlti1rLTJ5eXJGOW8TZCwlPJzffeN6AcIdUfxfzHn95OfzvfZEi8wDDwD338/LQ4awz7epREzKCUnVqpxNR2lpVSLNmoWKSEXVXk3Q8MFmmAi1Q/yFC8XHxinJyuL+TiXIzNDXWnHTJFz0S0KOEhNDU/2xuLeAZYKMjw+74sbHA5gyBQ/hT+xBN5CDgizv3gUaQJ9yyQ0WLUsZO5b7+jp12CKzn5M2FS3IBAfz9bz/PtCsGUY1PwqdznZ3RjFWKONQA/0LUBtXcTYrEoBnKcfLlQNS7rIgUx53EHAiushCZMbjqsh8uh29+MJPnza4MduKcOFWPWIeLsjUg77RTz9d9I6tW8vL06cbFvMuXTMbZ2op6encR6anc5/i48MWmbxAL4klsBBNkHETIsiuuAe2ONcyMfEx10kIofvoUdPfe5tFRiAEmfh4dZrrojC2OPj5cQBh//62/b7eRReNGsk16wDPF2TEgCk054Ki0labQlnkUAgyygm6OUHG2CIDsKbcx8d1hRxdTUQEu/AZZ8EVgoyY8JtzFU1JYa28JwoyliIKrnbtan7f3m92Ur1UH2I6/KB2Ijf3vBnvZ3jWidgpvWtXYNCgIo9zNvHxcgIX0Q+bsswIjyfAvGuZIDyctdmffcbr+9AJ0rlzDomTERYZql7Ds2biRSDSVI8fz5+FXAyVKIJrah5dD8D29P1CWaic9EvxNxCAPMTqA8E9adwICQHSsljT1AX/QsrNLdKyYCxLhIWxkH2KogAAVzafQ3Awp/W3FVGQVMR3AuCJj4c+c6rY5eImGD4+bMmMj2e/M32n3grHHBXGhhEjeI5y6RL/X319NUFGwwH89RcHvIpOtDhBxhLXMnOCjNCUGwemC7xVkClXju/PuXOWJcpRCjuJibLF21zmuKIQgkxwsNptzUP7VgNiwDSuVl2kdtIEYiIA8HNjXEPFnHuU0iIjJuSXL7M1pgRZu1WI6zR2gxKCjCTxuGbOE0O48JmqW+MtBAbyM7JhgwU7+/iwr8qECcCSJaiAO+iKvapdhJJGklSKzUKIPtAwR71wgTvIcePc6mtftWrRAoyyb6laleWu0FDO3GYJFSvyJb7+Oq+nQT+btldNDkCXcgf3YwPQto35nT2Ahx9mrb5I8lKse6Zi4A395G1E4aRZ9+WiEBYZlfXiiy8AAJcRibJlPavfCwkB0jN9UODrj97Qm6GKEGRMvTbVqgHHszgw+MrmcwA4UYytXOcQOXXiBg+3yHTEPrw/4F/zmrnmzfnFrlABuHQJBT5+eAh/IOtObvHHWYgySaFSkMkPKmWCjCRJgZIkHZAk6ZgkSackSXpPv72uJEn7JUm6IEnSSkmSSkb6Aydz//2celFI7eYEGXOuZZYKMkW5oHmrIKM0Mxv7mZtCKcgotdm2Zn5TCjJKa48nadZMIdpnLGxYI8goERYZJUIDWRTK51q4WJ46VXLdygD5+VNaDwoK+F6I/4mvr3mLjHjfPV1gNkelSqbdPU3SuTPw7bfAI4+AfH2x8pl/DNn2Ro/mmGLxPH/0URHnuHYNNV8YigpIkV18xOyqiEmaO1D2wxcvqrNlRUTw9d65Y7krqBBkBFmw0IfRAppFL0MEbkP6v/+z+1yuQmnJKlaQmT+fI/z1AeuPYaXdgozhfyssgQCi0c7j3uXgYH38rm8Z+EIHatCgkNvDiy8W/QxWrQpcSggGatdGvehfEYAcu2SOuDjZTRJ79nAu640bPXbACAwEDqAj4mp3tu7AkBDcbNobY/EDgj582/z+FqDMgBYcLLuW6YK8bMJnBkvUUDkA+hBRKwCtAQyQJKkTgI8BzCWihgBSABTjDKhhjHgHixpPiArXkVFiqSAjNLhFFQ0WnbO3pF8WdOwoL1trkVFia/ZBpSCjnMh7uiAjnhtjQcbWgSY0tPCAZm7AVz7XQpCJi/PYcckhiPurfF/Fuy+eGT8/8ywVcMMAACAASURBVIKMUEh4+nPmFEJDIbVtiyrndgE//ogKD/bEfR1SkZEB/P23mWPHjUP9k2swrOJOuQbP4sXsb9SkibNbbjHKfrh+fXW/JZS71mjvjePXsqF/ybKyEBvLFp+JE21oaEwMBux5Cyd9W3pUkL85lMJzsVbNRo04m0Tr1kC3bhgsbTD0+daSns7/V4P1IiEBSErCVw2+QDLCPE6QCQnh8S3fl7V8kjKWQ8+8eUVnWKxaVW/1euUV1Ew+gcFYb3MsKsBjQ7VqgF9GKtC9u5wtxFbtm5MRylFbxrOTT7AmJnjtTw5pi3K8uX4deO3Ga+iGvcgPsTC7j5dgVpDRF9gUUxN//R8B6APgd/32pQCGOqWFJRRzFpmCAhZm7A32F4NYdrbpApDeapGpW5ctKyEhliUcMhZkRGIaW7ODiEGtbFm1Zs8Ts0kpKcoiYy1iUA4JKRyvYU6QUVpkypeX/zclWZAx9b4bCzK+vuZdy0q1IAOw9eTAAeCpp4Bdu9BR4jzgq1ebOW7bNgBApWr6f0R0NGvblX6SHoCpfvj55203GhlnFFRaZFasYGvBggU2nPjnnxGUm4a3K35tW8PchHDVCw62who/cCDa0GH43LJtNp6WZqQo1PszJ0Swf6CnCTLCIuNToFfnW1ktu1o1fuayx00CAKzCI5DOn7O5PfHxeplF/w4bEMXMPAzxXNnSR2dHtcPL+AL+N69zZgo7uXtX1jPUqJyHUQmfAwBudS5Z03WLHIMlSfKVJOkogAQAWwDEAEglIjHsXgdQo6jjNQpjTpARJkF7g/3FIKbTmd7XWwUZSeL8/LdvWzYgGVsctm/nLEtpaeYn3qZQWmSUE/liA0g9ADFo5uaydrJBA2DwYOvP888/HEgYGqpWjAUEmHctM7Y0iiDckizIiPdYWbXdOOW0Ja5lpV6QqVVLZcYOT+BJ4eHD8i6pqXwfp0/X+9dfu2bIRhHqq+8Ep07lF6C4MvBuwJRl/Kuv1AVnraFiRSPhWWGREfUqLHbxUxIdjRuhjXCyQjfbGuYmxLVapXAaMAAAUCfWtrRlaWlG91gvyCRVjVK1yVMICeExMShX7xMntH4WIty+Y64F4FgAz6KjYtfZ3B6DK+WePbyhRw+OMfLQYo6iWbaMZ4GBwGbcxyu9e6sHDBu4e5dzmRw8CEwfdRUAcBthuH3vY3ad19OwSJAhogIiag2gJoAOAJqa2s3UsZIkTZAk6ZAkSYcSLU0vUwow51omnl9HBfsDpieYqanssuCNk8iICNtjXCpX5hoggG2ThKIEGUsyqLkTf3+5oy1blmsBrl1r/Xl69eIib5KkHogjI62zyAAsUAIsmJZUTCkujIUSa1zLPE2L6zIMfmFM6KqlkKBTJTNJSuI5z0cfsTVDWagxRMpkqeeff1iYsbSAkotwdPyysfuUwSKTlWUYGzIyzFimMzI47Zn+4YyPB1K2HcahgrZeJ1ALi4xVyTKiolAAH4QlnrXpN9PTTQgyFSsiP5yfZU8VZAzUrGnV8SJm9fx5YEjZrchAMGqnnbS+IVeuAFOnIuNOAd+jc+fYOrRzJ/Dyy9afz0WId8mWdyMwEDgLvavrpUsW+MwW3w6R3K19e6DMdS6c9zBWwT/ExomTh2JVqhYiSgWwA0AnABUkSRIpGWoCuFHEMQuJqD0RtY/w9FmeCzFnkRGCjL3B/sqij6YEmcuX2eXUk7KmOItvvlGngezRg+/vP/9Yf66igv2LEjw9CaH1FWl/7f3fK4+vU6d4QUbEfimfa6EUV8Y9lTSKs8hormVWYJQn3PfkMTyOn1UT8dRU2RIdkJvBRUT0k7EaeZdkX9TevV3RYqsQ7+a99zrmfMr4tRdfVFhksrMNqcB1OjOx/ytXctozfXq4pZ/dRsW0q9h91/sEGSE0WBUTWqYM4svURZXUczhzhutCWlPn0aRrWbNmCAmVVG3yFIKDoY4Hql7dquMbcsIynD0LXE8rh73oiiY5x6x3427cGPjkE1RKvsj36MIFywJi3YwoRWCLwYgVyhK2f6efuAmzqQ3k5bHuwaD0usg1ny6iQYkrc2BJ1rIISZIq6JeDAPQFcAbAdgCP6Hd7CoAdmcJLH/a6llkiyKSlceHCtm153ZQgExPj0bXMHMrEiXINC4DvYe3a7HliLZmZrEEPCCg6kYKnohRkHE1YWPGCjJgAKAU+EZBsUTpeL8USi4zmWmYBSkFm9WqgRg087vurapfUVI6nBoDhcV/wjOr99wEAbRI3yzu28by0weXLs2Vy0SLHnE85Sa5bV22RSUvV4QV8aT618KFD/KmXdvIPciavaLTzuudQvIfWjnlxIY1RPf0chg3jhGPKbHLmuHtXYWnT6YCTJ4GoKIM7t6e5dYeEGMXeWlNkDPzMRUQAx4/z5R4P7oKWOIbMq8UHZup0wIoV+t/OzpY7y8xMVAzOBWJjvUKQEQKbLfMC8T6lVKzHbg/Ll8v5p62kkPV++3Yk+VdFPKoVKszs7VhikakGYLskSccBHASwhYjWA5gK4FVJki4CCAOw2HnNLHmIDtVW1zJzrmkA+0USAX378ropQSY2tvQIMqaoVKno7CvF4cFp7M0iBk5nuCeFhBQfI1OUgF6hgnfXRjGHJRYZLWuZBShdy8LDgT590Jp4Yi36xNRUztAahLtomrwHaNGCc94DaHz7X95p3z6P9adt2tRxTVP2UeHhskWGsrLRLPVffImXcBItkJFejLrcSJCpeIETLBxBG69zcRR1NUS2REu5GNEZDe8eQ/34PVb/Zk6O4v95+DA/oN27GxRKnqYIKzSu2eBmEBYm1206GzkAPiDk/la8rnvvXn5Nn3oKqgreQRmJ6HTtN+4cO3Swui2uRggJtng6BMohbByEeusWpy+0ISuRKlV/Xh6wcSP+DR8Cgk/pE2SI6DgRtSGilkTUnIje12+PJaIORNSAiIYTUTEVUTSMEQ+sra5l5gQhANi0ifsgfaxioQlmRgb7k9eta1mbSyLG6UktRaVlA/D558D33zuuXc7EGRaZS5dYcC7kX22EOQG9pGKpRaY4l5WMDOCVV9THlDqMUwTWr4/quusIQA5q1+bNqalA1KbPkIZyaBG3CejUCfD3Ry70D1337iXbj1GBUtsfHi5bZHLvZKFJgRy3kHPJpGc4P7DHjhmWiYCWt7bgCFojFRU9zppgjkmTuNi6tSmnt7d8GVf96mJh2mOQoLMqFXN2tiKWc8kSnuH262cYPzzNrdsR/9OKFWVB5m7UPTiMNgj5/L0i9//iC+BrfQK8334DDyZ6/sAwDF8zirPK3H+//Y1zMkLmsKXGrkpBPWoUv7S5uWxVthKVIHPyJJCZiWOV2J221AkyGs7BUteyoiZ8fn78V5wgs2UL0K2bnFXKWJARPtKenjLYmTjKIvPKK+w77Q04wyITGckBhaGhxQcPm3OZLKmYssiI506pmS3OIrNwobxsa5ILr8fHR04NGB4O1KsHHxAicdkgyFBMLB49PA0JqIxjFXoCY8cCAHRiuPOgApjORtlHhYXJFpns1GxE4ZThO+nQQeNDmZMn5Zc2JweZlxPRSbcX2wP6A1C76noDtWtz/LS1Y15AWCi+CXgFNXADYUiyWpAJDAQPuAsWABMmAFWqGCa6tkx4nYnoj4ZF7MGLUTYEkIIFGeHeWaeuD1bhYQTEXy1ywvPKKxyKBejdIRWCTFnoNT7Ll8MbgjscYZHJzga/vP/9xxv2WG8JVAky+vs57KMOaNuWdTslCQ97hUoP9gb7A/zQFyXIFBTI8TGiYxIVhgVi3dOCDV2JcZ0FS/Fm1zJnxsiEhLAQU9RAX1otMuI9Tktj16FNm3hsCQqSs/yYcy1TCi+epsV1KUePcjFLvSADAPUQixo1gEApB0O/HQg/ykdH7Mf4BjuAzlxhOxD6zraUWGOAwq5lwiKzeU0W2uIwroW3RgaCEby7iOxIwq0MAHJykLdkGfyRj/DJT2LZMo/LXu00qlcHLt5ljWB13LBNkImJ4c6xPwuBYsLraYKMSBCxOrErzlXvY9c5AA5FS4beb1jh/pCfDwwcCOzYoT42KAjcOUZFGbbF9hgD9LGtLa7Gnv+rGJMN8zpR/PNGERbTYlAJMnv2AOHhiBpcF9HRnheXZS8e9gqVHszFuFiiuS5OkLl2jYWkxo3lTkVYYASaIMMWmdRU602t3izIODNGRmTqvHrV9Pel1SLj48PC25Ej7CXw+uvA/v1sxRJKRnOuZVZlWirJVKsGjBvHy/Xr8wdiUKEC8EjZDQhPPo8bAZG4jlqGPk9lIbQ2QMKLMXYtExaZmH9voiP243KzQViLIYjY9bvpTtBIkJEO7MdF1EdwhyiMGuV5k3BnUasWcAOcvasa4s1mC1WSk6NXQghfK70vt3gmPU0poYxVtHWME3OO8HBWpqZAv0GhNbx8ma1jIv2+ICA7jVMtK7IKUrj3ZLy15/+qipEBeNAIDwdu3rT6XAZBJlAHbN4M9OvneQ+bgygl3ZDnYalFpjjNdXGCzDl9Id3GjVnKL1OmcDV3TZDhTpuIi24ZExPDna0pvFmQcaZFRiSVOX/e9Pel1SIDsPAmauXUrs2ZjxRKR7OuZSXNr9khVKmCbN+yqIdYlC8P1A5kf5aXav4BQFYAqwTEyEjXttGNKPuo0FAg35dnStPwMfxQgLs9BmAL+qFMepLpl1apGc/JgW/seZxDY5XGvTRQsyYQD7bIvI2ZyMywLPiaSGGRMRJkRJKdli0d3Vr7UP5v7RVk7rlH7/VgwiIjCtcbx1Q2y4rmG9erl2GbVMV7BBkhf7Vubf2xfn6sHFDN66pWtUmQEVbDsLjjnDTgvvusb5CXoAkybkJopO11LSvq+N9/52ObN2chPCyssAuVJsjI2qdbtwp/16BB0YkQvFmQEZoaZ2QJEzUExJxo+XJgyhT5+9JqkQFYmSBSfYeHs2KhRg35e3OuZaqUqBqMJCG5fD3URwyqVweqBnAndzKfi8olJ/Mzp+onS5FpSxlSIEkw1C4RtHquM/ZB7zC/b5/64D//5EB/kS0mOxtB187jPBqVakGmK/5F4PnjFh2Xn88KiMBAsHUrLMxQlXPAAL7l1iYecDaOEGTEGNOuHWekNGWRiY1VH/Pzz8BrrwGtcvXxMT16GL7zq+o9gszjj3N8kC0erJLEz4qqr69WjavQWsmlS4Av8tF45kjeoAkyGo7Gx4cnc0WZqM0F+wNFW2REPvbRo+XJaqVKmkXGFPfcw587d1p3nDcLMqIurT4G2qGUL8+x2OfP85j15JPAnDny96XdIiMQAo2y1pw51zJNkDFNRMd66Fc3BuPHA1X8kpDtE4TkrCBIEit2ExOLVviUNpSuZomzvkXZUF+cQ2PkBgSrUt4iJ4dr79SrB7z1Fm87fBh+OXdLrUUmB4G4BvadDT+z26LjxPgcnnWN03EZdbodO3qet48jBBlReqZ/fxamc8oaWWTu3sXNU+oJSa1aLPQ0xwlQzZryQAXAv0Zl2xriJuyp/S4S5hiw0SJz5gzQP+Rf+J8/zQ+wlYVNvQlNkHEjtWuz69KJE4W/Ew9ycR1JUYLM9es86RGTdMB0di5NkAGaNOF4uq1brTvOmwWZd97h4tJNmzrn/I0asSATEyNvE37DpVmQUQbrC4uV0iJjzrVMCDJffOH4tnkz/m1boOzVs/BNT0W4lIRU3zBkZBjCZ3DzJs/L6+Ayfp550b2NdTPBwcBoLEVbRMNv0gQEBQEEHyRXaiBXeYyL45t39Cjw0ENyR7d0KXIDgvE7Hil1goy4BbVxFZcQiXbb55j2RzZCCNB1r+3iTlBfz8iT8feXr9fWMW7yZDZAdevG61RB/8BMn84pqIODMX6JOntgeDj3kY1xDgUNm6i+C2xQ07aGeCHlyhklZhLmeys5cwboXFk/CG/b5pjGeSiaIONGGjViJU3LloWt+qKP1FuhTVKUICMmScoiuMW5lpUiT4tCSBLHERUVnF4U3izIlC0LNGvmvPM3asRzImXfKwZ0YYH01ntnD0qLTFwcfyqVZOZcy+7eZWHnpZec0z6vZcAAvnGbN6MiJSMJYcjKkmMQhCBzFXWQW6u+e9vqZkJCgGUYjSNoi9BQnrT6+QGJFRrKgswHH7AZa+5c4M03VVqHuKrtkSyFl+IxQ8Lz+BoV7lwF/bkaixfzrSoKMT7XurKHb37z5q5ppp2I+Elbs1v5+7NbmcC3kn4iExcHPP00AKBa+nk0hlwfpXp1oEwAoTHOIb9eY9X5yrZ34oDlYYSGGpXKKFNGnbffQq5eBaKCYtn9p4THBWqCjBtprHhXxcRGIIQMawUZIuDZZ3lZxCsARbuWBQd7XmVhV1OtWuHshsXV5yHiSaW3VbV2FY0a8eRRmShBPM+igy6NEyFTtV+MLTLmXMuCgjzPFcXtdOrEL/GCBQjPuY5beezGIpKT3bolC9L/396ZB1lV3Xn884Puhu5mafalWQyCKCqoRRCNCwbcDWhKosQkZDQ6xlShiTOTEJNY0aRqrFmMM2N0HDWicUMlkZmKRgQSzAzigmERZBER2aSBhmZvsM/88bun732v3+t+r+1+C/f3qeq67557nx7euefc8z2/5cR2/50APzEtLw/jZyoqYHvXEepUX1+vGY6uvFKX1Xv21Acu+OFqyofQvXt8spWlYn7p5eypHMj+5/+H73wHTj9dPaZSuS8ePgwD2MqIxbP0Ny2Sl60XMm214NStR0fuOeWZ0G9t5kyOSCf+UDqFCSwEoGsXx6iVs+lOHUdOSBQynSoLf/+YtqKJRaasTPtlus3Z0lBXB9VHNqjrz3HuAhHj4Sj/RFfFk/PSe4tMc25fqYTMnj0aRNejR+JqbzrXsji7lXkGDtSJdzQrVHKq6iiHDumYEkerQib4lfB33w3L/MDsXSaPtzz2mZCc4KBTp0R/9Excy9oj01zRU1Kik+6FCxm8Yym70N0Oky0yEE8hM22axqpB2O+iCwnl5fBBry+piv6v/1Kf0GDvnUaCCfj2TkNj+8547jndS7Vff+Gv1VdSvuiPlFLPp5/q+zXVpvOHD8OF/JmS+kPwwx/mvtKtxMfWttU7rkcPmFM2DebM0YfwttuY3Hke/TrsZCFfZnGfybBwIZMevR6Aui9OBOCfrvk/zqta2TaVKBKaCBk/aPnA6Uw4coRT6pbQ78CGcCA8jjEhk0ciSTmaZM3au1f7e3MLOJ06NRUy3sz97/+euHLbvbu+zKOrRnv3mpABXcw9dizRYtWckImze1Qm9O6tx2g2V7PINN2Uuro6sY9mkrXMhEwa/v7vYckSFl9+D//CnYA+h6Wl+szFWcg88ww8+aR+9kImOu5XVMA7vS7TLB13362FyVt/B4Pe1pIhsVyEALjuOvjf/9Vxf0mfqyg5uI8LWNR4PVUYwpEjMJrlNJSUJuZaL3D8GN5W77iqqiDOf8IE2LePfd0H8dqh83lspsasja/5b3j4YQB20ov9gzWAc3nlOWzrWTy/W1vQxLXMr4Clci/zaRmTaPjB37GofjxDt71pQsZoX6KuXzt2JF7bu7d5tzJIbZHZuVOPfiDy+P9WND7RLDLKAM2qyfPPh2XNCRlvPTMhkxr/7Pm9jCAcmONskUnOOpacRCYT1zJzZ0yDCIwbx9rrfsqSIJ1wZaU+Z3EXMlFSCZnycjhwuKOmjdq1Sx/EaIBDhE0MieUiRJTKSnir6yTqS8qZwsupb3rpJZgyhfrd+/kaszkw+JSiyjnfS42aWRkBmqNRyAR4V+5ew3to5hmAF15g/4ARDOYTrp0qzJ8fzzlKStcyCAexhgZ45BHtq716qTU6iYY3/hKemJAx2hMRTQ/eo0dqi0xLHdiETNswZoweH3hA3f1eey1x0E3GhEzz+Gcv+kwnW2Ti+Nv57IQ+NWc0PgZadi07eNAsMi3RN5KltXfvMJWpHydNyOgxKkYqKgKDi++4p5+eVjF/1DA09kKmogJ2H65gxYBLmSbP0Y299GQXd3A/7rPAP/naa2HuXE68ZzonsoHtV9+a30pniR+bkl3eW0tVlY793n3bC5mBA0lY0a0ZdwWHKWfVKpg0KZ5zlLSuZd4i88orGgjtU8KlMAU2RN8j6TbDO44wIZNn+vfXfpy835FZZHLHiBFw882wfr2mLLzzzsRA9eRVcm+taal94opfzYNwP4FojExFRdHEvLYpt9yiE8lp0/TcB6N7zLXs8xN11z3nnNAi4+MD45Y2OJl0QubQIUKFffrpab//Yf3g2AuZykp9pu4ruYvebiff5gl+yV3czw/Y+9wr8JvfNN7b9y9zeJnJ1F5fYLtetsDdd8Ptt7ddtmj/zHiLfIKQKS2Fs84CoPaCqxO+V1sbPzfkbt1UszSGASS7li1erMcPgoxv5yamsQZoqI+Y0qLpa49T4pMKooAZNgyWLEks27u35Z3XsxEyVVV6jLpMmZAJ8RNu0Mnkn/4Unh8+nOgK5S0N0e8YIWVloZ/v+PHw8suJFpm4vZg8Dz8MDz6oyaH694fbbku8nmnWMiM9lZXw85/r3KhTJ+23+/eH+8n175/f+uUbL+SiFufy8uDcr8w0k6q15mAlJ8fQLTRKRYUmMnmXsdT0Gsnlu15hE0MA6HrTVDgS+JAOGcJHIy/j5nn3Mq/ILIHdu7ftflV+nrFvn35ukn7+jTdg6VIOyXkJ31u2DK5O1DbHPf79WFcXrC1EXcucg7lzE7/QuXPieUMDpVs2AnC0ohulZ57ZrvUtBFq0yIjIYBFZKCKrReR9Ebk9KO8pIvNEZF1wjPlaV+s56SS1ACQH4mdikfEZtDw1NVqe7BlgFpnmibqkrF4Ns2eH58li0YRMy3jfap/8yLuU7dsXz/gYUFfSkhK1AM6c2bR/t+RaVsx7F+WSn/1Mf18IBfW2bSpuWlocOt6ZOrVpWaNrmR/omlF7cV6I8ET74NErruZSeY1LOuiOyh29iJk8GTZs4E/T/pMa+sb+N/PzDL+gtXWrPkeNv0tFBZx3XkrXz+uuy0kVCwbvaffXvwYFUdeyiy5SH+X77gMvUJIDmTZtouOhA9zFL3hrzpZY5OvPxLXsGHCnc+4UYDzwPREZBfwImO+cGwHMD86NVnDSSSpGojuhZyJkunbV70WDiLdv1wl28rObLGScMyETJZUomagZIJsIGZ+YwXtiGE0ZO1aPEyboBD3qWhb3l3o6WnIty2RMMBKJWmT694/FO71ZhgxRl/pnngnLysuDd4j3+0nVQZ9/HmbPNiFDKGS+8Q0Y+PDPkIEDOaHhIwCWj70Rhg5lyy0/58XfdTRLYEAqIZOc7ARSx7BFkyLFgfPP10WXefOCAm+RWb8e/vxndf2cMQOWLtUOnZzNbKWmq17IRVT2i8eqYYtCxjm3zTm3NPi8D1gNVANTgFnBbbOAmBkA245TNNNgYzAwZDZpSR4cAD75RJPPJJMsZA4e1MA7EzJKsiveU0/B9On6OZVFpmfP436Pqc/F66+rMD/7bJ34mGtZy7TkWmZCJnu8RcYLGUMXdYcMCc8bLTIzZuigNmlS0y997WvUT5lKfb31X+/t0KtXcPLVrwLwbsk4Hhr7GGzcyMjrzmDqVLUEdutm2Qb9M+Mt89kImeS09cc7lZUq3jZsCAr8j+Jdyl54IXQnKy1tapEJTDnvc2ps+mpWwf4icgJwJrAE6Oec2wYqdoC+6b9pNMdpp+lz6eNk6ut18pyJRQYSc45v3gyDBjW9t1s3XY30vtF+YmlCRkkeVPv2DceKZCGzcaO5lbVEp05h1sdoFpaamsRkAEZI1LVs3LjEbSe8BdWETHZ4i8zWrSZk0tEoZM49V18+Ph99Et8N4tXjmKgjin8vNPbFWzUjWWVpfWPsh8/2tXFj2p8zVmRqkYnLxLslKivD/eoaLTJvv60/UDR4P5WQefVVdgw+izq6x2Z+l7GQEZEuwEvAHc65upbuj3zvFhF5R0TeqfG7NRoJlJWpK86bb+q5t5pka5FxToVMKotMx446OffZ0bygsYmRcvLJsG5dOGb065coZGprYf58tZr94Q9FtbdZ3unWLRTb6Z5PI9G17O23YdWq8NqBA3rN+mt2eGvghx/CiSfmuzaFiRd7PjVuOp59Vo9xH/t8H22MsR41CmbP5r7Tn26Stn/tWhPQkDhXcS69kBkwQA0Pze3jFgcqKiKpr/2kZOVKGD060T82Wcjs2weLF7Nm2BVAmOTpeCcjISMipaiIedo5Nyco/lREBgTXBwA7Un3XOfeIc26sc25sHwsqSMvw4eoWBpkLmWSLTE2NJgxIN1Gsrg6zhaxerccYZObLmOHDQ3fTqEXmyBH9nSZNgrfe0jIfTGy0jLfI7Nunz3Yqi6HRvGtZpmOCkUiXLtp/Dx60sS4dffro5DzV5HHPnnCeNGQIXHUVfOUrua1foeHfEQmuxVOnUjdoFLW1mpXQs26dCRlIzFq2e7f2yVRCBvT56t5ds136eUrcSLDIRP3tkjOQJQuZ5cuhoYE1Pc6msjI+7u+ZZC0T4DFgtXPuXyOX5gJBFAHTId0Wt0Ym9OmjQeTOZS9kvEXGd/rk/Sk8USGzbBl06KCLSUaI39C6T59QyOzcGaa1XrtWj2ZVyBwvZPyzZ0ImNamyli1apC80EzKt4+tfDz/HLWg4U3zGxh1JS5HHjmm65ltu0fNt22Kxt16L+HmjXyj3VFXpJvWTJ4dlDQ32roDEuUrCHjLNMHmyekrEkZQWGWiMx2okWcg8/DAAa8pGx8YaA5lZZL4EfBP4soj8Nfi7AvhH4GIRWQdcHJwbraRvX13p8avWkLlrmbfILFumxzPOSH1/VMgsXw4jR9q+FMm8+qpOHktKwjTB778fXl+7Vi27cd9YLxu8e8/mzXpuQiY13rWscSUOuPBCjU0wIdM6hg+HBx7Qz83s8xhrvJDxaeU9/n3yxBP6TNbVWbwHpLHIEL4TgqRRa4NugAAAD4FJREFUjZj407l4585q4fPPmVmq0pMyRgY0pVmUsrLwgVyzBn77WwA2HB1sQiaKc+4vzjlxzo12zp0R/P3BObfLOTfROTciOO7ORYWPV7zX3Y4drXctW75c/zvpBojqati1S1NtLl+u7pZGIr17h2OF33PCv9BBXQWqquKXSeXz4C0yJmSax7uWeeufZ+lSEzKfhxkzNM7NEnSkxv8uyRaZN97Q45AhYWylCZlwt/uLL04sjy5uRUVzM/uLxooePXSe7X83S/qSngSLTNS1LHniEbXIvPuuHl94gT17JVaLrVllLTPaDy9kampaH+z/4YfqPpFurwQ/gXzpJfXjNSHTPKmEzNq1TVM1G83jg/29kGnJpSCueNeyZCHToQNs2qSf4/RyaktSpXU1lHQWGW+J/uwzePpp/TxyZO7qVaicf766gCe7cEc3yoy6l5lFRunZk8Z9dfy5kZq0FplkokJm2TI9nzyZ2tr4BPqDCZmCwb9Mamr0D1qeMFdU6CQn04xQ1dV6/OY39ThmTOvrGwe6dNEFkDVrwrKjR03IZEv37vqMfvxxYuyRkYh3LUsWMiLw6KOaLcomkkZb07u3vkd+9ztNF+xZt06PW7bA3XfDtGlwzjl5qWJR4PcT/fa34eabw3KzyCjJwsUsMunxKdEbGshOyIwaBWVl7NljQsbIA75T79ypK2OdO7ecU10kjD9oLvWyxwsZ0Diayy///PU+nonGwkTdoWwAzg7/AluxwtzKmiOda5mIpmK+5BLbmd5oezp2hKFDYcGCxKRI69cn7hnzyCO5r1sx4V26J05MfHdbHKqSbE1ubn4ed7x17/BhmjcnJwuZMWM4elQXw03IGDnHD3wHDqivcr9+mU1a/M7Vu3bpQ9/cRDEqZH71K12FM5rHu+1973vhwBvdFdtoGS9kli41IdMc3rUsebutLVt0XLAYD6O98FZSn4J51y597n7wg/Aen/zESM2NN8Irr8ANN9jGjqmw+UbmVFTo8cABmld8ZWUqZN54Q/32xoxhzhy1Dl56aU6qWhDYo1UgeAW+f79aZDKdtPhAar8HTXMTxegur+PHt66ecePIET1eeWU4uNh+FNnhLVhHj8IFF+S3LoWMj+P0QdcXXaRHL2xMyBjtxeHDevTP4OzZevz61+H3v9cNWo3m6dABLrtMFyBLS9W69dBD+a5V4bDb0kFljJ8P7tgBo8cGQubHP256Y2mpZi279VY9P/dcVq7UZ/CKK3JT10LAci8VCJ0760B44IAKmUxX/b1FJpMNLkXUEjNqlAW/Zsupp6plYc8eEzLZEnXFu+22/NWj0PFuPNu3ayzRggWaOviOO7Tcx9EZRlvzk5/ATTeFrsmLF+ui2Jgx6dP5G82zdGm+a1BYmJDJHL9ounAhrHi/Axec71j0yxQ3lpbqKsTWrTBtGvcvHs8vfhHOJ+NCjP6phY2IqnAvZDKdtPiMUCtX6mpaS8HAt9/eNG2kkZ7XX4cnn9RBYdo0LbOA6+yIBnn6Adpoivel37QpTCgRXdAwi4zRXtx4I0yfrtks//Zv4amn1BXZYrKMtuKee9Q98fXXw0zBRmq898yHHyaeN6G0VHOjHzkCkyY1uoLGbaHaLDIFRGWlipKamswnLV27qhhfsUJ3wbUAurZl4sTw8z33aMa3YcPyV59ixJIjZIaPQdi4MdyvY+jQ8LpZZIz2xAcH+6B+27DQaEuuuSbMsGo0j0/ZvWCBHtPGXEV3ZR01qvFj3ISMWWQKiC5dNNbls8+yi5HZt09FuQWhty8dOpg1pjX4/ZDiZOpuDVEhk8oiY5sRGu1J8r5ltgBhGPnBC5nly/VYX5/mRi9kzj4bxo1rLDYhY+SNysrQlJiNRaauTl3SLKuMUYiIwKxZ6v5opMevuh07FqYq7dVLA4h/+9ummzobRlviN+Dz+574RCeGYeSWZCGSnJK/ER9YecklCSuFcVs0tFdjAdGlS6jAs7XIdOmSuLOwYRQS3/pWvmtQ+EQXIryQEdGUrobR3vhU86eeqlZBv8GjYRi558ILNczg5JMTN+VOIPDV23ioLydEip1r79oVFjHTbYVNZWX4AGbqD9+1a7gbuAkZwyheokImTpuZGYXBT3+qQf+PPaZeKvfem+8aGUZ8WbhQF7Z79266t1gjtbUAzPznXrz3XlgcNyFjFpkCIipEMrXI+AnP4cPmWmYYxYwJGSOfDBqkIgZgyZL81sUw4o6Ieo4NHar7ydTVpcheFgiZWnqwYUPu61gomEWmgPATmbKy0LWkJaLBwGaRMYziJZVrmWEYhhFfRo/W44oVKS5GhMyqVWFx3CwyJmQKCD+Rqa7OPFjLB2aCCRnDKGaiKTbNImMYhmGMGaPHZcuaXms49XQAtlDN66+H5SZkkhCRx0Vkh4isjJT1FJF5IrIuONr6YRvg94zIxkUsapEx1zLDKF6iCxEmZAzDMIxBg/TdsG5d02sb7nyQL/IWWxjEokVhuQmZpjwBXJZU9iNgvnNuBDA/ODc+J6ecosdsNo2KTn7MImMYxUt0bzNzLTMMwzBEdMH644+bXlu7uYJ3+GLuK1VgtChknHOLgN1JxVOAWcHnWcDVbVyvWNIaIQNhznETMoZR3IwapdYYvyGaYRiGEW+GDtUUzK+9lmhtWb9ej37u9+yzua9bIdDaGJl+zrltAMExbbJgEblFRN4RkXdq0uaQMwCGDdN9jZ57LrvvDR6sR3MtM4ziZulSTbVZXp7vmhiGYRiFwNChsGoVXHopPP64li1Zovs9lZTABx/A9u1w0kl6LW6uZe2eftk59wjwCMDYsWNj9vNmR8eO8Mc/Zv+9QYNUmX/2WdvXyTCM3JG8o7NhGIYRb6JJnZ55RuNl7rtPzwcM0DkgwK5deoybkGmtReZTERkAEBx3tF2VjGyZMUOPJ56Y33oYhmEYhmEYbYe3tAAsWBCKGEjcPL1z59zVqZBorZCZC0wPPk8HXm6b6hit4Zpr1BpjQsYwDMMwDOP44eST019LJWTMIpOEiDwLLAZGishmEbkJ+EfgYhFZB1wcnBt5JNN9ZwzDMAzDMIziYPjwxPPTToP779fP+/eH5SVBsEjchEyLMTLOuWlpLk1s47oYhmEYhmEYhhFQVgY/+YnGyrz0Ejz4oJZ///tw8GB4n0/4dMMNOa9iXmn3YH/DMAzDMAzDMFrHvffq8aabwrJf/xomTAjPKyqgtha6ds1p1fKOCRnDMAzDMAzDKCK++92mZVVVua9HvrHICsMwDMMwDMMwig4TMoZhGIZhGIZhFB0mZAzDMAzDMAzDKDrE5TBPm4jUAB/n7H/YPL2BnfmuhJE11m7FibVb8WFtVpxYuxUn1m7Fh7VZ+zLUOdenpZtyKmQKCRF5xzk3Nt/1MLLD2q04sXYrPqzNihNrt+LE2q34sDYrDMy1zDAMwzAMwzCMosOEjGEYhmEYhmEYRUechcwj+a6A0Sqs3YoTa7fiw9qsOLF2K06s3YoPa7MCILYxMoZhGIZhGIZhFC9xtsgYhmEYhmEYhlGkxFLIiMhlIrJGRNaLyI/yXR9DEZHBIrJQRFaLyPsicntQ3lNE5onIuuDYIygXEfm3oB2Xi8hZ+f0XxBsR6Sgi74nI/wTnXxCRJUG7PS8iZUF5p+B8fXD9hHzWO86ISJWIvCgiHwT97hzrb4WNiHw/GB9XisizItLZ+lrhISKPi8gOEVkZKcu6b4nI9OD+dSIyPR//ljiRpt3+KRgjl4vI70SkKnJtZtBua0Tk0ki5zTNzROyEjIh0BB4ELgdGAdNEZFR+a2UEHAPudM6dAowHvhe0zY+A+c65EcD84By0DUcEf7cAD+W+ykaE24HVkfP7gPuDdqsFbgrKbwJqnXPDgfuD+4z88ADwqnPuZGAM2n7W3woUEakGZgBjnXOnAR2B67G+Vog8AVyWVJZV3xKRnsDdwNnAOOBuL36MduMJmrbbPOA059xoYC0wEyCYn1wPnBp859fBgp7NM3NI7IQMOhisd85tcM7VA88BU/JcJwNwzm1zzi0NPu9DJ1XVaPvMCm6bBVwdfJ4CPOmUN4EqERmQ42obgIgMAq4EHg3OBfgy8GJwS3K7+fZ8EZgY3G/kEBHpBlwAPAbgnKt3zu3B+luhUwKUi0gJUAFsw/paweGcWwTsTirOtm9dCsxzzu12ztWiE+rkSbbRhqRqN+fca865Y8Hpm8Cg4PMU4Dnn3BHn3EfAenSOafPMHBJHIVMNfBI53xyUGQVE4AJxJrAE6Oec2wYqdoC+wW3WloXDr4B/ABqC817AnsjgH22bxnYLru8N7jdyyzCgBvhN4BL4qIhUYv2tYHHObQH+GdiECpi9wLtYXysWsu1b1ucKjxuBV4LP1m4FQByFTKrVKEvdVkCISBfgJeAO51xdc7emKLO2zDEichWwwzn3brQ4xa0ug2tG7igBzgIecs6dCRwgdHVJhbVbngnciqYAXwAGApWo+0oy1teKi3TtZO1XQIjIXagL/NO+KMVt1m45Jo5CZjMwOHI+CNiap7oYSYhIKSpinnbOzQmKP/UuLMFxR1BubVkYfAmYLCIbURP6l1ELTVXg/gKJbdPYbsH17jR1wTDan83AZufckuD8RVTYWH8rXCYBHznnapxzR4E5wLlYXysWsu1b1ucKhCDRwlXADS7ct8TarQCIo5B5GxgRZHkpQwO15ua5TgaNcRWPAaudc/8auTQX8NlapgMvR8q/FWR8GQ/s9WZ7I3c452Y65wY5505A+9MC59wNwELg2uC25Hbz7XltcL+tVuUY59x24BMRGRkUTQRWYf2tkNkEjBeRimC89G1mfa04yLZv/RG4RER6BNa4S4IyI4eIyGXAD4HJzrmDkUtzgeuD7IBfQJM1vIXNM3NKLDfEFJEr0BXjjsDjzrlf5rlKBiAi5wFvACsIYy1+jMbJzAaGoC/yqc653cGL/D/Q4MeDwN84597JecWNRkRkAvB3zrmrRGQYaqHpCbwHfMM5d0REOgNPoTFQu4HrnXMb8lXnOCMiZ6AJGsqADcDfoAtc1t8KFBH5OXAd6uLyHvAd1P/e+loBISLPAhOA3sCnaPax35Nl3xKRG9H3IMAvnXO/yeW/I26kabeZQCdgV3Dbm865W4P770LjZo6h7vCvBOU2z8wRsRQyhmEYhmEYhmEUN3F0LTMMwzAMwzAMo8gxIWMYhmEYhmEYRtFhQsYwDMMwDMMwjKLDhIxhGIZhGIZhGEWHCRnDMAzDMAzDMIoOEzKGYRiGYRiGYRQdJmQMwzAMwzAMwyg6TMgYhmEYhmEYhlF0/D96Sm3nqdMPLAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "Prediction MAPE: 0.129854\n", "Prediction RMSE: 4.22777\n", "Running time: 92 seconds\n" ] } ], "source": [ "import time\n", "start = time.time()\n", "steps = 5 * 144 / 12\n", "alpha = - 0.005\n", "window_length = 144\n", "forecast_horizon = 12\n", "rho = 0.05\n", "maxiter = 200\n", "\n", "data_hat = np.zeros(np.int(steps * forecast_horizon))\n", "for i in range(np.int(steps)):\n", " y = data[: 4 * 144 + 12 * i]\n", " y_tilde = scale_to_data(y, alpha)\n", " x = SLRMC(y_tilde, window_length, forecast_horizon, rho, maxiter)\n", " y_hat = scale_back_data(x, alpha)\n", " data_hat[forecast_horizon * i : forecast_horizon * (i + 1)] = y_hat[- forecast_horizon :]\n", " \n", " if (i + 1) % 60 == 0:\n", " fig = plt.figure(figsize = (14, 2))\n", " plt.plot(data, color = \"b\")\n", " plt.plot(np.arange(4 * 144, 9 * 144), data_hat, color = 'r')\n", " plt.show()\n", "\n", "mape = compute_mape(data[4 * 144 :], data_hat)\n", "rmse = compute_rmse(data[4 * 144 :], data_hat)\n", "print('Prediction MAPE: {:.6}'.format(mape))\n", "print('Prediction RMSE: {:.6}'.format(rmse))\n", "end = time.time()\n", "print('Running time: %d seconds'%(end - start))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This is a case of single time series forecasting. How about multivariate time series forecasting?" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.3" } }, "nbformat": 4, "nbformat_minor": 2 }