Stepwise

Stepwise.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Stepwise Selection \n",
    "- [Statistics - Forward and Backward Stepwise](https://gerardnico.com/data_mining/stepwise_regression)\n",
    "- [変数選択法の種類](http://www.snap-tck.com/room04/c01/stat/stat07/stat0703.html)\n",
    "\n",
    "[AIC](https://gist.github.com/shotahorii/99342b1b4d63b2e743301b7ab8e56f19), [BIC](https://gist.github.com/shotahorii/4d05acadcd5e2a76568532b678f6f4d9)と、モデル選択基準について続けて書いたので、これらが実際に使われるユースケースの一つとして、ステップワイズ法による重回帰分析の変数選択を実装する。  \n",
    "  \n",
    "実際、AICはネストしているモデル間の比較で有用なモデル選択基準である(ただし一般的に「（ネストされてなくても）たぶん問題無いだろうという理由で使われているのが現状」(久保緑本)とのこと。)ので、ステップワイズ法による変数選択でAICを用いるのは理にかなっていると思われる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "基本、パラメータを最尤推定しているモデルであればAICは計算できるはずなので、例えばGLMでどのような確率分布を用いて回帰を行っても(e.g. ポアソン回帰, ロジスティック回帰, etc.)良い。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 179,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import statsmodels.api as sm\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "from functools import reduce"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 180,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# return aic or bic\n",
    "def xic(result, criterion):\n",
    "    if criterion == 'aic':\n",
    "        return result.aic\n",
    "    elif criterion == 'bic':\n",
    "        return result.bic\n",
    "    else:\n",
    "        print('criterion must be aic or bic')\n",
    "        return None"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Forward selection method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 212,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def forward(y,X_full,c='aic',dist=sm.families.Gaussian()):\n",
    "    # initialize as null model\n",
    "    X_curr = pd.DataFrame(np.ones(len(y)),columns=['intercept'])\n",
    "    res = sm.GLM(y, X_curr, family=dist).fit()\n",
    "    xic_curr = xic(res,c)\n",
    "\n",
    "    remaining = list(X_full.columns)\n",
    "    selected = []\n",
    "\n",
    "    # start forward\n",
    "    while len(remaining)>0:\n",
    "        scores = []\n",
    "        for f in remaining:\n",
    "            X = pd.concat([ X_curr, X_full[[f]] ],axis=1)\n",
    "            # modeling\n",
    "            res = sm.GLM(y, X, family=dist).fit()\n",
    "            scores.append( ( f, xic(res,c) ) ) # store the score (feature_name, score)\n",
    "        # select the best feature\n",
    "        best = min(scores,key=lambda tup:tup[1])\n",
    "        best_f, best_score = best[0], best[1]\n",
    "        # break if score doesn't improve \n",
    "        if best_score > xic_curr: break\n",
    "        # else go on\n",
    "        xic_curr = best_score\n",
    "        X_curr = pd.concat([ X_curr, X_full[[best_f]] ],axis=1)\n",
    "        remaining = [f for f in remaining if f != best_f]\n",
    "    # list of selected features\n",
    "    final_features = [f for f in X_curr.columns if f != 'intercept']\n",
    "    return {'features':final_features,'score':xic_curr}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Backward selection method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 215,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def backward(y,X_full,c='aic',dist=sm.families.Gaussian()):\n",
    "    # initialize as all model\n",
    "    intercept = pd.DataFrame(np.ones(len(y)),columns=['intercept'])\n",
    "    X_curr = pd.concat([intercept,X_full],axis=1)\n",
    "    res = sm.GLM(y, X_curr, family=dist).fit()\n",
    "    xic_curr = xic(res,c)\n",
    "\n",
    "    remaining = list(X_full.columns)\n",
    "    removed = []\n",
    "\n",
    "    # start backward\n",
    "    while len(remaining)>0:\n",
    "        scores = []\n",
    "        for f in remaining:\n",
    "            X = X_curr.drop(f, axis=1)\n",
    "            # modeling\n",
    "            res = sm.GLM(y, X, family=dist).fit()\n",
    "            scores.append((f, xic(res,c))) # store the score (feature_name, score)\n",
    "        # select the best feature\n",
    "        best = min(scores,key=lambda tup:tup[1])\n",
    "        best_f, best_score = best[0], best[1]\n",
    "        # break if score doesn't improve \n",
    "        if best_score > xic_curr: break\n",
    "        # else go on\n",
    "        xic_curr = best_score\n",
    "        X_curr = X_curr.drop(best_f, axis=1)\n",
    "        removed.append(best_f)\n",
    "        remaining = [f for f in remaining if f != best_f]\n",
    "    # list of selected features\n",
    "    final_features = [f for f in X_curr.columns if f != 'intercept']\n",
    "    return {'features':final_features,'score':xic_curr}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Stepwise forward selection method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 219,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def stepwise_forward(y,X_full,c='aic',dist=sm.families.Gaussian()):\n",
    "    # initialize as null model\n",
    "    X_curr = pd.DataFrame(np.ones(len(y)),columns=['intercept'])\n",
    "    res = sm.GLM(y, X_curr, family=dist).fit()\n",
    "    xic_curr = xic(res,c)\n",
    "\n",
    "    remaining = list(X_full.columns)\n",
    "    selected = []\n",
    "\n",
    "    # start stepwise\n",
    "    while len(remaining)>0:\n",
    "        xic_prev = xic_curr\n",
    "        \n",
    "        # forward\n",
    "        scores = []\n",
    "        for f in remaining:\n",
    "            X = pd.concat([ X_curr, X_full[[f]] ],axis=1)\n",
    "            # modeling\n",
    "            res = sm.GLM(y, X, family=dist).fit()\n",
    "            scores.append((f, xic(res,c))) # store the score (feature_name, score)\n",
    "        # select the best feature\n",
    "        best = min(scores,key=lambda tup:tup[1])\n",
    "        best_f, best_score = best[0], best[1]\n",
    "        if best_score <= xic_curr: \n",
    "            xic_curr = best_score\n",
    "            X_curr = pd.concat([ X_curr, X_full[[best_f]] ],axis=1)\n",
    "            selected.append(best_f)\n",
    "            remaining = [f for f in remaining if f != best_f]\n",
    "\n",
    "        # backward\n",
    "        scores = []\n",
    "        for f in selected:\n",
    "            X = X_curr.drop(f, axis=1)\n",
    "            # modeling\n",
    "            res = sm.GLM(y, X, family=dist).fit()\n",
    "            scores.append((f, xic(res,c))) # store the score (feature_name, score)\n",
    "        # select the best feature\n",
    "        best = min(scores,key=lambda tup:tup[1])\n",
    "        best_f, best_score = best[0], best[1]\n",
    "        if best_score <= xic_curr: \n",
    "            xic_curr = best_score\n",
    "            X_curr = X_curr.drop(best_f, axis=1)\n",
    "            selected = [f for f in selected if f != best_f]\n",
    "            remaining.append(best_f)\n",
    "\n",
    "        # if forward and barkward didn't improve, break\n",
    "        if xic_curr == xic_prev:\n",
    "            break\n",
    "\n",
    "    # list of selected features\n",
    "    final_features = [f for f in X_curr.columns if f != 'intercept']\n",
    "    return {'features':final_features,'score':xic_curr}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Stepwise backward selection method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 220,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def stepwise_backward(y,X_full,c='aic',dist=sm.families.Gaussian()):\n",
    "    # initialize as all model\n",
    "    intercept = pd.DataFrame(np.ones(len(y)),columns=['intercept'])\n",
    "    X_curr = pd.concat([intercept,X_full],axis=1)\n",
    "    res = sm.GLM(y, X_curr, family=sm.families.Gaussian()).fit()\n",
    "    xic_curr = xic(res,c)\n",
    "\n",
    "    remaining = list(X_full.columns)\n",
    "    removed = []\n",
    "\n",
    "    # start stepwise\n",
    "    while len(remaining)>0:\n",
    "        xic_prev = xic_curr\n",
    "        # backward\n",
    "        scores = []\n",
    "        for f in remaining:\n",
    "            X = X_curr.drop(f, axis=1)\n",
    "            # modeling\n",
    "            res = sm.GLM(y, X, family=sm.families.Gaussian()).fit()\n",
    "            scores.append((f, xic(res,c))) # store the score (feature_name, score)\n",
    "        # select the best feature\n",
    "        best = min(scores,key=lambda tup:tup[1])\n",
    "        best_f, best_score = best[0], best[1]\n",
    "        if best_score <= xic_curr: \n",
    "            xic_curr = best_score\n",
    "            X_curr = X_curr.drop(best_f, axis=1)\n",
    "            remaining = [f for f in remaining if f != best_f]\n",
    "            removed.append(best_f)\n",
    "\n",
    "        # forward\n",
    "        scores = []\n",
    "        for f in removed:\n",
    "            X = pd.concat([ X_curr, X_full[[f]] ],axis=1)\n",
    "            # modeling\n",
    "            res = sm.GLM(y, X, family=sm.families.Gaussian()).fit()\n",
    "            scores.append((f, xic(res,c))) # store the score (feature_name, score)\n",
    "        # select least RSS feature\n",
    "        best = min(scores,key=lambda tup:tup[1])\n",
    "        best_f, best_score = best[0], best[1]\n",
    "        if best_score <= xic_curr: \n",
    "            xic_curr = best_score\n",
    "            X_curr = pd.concat([ X_curr, X_full[[best_f]] ],axis=1)\n",
    "            remaining.append(best_f)\n",
    "            removed = [f for f in removed if f != best_f]\n",
    "\n",
    "        # if forward and barkward didn't improve, break\n",
    "        if xic_curr == xic_prev:\n",
    "            break\n",
    "\n",
    "    # list of selected features\n",
    "    final_features = [f for f in X_curr.columns if f != 'intercept']\n",
    "    return {'features':final_features,'score':xic_curr}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Bostonデータで試してみる。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 213,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# data (using boston data)\n",
    "from sklearn.datasets import load_boston\n",
    "data = load_boston()\n",
    "y = data.target\n",
    "X_full = pd.DataFrame(data.data, columns=data.feature_names)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 218,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'features': ['LSTAT',\n",
       "  'RM',\n",
       "  'PTRATIO',\n",
       "  'DIS',\n",
       "  'NOX',\n",
       "  'CHAS',\n",
       "  'B',\n",
       "  'ZN',\n",
       "  'CRIM',\n",
       "  'RAD',\n",
       "  'TAX'],\n",
       " 'score': 3021.7263878250615}"
      ]
     },
     "execution_count": 218,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "forward(y,X_full,c='aic',dist=sm.families.Gaussian())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 216,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'features': ['CRIM',\n",
       "  'ZN',\n",
       "  'CHAS',\n",
       "  'NOX',\n",
       "  'RM',\n",
       "  'DIS',\n",
       "  'RAD',\n",
       "  'TAX',\n",
       "  'PTRATIO',\n",
       "  'B',\n",
       "  'LSTAT'],\n",
       " 'score': 3021.7263878250615}"
      ]
     },
     "execution_count": 216,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "backward(y,X_full,c='aic',dist=sm.families.Gaussian())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 221,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'features': ['LSTAT',\n",
       "  'RM',\n",
       "  'PTRATIO',\n",
       "  'DIS',\n",
       "  'NOX',\n",
       "  'CHAS',\n",
       "  'B',\n",
       "  'ZN',\n",
       "  'CRIM',\n",
       "  'RAD',\n",
       "  'TAX'],\n",
       " 'score': 3021.7263878250615}"
      ]
     },
     "execution_count": 221,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stepwise_forward(y,X_full,c='aic',dist=sm.families.Gaussian())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 222,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'features': ['CRIM',\n",
       "  'ZN',\n",
       "  'CHAS',\n",
       "  'NOX',\n",
       "  'RM',\n",
       "  'DIS',\n",
       "  'RAD',\n",
       "  'TAX',\n",
       "  'PTRATIO',\n",
       "  'B',\n",
       "  'LSTAT'],\n",
       " 'score': 3021.7263878250615}"
      ]
     },
     "execution_count": 222,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stepwise_backward(y,X_full,c='aic',dist=sm.families.Gaussian())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[ここ](http://localhost:8888/notebooks/development/gists/BIC.ipynb)と同じデータで試してみる"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 223,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# create matrix of X \n",
    "# for example, if fitting with 3 degrees, return [[1, x, x**2, x**3],[1, x, x**2, x**3],...]\n",
    "def create_matrix(x,d):\n",
    "    \"\"\"\n",
    "    x: predictor variable data(array)\n",
    "    d: dimension (int>0) including intercept\n",
    "    \"\"\"\n",
    "    X = np.zeros((len(x), d), float) # make a matrix size of (data size - dimension)\n",
    "    for i in range(d):\n",
    "        X[:,i] = x**i\n",
    "    return X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 224,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# generate random dataset\n",
    "x = np.linspace(0, 1, 100)\n",
    "y = 0.5 + 0.4*np.sin(2*np.pi*x) + np.random.normal(0.0, 0.1, len(x))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 225,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "d = 10\n",
    "X_full = pd.DataFrame(create_matrix(x,d), columns=['x'+str(i+1) for i in range(d)])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 230,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'features': ['x2', 'x10', 'x3', 'x6', 'x1'], 'score': -169.97748050293859}\n"
     ]
    }
   ],
   "source": [
    "fw = forward(y,X_full,c='aic',dist=sm.families.Gaussian())\n",
    "print(fw)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 231,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'features': ['x1', 'x3', 'x4', 'x5', 'x6', 'x7'], 'score': -176.64499910783337}\n"
     ]
    }
   ],
   "source": [
    "bw = backward(y,X_full,c='aic',dist=sm.families.Gaussian())\n",
    "print(bw)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 232,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'features': ['x2', 'x10', 'x3', 'x6', 'x1'], 'score': -169.97748050293859}\n"
     ]
    }
   ],
   "source": [
    "sfw = stepwise_forward(y,X_full,c='aic',dist=sm.families.Gaussian())\n",
    "print(sfw)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 233,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'features': ['x1', 'x3', 'x4', 'x5', 'x6', 'x7'], 'score': -176.64499910783337}\n"
     ]
    }
   ],
   "source": [
    "sbw = stepwise_backward(y,X_full,c='aic',dist=sm.families.Gaussian())\n",
    "print(sbw)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 235,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# modeling\n",
    "model_sfw = sm.GLM(y, X_full[sfw['features']], family=sm.families.Gaussian()).fit()\n",
    "model_sbw = sm.GLM(y, X_full[sbw['features']], family=sm.families.Gaussian()).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 237,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# visualizing data\n",
    "x_vis = np.linspace(0, 1, 1000) \n",
    "y_vis = 0.5 + 0.4*np.sin(2*np.pi*x_vis)\n",
    "# make X matrix for visualizing data\n",
    "d = 10\n",
    "X_full_vis = pd.DataFrame(create_matrix(x_vis,d), columns=['x'+str(i+1) for i in range(d)])\n",
    "# model prediction\n",
    "y_sfw = model_sfw.predict(X_full_vis[sfw['features']])\n",
    "y_sbw = model_sbw.predict(X_full_vis[sbw['features']])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 239,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 244,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmAAAAGoCAYAAADhOFmUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XdclfX7x/HXjYBiimliwxIsLXMvGmaJK0cWlvZTNCe5\nc5ZfNU1tWKaWlXugliVamis1y4HmStzbHIEr46iFCxTh/v1x1BBQGeecm/F+Ph48HnLu+9z3deAI\nF5/P9bk+hmmaiIiIiIjruFkdgIiIiEhOowRMRERExMWUgImIiIi4mBIwERERERdTAiYiIiLiYkrA\nRERERFzMIQmYYRghhmH8bRjGrtscb2EYxs7rH+sMwyjniPuKiIiIZEWOGgGbDtS7w/GjwAumaVYA\nPgKmOOi+IiIiIlmOuyMuYprmOsMwfO9wfFOiTzcBRR1xXxEREZGsyIoasDeBZRbcV0RERCRTcMgI\nWGoZhlETaAdUd+V9RURERDITlyVghmGUByYD9U3T/OcO52lzShEREckyTNM00vocR05BGtc/kh8w\njGLAPKCVaZpH7nYh0zT1kY6PIUOGWB5DVv7Q109fP33tsuaHvn76+ln5kV4OGQEzDGMWEADcZxjG\nMWAI4GnPpczJwHtAIWC8YRgGEGea5lOOuLeIiIhIVuOoVZAt7nK8A9DBEfcSERERyerUCT8bCQgI\nsDqELE1fv4zR1y/99LXLGH39MkZfP2sYGZm/dAbDMMzMFpOIiIhISgzDwLS4CF9EREREUkEJmIiI\niIiLKQETERERcTElYCIiIiIupgRMRERExMWUgImIiIi4mBIwERERERdTAiYiIiLiYkrARERERFxM\nCZiIiIiIiykBExEREXExJWAiIiIiLqYETERERMTFlICJiIiIuJgSMBEREREXUwImIiIi4mJKwERE\nRERcTAmYiIiIiIspARMRERFxMSVgIiIiIi6mBExERETExZSAiYiIiLiYEjARERERF1MCJiIiIuJi\nSsBEREREXEwJmIiIiIiLKQETERERcTElYCIiIiIupgRMRERExMWUgImIiIi4mBIwERERERdTAiYi\nIiLiYkrARERERFxMCZiIiIiIiykBExEREXExJWAiIiIiLqYETERERMTFlICJiIiIuJgSMBEREREX\nUwImIiIi4mJKwERERERcTAmYiIiIiIspARMRERFxMSVgIiIiIi6mBExERETExRySgBmGEWIYxt+G\nYey6wzlfGYZxyDCMHYZhVHTEfUVERESyIkeNgE0H6t3uoGEYDYDHTNMsCXQCJjroviIiIiJZjkMS\nMNM01wH/3OGUQOCb6+f+DhQwDON+R9xbREREJKtxVQ1YUeB4os9PXn9MJNuw2WyEh4djs9msDkVE\nRDI5FeGLOEBo6Bx8fUtRt25nfH1LERo6x+qQREQkE3N30X1OAo8k+vzh64+laOjQoTf/HRAQQEBA\ngLPiEskwm81GcHBXYmJWExNTHthFcHBN6tSphY+PT4auGxERgZ+fX4auIyIijhMWFkZYWFiGr2OY\nppnxaADDMPyAxaZplkvhWEOgm2maLxmG8QzwhWmaz9zmOqajYhJxhfDwcOrW7Ux09Nabj3l7V2bF\nikn4+/un65qhoXMIDu6Kp6cfV69GEBIynqCgZo4KWUREHMQwDEzTNNL8PEckO4ZhzAICgPuAv4Eh\ngCdgmqY5+fo5Y4H6wCWgnWma225zLSVgkqXYbDZ8fUsRE7MasI+AeXnVJDLyQLpGrhx9PRERcZ70\nJmAOmYI0TbNFKs55yxH3EslsfHx8CAkZT3BwTTw8fImLiyQkZHy6k6WIiAg8Pf2uT2cClMfDw5eI\niAglYCIi2YTDpiAdRSNgklU5qmZLI2AiIlmHpSNgImIfCXNEguToETUREcl8NAImkklpFaSISOZn\naRG+IykBk5xKCZeISNaT3gRMjVhFnCQtnfHVyFVEJGfRCJiIE6Slj5eK7kVEsi6NgIlkEok740dH\nbyUmZjXBwV1vOxJ2o+2EPfmCxG0nREQke1ICJuJgaU2o/Pzso2Sw6/oju4iLi8TPz8/JkYqIiFWU\ngIk4WFoTqhttJ7y8auLtXRkvr5pqOyEiks2pBkzECW7UgCXu43W3vRy1ClJEJOtRGwoRB8toQpT0\n+UqwRESyHxXhiziQI9pC+Pj44O/vj4+Pj9pMiIjILTQCJpKEo9tCqM2EiEj2pREwEQdxdFsItZkQ\nEZGklICJJOHothBqMyEiIkkpARNJwtFtIdRmQkREklINmMhtOHrVolZBiohkP2pDISIiIuJiKsIX\nERERySKUgEm2ZbPZCA8Pv+0m2CIiIlZRAibZkhqfKgEVEcnMVAMm2Y4an/63F6Wnp70FRmr2ohQR\nkbRLbw2YuzOCEbHSjcanMTHJG586OwG7HHeZzSc3sydqD3+c/YPI6EjOxZwjOjYaE5NcRi7y587P\nA/keoGj+opTxKUP5+8tT7v5y5PXI65AYbDYbwcFdiYlZff1rsIvg4JrUqVMrxySgIiKZnRIwyXZu\nbXxqT0Cc1fg0NhZW7NrLlss/sPzIcnb9vYtyRcpR8YGKFPd+nD9X16RkvkI8cO+9PHC/wcPFrlHA\n5wJnYk9z/Pxx1h9fz8StEzl45iD+Rf2p+2hdXnniFcoWKZvumByVgKpthoiI8ygBk2znRuPT4OCa\neHj4EhcX6bDGpydOwMqVsHL9Pyy3TcP2cAjueS/SrWZThtUaxtNFn+Yez3sAuHoVPLbA2bNwdh/s\n+BkOHrQ/fuwYGIkGrC9evciaiDX8cuQXGnzXgPu87qNV+Va0qdiGwnkLpylGRySgmsIUEXEu1YBJ\ntuXoERzThNLV/uTq08P5q9D31CzaiH61ulDd7xncjNSvZ7l2Ddzdk8cXHW0fUfMpksCaiDV8vfNr\nFh5cSPMyzen9bG8ev+/xVN/jRgKVOAFNbQKlGjoRkdRTI1YRJzp5/iQfrf2I7/d9T9eqXen2VDce\nyPdAuq+X0ghToULNaNYMAgKgQwdo0ACiLp9m3OZxTNw6kYYlG/JBwAf43uubqnukNwENDw+nbt3O\nREdvvfmYt3dlVqyYhL+/f1pfqohItqYETMQBLl6EqVOhUCFo3Rquxl9l9MbRjNgwgvYV29Over80\nTwkmdacRpjx5fPjhBxg71j5V2bcvBAVBTEI0ozaMYvyW8bSr2I4hNYaQP3d+h7zmtMSnETARkVup\nE75IGiXukxUbCyNHwqOPwvr1UK4crDu2jkqTKrEmcg2b39zMyBdH3pp8mSacOweHD8OWLbBqFaxe\nDWvXwsaN9oKvf/+1n5fIjSJ5e3IDiYvk8+eH9u1h61b44gv49lvYsQMK5CnAh7U+ZG/XvZyLOUeZ\n8WVYcGDBbV9PRmjzcBER59MImOQYiafkVqxYdb1Gyo+YmCoUKDCaatXuYfhweLTkFQavHszMXTMZ\n02AMrz35GsbFi7B5sz2x2rnTnnQdPgy5ctmHywoWBG9ve7IVHw9xcfbq+7//tg9lPfoolC4NpUsT\n/dhjlO3cixOxa0jvCNOaiDV0XtKZ0j6lmdRoEr8uXOnwonmtghQRuTtNQYrcQeKaqytXjpKQYHL1\n6lrsCdBZPD3bcuLENM4aZwmaF0Qx72JMKd6DIr+sh8WLYf9+qFgRnn0WKleGkiWhRAm499673/zy\nZXuytm8f7N0L4eFcXbuWw7FXWedekEXGFVpNnUCzVi3T9JquXLvCoFWD+G7Xd5yddpGr+9ehKUMR\nEddSAiZyG8lrmmYBQ4BDN8/x9q7MwG+bM3Lvp3wU+ywdv96DkcsdAgPh5ZehWjXw9HRcUNeu8U9Y\nGJfmz+f+33/H4+hReOklaNMGatUCtztXB/TrB8WKQZcuMGH5OLqv7I25uwes+AQSPFQ0LyLiIkrA\nRG4j+ao+G+AHbATKg7GNx+s8Q0wVN+Yuys1TtdtAu3b2ES8jzf+n0ufUKZg7F6ZNg+hoeyFYx45w\n//0pnr5vn32lZO7cMHr0WZ6pWYLYhmXB3QN+GIyX2UQjYCIiLqAETOQ2bDYbxYqVIjZ2J/AwsAsP\nj+fwzOXBq575OP7Sca555+bHZz6jSPNgyJPHumBNE7Zvh0mT4Pvv4f/+D95+Gx5P3gMsPh6GD4ev\nvoIWLdYzcfLLJNTwJK6UjWHlhjGgXX8LXoCISM6iBEzkNv79F+rUOcm2bVHcc89g4q+tZ2n3jpRa\nPo+XA05Rwq8K33RfgYe7A6cYHSEqCsaNgwkToE4d+OADe91ZEps3Q8uW8NZbF6hW7QC7ru2i/7r+\nTH15KoGlAi0IXEQk51ACJpKC/fvtZVz16kH//jb+CfuFxydO5GTsX9RvEktT/9Z8VGsYhqumGtPj\nwgX48kt7X4qmTWHIEHjwwVtOOX/eXjaWL5/98/CT4QTODmRIjSF0qtrJgqBFRHIG9QETSWLhQqhR\nA959F8Z8epmiX42ibO/eHG1Wm+eDYulRawDDan+cuZMvgPz5YdAge1+x/PntTco++8ze6uI6b+//\nki8A/6L+/NbuN0ZtHMV7q95Df9SIiGQuSsAk2zp5EpYsgbaPrLQnLcePs3/tPGpfncKwWsPo9lQ3\nq0NMm/vus3eL3bABfvnFvkggLOy2pz9W6DHWt1/Pz0d+5q2lb5FgJrguVhERuSNNQUr2deWKffhr\nzhyYMoW9VYpRd2ZdRtQdwRvl38jw5V3dqPSW+xUuDPPnQ69eUL++fUQs/61bE128CD/8AK82j6bh\nrAaUK1KOCY0mpGnjcBERuTNNQUqW46itc1K0fz888wwcPQo7d3LAvzh1Z9Zl1IujHJJ8hYbOwde3\nFHXrdsbXtxShoXMcEHQa7jf7e3jtNdizBxISoHx5+zZIiVy8aC8bG9CnAEuaL2ffmX10WNSB+IR4\np8YqIiKpYJpmpvqwhyTZ3axZs00vr0JmgQKVTS+vQuasWbMdd/HvvzfNwoVNc9Ik00xIMI/9e8ws\nNrqYOX37dIdcPioqyvTyKmTCTtPeN2Kn6eVVyIyKinLI9dN1vyVLTLNoUdPs2dM0Y2NvPhwdbZp1\n6phm06amefbCBTNgRoDZZn4bMz4h3imxiojkNNfzljTnOxoBE5ez2WwEB3clJmY10dFbiYlZTXBw\n13SPhJ09Cw0awF8n4qF/f/jf/+w1Uh07cibmLC9++yI9n+5J24ptMxx3eHg427dvv+1m2s5wp827\nb2rYEHbvhmPH4Lnn7CN/2Ivzf/rJvh1lm6B8zHttCUf/OUrPZT1VmC8iYiElYOJyqUooUikqCmrW\nhKdKRfPAmy/Bli0QHg6VKnHhygUafNeA10q9Rp9n+2Qo5sRTgIGBzYiJOQLsun50F3Fxkfj5+WXo\nHrfj52ffYPuu9ytYEObNg9at7dOvc+cC9m75c+fCPffA1Al5WRy0mPXH1zMkbIhT4hURkbtTEb64\nXPK9GdO3efRff0Ht2hBc7wR9VjbEeOEFe9GTuztx8XG8NOslit9bnImNJmao1URK8Xp4PIe7uyce\nHr7ExUUSEjKeoKBm6b7H3dzYTDzV9wsPh2bN7HVin34KuXIRH28vF/PwgKhLUTw//Xk6V+lM72d7\nOy1uEZHsTo1YJUtJc0KRxJkz8MIL0LPWbjotegm6d4d33gH7fwS6LOnC8fPHWdh8Ie5u7hmKNfle\nkvbNu3/4YTgFCxa0ZhVkau539iy8/rp9a6XQUChQ4JbDx6KP8fz05/mw5oe0rtDaSVE7jqtXnYqI\npIYSMMlyMvILdfZsuLBsHR2WvWbfDLF585vHPt/4OTN2zGB9+/Xkz53/DldJfZyOGLGzRFwc9O4N\nK1fCokVQsuQth/fb9hPwdQCzm8ymZvGaFgV5dzcSdk9P+3Sss0ccRURSSwmYZHoOHcFYtco+xTZr\nFtSte/PhhQcW0nVpVzYGb6RYgWIZjPg/GR2xs9ykSTB4sH2D7xo1bj584QJs+nsVbywIYk3bNZQq\nXMrCIFOWpRNgEcn2LO0DZhhGfcMwDhiG8YdhGP1SOO5tGMYiwzB2GIax2zCMto64r2QdjuibdWMV\n4r9z5thHvH744Zbka+fpnby5+E0WNFvg0OQLICioGZGRB1ixYhKRkQeyVvIF0KmTPVl9/fWbxflg\n31Zy1rBafFJ7OC/NeomoS1EWBpkyRy7aEBHJLDI8AmYYhhvwB1AbOAWEA81N0zyQ6JwBgLdpmgMM\nwygMHATuN03zWgrX0whYNuOIEYwbI1AvGwUZc/koO4cMpe7QwTePn4s5R9XJVRlWaxhB5YKc80Ky\ngx074KWXYMAAeOstLl60ryJt0AASAgax8s+VrGq9Ci8PL6sjvUkjYCKSmVk5AvYUcMg0zUjTNOOA\n2UBgknNM4EYxTn7gbErJl2RPGR3BsNlstG8/hGoxnzDm8nkaMZPAEV/e7BsWnxBPi3ktaFyqcY5K\nvtK1k0DFirBuHYwZA+++S757TJYssQ+OPfzHh/jd68ebi9/MVD3CfHx8CAkZj5dXTby9K+PlVZOQ\nkPFKvkQkS3NEAlYUOJ7o8xPXH0tsLFDaMIxTwE6gpwPuK1lEqvtY3cbIkZepcnUgoQykCfMIp+Ut\nCdzQsKHEXovl0zqfOiH6zClDU7rFi9uTsBUroFs3ihRO4Oef4f33DZq4h7DPto8vNn3hvODTIctP\nAYuIJOGqRqz1gO2maT4EVALGGYaRz0X3FotlZARj0SIIDznHvIT2tOZ91vE8iRO4BQcW8PXOr5nT\ndA4euTyc/2IyAYfsJODjY0/Adu+Gjh0pUTyehQth0295md9sPp+u/5RVf65y3otIBx8fH/z9/TXy\nJSLZQsYaJNmdBBJXPD98/bHE2gGfAJimecQwjD+BUsCWlC44dOjQm/8OCAggICDAAWGKlYKCmlGn\nTq00rYIMD4cP2x5hvXsjwrt3Zc3U9/D2mHpzFeK/bv/ScXFHFgct5v5897vgVWQON6Z0Y2KST+mm\nKTnx9oZly+CVV6BdO56aNo2nnnIH/JjVZBYt5rXg9zd/x/deX2e8DBGRLCksLIywsLAMX8cRRfi5\nsBfV1wb+AjYDQaZp7k90zjggyjTN9w3DuB974lXBNM1zKVxPRfjCqVPwYpWzbHSrRv73ekPnzre0\nsfAu6M2zIc8SXCmYbk91szpcl3J4Ufrly9C4MRQqBDNn2lvlA6M3jmbmrpmsa7+OvB55HfoaRESy\nC0v7gBmGUR/4EvuUZohpmsMNw+iEfYfwyYZhPAjMAB68/pRPTNMMvc21lIAJi+deocL/XqRY06dh\nxIhkx3ss68GJ8yeY93/zMrTNUFbl8L5ksbHQpIm9W/7MmZArF6Zp0mp+KzxyeTA9cLrjghcRyUbU\niFWyD9O0bygdE2NvHOp2a6niggML6L28N9s6bqOgV0GLgrReRhvbJnt+TAw0amQv0p88GdzcOHfx\nIs997c+A6gOyxHZFIiKuZmkjVhFHuNFW4VK/fvDHH/DNN8mSr8h/I+n0UydCm4Tm6OQLMlaUnuIq\nSi8vWLgQ9u+3b19kmjR/LR8d7v2et395m/22/Xe/sIiIpIpGwCRTuDGl1tQowEcxkWwdP55XO3e6\n5Zy4+DhqzKhB41KN+d9z/7Mo0qzrxohXvnz5qFKlerIasq1b13Hx4kWKFyxI4WbNoF49tjb5mPr1\nodfXIYRGjGZzh82qBxMRSURTkJJl2Ww2ihUrRfHYqYTRiQZ8xX6vbsmKygeuHMi209tY0mIJbobj\nB28duldlJpN4M+vY2MO4ufkSE7Pr5vE8eR7FNP8hT55HuXo1gpmjh9Pkqy/gjTeYUngAo78wKT+k\nFffkzk1IYIiFr0REJHPRFKRkWYcORZI3bg7z6U9fRrKN5sk65a87to5pO6YxI3CGU5IvR+xVmVkl\n7Rt25cpCYmIO819j3DBiY09z5cqam33FWvXuz9nZs2HyZN50n8EzTxvEL5zI+uPrmblzppUvR0Qk\nW1ACJpb7ZkYpvo4fyK9U4hvakLRT/oUrF2g9vzUTX5rolH5fDmlsmokl3woqgDx57id37hp4e1cm\nd+5AvLxKkHSrqKOxsbBsGcaA/kx4ZRkHduWmifkpvX/uzZFzRyx5LSIi2YUjGrGKpNvcufD4D6Op\nWuICLU/8grdn5ZttFW5MA/Ze3ptaxWsRWCrpFqOO4bDGppnUrVtB2Wu+DOM827Zt4OLFizdrwhIf\nv5kA+/jAjz9i1q9PgStejHvvYS6Xj6E+Ddjfdx/ubun/EZKdp3xFRO7KNM1M9WEPSXKCQ4dM85UC\nYeaV+x4wzZMnzaioKHPz5s1mVFSUaZqmGRUVZY5cPNL0/dzXPB973mlxREVFmV5ehUzYadp7YOw0\nvbwK3YwjO5g1a7bp5VXI9PauZHp5FTJnzZqd6uNRUVHm6575zJP4mI9y2MTYbrq1cTf/t+R/GY6n\nQIHKKcYjIpJVXM9b0p7vpOdJzvxQApYzJCSYZr0qNvPCvQ+b5rJlyY7PmjXbzHPfvabR1930LJHf\n6b+g75agZAdJE9zUHt+8ebNZoEBlsyMTzT8oYd6Hzcz3YFmz0MeFzE3HN6Urjuye8IpIzpHeBEyr\nIMUapsml2i+Tt0ppjJG3drq32WwU832C2Mbl4O/qsKpZxrbaSSVNiaUs8dZHnzCLZ1nBK3mOMnr1\nKD4J/4TtnbaTzzNfqq8XHh5O3bqdiY7eevMxb+/KrFgxCX9/f2e8BBERp9EqSMlaRo/mnstnMD4e\nluxQREQERuX84H0B1gwhcU2WM2WksWl25uPjQ0jIeLy8ajI8/3LOu+1i2zMVqOTZjjL5qtP7595p\nut6tNWmQdNGFiEhOoARMXG/bNhg+HEJDb278nJhbITdinjsOPw6CeE/0C9p6QUHNiIw8wK8rJ/PM\n4YM89s8/xI/6nH2ff8XKo6tYcGBBsufc2Nkg6WrSxAmdt3dlvLxq3rLoQkQkJ9AUpLhWbCxUqQID\nBsAbbyQ7bJomL377IoXP+7Cw73LHbTYtjnXsGDz7LJ89Pomw0gXZ8mhTdnfZTeG8hYFbG79evRqR\n4vdPU74ikh2oE75keocOQcnJfeHPP+GHH8BI/n6dsnUKk7dNZmPwRv45+49+QWdmv/9OQqOXecVr\nBfT9Gu+H/mJWk1m31Iwl3urI2TV8IiJWSG8Cpj5g4hLbt8OQmmuZn+c7cu3emWLydSz6GO+uepfV\nbVbj7uaOj4+PfmFnZk8/jdtXXzLv7UDKDltLfL/azN8/n4cvPpyt+6qJiDiCasDE6S5fhjebXWCW\nZ1tyTZlkb+6ZhGmadFzckZ5P96RskbKpuu7taozEhYKCyP3G6/xaIJgej0yl29JueD/grSJ7EZG7\nUAImTte/P4ww3yHfyzXh5ZdTPGf6julEXYqi33P9UnXN7Lx3Y5bz8cf4PZJAr63LaF62Oe9vfj9d\nRfZKqEUkJ1ENmDjV6tUw+fVf+TZvB3Lt2QXe3snOOXH+BJUmVWJFqxVUeKDCXa+pGqNMyGYDf38u\nj/yYiqeH8mmdT6leuHqqa/hSU7QvIpIZqQ+YZEqjP7pEiHsnck2emGLyZZomnX7qxFv+b6Uq+YKU\nNpd2TZ8wuQMfH5g7l7xdezKtwnt0W9oNt3vcUtVXLbtvhi4ikhIlYOJU8ysMIW+dalC/forHZ+6a\nyYnzJxjw/IBUX1ONPDOpqlVhxAiqd/iIZiUb0/Pnnql6mhJqEcmJlICJ82zdSq7vZsLo0SkePnXh\nFO/88g4zAmfgmcsz1ZdVI89MrF07qF2bjh+dYM2RjSw7tOyuT1FCLSI5kWrAxDni4uCpp6BPH2jV\nKtlh0zQJnB1Ihfsr8GGtD9N1CzXyzKSuXOHsE9X4X5GnWdViGXu67uEez3vu+JQbNWBqvCsiWY0a\nsUrm8umnsHIlLF+eYs+vWbtn8fFvH7O141Zyu+e2IEBxJvPwEc6Xfob6vf15rnZpRr046q7PUUIt\nIlmREjDJFHbuhJK5jpI34CnYvBkefTTZOVGXoig/oTyLgxbjX9TfgijFFc5MmkfkO32o924My9ss\no8pDVawOSUTE4ZSAieWio6FMGdjp+wr3NXrWvt9jCoLmBfFw/ocZ+eJIF0corra7Vg++u7KeX1qZ\nbO64GXc3bb4hItmL2lCI5f73P+hXbgn32Q7Ya79SsOjgIrac2sL7Nd93cXRihTJLRtL1MBSIusSX\nm760OhwRkUxDI2DiEGvWQPsWsfyRuyy5JoyDevWSnfNv7L+UHV+W7177jhp+NSyIUixx9CiH6/nz\nTLt4wrtsp3jB4lZHJCLiMJqCFMvExED58rDI/0OejN0OP/6Y4nkdFnXA3c2dCY0muDhCsdy8eQz/\nugNhr1VmWZtfMVJYmCEikhVpClIss2ABvPh4BE/+8uVte36tPLqS5UeW82ndT10cnWQKTZrw9oNN\n+OuPbYTuCbU6GhERy2kETBwiIfBV3PyrwKBByY5dunqJchPKMbbhWBqWbGhBdJIZ2CIi2Pv6s7R4\nJYb9fSMpkKeA1SGJiGSYRsDEqWw2G+Hh4Snvz7dyJW57dsE776T43EGrBvFcseeUfOVgoaFz8C1d\nhcH7H6Lmlou88UWg1SGJiFhKCZjcVWjoHHx9S1G3bmd8fUsRGjrnv4Px8fYVjyNHQp48yZ676cQm\nZu+dzeh6KU9NOtIdk0SxTOLNtn+7tJV7V/Rn/dm1rNr9q9WhiYhYRglYNuGs5CPxL8/o6K3ExKwm\nOLjrf/eZMQMKFIBXX0323CvXrtB+YXu+rP8lhfMWdmhcSd0xSRRLJd1se8LlD2i76jHeCX2DBDPB\n0thERKyiBCwbcGbykfSXJ5TH3d2XDRtOwcWLMHgwfPZZitsNDfttGE8UfoLXS7/usHhSctckUSyV\ndLNtkz38sKMn7mfOMW3uQEtjExGxihKwLM7ZyUfSX56wiytXyjNoUGnMT0dAzZrgn3w7oZ2ndzJx\ny0TGNRzn9JYDKSWJHh6+REREOPW+OV1qR119fHwICRmPl1dNvL0r4+VVk0Hji/LA6hEM3DKCM3//\n6aKIRUQyD+0LksXdSD5iYpInH47Y0PjGL8/g4Jp4ePhy9eo58ubdx9QhpzE6jYPt25M951rCNYIX\nBTO8znC05yYkAAAgAElEQVQeyv9QhmO4m1uTxPLALuLiIvHz83P6vXOq0NA5BAd3xdPT/rUPCRlP\nUFCz254fFNSMOnVq3bLZtpcXeIeFMGBkfaaMOujC6EVErKc2FFmczWbD17cUMTGruZF8eHnVJDLy\ngEMSsMT3iYiIYOrUMly7lpeQuNbwyCMwbFiyc0esH8GvR3/llzd+cVnDzRsJgYeHL3FxkXdNCCT9\nbvee27p1HRcvXryZYKVG9JkTlB7px9yKn/BsUF+nxi0i4gzqhJ+DuSr52LYNGjSAA7N3ULBFA/jj\nD8if/5ZzDpw5wPPTn2fzm5tdvuXMjSQxLQmApF14eDh163YmOnrrzcfy5HkU0/yHPHkeTdWIWGKz\nfhjMyLWfEP7ecdyLPOCssEVEnEIJWA7niuSjdm1o2RLaz21oz8S6d7/l+LWEazw37TnaVmhLF/8u\nTolBrJd8BCwMaAhsIj2jsKZpUnuwL4G2QvScsD3FBR0iIpmVGrHmcD4+Pvj7+zt15Ofbb6HtY7/B\n/v3QsWOy4yPXj8Q7tzedq3Z2WgxivaRF9blzB+LlVYL0LoIwDINxnRbxUcE9RM2c6KywRUQyFY2A\nSeqZJjz/vD35at36lkO7/95NrW9qsbXjVooVKGZRgOJKN0Zd8+XLR5Uq1TNch9j6s1a4b57LtM8P\nQ9GiTotbRMSRNAImzrd0Kfz7r30eMpG4+DjaLGjD8NrDlXzlIDdGXZ988slkbSZCQsaneTT2n43j\nWPSYO5t7NrUn+yIi2ZhGwCRFyWrKEhKgUiV4/31o3PiWc98Pe5/fT/7OkhZLXLbqUTKfjNYhRkZC\n+ZYhPFGlB5ueHIVbZ+vrCLWwQ0TuRiNg4jA3OuvXqdOVYsWud9afPRu8vCDw1k2Ut/+1nXHh45jy\n8hQlXzlcRusQfX2hb4N2HHMvwTff9oUjRxwcYdpoeysRcSaNgMktkq9wiyJ/nlL880ABcoWEQK1a\nN8+9cu0K/lP86VutL/Ufqq+RAsmwK1egZI1wYl6szeENFSjwyxpwc/3fia7qryciWZ9GwMQhkm/r\nU4T25ONikSK3JF8AA1cNpEShEuTa66GRAnGI3Llhyvv+5DrdlA+KHobJky2JQ9tbiYizaQRMbvHf\nX/6HgEJ4sJXDxlPkX/ITBRs0uHner0d+pf2i9qxouoJKT1TTSIE4VMSZv/Gf8iRrpyXw5MpdUMy1\nizs0AiYiqaURMHEIHx8f+vWbA8SSP391Onq8gEfZsrckX2cun6HdwnbMCJzB+dPnNVIgDudX+H4G\n1hxMr1Y+mB07uHxVZEobiKdnZaeIyO04ZATMMIz6wBfYE7oQ0zQ/TeGcAGA04AHYTNOseZtraQTM\nQvHx4O8PHTuep2r5PVRq3pxcs2dDtWqAvWv5q3NepWShkox8caRGCsRp4uLjqDChPB8vjaXx/w2B\ntm1dHoNWQYrI3Vg2AmYYhhswFqgHlAGCDMMoleScAsA4oJFpmmWB1zN6X3EO04S334ZOnbypum8f\nuZ544mbyBTBl2xSORR/jo1ofARopEOfxyOXBVw3H0KdmHDED+sKpUy6PwRU7TIhIzpThETDDMJ4B\nhpim2eD65/0BM/EomGEYXYAHTdMcnIrraQQsM4iLg8cfh5kzoXp1APZG7SXg6wDWtl3Lkz5P3nK6\nRgrEWV6d3YRKe88xeJc3LFigvSJFJFOxsgasKHA80ecnrj+W2ONAIcMwVhuGEW4YRisH3Fec6Ztv\n4NFHbyZfF65coMn3TRhVd1Sy5As0UiDOk2/9Z3zqvotjp/bDHGtX2NpsNsLDw7HZbJbGISJZn6uK\n8N2BykADoD7wnmEYJVx0b0mruDgYNgyGDAHsdV8df+rI88Wep03FNhYHJznN4B5+mJvfokcrX+jV\nC6KiLIlDjVlFxJHcHXCNk0DiNeIPX38ssRPAGdM0Y4FYwzDWAhWAwyldcOjQoTf/HRAQQEBAgAPC\nlFT79lt7W/IXXgBgfPh4Dpw5wIb2GywOTHKikiWhc9l+TLKV4re2tXi+e3eXj4TZbDaCg7sSE7Oa\nmBj7YpPg4JrUqVNLo74iOUxYWBhhYWEZvo4jasByAQeB2sBfwGYgyDTN/YnOKQWMwT76lRv4HWhm\nmua+FK6nGjAXmzsX8uWD+vWxL4MsVQqmTIGAADaf3EyjWY3YGLyRxwo9ZnWokkNduADFXprF/YEj\n2TfpIm4jRyXbFsuZwsPDqVu3M9HRW28+5u1dmRUrJuHv7++yOEQk87GsBsw0zXjgLeAXYC8w2zTN\n/YZhdDIMo+P1cw4Ay4FdwCZgckrJl7je+fPQowcULHj9gXnzoEgRqFGDvy78RdPvmzL55clKvsRS\n+fPD6PZBnDjuxfShr8Fbb9nfvC7i5+fH1asR2H+EAewiLi4SPz8/l8UgItmLOuHncP37w+nTMGMG\n9h4UVarA++8TU78ONWbUIPCJQAa+MNDqMEVISICv5m1mxLHGHDxQh/x5CsCYMS67f2joHIKDu+Lh\n4UtcXCQhIeMJCmrmsvuLSOaU3hEwJWA52JEj8PTTsGsXPPQQ8Msv0KcP5s6dBM1vSS63XHz76rcY\n6Vz2r9YU4gyt57emqEchPukwB378EZ591mX31ntaRJJSAia3SM0viiZN7ANe7757/YFataBdO95/\n5CjLDi8jrG0YedzzpOv+N0YLPD3tUzcaLRBHOXn+JOUnlmfLA0MoPnIKbN0Knp5WhyUiOZT2gpSb\nUrNc/sIFuHoVeve+/sDvv8ORI0x+4iJf7/yaBc0XpDv5SrxiLDp6KzExqwkO7qreSeIQRb2L0vuZ\n3vTNs9a+SfeoUVaHJCKSZhoBy2bSvTfja6/xw3P30stYzpq2ayhRKP1t2rRiTJwtJi6GJ8c9ydfP\nfkqNl7rBhg32nRtERFxMI2ACQEREBJ6eftiTL4DyeHj4sn379tt38N6/n1+OreKthJ9Y2mJphpIv\n0IoxcT4vDy8auH9K+7WfED/wXejUyb6IJIPU6V5EXEUJWDaTUvITE3OYxo2DbjsluXRMD954+Rrz\nmv1IhQcqZDgGbdAtrtDQ9//4KyI/U57NBxcvwvTpGbqeOt2LiCtpCjIbSrxc/urVP0lIMLl6dS0p\nTUnOXzeFTos7sajNzzxT+kWHxqEVY+JMpglPNd7EvoqvsPe5EPzeCIbdu+H++9N8rXRP3YtIjqdV\nkHKLG8nPP//8w//934AU67HCCefDn/7HkkuBVB71nYXRiqRdaOgc2rX7nCv1S5DryjJ2eVSndL58\nMGtWmq+lukURSS8lYJKi5H/ZnyJP3rK0nPkaG06uY/Hov3ls5TYoXtzqUEVS7Zb3dT4f6FqGAt/E\nczZ3QXJNmAANGqT/ehoBE5E0UBG+pOhGPZan51AM4zi5i1TEb8gDnLh8go0J7XmsSh0lX5Ll3LLY\n5OKDsKEvl2q7cfidd6BLF3tNWBqoblFEXE0jYDnAtWtQrvxVyr0xmtXuo3in2ju883RvcpV6EmbO\nhGrVrA5RJE2SjVi5h2O89QxzW//Aa5/PhwcfhBEj0nVd1S2KSFpoClJSZJomb3/1GyGR/Sj1ZALj\nXxpPlYeqwIIFMHw4bNwI6dxqSMRKSfdm7PRFG1bGr2B745/JVaEirFwJ5cpZHaaIZHNKwOQWZy+f\nZeHBhUwIn8T2/ecYXHMgg15ujZtxfdb5hRegWzdopu2BJOtKPGJVuHBhAr4OoEXZFnTaCnzzDfz2\nG7ip0kJEnEcJWA735z9/8tHajzgXe44DZw5w4vwJ6j1Wj5blWlLW4xVKlsj138lbttg3gjxyBNzd\nrQtaxMG2/7WdBt814EDXfdxb+yUIDoY337Q6LBHJxpSA5XBRl6JYeGAhBb0KUrJQSUr7lMYjl0fK\nJ7dsCZUqwTvvuDZIERcICu1AwbzejPdrDXXrwt69oHouEXESJWCSOidOQPnycPQo3Huv1dGIOFzH\nPn8zM38ZdvbYwOMfT4KzZ2HGDKvDEsmU2s9rT7k85XjjqTe08CSdlIBJ6vTvD7Gx8MUXVkci4hRn\nz4JvixH4N1nH6qDvoEwZez1YQIDVoYlkKjO++4Z2e9qS/5uyXPvnJCEh4wkKUl1wWikBk7uLiQFf\nX9iwAUpkbMNtkczs8y+vMOivMix4cwIv7r4I774LO3eCp6fVoYlkCjabjaL1ixP3ZEX4bh1qPpx+\nasQqN02caP9IZs4cqFpVyZdke9275qZg+Cg6zOvFtVca2d/zo0ZZHZbL2Ww2wsPDsdlsVocimUxE\nRARmuVywq+v1R8rj4eFLRESElWHlKErAsplz52DwYKhePckB04SxY+GttyyJS8SVPDxgcp9AzkY+\nyIQtk2DMGPj8c3vtYw4RGjoHX99S1K3bGV/fUoSGzrE6JMlEvIp4ce3BC3DgseuP7CIuLhI/Pz8r\nw8pRNAWZzfToYe98P358kgO//w4tWsChQ+qLJDmCacLK3XtosbwW+7vt574xUyEsDJYuzfbNh7W3\npdzNyPUjWbJ5KZsH7rrZzFg1YOmT3ilINYHKRvbtg9BQ2L8/hYNjx0LXrkq+JMcwDKhTviyvH3+d\noWFDGdPnc/vWW/PmQdOmVofnVDf2yoyJKX/9kf+ml5SAiWmaTNsxjSlNpvDEG09o+y2LaAQsmzBN\naNAA6teHXr2SHIyKgieesDdeLVTIkvhEXCGlvRzPXj7Lk+OeZFWbVZQ9+A8EBdn/WvH2tjha59EI\nmNzJphObaD2/NQffOoiRzUeDXUFF+Dnc1atQoYJ9kCuZqVPtne+VfEk2druap/vy3sd7L7xH7+W9\nMatXh3r17IWS2ZiPjw8hIePx8qqJt3dlvLxqEhIyXsmXADBt+zTaV2pvT75OnoS4OKtDypE0Apbd\nXbsGxYvDokX27ve3kdLIgUhWcbsRn61b13Hx4kWKPlKUOvPqMLzOcF65r5q9N9iyZVC5stWhO5X+\nX0vS90DEqQgqzKjAuhbrKPPAE7g97Q9Dh8Krr1odapalETBJ2aJF9t5fd0i+tFpKsrobNU/25Aug\nPKZZgEqVqlG3bmdKPFqOe9b1ofP8Ply5Nz8MHw6dO0N8vIVRO5+Pjw/+/v5KvnKopD/bu3fvxeOB\nZbl8EJ56sg5fPTiYSwUfhsaNrQ41R9IIWHZXqxZ06GCve0mBakUkO0j+Pg4DGgKbuPG+9vT8GI/W\nlxjY6nkGVH/H3hm/WTPo1s26wEWc5Lb/J9qWgd/7U3Z/AVbyOm471lC4Qvk7X0zuSCNgkty+ffYl\nkU2a3PaUlEYO1IxPspqkNU+5cwfi5VWCxO/r3Ln/oMSR9xkWNoLTl6NgwgT71Mvp0xZGLuIcyX+2\n3wMF7wefSPjjJSbzHh95PsqfV69YGGXOpgQsO5s4Ed58847br/j5+XH1agSw6/ojasYnWVNQUDMi\nIw+wYsUktm/fBJwk8fv62rVIvhrqR8K2dryzbKC9Diw4GN5+28KoRZwj+c/2S1DxBOyqR+f46SSQ\nQIjbn/pZbyElYNnV5cvw3Xf2BOwOtFpKspMbNU9PPvlkiu/rF14oRFDRQczfs5Qtp7Zg69yZK6tX\n8+/cuVaHLuJQSX+258n7Gvmq5+HBPVv4gPfo4fk3U6dN0M96C6kGLLv65huYPdve9TsVtFpKsqOU\n3tc2G3ScEMK+3CM59n4UgUZBPoyJYNv0GTRr08riiEUc68b/gSPGET7b/hnLvy/KwYQHKTH7A/2s\nd5D01oApAcuuqleHvn0hMNDqSETSxZl/FJz++zQPDX0Yc90w2NOPH6nFHvdNdD4VqV9Kki0Fzg7k\n5St+vPnhEti9G7y8rA4p21ARvvxn717480946SWrIxFJF2e3Rjl+7Dh51z4GdceBx2V6MoOu8Vf5\n67ffHHofkczgePRxfov8jaCPFtg3ClbylSkoAcuOpkyB9u3BXVt9StZjs9kIDu5KTMxqoqO3EhOz\nmuDgrthsNofdw8/Pj4Q/z8Dx0vDcCI7zL6NyeVJqzBj7vl4i2cjUbVNpeekx7qlaDV580epw5Dol\nYNlNTAx8+619ddcd2Gw2wsPDHfpLTcQRXNEa5UaBcu61mzCe/ojcRV6gwrQpeJ49C99/77D7iFgt\nLj6OqeGT6DT7EIwaZXU4kogSsOxm3jzw94c7LC1W53vJzFzVGiUoqBnH9xzixXs7U+fTWjRv1dLe\nG6xPH4iOdui9RKyy+I/FFD0eR7Gm/aFoUavDkUSUgGU3kydDx463PeyK6R2RjHBlaxQfHx+KRn7K\nmj/D+S3yN3juOWjYEN57z+H3ErHCFws+osVvHsS91dvqUCQJJWDZyf79cPgwNGp021PU+V6ygsRN\nVSMjDxAU1Mxp9/rk/XvglxF0XtiT+IR4+z6R338PW7c67Z4irvDHX3vZ988uHn5qIvc9lNvqcCQJ\nJWDZyZQp0K4deHjc9hR1vpeswlUbSRcpAkOaNOfUMS9m7JgB992XYzbrluxt+Miu1NnrS+A4bbad\nGSkByy5iY2HmTHW+F0mHHj0MvNd/Sd+fBxEdGw1t2kDevPbtvESyoPNHjrLY/TdeeWXynf4mFwup\nEWt2sX8/jB5trwFLhYw2uVTnfMluli6FDj8F0zzwXj6r95m9n15AgL1p5QMPWB2eSJqMaVWNkHtO\ns2PiUatDyfbUCV9cJjR0DsHBXfH0tE9nhoSMd2qNjoirHD8XReWQMoS1CaNMkTIwYABERsKsWVaH\nJpJq5saNVJlVg486zKZh+desDifbUwImLmGz2fD1LUVMzGrshfy78PKqSWTkAY2ESbYwdvNY5u2f\nx6rWqzBiYqB0aZg6FerUsTo0kWSSzUaYJmsblaNDtTPsf/cUboYqjZxNWxGJS2gVpWR3nat25lzM\nOebsnWOvAxs7Frp2tddZimQiKfZ0nDuXL4ueoGedQUq+MjmNgEmaaARMcoJ1x9bRfG5zDrx1gHye\n+eDVV6FSJRg82OrQRICUfxYXyBPA1sfy8nTQBSL6nrS/d8XpNAImLqFVlJITVC9WnRrFavHBmg/t\nD3z5JXz1FRw6ZG1gItelNBvRNaEQX/ibtPPvqOQrC3DICJhhGPWBL7AndCGmaX56m/P8gQ1AM9M0\nf7zNORoBywK0ClKyu0bNT7OmTFnCu6yjVOFS9n30fvkFli8HI81/7Io4VNIRsEL8xhbPQCoNSmBn\n95343utrdYg5hmUjYIZhuAFjgXpAGSDIMIxStzlvOLA8o/cU67mqSaaII6VlE/qBPR8g1/qBdF3c\nHdM0oWdPOH3aoZt1pyUekcSSzkYMNfrx0XMlqf1EbSVfWYQjpiCfAg6ZphlpmmYcMBsITOG87sBc\nIMoB9xQRSZO0bkL/7LPwks9b7I74ix/3/2jfYWLiRIdt1p3WeESSurFl14IRQ2nOQcIanqXX072s\nDktSKcNTkIZhNAHqmabZ8frnbwBPmabZI9E5DwHfmaZZ0zCM6cBiTUGKiKukd/HIyZNQusEa8rVu\nzR899nGP5z3QoQN4edlrwlwcj0hKNj/8KsurFWBBnd1s6bAFQ1PkLpXZi/C/APol+lzvDhFxmfS2\nTylaFPo1r4Hbief4+LeP7Q86YLNutXMRR9n51RoeOr2dn+scom+1vkq+shB3B1zjJFAs0ecPX38s\nsarAbMP+zigMNDAMI840zUUpXXDo0KE3/x0QEEBAQIADwhSRnOrWTejtI06p3YS+Tx8otnwUvbaW\np23FtpS8r+TNzbptixcTcfx4mhejZCQekZsSEig3vQ8LP2rHXzFf07R0U6sjyhHCwsIICwvL8HUc\nMQWZCzgI1Ab+AjYDQaZp7r/N+ZqCFBGXu7GFloeHL3FxkWneQmvk+pGsiljF0hZLMYCo0mX45Egk\n0/OWSteWXBmNR4SZM2HsWAJ7PUC9EvXo6t/V6ohyJEu3IrrehuJL/mtDMdwwjE6AaZrm5CTnTgN+\nUgImIq6WkfYpV+OvUnFiRYbVGkb1wtWp/0gJfr6Si/Ls5TS2dNVwqZ2LpNvly1CqFPumfEyt3e/w\nZ88/8fLwsjqqHEl7QYqIOFlYRBit5rdi5jMzadzgbfpFv4gvkbRkFt7elVmxYhL+/v5Whyk5wbBh\nsH077Vrl57GCjzHohUFWR5RjKQETEXGBtgvaksfMwzctf4CYpeylOR3pz3qvd7WKUVzj9GkoU4YT\nqxdSftkrHO5xmEJehayOKsfK7KsgRUSyhZJ/juL7A/N5b1xf8GpIPy83xhtdmT7hCyVf4nQJCbDi\nuSGcf7UNo/+aT9uKbZV8ZVGOWAUpIpJjPPZgYfKFfsKPL03k6J97OX7sOMWGDqVkxFGrQ5McYOGw\nPdQ4Np+4939n+rdV2Nl5p9UhSTppClJEJA1ME2rVTuBY7Rr0qducbk91g+PHoVIl2LgRSpa0OkTJ\npqKjYWuRBpR4qz4hjc5x/PxxpgVOszqsHE9TkCIiLmAYMG6sG+e+mcjg1UM5deEUPPIIDBgA3brZ\nMzQRJ/iuzS+Uzn2Iw688wdjNY1V4n8UpARMRSaPSpaFDYBke+qsDvZf3tj/Yo4e9OHqO9nSU9Lvd\nBu3bwuN5YfHb9Iq9SMNP2hEdfpHffw63KEpxBCVgIiLpMHgwVL44iC0nt/Dz4Z//26z77bfh33+t\nDk+yoDtt0O729Tj+MW3MyfUDVyrEEb96HsHBXZMlardL4CTzUQ2YiEgG/Hz4Z7ot7caeLnvsjTC7\ndLEvVZs0yerQJAu54wbtXl5cLV6cejH3EVY1CAr+CQtmJOs9d2N3BU9Pv3TtziDpoxowEREL1C9R\nnyoPVuHDtR/aHxg+HJYsgTVrrA1MspQ7btA+ciQJL7zAJo+/4ekvYO1Aku4farPZCA7uSkzMaqKj\ntxITszrFETLJPJSAiYhk0Jf1v2TqtqnsOL0DChSAceOgQweIjbU6NMkibt2gHW4kWI/mzg1jx5Ln\n889p9GFdch25jPe1Znh51SQkZPzN3nN3TOAkU9IUpIiIA0zfPp0xm8ewucNm3N3c4fXX4fHH7VvG\niKRCihu0/7ocihQheugASowpwU+Nf8LtX7dk+4fecQpTDYKdSlsRiYhY6PBhk5Y/v8hrFevQr3o/\n+4rI8uXh11+hQgWrw5Ms4sYG7UWL+rHq85O0/LY+xsGDDN72GceijzGj8YzbPjfFBE41YE6nBExE\nxEI//QTdBv3JpRb+bAjewOP3PQ4hIfZi/I0bIVcuq0OULGTkCJM6I+pS8f3XiGrThNLjS7O141b8\n7vW74/NuJHBJR8jEeZSAiYhY7PXX4d9SX3L1sR9Z3WY1bhhQqxYEBkKvXlaHJ1nEsWPQr+wSvvZ5\nG88Du+mx4m3cDDe+qP+F1aFJCpSAiYhY7NQpKF8xnoffq07nZ9vQuWpnOHQInn0WtmyB6yvWRO6k\nSeA1Jm4sj0/Ip/xZvSxVp1Rlf7f9FLmniNWhSQrUhkJExGIPPQQfvp8LY3EI761+j+PRx+17Q77z\nDnTurG2K5K5+/BHK/j6VQqUfgEaNGBI2hO5PdVfylQ1pBExExIESEqBOHSjV6UOOxf/O4qDFGNeu\ngb+/PRF74w2rQ5RMrFf784xY+ASeK5ay+yF36sysw6Huh/DO7W11aHIbmoIUEckkTBPiEq5SdXJV\n+j3Xj5blW9qnIBs1gt27QcXRcjsDB8KJE/D117wc+jK1i9em1zOqH8zMlICJiGQyW09tpcF3DdjR\neQcP5X/IPgJ26hTMmmV1aJIZHT8OFSvCzp2svHqQDos7sL/bfnK757Y6MrkD1YCJiGQyVR6qQpeq\nXXhz0ZuYpgkffGAfCVuwwOrQ0kUbPTtZ//7QtSvxDz1In1/6MKLuCCVf2ZgSMBERJxr0wiD+vvQ3\nIdtDIG9emDYNunaFs2etDi1NQkPn4Otbirp1O+PrW4rQ0DlWh5S9rF8Pa9dC//5M3zGdArkL0OTJ\nJlZHJU6kKUgREScyTZi3di9dNgew+c3NFC9YHHr3hqgo+O47q8NLFW1z4zwHD4IbCZRs+RT07s2F\npq/wxNgnWBS0iKoPVbU6PEkFTUGKiGRCly5BzxZlaPrg/2i3sB0JZoJ9f8jNm7PMVKQ2enaO+Hho\n0wZOfvw1eHpCixYMXzecuo/VVfKVAygBExFxonz5YPx4+PWDPsRdi+er37+yT0VOn55lpiL9/Py4\nejUC2HX9kV3ExUXip8ayySStk7tT3dzYsXCv23lq/DoQvvySyOhjTNw6kWG1tIF7TqAETETEyQID\n4amquXjiwAyG/TaMA2cOQPXq0KwZ9OhhdXh35ePjQ0jIeLy8auLtXRkvr5qEhIzX9GMSSevkunfv\nddu6uUOH4MMPIbTsMIx69cDfnz6/9KHn0z152PthC1+FuIpqwEREXODMGShXDt74YgJrzk9nQ/AG\n3GOvQoUKMHIkNG5sdYh3pY2eby95nVwY0BDYRNK6uUKFfKhRAzrUPEyb8U/Dnj38fGkn3ZZ2Y2/X\nveRxz2PhK5G0Ug2YiEgmVrgwjBkDB2d1ppBXIT5Y84HLV0VmtI2Ej48P/v7+Sr5SkLxO7h7gEVKq\nm9u50z413XrXO9C3L1d8CtF9WXfGNBiTLPlS64/sSwmYiIiLNG0Kc38wmNF4BlO2TWFt5Fp4/nn4\nv/9z+lSk2kg4V/I6uUvAcVKqm6tcGZb1+RVjz27o1YuRG0ZStkhZGpZseMs19T3L3jQFKSJigaWH\nltJlSRd2dNpBQTO3vQP6sGHw+usOv5faSLhGaOgcgoO74uHhS1xcJMHBbxAS8u3Nz0NCxhMU1Azi\n4m5+vyMCKlJ1clW2dNyC371+N6+l71nWoSlIEZEspGHJhrxa6lU6LO6A6eUFM2fCW2/ZtypyMLWR\ncJ7EU4RBQc2IjDzAihWTiIw8wJgxX97yeVBQM/uTvvoKihaFwEB6/tyT3s/0viX5An3PcgIlYCIi\nFhleZziHzx22d8l/+mno0gXatbN3b3UgtZFwjpSmCJPWySWrmzt5Ej75BMaO5ccD8zl45iDvVHsn\n2WAEWvQAAB2NSURBVLX1Pcv+lICJiFjkWmwehpYNZcDKAfbWFAMHwr//wrhxDr2P2kg4ns1mIzi4\nKzExq4mO3kpMzGqCg7umWCwfFwc//XT9kz59oHNn/nnEh+7LujP1lakp7veo71n2pxowERGLhIVB\nq1bQ57vJfLN/ApuCN5H7aCRUqwbr1kGpUg69n9pIOE54eDh163YmOnrrzce8vSuzYsUk/P39bzl3\n0CDYvh1+6rWChDeD2fHtt4z+axIF7inAuJfunGzre5b5pbcGTAmYiIiFeveGk6dM4l9ryoP5H2Rs\nw7EwcSJMnQobNti3qJFMJ7VF8uvW2ddV7Pj9Cl7PFCf47AWWPf4Al188ypRKU3mlXiMlWFmcivBF\nRLKgjz+GvXsM6l0J4efDPzNnzxzo1Anuv9/eKl0ypdRMEUZH20c4J0+GfJM+YN3fZ5hrruBSrQTM\nxV/SpX1PtZnIwTQCJiJise3boV49mLF8O21WvMhv7X6j1LV77a0KfvzRPiUpmdKdpghbt7b32p3Y\nP4JrFSpQxXyEXdVehHynYd6XgB+wEbWZyNo0AiYikkVVqmSfijywuhKf1P6Ept835VKh/DBpErRs\naS/Ml0zpdrsDREfbPz4bZULPnlzp0oUDDx6DMt/C0jHAr8BDqM1EzqURMBGRTMA0wTDANE3aLWxH\nvBnPN42/wejeHaKiYM4c+wmStcydC++9R/SmNZQYW4b/b+/Ow6Oq8jSOfw8ESBQEIiCoQETENNBo\nSGNDuxAUGmhHUVpAaLDFqCBLa/vYiqPO0O3QgOLGIAqaaRaVTcWdHSKCIAFlDyDIoqikWIwdTCDL\nmT9uASFkqVRVblUl7+d56nmqKrfuPXVMvC/3nPs7P8/KIfrbKzh5ci8FBZaTJ1eiK2CRTZPwRUQq\niV9yf6Hj6x0Zcc0I7m9zl1MjbPhwuP/+UDdNyuPYMWjbFubO5W7Pa0RHRfP0b58+PWS5dOnysyrn\nn66ULxFFAUxEpBLZdWQX16Zcy4SECdxy/hXE9uoFK1Y4J3SJDPfdBzVrMn9EV/625G9sHLqR2jVr\nn7WJykxEPn8DWFRFNEZERAKzYfFX/PxGDoMz7qPW9PNZ8cc76divH6SlOTO7JbytWAELF/LDF8t4\n4K0bmN9v/jnhC5w5ZApeVZOugImIhBmPx0PTprdy4sRcuGE6XDGP6NkHONbjJqJjY526BhJWvvsO\nXngBJkwAk5MN7dqRP+EZfp85ieubXc/opNGhbqJUEN0FKSJSSTgLMccCTWHlE/BzK/J7VmPr8OHO\nlZU5qhcVTvLy4M474cILvfdJ/OMfkJDAmHpbKLAFPHXDU6FuooQhBTARkTATFxdHXt5aYB9g4L1H\nyWuYyaL8z2H2bBgxAnbuDHEr5ZQnn4Tzz4dRo4Avv4SUFFJH3ckr61/hzd5vUr1a9VA3UcKQApiI\nSJg5U2U9kZo1U6iWn8NzHV5k0qZJLKl3FMaMgd69ISsr1E2t8t55B2bNgpkzodrJHLjrLjwT/s7A\n1L8wrdc0Lq5zcaibKGFKc8BERMKUx+Nhz559PP74VSQm1uTWkSvpM68PywYtpe2o5+HECXjzTdUH\nC5H0dLjhBliwAH7zG+Cxx8jf/TV/uO04CU0SGNd1XKibKC5QGQoRkUrqyBFYuhT69YNZW2Yxatko\n1gxYzsXd74B77oGRI0PdxCopK8u5KbVLF2D1arjjDkZNuYO0n7azaOAioqqp0EBVENIAZozpAbyI\nM6SZYq0dX+TnA4DHvC//DTxgrd1Swr4UwERESjH2s7HM3T6XlUnTqXNDN5g/X+tFhtLx43DVVcx9\n4jYe+/kd0u5Lo8F5DULdKnFJyAKYMaYasAu4CfgeSAPutNbuKLRNRyDdWpvpDWujrbUdS9ifApiI\nSCmstQz5aAjf/vwtH9Z9gKgHhsP69XDRRaFuWtU0fDibTuyn65VfsHjgYhKaJIS6ReKiUJahuAb4\n2lq731qbC8wGehXewFq71lqb6X25FrgkCMcVEamSjDFMvnkyAMPsR9h7BkOfPnDyZFCP4/F4SEtL\nw+PxBHW/lcrChRxZNJ/b22zlpR4vKXyJz4IRwC4Bvi30+jtKD1j3AguCcFwRkSrrwL4o5vxxLl/+\n8CVPXp/nFKEaNsxZ1TsIZs2aQ/Pm8XTrNpTmzeOZNUu1x7Ztc0YbT/vxR04k383tIxpyR9u+DPj1\ngFA1TSKQq2UojDFdgMGcmQ8mIiLlZC3cdRe8OrEOCwcuZP7O93hmZHtnRvjEiQHv3+PxkJw8jOzs\nFWRmbiA7ewXJycOq9JWwffugWzdYt877RkEBBXcNYvC9DWh08RW641HKLRi3aBwEmhV6fan3vbMY\nY9oBU4Ee1tpjpe1w9OjRp58nJSWRlJQUhGaKiFQOxjj1WDt2hCuuaMCSQUu4/l/XU/fZexny53EQ\nHw/du/u9f6cSfxzZ2e2877SjRo3m7Nu3L6zXLayoha0zM+Hmm51Cq126eN987jmeapLOvqZNWXb7\nTKoZldWsKlJTU0lNTQ14P8GYhF8d2IkzCf8HYB3Q31qbXmibZsAyYJC1dm0Z+9MkfBERH2zYAD16\nwAcfQKMr99B5WmeebXYv/YdNhpUrnSDmB4/HQ/Pm8WRnrwDaAZuJienC/v07wjaAzZo1h+TkYdSs\nGcfJk/tISZlM//79At5vbq4Tvq68Ev73f71vpqXx+iM3Mq7XhawZkkbD88OzT8Qd4VCG4iXOlKEY\nZ4wZAlhr7VRjzGtAb2A/YIBca+01JexLAUxExEcLFsDgwbB8ORQ02ErXGV2ZGH07fV9aCp9/Dn4G\nplOBpkaN5uTm7g9aoKkIFRUYCwqcod7MTKfSR1QUcPQo7/b+FcO7neTToV/Q6sJWwfoaEqFUiFVE\npJIraYjtzTehXj3nSs3mQ5vp/kZ3nj92Df2XHnKS2XnnBfV44SYtLY1u3YaSmbnh9HsXXNCepUun\n0KFDB7/3m5sL48fDww97u7CggEV/+i2DWm1j4dBVtG/SPgitl0inACYiUomVZ4hta8ZWfj/z9zyz\npwUDv4uFd9/1Xr6JLL4GQLeGTD8bPZjeuW/wXvJSrm3ROWj7lcgWyjpgIiJSgcp7V2LbRm1ZetdS\nHrt8L9PrfAPDhwetPIVbylMG48zi5V244IL2xMR0ISVlclDD1/q5L/LHnBm8dfsbCl8SFLoCJiIS\n5vwdYttxeAfdpnflkc/hwYSh8OSTbjQ3YP5e0aqoIdM1a+fRa34/Xu80lltvUxUlOZuugImIVFJx\ncc6wI2z2vrOZ3Nz9xMXFlfiZNWugaUw8q5JX88r10Tz+1QTs1KkutDZwp8pgOOELCpfBKE3Dhg3p\n0KGD3+Hrww+dCfeFpW75kF7v92dGk2EKXxJUCmAiImHOnyG2GTPgllugYc3mrLp/LSu6xHFP6kPk\nzpzmWrv9VVLgrF27doUtjTRjBgwdChkZZ95btPNj+s7+I7Nze9HjoUlBP6ZUbRqCFBGJEOUZYsvP\nd8pTHDjg1AmrHn2cPv/qifliHbO7v0adfoNcarV/ipbBSE4eRErKzKDX+QKYNg2eeAKWLTtTOm3O\n1tmMnDuY+bsSuPbNz6B69aAcSyof3QUpIiJnyc+HESOcFYoWLIB6sbmMmNmf1V++xwdJU2jROznU\nTSzVqcBZu3ZtEhOvq5C7HF94wXksXuyEL2stz615jolLx/DRkoa0+2QD1KkTlO8jlZPmgImIyFmq\nV4fJk51Vibp3h2rU4NU/z2Po9X/ld1/cx6fzng11E0t1ak5XVlaWX3PCivJ4PGcNYb73HkyZAqtW\nOeErvyCfkQtGMuPTiXz+ZjTt3lqm8CUVRlfARESqgG++gRYtzrxe+uFL/GnVX/mvK+5jWPKrGFPu\nf8C7Jhh1voqro9anTz9+/hliY+HIL0cY8O4ACjIO8faL31N3YSq0bl2RX0sqCQ1BiohIuexe+R59\n5/bh8ssSeX34IupG1w11k0oUyNJIZQW4jT9upPec3vSu3YFx/7mcqA8/hmuKXS1P5BwKYCIiYSYS\nlvLJ2bSBR8bcwIJ25zP3ngUkXpwY6iaVyN/+LKmO2pIlr5JeM51HljzCxJYj6T9kEsyaBTfdVBHN\nl0pKc8BERMJIeSq5h8rq1VDtV4lM+p+vGL/E0vO1zvxz5RjyCvJC3bRi+VvnKy4ujuzsxsA27zub\nOVltL+O+Hsf41eNZ1vZZJ3xNn35W+Co6Z0wkmHQFTEQkyNxamzAQ1kLv3jns3p3HjBknSGh8kgN3\ndOPeG37i2GWNmXbbdNo0alOufYbjFb+CAnjmGRg7NocTJ3pRq5aHnIt2U/fP0fS5ug/PxtzKeXcO\nclY079bt9OfKs/amVG26AiYiEib8reReXkWv0JTnis3s2XNYuPASdu9+nfbt83n0ha9ptnANiza2\n5f7VJ0ia1pmnP32anLwcn9oSjlf8PB64+Wanwv3WrdFs2z2Zri9eRuz955NyRwovV/sPJ3zNmXNW\n+Crv2psi/lAAExEJMn+WDiqvooFn5MiHfA5ApwJGTs4KcnIeAn5mwoQW3PNgFL/M+pD7anViw9sX\n8uXulbSZ3IYPdn5AaSMT4RhYFi2C9u3hqqtgxQpLWta7dJ7bmQYNGpA+Ip1bPjvkVKp9/33o0uWs\nz5YWoDUsKUFjrQ2rh9MkEZHI9tZbs21MTKy94IIEGxMTa996a3bQ9p2RkWFjYmItbLLOYOIKCzGF\nXm+yMTGxNiMjo9jPr1u3ztat2967rfOoU+d626tXht21y1pbUGDtK69Y27ChXfTm3238pHjbfWZ3\nu/nHzT7v74ILEuy6deuC9p3L67nnrF2yxNr1B9fbpGlJtvXLre2n+z51vttTT1nbooW1O3cW+9lz\n+9fpz1dfnWpjYmJt3brtg/7fVCKXN7eUP+/486GKfCiAiUhlkZGRYdetW1diEPLXuYFnnYVWPgeg\nkgLGOe38/HNrL73UnvivJ+zzqyfYRs82sn3n9bXbM7af9f22b9/u2/5ctPPwTjvw3YG28YTGdsr6\nKTY3P9farCxr+/e3tkMHaw8dKvXzRQP0qfAVTt9RwoO/AUxDkCIiFcTfu/bKcu4Q53HgW3wd8vR5\nce9OnSAtjZqffc5fR73Pnl7LSWicQOdpnen0XCeadmpJt25DSUy8juTkgeVaLLyibPpxE/3e7se1\n/3ctLeu3ZNeIXdyfeD9Re/ZCx45QsyZ8+ik0alTqfvr378f+/TtYunQK+/fvoH37q12Z1ydVh+6C\nFBGJQOcuVj2QlJQ3ylWo1Ne7Fh8amU9y5vO0XfAM5oUX+CapE/H925HbPhZ+uQTW3UL0vuf4cu1q\nsrKyXLsLcv9++Oc/4ZbbT5B16XymbpjKjsM7eLjTwwxJHEKdWt5lhObNg+HD4e9/h6FDwY+q/5Fw\nZ6uEhgqxiohUMUUDVEWVgfj4Y3jqKWj1y0amHv8TeXEX0nnTUbZmbYJWH0HiVGi2iJ6X/54Hkx4k\nKS6JWlG1gnb8ovbuhX+Oy2POmlW0uu0d9teeQ7vGv2ZI4hB6XdnrzLF/+unMauRvvAEdOgR03ECq\n8UvlpQAmIhLmygpI4VhH6xRrYf58mPjsCW7ePp7BP49lPA/wIuPJI53oCzvz5OxH+Hj/x2zzbKNL\nXBd6tuzJ75r+jtYNW1O9WvUAj2/5fMc3DH9mJTuyP6V6/Ce0bNSUvm17069tP1rGtjy7se+/D3/5\nC9x6q1MI7LzzAuwBRzj/N5LQUAATEQmx0k7OZRX2jKTCn+vXw6KX53DNjEHEUY3/jqrBLf96jf4D\n7gTg8C+HWbxnMQt3L2Ttd2v5IesHEpsk0rphay6vfzktY1vSpE4T6kXXo150PaKjoskryCO/IJ/s\nvGwyjmdwKOsQB/99kHRPOts829iSsQWsoVlBZ/pfewO9WvfgsvqXndu4PXtg5EjnMtnLL8ONN7rc\nO1LVKICJiIRQaQGqrPlDkTq/yJORwdE5c2jx2mvUqFULnn6alTHdycs3JCRA/frOdkezj5J2MI1d\nR3ax++hu9hzbw6Hjh8jMyeSnnJ/IzssmqloU1U11qhNNdH4jWjVpTLP6TYhvEE+bRm1o07ANzeo2\nw5Q0f+vgQRgzximq+thj8NBDzoR7kQqmACYiEiJlBaiSFoNeunQKHTp0KPPnYa+gwJnoPmYMRw5b\nptZ+mAnfD6BOg1q0aAHNm8Mjj0CbYlY2ev55WLkSDhyAnTvh0kshIQHGjoXLirnAdY70dJg0CWbP\nhuRkePRRaNAg6F9RpCT+BrCoimiMiEhVcqpyenb2uSUKGjZsWKRshBPQCpeJKOvnYa9aNejXD/r2\n5cIlS3j8+ecZdfRvZP6uD7t+exebz+tInTrFn59atYK4OGjaFOLjoU4dH47373/Dhx9y8tVXMdu2\ncfLuuzl/2zZo3DioX0ukIukKmIhIgHwZQizrDrpKd4fdgQPOAtczZ8KxY9Cjh/Po2BGaNStfKYiC\nAti6FT77DJYsgeXLOXh5S57YtotPYq4gK/dA5PeXRCwNQYqIhJAvASqS74IMyJ49sGABLF7slITI\ny4Orr3bGJps2dYqiRkdDrVpw8iRkZjqhbe9eZ1xyxw5nm+uvhy5dONyxI83a/Tbi5sxJ5aQAJiIS\nYpU2QAWTtfD997Bli3OV7MABOHwYTpyAnBxn4nzdulCvnjM22aoVXHklFOrPiJ8zJ5WKApiIiISV\nigqkkXrXqFRO/gYwrQUpIiJBN2vWHJo3j6dbt6E0bx7PrFlzztnG4/GQlpaGx+Mp1759XcuyvPv3\ntz0i/tAVMBERCary3JQQSOHZQArfFhVJhXAlvGgIUkREwkJZc7QqegixvPvXkKYEQkOQIiISFs6u\nawan6prVrl2btLQ0vvrqK2rWjMMJO1C4blownKrL5uv+y7u9SDCoEKuIiATVqTlaycldTpflSE4e\nRGLiddSsGceJE99QUGCpqMKzvha2PTWEWbt27cguhCsRSQFMRESCrn//fnTteuPpgJOYeB3Z2Su8\nqwVspkaNa4mJ6XJW3bRgDfcVFwCL7r/onK/k5IGkpFRMe0SKozlgIiJSoUqaEzZv3jjq169fYXXT\nSpqkX9Kcrw0bVpGVlaU6blIuWgtSRETCUklDggkJCRUadBo2bFjs/ktauzMrK0uFXMU1moQvIiIV\nyte6XW4p6SYBzfkSN2kIUkREXBFOSzVVusXPJWRUB0xERKQcwikQSuRSABMRkYAokIiUnwqxioiI\n33xZu1FEgkdXwEREqjgtxSPiP10BExERv2gpHhH3KYCJiFRxKssg4j4FMBGRKi7c6nSJVAWaAyYi\nIoDughTxR0jLUBhjegAv4lxRS7HWji9mm4lAT+A4cLe1dmMJ+1IAExERkYgQskn4xphqwCSgO9AG\n6G+MiS+yTU/gcmvtFcAQ4NVAjysiIiISqYIxB+wa4Gtr7X5rbS4wG+hVZJtewAwAa+0XQF1jzEVB\nOLaIiFQSHo+HtLQ0PB5PqJsiUuGCEcAuAb4t9Po773ulbXOwmG1ERKSKUiFYqWp0F6SIiISUx+Mh\nOXkY2dkryMzcQHb2CpKTh+lKmFRqUUHYx0GgWaHXl3rfK7pN0zK2OW306NGnnyclJZGUlBRoG0VE\nJEydKgSbnX1uIVjdjSnhJjU1ldTU1ID3E/BdkMaY6sBO4CbgB2Ad0N9am15omz8Aw621NxtjOgIv\nWms7lrA/3QUpIlKFaCkkiWQhuwvSWpsPjAAWA9uA2dbadGPMEGPM/d5tPgH2GmN2A1OAYYEeV0RE\nKgcVgpWqSIVYRUQkLKgQrESikBZiDSYFMBEREYkUIRuCFBEREZHyUQATERERcZkCmIiIiIjLFMBE\nREREXKYAJiIiIuIyBTARERERlymAiYiIiLhMAUxERETEZQpgIiIiIi5TABMRERFxmQKYiIiIiMsU\nwERERERcpgAmIiIi4jIFMBERERGXKYCJiIiIuEwBTERERMRlCmAiIiIiLlMAExEREXGZApiIiIiI\nyxTARERERFymACYiIiLiMgUwEREREZcpgImIiIi4TAFMRERExGUKYCIiIiIuUwATERERcZkCmIiI\niIjLFMBEREREXKYAJiIiIuIyBTARERERlymAiYiIiLhMAUxERETEZQpgIiIiIi5TABMRERFxmQKY\niIiIiMsUwERERERcpgAmIiIi4jIFMBERERGXKYCJiIiIuEwBTERERMRlCmAiIiIiLlMAExEREXGZ\nApiIiIiIyxTARERERFymACYiIiLiMgUwEREREZcFFMCMMfWNMYuNMTuNMYuMMXWL2eZSY8xyY8w2\nY8wWY8xfAjmmiIiISKQL9ArYKGCptfZKYDnweDHb5AEPW2vbAJ2A4caY+ACPK8VITU0NdRMimvov\nMOo//6nvAqP+C4z6LzQCDWC9gOne59OB24puYK390Vq70fs8C0gHLgnwuFIM/REFRv0XGPWf/9R3\ngVH/BUb9FxqBBrBG1tpD4AQtoFFpGxtj4oCrgS8CPK6IiIhIxIoqawNjzBLgosJvARZ4spjNbSn7\nqQ28DTzovRImIiIiUiUZa0vMTGV/2Jh0IMlae8gY0xhYYa39VTHbRQEfAQustS+VsU//GyQiIiLi\nMmutKe9nyrwCVoYPgLuB8cCfgfdL2O7/gO1lhS/w70uIiIiIRJJAr4DFAnOBpsB+oK+19idjTBPg\nNWvtfxhjrgVWAltwhigt8J/W2oUBt15EREQkAgUUwERERESk/EJaCV+FXP1jjOlhjNlhjNlljHms\nhG0mGmO+NsZsNMZc7XYbw1lZ/WeMGWCM2eR9rDLG/DoU7QxHvvzuebfrYIzJNcb0drN94c7Hv90k\nY8xXxpitxpgVbrcxnPnwt3uBMeYD7//3thhj7g5BM8OSMSbFGHPIGLO5lG103ihBWf3n13nDWhuy\nB87csUe9zx8DxhWzTWPgau/z2sBOID6U7Q5xn1UDdgPNgRrAxqL9AfQEPvY+/y2wNtTtDpeHj/3X\nEajrfd5D/ed73xXabhnOjTe9Q93ucHn4+LtXF9gGXOJ93SDU7Q6Xh4/99zgw9lTfAUeAqFC3PRwe\nwHU4ZaA2l/BznTcC679ynzdCvRakCrmW3zXA19ba/dbaXGA2Tj8W1guYAWCt/QKoa4y5CAEf+s9a\nu9Zam+l9uZaq/ftWmC+/ewAjcUrOZLjZuAjgS/8NAN6x1h4EsNYedrmN4cyX/rNAHe/zOsARa22e\ni20MW9baVcCxUjbReaMUZfWfP+eNUAcwFXItv0uAbwu9/o5z/0MX3eZgMdtUVb70X2H3AgsqtEWR\no8y+M8ZcDNxmrX0Fp2agnOHL714rINYYs8IYk2aMGeRa68KfL/03CWhtjPke2AQ86FLbKgOdN4LH\np/NGoGUoyqRCrhKpjDFdgME4l57FNy/iTCc4RSGsfKKA9sCNwPnAGmPMGmvt7tA2K2J0B76y1t5o\njLkcWGKMaadzhrilPOeNCg9g1tpuJf3MO6HtInumkGuxQxbeQq5vAzOttSXVGqsqDgLNCr2+1Pte\n0W2alrFNVeVL/2GMaQdMBXpYa0u7bF+V+NJ3vwFmG2MMzhycnsaYXGvtBy61MZz50n/fAYettTlA\njjFmJXAVztynqs6X/hsMjAWw1u4xxuwF4oH1rrQwsum8EaDynjdCPQR5qpArBKmQaxWQBrQ0xjQ3\nxtQE7sTpx8I+AO4CMMZ0BH46NdQrZfefMaYZ8A4wyFq7JwRtDFdl9p21toX3cRnOP5qGKXyd5svf\n7vvAdcaY6saY83AmQ6e73M5w5Uv/7Qe6AnjnL7UCvnG1leHNUPJVaZ03ylZi//lz3qjwK2BlGA/M\nNcbcg7eQK0AxhVz/BGwxxnxFFS/kaq3NN8aMABbjBOgUa226MWaI82M71Vr7iTHmD8aY3cBxnH8V\nCr71H/AUEAtM9l7JybXWXhO6VocHH/vurI+43sgw5uPf7g5jzCJgM5APTLXWbg9hs8OGj79//wNM\nK1Qq4FFr7dEQNTmsGGPeApKAC40xB4D/Bmqi84ZPyuo//DhvqBCriIiIiMtCPQQpIiIiUuUogImI\niIi4TAFMRERExGUKYCIiIiIuUwATERERcZkCmIiIiIjLFMBEREREXKYAJiIiIuKy/wfzFZoCQMz6\nWgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1131547f0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# visualise\n",
    "fig, ax = plt.subplots(1,1,figsize=(10,7))\n",
    "vis_samples = ax.scatter(x,y)\n",
    "vis_true = ax.plot(x_vis,y_vis,'b--')\n",
    "vis_sfw = ax.plot(x_vis,y_sfw,'r')\n",
    "vis_sbw = ax.plot(x_vis,y_sbw,'g')"
   ]
  }
 ],
 "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.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}