{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Dimensionality reduction using Gordon's theorem\n", "\n", "The Johnson-Lindenstrauss lemma guarantees that with high probability a random projection $\\Phi$ will preserve the inter-point distances and angles in a data set $X\\subset \\mathbb{R}^d$ (of $n$ samples) up to a distortion parameter $\\varepsilon$ in the sense that\n", "\n", "$$\n", " (1 - \\varepsilon) \\| x - y \\|^2 < \\| \\Phi(x) - \\Phi(y) \\|^2 < (1 + \\varepsilon) \\| x - y \\|^2, \\quad \\forall x,y\\in X\n", "$$\n", "\n", "as long as target dimension $m$ satisfies\n", "\n", "$$\n", " m \\geq \\frac{4\\log n}{\\varepsilon^2 / 2 - \\varepsilon^3 / 3} .\n", "$$\n", "\n", "Note that this lower bound is a worst case scenario estimate as it doesn't depend on $X$. This has been implemented in skicit-learn module `random_projection.py` as function `johnson_lindenstrauss_min_dim`.\n", " \n", "Gordon's theorem [1] improves the lower bound on target dimension by taking into account the geometric complexity of $X$. More precisely, there is a universal constant $c > 0$ for which the lower bound\n", "\n", "$$\n", " m \\geq c \\, \\frac{g(T)^2 + 1}{\\varepsilon^2}\n", "$$\n", "\n", "guarantees that, on average, the distortion doesn't exceed $\\varepsilon$ in the sense that\n", "\n", "$$\n", " \\mathbb{E} \\max_{x \\in T} \\big| \\, \\| \\Phi(x) \\|^2 - 1 \\, \\big| < \\varepsilon .\n", "$$\n", "\n", "Here the geometric complexity of $X$ is measured by\n", "\n", "$$\n", " g(T) = \\mathbb{E} \\max_{x \\in T} \\big| \\, \\langle \\gamma , x \\rangle \\, \\big| ,\n", "$$\n", "\n", "where $T$ is the set of normalized difference vectors of $X$ and $\\gamma$ is a standard $d$-dimensional normal variable.\n", "\n", "This presentation follows [2], which contains numerous ideas for further improvements.\n", "\n", "References:\n", "\n", "[1] Y. Gordon: On Milman's inequality and random subspaces which escape through a mesh in R^n.\n", " _Geometric Aspects of Functional Analysis_, 84--106, 1988.\n", "\n", "[2] J. Bourgain, S. Dirksen, J. Nelson: Toward a unified theory of sparse dimensionality reduction in Euclidean space.\n", " _Geometric and Functional Analysis_, 25(4): 1009--1088, 2015." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "# Define a function to compute the set T of normalized difference vectors\n", "def diffvect(X):\n", " \"\"\"\n", " Input: Data matrix X as a numpy array of shape [n_samples, n_features]\n", " Returns: Numpy array T of normalized difference vectors with shape [n_diffs, n_features], where\n", " n_diffs = n_samples * (n_samples - 1) / 2\n", " \"\"\"\n", " n_samples, n_features = X.shape\n", " meps = np.finfo(float).eps\n", " T = np.array([(X[i] - X[j]) / max(np.linalg.norm(X[i]-X[j]), meps)\n", " for i in range(n_samples) for j in range(i)])\n", " return T" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "# Define a function to compute the geometric complexity of the dataset X\n", "def geom_complexity(X):\n", " \"\"\"\n", " Input: Data matrix X as a numpy array of shape [n_samples, n_features]\n", " Returns: Geometric complexity of X\n", " \"\"\"\n", " n_samples, n_features = X.shape\n", " T = diffvect(X)\n", " N = 1000 # number of samples from the normal distribution, increase for more reliability\n", " G = np.random.randn(N, n_features)\n", " # G.shape = (N, n_features)\n", " # T.shape = (n_diffs, n_features)\n", " # np.dot(G,T.T).shape = (N, n_diffs)\n", " # np.amax(np.abs(_), axis=1).shape = (N,)\n", " g_T = np.mean(np.amax(np.abs(np.dot(G,T.T)), axis=1))\n", " return g_T" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "def gordon_min_dim(g_T, eps=0.1):\n", " \"\"\"\n", " Inputs:\n", " g_T: A numerical value describing the geometric complexity of the dataset\n", " eps: Allowed amount of distortion on average (default is 0.1)\n", " Returns: A suggestion for the target dimension m\n", " \"\"\"\n", " c = 0.7 # experimentally adjusted\n", " m = c * ((g_T**2 + 1)/eps**2).astype(np.int)\n", " return m" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "from sklearn import random_projection\n", "\n", "# Define a function to measure the distortion induced on the dataset by a random projection\n", "def test(X, m):\n", " \"\"\"Measure the distortion on X from random projection with target dimension m.\n", " Note that Gordon's theorem states that the expectation of the distortion doesn't\n", " exceed eps, whereas here we consider individual random projections.\n", " \"\"\"\n", " n_samples, n_features = X.shape\n", " Phi = random_projection.gaussian_random_matrix(m, n_features)\n", " T = diffvect(X)\n", " # T.shape = (n_diffs, n_features)\n", " # Phi.shape = (m, n_features)\n", " # np.dot(T, Phi.T).shape = (n_diffs, m)\n", " # np.linalg.norm(_, axis=1).shape = (n_diffs,)\n", " norms = np.linalg.norm(np.dot(T, Phi.T), axis=1)\n", " # np.abs(norms - 1).shape = (n_diffs,)\n", " return np.max(np.abs(norms - 1))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will then compare Johnson-Lindenstrauss lemma with Gordon's theorem in dimensionality reduction of two datasets of $n$ samples, namely the Olivetti faces dataset and a random dataset of the same size." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "from sklearn.datasets import fetch_olivetti_faces\n", "\n", "n = 150 # number of samples (max 400)\n", "eps = 0.1 # allowed amount of distortion\n", "dataset = fetch_olivetti_faces()\n", "X_faces = dataset.data[:n,:]" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD8CAYAAABXXhlaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztfVmsZNd13dp1a6439euJze6WmoxompRgURIj01FgTZYgO4b1ETuwrQR0QIAwYAcy4sCSkiCwAweQf2z5I3BAWLL5YVuDJwqCIZtgxARBHEqUREkcxEFks9nsud9c83DyUfXqrL3fu7er2d31KNZewMO7VffWueeee0/V3mftvbaEEOBwOGYLub3ugMPhmD584jscMwif+A7HDMInvsMxg/CJ73DMIHziOxwzCJ/4DscM4pomvoh8RESeFZEXROST16tTDofjxkJeawCPiCQAngPwIQCnAXwDwC+FEJ6+ft1zOBw3Avlr+Oy7AbwQQngRAETk8wA+CiB14hfKtVCsLQMAgjnzIKEXoveFQvxyyhf6sb2kr47LyyA2IfoLLYTYaELHBXuyDHCLkvK+3Sc79vLn9LlzfKzw+wN1HH+qIHoMcnRtuYxz9wIPuD6Ox7EXolGYM2M6oH19cy1p42NHm9u34PHJGse08wL6OejTua7mvvBr+zvJRw7o1SDoNvg6++aau714L0Kfrrlnno8unVff9vEgd+or6LbqV3yor2XiHwXwCr0+DeDHsz5QrC3jbT/9GwCA1rK++PZi3O5X9Oh2jsQrPnjT+nj72PyaOm6x2BpvV5Ku2telb5aFfDO+ryaAfph5EgF6svCkGpgHJU93xU5M1Sdz7lKuR+eO7Zdz+lpKEo87XFhX++aTeG3FjHNf7s2NtxNzncvJ1nh7pR+PK4vuR31QGm9vDspqX0Ljw18KiZlwl6gfFjxBbB8ZfC/seBdoTNd71dTj+B7umJh0n3oDvY/vU3sQp9NmT4/HeqdM2xW175VLS/FcW8XxdvF8QR1XOR/7WNg0X9at4esnv/oZTIJr8fF3+1bZ8TUqIveLyOMi8nivVb+G0zkcjuuFa/nFPw3gOL0+BuCMPSiE8ACABwCgevB46FaH3xdd80U/KNEvRFl/fyTl+O1czsdv8HLSU8flM34VGO1B/CZ9rb/WbBnUkrbaZ01iBreZD8aioF+TatIZb/MvPKB/xewvYVH1Px5nf8X25+OvesdYHvzLfjC/Md7mX3hAWyKtoH+dBikmfMv6eBnI+pVX56J7aK0o9pK4PXtv29Qve14+divoMbAm/bgfg2TX94Gdz8d8LVqqK2tkKWR4N3Z4Mzym3ftwdYcrfAPAbSJyi4gUAfwigC9fQ3sOh2NKeM2/+CGEnoj8OoC/B5AA+FwI4anr1jOHw3HDcC2mPkIIfwfg765TXxwOx5RwTRP/WmBdIFoQRX9O+1jlYvRV54rRny4l1veNvpilsvSKeWx/YHxC9uesL8a+dppvdzWwq/V9iZ4Xn4t9+h1tmJV29usLiNfSzbjVR/Or6nWX/XAaAnuuBvn8+4kJAICLvQVqL45xmu9/Lcgaqy49WHycXQvIYgbY/+e1F4szrbg6z8/i6J2MfRGl5cjK4Py82pc1dNskVsbykoKH7DocMwif+A7HDGKqpr4MgEI9jLa1qTwoko2St6Z+NDGLZCZlUSZZ1F4pw3Rm8z6Llkv7DKDdALtvB91EqOaiGZkodySXelzffHezeX+IzO9WaKnjLvSjGdk1FJsKuKF+lKFN/WpO05gMdqe6/dh+25yr0Y8BK3OGFrVU67h/NsAGcUzbfd2+Giv6XBZVa/vILoLdp9qkZ7PXT/9Ntc9EjZ7v1ZXIc+dqhqpdj23ax6iwNTy3DCZ8Zic6yuFwvKHgE9/hmEH4xHc4ZhDTpfNyQK8y9Nt2RG6SO1djSgPA4fnN8XY5H/2hHZld5BLlcsbv5sypfqRWLCWYRdOl+f+WlmNYf5QTOawLm+b/V43vW8o4X4NCStmvt9lzh5LN1H0douku9iMtN5/T94XDg1sDHbKb5kNbmpXH215/WshuVkhtlg8+KWyIdNa6DKNC6wl2/YmfF/uMMUWdy9NxekjV8zIo6jbqN+VHn5mMZvZffIdjBuET3+GYQUzX1B8ASXtoyiRNY5LEoCfMlbVpy1l4rV60f9jsB4BWP/1yeJ/N6mOwKZplUiuz35ivyvQ0Q5xFEabRV0xJWVgzlCP3bMacOo6ouI1+Ve1byLXs4cNzQZ+LzXtLKzYGkaab1FRe6+p+1PK704XWVM4l6WPK5+b7Yt0Fdgmy+mvvdWvCa8sCuwVVytSrX9CZgHxqm4/fWRy50BP+lPsvvsMxg/CJ73DMIKaepLNt6Vorur8c37iptqn2cRRYi1bkYSz2TlYkXy62YWW50tA2y6qVlASNDSOzNEcmql1VzjI3GbzP9oNhzWGO3OtnyEn1jVTWJOC2LRKjC8juyTqiCb/V1+YrX5tiPAAMetcucsFuAUcGTup+ANn3jN2zVXKZriaJiwVYEl79rxitRTrODBVKq6OI2PRbpOC/+A7HDMInvsMxg/CJ73DMIKbq4+f6AaW1oRPS2qe/c/KldOeE/XrOvrK+XpbAAfv1q0bemLFYiHRKVlQfZ/g1DW3GNNcO6knSBR/4NfvkVqCykBExtybRz8yS1+b1ECuwwePNlKBtj4U4OqaPq73aeDuNprwS+F5zG5ZS47WBrMxLphh3RM+R/28jJbf6cT2kKnqdR4mATOhg78jY7O++3iAVey3xOu0SRbE+XA9wH9/hcKTCJ77DMYOYeuRerjM0SXJdbe6UK9GEqveK+mOcyEFmUZLTdEctH9vIm31cyeTl9X3j7cNzWituvpAuLtGcMBqNq6jkbcUWptXM6DOtw5RgwfCWXBHGtmEj6LZh6TZGN6cbWTORfON+WC26DJpR9SkjnIyjI61LoFw5elyySnnVDV3IbgHTws2+fsbYXWCXAAAWqTrRJVMQgjX4+JmzbiJrmExK9YlJNONLS0xwZdIamfoT1sL0X3yHYwbhE9/hmEH4xHc4ZhBTF+IYFIffNduim9sYkH9kK5Km0R0l48er0N6e9j/P1aO4ZLMdfbh2xdBQ7ejfXjS+aZGoog75hNZnYxrK+npVWocYhAW1j6v9MrKqt9qw5YZER3A+SRfiYH+90dN+cVp9u7IRqOA2bQZhWhs2LDeLAmsyrUhj2jM0Lu+z/j+Ls+ayBEzz6RmQFxJ6dvo2jDs67/P51q7bALCWixTyat+Khezu8w96uo+89GBD3ovrwzekf518fBH5nIhcEJEn6b1lEXlYRJ4f/d+X1YbD4Xh9YRJT/08BfMS890kAj4QQbgPwyOi1w+H4IcEVTf0Qwv8WkRPm7Y8CeN9o+0EAjwL4xBXbEkG/vPt3jZCJ1urpbgUyhfg4mynVIBrQmsdsjrNOf6NryzvH42w/NpuRpuuRGdbtpA9jkhjKsRppOu4HAKxTRGFuLl5nc6BNbNZ22zBZdrmU8cmZsWIqbl++rvb1qeZBklGrWQmOGGqP3YJ1et+a2CtEfZ7a0obj+Q0ysRuaYmOEPkX15XV/c+QO5guxv6WCHlN+rvbXGmpfFv3G+9hVsyIx7L5a92+1Ee97IU+ZgEXjylKTlRWzr9Hd7lBqX9XxEx21E4dDCGcBYPT/0Gtsx+Fw7AFu+OKeiNwP4H4AKFaWrnC0w+GYBl7rxD8vIkdCCGdF5AiAC2kHhhAeAPAAAMwtHw/bC7J9Iw/MplB+kG6IFMl0W2vpZBs21zo9vfLLZYoKlBC01tSmMn9uYPrR2IjHhnbcJ2b1NbCYQkmbZBtUWqk/Z9iAbrwdlXw0cxcKeoW4naTftqVCNFOzIuZUuS4T1cfm/Qol29gVeJvMwuDoQitUwrjYipFwFzZ1VFx9PX5OtuI1J4306xqUtKnbI9OfdT2aVW1u5wrprFKVnh1r9nOprGY3PTmrVozjXStoBiFP7iC3Nz+n5cy3luJ4FOr6XiTrQ3dN+ukRmozXaup/GcC9o+17ATz0GttxOBx7gEnovL8A8I8AbheR0yJyH4BPA/iQiDwP4EOj1w6H44cEk6zq/1LKrg9e5744HI4pYepim9sYpLMzKCbaf+HIpg5F8bW6uvssVJjLpfs6+0rRD7a+GPtpkjf9WIqf6/XSs/M67cmy1mwbSWH3MktWRJSz2GpGAJSj+ric1JzRyldCnMbws9TfNmzWWlaGIu/jbMXNro4SZH+6b0pLc3ZaqFIJ6oKh7JokQlk2ApVlitwjP94KeQq9rBqa1T6PjG6OBFM6JBjTTZ9aliZm+pTXqfKGCu5X4r5+yRjr21l5k7F5HqvvcMwifOI7HDOIqZr6IQd0q8PvmkGiTS2mzmzSQhrFFsxxSxRxZYU4Dpaj4MbRytp4+0fnz6vjOAnDlnRa70aTlaMEl4qaduF9LAAC6IQjm1DCVVMXKAosb0xvpveWCzrqjqP1OCJv3pj6rKVnTXim85jCs5r4/LMxb5UhqMucsGKTXLic2YF5fS31UnRj7L1msLtQNhF55Xx8XaLtdi/90T9Q0eIs+8k1tLUV6pTg9Gpjcbx9sV5TxzXJDbAuXp5cygK5FTb6L5CrUj9s6OoXR89ZzqvlOhyOFPjEdzhmED7xHY4ZxJ7ReT0dnal8G+vPJUp0MPo5h0wWVYV8oruWTqt9t5QujrdrVCL6cl935FI3hsrurJ0X2z/diHkHHHYKAOcoq6y+rsOKOZMsKWuaaHlJ+5bbKJtaf+zjW996jl4fLayOtwtGRONib4G259W+RHbXpt8hOEKiGonodQIWuVwgH9/q3l/qxLGzZc7r5Bc3WiR0arMh6fHYMpQga08WyvHcSwv62eGQ2o4tTEfYb9ZUbipujLc5686u35zbpOeqpZ8r9vGz1jLyS7GP3Zpef+ocGq4phJcnqwnov/gOxwzCJ77DMYOYrqkfgFxvaAK1l/UuzrqzVFylEE1dpquWS9pcY1rtLWVN0y0n0Yx+snl8vP2t9ePquFe3IiVz4bLWxJPzkbqZezl+Z3YW1WGoXIhmXm1Lm3yt/URbmujFteVI/V26OZ57335dNpxLXFWM1t1cJZrVXF6rBW1ecskrK47RSnksdkTq0aVZ6pNNetbOt20w3XlqRQtxNC/FNpPN+LmCZk/R2R+fieop3X75Eo9/HPDimqZZG7U4BmsL2tx++g4S8Dikn7nbD8XE1DfVVsbbluI9L9Gl4RoSFhy51zVZgoeXo1uxWtXjLdsCHK6r73A40uAT3+GYQezZqn53Sa9os1BG1QgVWJN+G0fK6+r1TaX42q5Uf/Xy28bbX3/hRNyxoU3gZH9c8c+/rM3BuZfjdr4ZTarNW7V5VdhgkQ7d5/U76Y2CcWlORlM093J0K1Ya2nyt3JpeEZZZCTbhrbw2m9xWsINdhLLVcSZw+zYijysSFzOqGJ8h16p1Vke7JS1KvqlSgs1xbUa/98RL4+3/hTvVvpCL19k4HsftLX9u5MBvj/3d95zet35HbEOe0M/V90vRhH/57dFVuWXfZXXcErlgm20dAWlLwW2jYN5n5quzqJ856Y6O9SQdh8ORBp/4DscMwie+wzGDmK6PL8AgP/Q1c3Pad2T/xZYzYmokK9OL8ciFH1Wvn3vm2Hi7eJl88Dt0tNy/vO2J8fYLbz6o9n33q7HNyvnoTB3+R33u+k1xu18xpasW4xrCW266qPY90419LJ+Jt0b6uo3Lm9EXXlvStM7ZVvSZD5Tite0r6HUSFttc7ek2uPQzU3ELpizUq60Yvdgz6wQXGtoXTsOl1XhcMAIb6MTXxYOx/4cW9T07VIp0Z36/9v+TH0QffO5kHNPuvD7Xxk/F9uvf12sN73/Xd8fb//fQLWpf8R9j/1vfihz1c2/X6zIH52OfrU9fLaSLeTIKtFbSM2tkjSPD9ajB05P9lvsvvsMxg/CJ73DMIKZq6g8SoL00NGUW5rVJxgkOtmpsjfTbl7nck+n997eOjLfPbOiou0Ba+u1D0cwzMuz41eVotx87pJNvPv+vvz3e/t0/iRqkvao2z1ivon5An2CRtNJPr+uQP3YDWmQNWlO/dTlSTyf36RBILgnWqsUBKtQMpUZjt9LRpi0FuOFoMSb62GhIdhdsBOQWVSRm8ZFuV5vA/U3SODTXyd5De4Pozbw2lV9pRBrtw//kWbXv79pvHW8XTsc2No/rh+fOm0+Nt3tH9O8hJyOJSb7ZenPsS2Er9r99Vj87fUoomytqV3aekq5YxKVrtBb59dxh7e70ysPnPcNTUPBffIdjBuET3+GYQfjEdzhmENMV25wboP2eIfVSyKiPZ8UlWdSABSRvyq+p476zEekwK2gwt59qypFYQ/sV7Yv93BP3jbd//61fUvuebcU1hPotkYKRrhF/oHp5UtS+dZN833bThAuT75osxvb7nfSxOrOq1wlY1KFRiuc639ZrHhu92H7bCE/UaZGiSDHHNrT3zsqr421LCZ7ejFTf5lZck+g1TC2EeappWNShyAUSzqxv7V6iHAC+e/7m8faxJf1M7FuOvvAK3af1su7H2qnYxkGTDfnChQPj7b6tp0D3qVOi0N627iOLihyq6vbnC7vXIGwZurrZi69vWtBtnDs6HO+sehWMSUpoHReRr4nIMyLylIh8fPT+sog8LCLPj/7vu1JbDofj9YFJTP0egN8MIdwB4B4AvyYidwL4JIBHQgi3AXhk9NrhcPwQYJLaeWcBnB1tb4rIMwCOAvgogPeNDnsQwKMAPpHVVpIMsDSis2wZazbnCxnZXFzS+VxvSe3bIN37iimDtFiOlMkiRQK+WNyv29iMJuuvfvNjqf3I1aIZmphSW0tEVfaMBhyb+oWSNm3zZOpXSVN+o25KeTdiG52WvoX5OSoZxeW6bXlnMtvtvjxp6ZcljuOLnUPqONYunDMls+dLpGsIcqeM2zJoxecgHNLjwRltXSpJ1Temfo5chNNr+pnotOPn8lQau2dprw5l8XW0iV2rxGuxpdNbBdLLL8Y2rEvTydDxZ33CcoZuX7vP4in6ArZ+ZPi82DLh6ee8CojICQDvAPAYgMOjL4XtL4dD6Z90OByvJ0w88UVkDsBfAfiNEMLGlY6nz90vIo+LyOO99d3z6h0Ox3Qx0cQXkQKGk/7PQgh/PXr7vIgcGe0/AuDCbp8NITwQQrg7hHB3frG62yEOh2PKuKKPLyIC4LMAngkh/D7t+jKAewF8evT/oas6sSkBPEeUhhWQZLzUjhlzK10dasp+z81z2ijh0ticzfXe5efVcac7kZzgUFBAryGwv9Wzfh/5c03jL5Zp7cGGfzIF2aC1AEshcbnnnCnhXKT6cKx8Y/3FLLCyDivwJIZmZW3+gVH4qVGG5VyN1GeMbzpYI3pzRdcgOEvXXSqnKwG1SKc+l9PXyZ/jdZ/ish43FnQ9Mbei9rFw6LMb2qPlsupc8zEs2TLck41/ljY/r33ZcN6kMtqXm+w8k/D47wHwbwB8T0S2c1b/I4YT/osich+AUwB+YaIzOhyOPcckq/r/B0Ba6P8Hr293HA7HNDDVyD1BQHFkmlo6IpehEsim9EovmvebPS1ayKWmqqbEMAt4HKOMszvLutTWOytRuPFbRS26cLIVqT+OhFtp67WLPvV3rqxprjbRUlVDOa43yZUg89XShUWir2wbTKOxoImlSNlU3EEbUSTfhZ6O+GNcojpo9v6x67ZUJXrTuEWco9nf0o/jYBCfkRyJV+SMkEU+H4+zJjWXzebr5JLkAHCkGl3DJSNawp87XNERc1t5Uzo8BTlyk2yknhp/Lvll7hnXm9jq6PNuRz2KZ+c5HI40+MR3OGYQ0zX1JeqN2ZqevCpsq9SWMrTdGUVa4bauBK9Ub/ajSX2yo3X1+tj9OCA9Em7H6iv1o2BW3Xnlvtk1UXdkynGZpVpJsxxc2ZVXgQFt3vO+eaOXVyc3ySbpcFTf883DsT1T6ZY/Z/dxZWGuYmzLo1Vpxb9b1P3gqrgcuVc17hOPqWWL0jTrG12dzbLZjeNhV8yzIknztI+vzbo0WbUF0p4l6z4VaYw7fd3H7eSsSdkD/8V3OGYQPvEdjhmET3yHYwYxVR+/18/h0taQjjswV7/C0RHs97Av2RS9FmC13Rm8hsCiEet9HS1W76XTM6zjz+0VjX/Lo5olKsJRXxYcSWZ9zAUSI2Vf2r7mLLueKU9tI+0Y3OcsnXcW6bQUWFp/m1V9zfUOZSva9RBa22hRBGTX+LdM2Vmf3tafG7dhfPDLrXgtds2G1yjmUkQzAD1uOyIIab3FriHwM83jbddesqIvrfDMleC/+A7HDMInvsMxg5iu5h603h3D0hNp0NSeNtOtWZ2GtW409SuJpsrYvLImNu/jc1nRD07gafX1EHObpbJ2Edhl4MhDa24rqixJpzqZwrTmZVqfAG16cvsHCjpqjZFVziyN2gOsqW/KhpO7U8ynU7UMa9pzyXW+L9at4DGw2oJ8PvucsnAGH2dpVn5eLOWYlXzD6GVEW05K440/f1VHOxyONwR84jscMwif+A7HDGKqPn5OwliIwlIy+rh0Cow13xcKOgy13U+/HPad2Bez9F2Wz8zn43NZ/5apsoHxF7OoOEaN1h4s9cahnHasCpIeGqqOozazMiMZXUMJLuaZwtMZihx2zedaMHURV1t6nYbB953DlK1/y+Gx/Yx6Dfw5FlwFtIZ9Ppfun9v7mUY1W2TtU/eQmq/3dFgx09VWrHY7XNiz8xwORyp84jscM4ipZ+eVCrubPDrzTX8f2UynbdiMszxF8tk2mAbc6EXzciGvTT7eZ2ku1txr9NLpK5VhZd0WMqst/VhQGnlEF2aYiVVDR26lRB7afiAjypHBFGZ3kP64WD0+pkmb/Wiy1kx/95Xj+F9qaA1FpsSystv42mzZKUXhURudjKi4ndGWg12PA9Iz66xpz8+fdZnYPbPPLSNtHrwW+C++wzGD8InvcMwgpmrqJzLAYmlonq+Z1Vwr0MBIi1SrGhnurBVt3pfGEgDa9OcIP0CbivtLsR929ZXBMtMW1s1gzOW5PJV2adgcZA1CYKeIyfgzhhlY66SvpnOfby7F6rPWnGfzNRnofTb6bRvW5WCdRCuJ3smIYtNtkEtgIvLYtWrQfbL94M9lRYBmCY7wGO9wrQhZzyl/zkb/sdtiWbGs+bPrea7qaIfD8YaAT3yHYwbhE9/hmEFM1ccv5AY4UB4KcKw0tf+cRVVwhBj7t1kinNYf7aeUha4l6cIKWeISTMlYX533WX+O99nIPT6WM+H2J1vquGdaN4+3n904rPalUU+WDjvfnI/9yKeLeczT+oL1W/la7L1oDIq0j8RT+no9ZJ3WGm5f1OUXOdqSo9ZYGBMwvrWJQkz1wTOOs751FpRPnrKOBGQLkzJsGW4G+/VWeGN7/oQJk/Su+IsvImUR+bqIfEdEnhKR3xm9f4uIPCYiz4vIF0QkfYXL4XC8rjCJqd8G8IEQwtsB3AXgIyJyD4DfA/AHIYTbAKwCuO/GddPhcFxPTFI7LwDYtjULo78A4AMAfnn0/oMAfhvAH2W11e4nOLU5rEBrtdFYZMCajYtJNKXZvLQRUNa8T4OOrEuPxGJKzaKVQpsBQC6kt8/7LB3JLs3BfKS2GgNt2n5vPZr65xtzal9C52MhC5uUwq8tBVbLuG4GU4d2PBr93Q1A6y7weFsRipvLkUpc7aaXWM8SAUmLyLPm9qTXbCPr2EXIcld1dJ420/k5Xuum06xcjZerIgNAc1RdeVLtvYkW90QkGVXKvQDgYQA/ALAWQtg++2kARyc6o8Ph2HNMNPFDCP0Qwl0AjgF4N4A7djtst8+KyP0i8riIPN5bTw9YcTgc08NV0XkhhDUAjwK4B8CSiGy7CscAnEn5zAMhhLtDCHfnF9PNGIfDMT1c0ccXkYMAuiGENRGpAPgpDBf2vgbg5wF8HsC9AB66Ulu9QQ6XNochpstz6VSZRYLom2X68fQ1tt7TPmGZ1g3miMLbsU6QY3omg9bJ8aY2drYoDLhs1iu4zaqhEsvCGVzx1rzQ0pQdU2BWkKFC2Y9MS9mQTu6zpRWPl1ficTTelprkfVk6/UzhWR//aG2djtO+OvdxgTIxbT9WUNv1M4D2p1lIJW/ayAqjzaIt02Bpy6zx4T5zf61Qa5b/3usmo2Mm6t5EPP4RAA+KSILh4/7FEMJXRORpAJ8Xkd8F8G0An53slA6HY68xyar+dwG8Y5f3X8TQ33c4HD9kmGrkXinp48T+oRlptdY4E8ua30xnKZrLiE4oemzCK6tKOo1jqRumrJiustQQi2NYyo5N3cSYpa0Q29/sRdEPmyXI5uBSRWfuzVOJJy6ZbWlF7rONVFOCG9R/W26Mj8vS7eNrtqW8skpQb9AYHCrGSMYsd89eJ2dfsnmfFT23A9SkdZmYwrPPLaNC49gcaDeAXRx2uy62NFWrqE+TnZeMy2SndkHBY/UdjhmET3yHYwYxVVOfkdiItgzNvRatcG8MormZtcI6n9MmcJ++49j8tsIVbJY2BunpB1mVS9V5MzTU+sYua9H5GmSi2tVuTqqpGqGPmytxlTxNfMS+ZrcC0BGES0nctkxJlsmdVgGWK+wCOvnGCmBobcR4P+2YsvltWRQ29bNW7rm/9rp4rKzrthXI3aQhts+mjiA0Wo4DPf5p6A/iONoV/nF13gmfS//FdzhmED7xHY4ZhE98h2MGMfUSWts+adPo0vPrrGwr9smtH1XLRSrLZost5GKeQJ3owflk8vwBXnvgqLtGP10YwoL7leX/MzVk/fMjlUhpHihpkY6sUtaM5aQ+3j7dWVb7jhcuj7dPFGIU30pP00u89mLHm6k4HjdL37W6sY1OBtXH0ZBzJuKRfXIrTApapuHxzqLeLDW5nI9jlSVGotrI8LUtpcnga7ZZfOzX7yyLPSGPt92/qzra4XC8IeAT3+GYQUzV1F/IN/HhA08DAL505l1qH0dAWVOfI8TKE+qhWVpnmXTr2FzjxBgAWOlrc5aRVgbJnotNynbIp+7LqoKrdPuMeXystDreflvlFbXvXG9pvP1qZ99425rit5fOjrdvKqypfS+2Y1Im3NE1AAAcqUlEQVTQJtGnfdNfdhfOtRfVvtON2I9FqpC7XKir45jCs5Fqq+1IH3K5tH15neDFY2XdLl1ZOI6jpeUYls7LOnYrRdzDugCq1JZxJfj+plV1BrSp37PJWaVhHyelmf0X3+GYQfjEdzhmED7xHY4ZxHR19aWPm/JDf9L6UQNJFyDg0Fmmq1igw2Ip0b4kn49pv67xwau0z9JtTEslJKjRyhjGvhFgzKKRkCJ6cXJL021vrb063rZrEp976Z+Nt9e2on9+19FX1XGLyfHx9nvnnlH7lvNxPYTXBmx484vNg+Pt5zYOqX2nVuL6QnM9Unsfe9dj6jj2SV+pL6l97O/yeEwqqmrBY29FUDhT0tLEfN38fADAAbr1m/14nXZtJ018FEivH2Az8HgdLEn0GByaG96zFyesoee/+A7HDMInvsMxg5iqqf/qxf34z//jVwAA7X2adjh49/nxdq6k97GJxmaXjbqrZZjpnA1YNzr1DDbna7kMXX0SzaiKpnvYNbFmaVZ5LaYFmXI8tzmvjvvD8+8fb9966LLaxxp8R5djpp6lkF5qHoj9T25V+zgjj7MaL3W1W3GqHl2Qi3WTddejrLuLcawePXebOu5AhSII17Spf2QxRihyjQNLr7GJneUGZGcTpmd6KvEUW5qNKM4Ci3vYU5HVzuW0hp+Lz0GjQxGs3XTNvcFAU30Xtob3Jkvbn+G/+A7HDMInvsMxg5iqqV/c6OPow8OkD2loM3rrn6Z/jk1iNqN59RkAchky3JwAwtF6A/Pdl5PJ5JNriP1ncRDbX2uWZkl2p0WIvevwafX6/51583jbVh3eX4tm+uFKZECshPbZ5kLsr4kM/JHqufE2uyZnWzo671w9uiAFs8o8V43js/mmOMZnL+o2Lpeji1Au6j6emIsJQnMZVXv5Xtt9fC/4WqzeIbuQNkIxK3Ivk6UhZJVcaxOLtdmNbosV21CJOd/W41h7eniduUuTTWn/xXc4ZhA+8R2OGYRPfIdjBjFdsc0QIP2hn9J/4SW16/S5d463D92qxSQ4W2/eCi1MiC7xKQskxGn9uSzwegD7bJmiDsbnZNrIildeIqELXgvYNPUDDsxFCmx/WUcosq4++/U1E6lWIe1/Syuy//tqJ64FXGppyo7LcNcK2g/eaO/uqxby+lyL5XgvbCnvWyqX4j6ibm3E5iAX78vOMtz026Z09fV6ArfZN7+HWbQu+/+qhLuZWkyz7SxFtnsJrbxZN9nYjM/LkWd0G/PfH66HJK3Jslcn/sUflcr+toh8ZfT6FhF5TESeF5EviEh6TKLD4Xhd4WpM/Y8D4KDu3wPwByGE2wCsArjvenbM4XDcOExk6ovIMQD/AsB/A/DvRUQAfADAL48OeRDAbwP4oys21huaKMnBg+rtw38fDYbyr2tzJUtnb1Is5SLNxabs1Zj6LYn9qGMyvXbrBmRF7q2S6c/n6vTTXYmtbnoUoqqIa2lFijI7UNC0KINNzwPGreCEEhsxVs/F+8nmfa2o+1GjugCLBe3GpSdW6fHgCM4dlXTJfdLJNrof6lymrBrTtUVDxypqOEPTj7X/rcAGj2Ojm077hYvxXtde0fcC3VG/JiyXO+kv/mcA/BZiIOJ+AGshhO1ROA3g6IRtORyOPcYVJ76I/CyACyGEb/Lbuxy661eNiNwvIo+LyOOdXmO3QxwOx5Qxian/HgA/JyI/A6AMYAFDC2BJRPKjX/1jAM7s9uEQwgMAHgCAxcqRyewQh8NxQ3HFiR9C+BSATwGAiLwPwH8IIXxMRL4E4OcBfB7AvQAeuuLZRBAKo1Me3q927ftmpG6sEEczx+WpyafNa9qP/bQC0ikqpnysycMlubP8/wL5ekXJ0Ek3PidTeFk0YL2X7rszWMQBABq96FsztWdLM7P/vyKapuOnIkv0c41KnVsf34pIjJs2QhFz1MfXWo8wm25L38fokHCGFWDltYEkQ+iiFdKnE2debpl72+rt7td3jaDm3EmiLTdMPYgxpTnZmtW1BPB8AsOFvhcw9Pk/ew1tORyOKeKqAnhCCI8CeHS0/SKAd1//LjkcjhuNPYjcG5qLoWTMm7MxE+upc5oguOf4yfE2Z+dlmcpFK3KRkRXHGJCZx2Ibme0bc/JiJ2atbfV1CWTWXmsaHTYbdTY+V6KvhaPkbGnpfG53atHWKmBTf6mgF13ZlM6iKtm8t+Y3m+1M4ZVMXQSmtizNxfeaxVMspZYFNrGTjHoKF3sxQhE5bUbbLFDGWp9cNxLYyKqZwOW/AaBO7hlHOW7V9bNz8w/idUtH9z+UR21MyE57rL7DMYPwie9wzCCma+oDwGBk8vS0iSpzcWW5/7wpYxWVoJXJutnXJhOXyepLus3Toe+7sjFl50k/r290zQrkWnBpKesScASX1bqbdOWak2iWi6Z0VTGa5m3DgLCJmVWOqUP7rPndTmKb7E5Zt4LLmbG5anGkGrXzGuY4fm3N49VuNKM5aclWBM6qOsxgU98m4vDqv72faWIe9rUSfzG3ma9tR+Reik5eOKNN/cqZeN0hMZ8ZXB1T7r/4DscMwie+wzGD8InvcMwgpu/jj8+sv3NCKfp6i8/rQ5vvJZ+LNm1knRXOZLRUZFb0F63o4qRg6q1uSzMzbWSyCdlftPvYj01Syi9bbIj2Aznij7O+Ljb1ukm7F8fj9KbWs08T2OAoO4t5k1lXonWCJRLYsGscDYpas4KgTHeum/UcBuvZ26w79sHnSYAlq6aBFfrg9YDOjtJYcbyV6CfSa0Osd/U94yjHNmnpz79kakOsxDWscQTs+ISTZ5kO++dwOGYOPvEdjhnEdE19StKRZrpW+fwret9zl6Nox6GbI6VhqZVORiQfuwVM51mTjE1/6wawbh8nhlw052Jz00aIsTnIpaoAHQXGJmWW+Igdg54qN0aRZIZC4vJMzbam2DYQTVHWut9X1RFtbLZbvbwl83r8fkG/z+b9Ql7vSxPOsJVo2TLvZugHsslux62fUTqNowZtxWBV3o36ZRNxOEmqZxKreBw3NqNLc/wlc997dG2GzpPtfddZiMPhcLyB4BPf4ZhB+MR3OGYQ08/O2xYFzBt/fED01QUtJLj2CtFNN6c3zxlzO0php3zF2ZBdXgvIEuJg331/Ymr4kR9vswQZ8yYLjMU9srLReJ+lBDehqaJtHK5u7vr+lVAjH3/e0HlMFxZNVqCl5uL7VmwztmmvpRLi62qSTiUOMkJ2+R5e6sasSRv2y2Nsw3mzskDT6jratYC1TvTdbRYmi6nmXqX1lbMb6jieM+N5NEIoj56JjFB1hv/iOxwzCJ/4DscMYvqRe9vmykCbOyEhc+eyNnHmn1seb6+8I2bxLeR1tNiruX3j7beWXlX7FNWXYQ21uAySoY0sNTfuX6JNdjbFmaIDdEah/RzTeywM8WJT1yBYo6y1y22tl7fe2d3Ut9hoTabpt0VUX6uX/rhs5PV5qxTxd7y2Nt7O5TXdxDSdNefZhOcSWtaM5vuZVUKLTfbVnh63w4X18bZ1E7mNLF29bgqVOuxXvO6WEUXpkLbewovx/WTFmPrXEf6L73DMIHziOxwziL1L0smZyCOKOAoVbWotPxtN7KcvHR5vHzm2ro7jiKuNgTY9uUIuR3PZiLYihYENjGlfD1Tmi/axKAcAnOlGl+O5xk1q35MrR5CGhVLs49nNuAK9vqHN0kE3jl1om3FskUuTkFk9bzTaqA2YMUCekk0KcTsx1VuRISoyoJXqk4Xoqi3UtHt2bD66AbfULqt9+wuR3eEV+SSjArH9JWPTn10Jq4XIkZJ2FZ/bGBjRDI7Wy5JEZ8l4riQMaG29wxfo2nqGEaJEtmCScuQGldByOBxvIPjEdzhmED7xHY4ZxN75+MZ/YWEBKelsseoTp8bbp168Zbzdvll3nzOsrCgH03llpZ1v2iBuqGZ8/GqO/frYx+fbh9Vxj63GPj51Tvv07QuRisu1tJ/2SiX6d9IXOk5fCzFlkJ5ugxPcOAlssGKzuah9q/NBbmKvFl/skP2ntQEboMheciAKby2vBUFWSrGU2hMHjqt9S/uij1+lCMJSXnf4TXOr4+0TFb1OsJynNog6tNF+7NfbfezjtwyVyK9ZCNYKh3JG3nrTCHFsxGcp36RBNutg2n/X7Vuf/0qYaOKLyEkAmwD6AHohhLtFZBnAFwCcAHASwL8KIaymteFwOF4/uBpT//0hhLtCCHePXn8SwCMhhNsAPDJ67XA4fghwLab+RwG8b7T9IIY19T6R+QkRhGRokkhX0w6KzjOJBmz6v/kr0eQ79Y5lddxtlQvj7f1IL3vE5rzVXuNvQqsPxzQga69baojNupuWdPTVRjlSSpumRNKgHW9H6FF5qqpJJOJudYzJSvQem9i5tinpRNZysPs6vJ3uEiTMzBlLk/Nt+uW40zKH7NL01g3NVY302HqJKNgF3ZG1Q/FzzWVtiv/YfIzgnE80lchg894mZ7UG6dr/rInP0Xq2svD5RnRxGnVN+wlRq0mLBs6a7xkJONIfuQgTSkhO+osfAPyDiHxTRO4fvXc4hHAWAEb/D03YlsPh2GNM+ov/nhDCGRE5BOBhEfn+pCcYfVHcDwDlwsIVjnY4HNPARL/4IYQzo/8XAPwNhuWxz4vIEQAY/b+Q8tkHQgh3hxDuLuZrux3icDimjCv+4otIDUAuhLA52v4wgP8K4MsA7gXw6dH/hyY64zZbYWt/8TlN2CFr7pe/HkX3v/etO9Vxb3v/mfG2zXwr5HcXubDCikr33vj/rRRBhn15LRxyrBrDUN9U0xzYPipJ3TBlsn+wdWC8zVr3rW5GRpjZ1+vu3sdBX/uH3TYdZ9YJhGoGhgqNgaklmNsiYQhD5+WILhyk61goinBQ1Pd9wPRmme5LTYt5LFWi727r0LEev65VkC5g2jB1Epjqs/eMMyU5/Hu1rbMy1+uxH8HcC9Br++ynYYcQxzYdPiGrN4mpfxjA38hwYSEP4M9DCF8VkW8A+KKI3AfgFIBfmOyUDodjr3HFiR9CeBHA23d5/zKAD96ITjkcjhuL6Wvu9VM06AY2LGx3yIFI4d3+xzpe6NWfiObx4YKm0ZSGPZlDtiRylq5+2r5bi3p5g0s1bZosQY4es27GXbUYofjMQhQXfKFuhDja0WzsGju60Y3XwwIPNquMwSWcAKDVjOZsLpd+X/oVosC6uv0+0ZFSjPfcCnHwq1JJm98L1TiO86U4VrZc1z7S8E/T+gN0ZOcOHb0MbUTW0rPl0ur93cuDn9nSC9ntVrwvdgzUCDNXu6P0dezjjhJaVwmP1Xc4ZhA+8R2OGYRPfIdjBjH97LyR32J9FEVjWOUROlaF857RvvXjZ2N216GS1k1neq9KvrUNrWxRHe4+tC9ZIDpomdYMrAjnoSSe+1RPhxVvZpR7XkoiLfju2g/G228q6YyzVztR4cdmi11oR6Wa88243eyZmm/9dI6tQdQn19Vjv92iUNb0Uo7UemqVON6Vgj6OS3LvL2tadLkYx7iWoavPZcSzwnLnaB/77YCm8GxJ7s1eXKdhfXwAaNG4dmi9ZW0z/T7bdZNBwhmQFLre0tcs8xQHYzP3JlwjG3/8qo52OBxvCPjEdzhmENMvk10ZmVhWV585NjFuAEUphRIJdpjmG1vRXHtqXQtgsC47a6hnlbHuGkEGNvU5qq8l2mzmfUfzhnKkbavbz6+5nNaCKbUFslITTQbhSDFGDa5Woml4uaMFMDiTbNWYr71KvJ56j8s76/Fgnf2yEcfgzMZS0tv1fUCX07YltPcXY4YlU3FMiQL6ntl9m5Q5mSW2wWga3fsLLcqs62kXoU0imquNOI59E0GZFOIzkSR6DFhLhU19MWXmWGwjNXJvwtA9/8V3OGYQPvEdjhnEnmnuhYI2Y0RFLJkVSl7xp9VMWdLRUUJfY9YsfXormv7JPGnFGyE51t+v2xJaFDnVyjAbeV3Z6vYt5eJKtY0a5Eg+Nm1tFVn+XN+cmyMKD+SjqZwzUYi8ql0yChtKp55NdnMcuwu2PkGT2t/oRnPbsijLpJ2fN6vdJXJ3uqTiZ1fu+Vo2jSgKr9ArU9+YxFukib/S0VmkXPLKsiNcimx9g/QUjTmfy8XX+bxhkpKUxJx8xvTMYL4mgf/iOxwzCJ/4DscMwie+wzGD2IPIvZEflxiqIiGqwnwfhTxlRA3Sa4v1t+LlWN93hYQRVojmsoIMa0k8rmMyuDjDj7PuLvc1Vca++v6cjkYrZmSB1SRSURxBCJioNeqWjdzjfu0oJ01gysrSV+y7L5BQv10P2ZcnURETCcdrA3m6ZituakudMwq0ppAEqudn+sH32kbd8drDFkXn2TLWHPG41TUZeF1arzDlxVkwlWsaVhYzriuxRQhSMvLS6uMBCCVzb8fzwmvnORyOFPjEdzhmENM39bfpOJtk0IvmG5v9AJRSQY5ovx3mDjVpaSM2+V5pxCQXa3py5JeNmLOm7jYaA23+sTlvNdr5dRfaleAy3Eo0whxXFUoySjTFxlTfxRDNV0tf8XXbaDoeO3YX8pZmpfG2Y8NtsjlvxzsL3cHuj+d6T+vZXe5G1826N9aNSQOLm1g3kc37hjH1u1tEiy7E+1It6QjCZif2o2QSlVTwJZv3lrIjiJ0/GZr7u8F/8R2OGYRPfIdjBuET3+GYQUzVxw9CQhrGX2QRzoAMIXaG9XMI1bz2sbieHYeQnmzsV8exD7qeGNEF8p+Z2qsa0cwcJhNFKMCEblL7NnOPsTlIF3l4LbDhtoys0s9NIb/VhPNym0xf2TY4lNiGJqdn1pk2VFixboPFMFmY9FJLh+Vylh1nHQJatKTTMveFqLhqOT5ziQk/5nDkasEIgtJjnGvQc2vXVMq0vtAxbRRH9+I6185zOBxvIPjEdzhmEFOm8wTIj75repYaou+gfDrVNyiSEIelMKjJ+YI2v9c70Wycy0czqWGyrZ6vx6K/J6pa666QEoHGohmApvcszVVQGX5G05+OZXfB0oX8ubW+prasEEUadujKE9gVYhOe3wc07Vc3NiaLbzDFVkl0/3SUoI52S6PisvreNDr33MYGPQNbHT2mnM3Z6OjzNjg6r63PzRRekbLugnFHysX4zFmqWZXQ6tK+DFcWeVsXYNTGhKzeRL/4IrIkIn8pIt8XkWdE5CdEZFlEHhaR50f/9125JYfD8XrApKb+HwL4agjhRzEsp/UMgE8CeCSEcBuAR0avHQ7HDwEmqZa7AOAnAfwKAIQQOgA6IvJRAO8bHfYggEcBfGLiM5uvnMAVWk3kXq7NpYOoQqvRHQOVJrLCE4vF3ZMmDhR1Es2zG9HU75kIrm41npvNzYN5LeXNboCNArvYj+IhVu8vDXVj6mcl4rBMNCfOcFVXQK9284o2oGWiW7TP6s0t0Jja8WYTnsta2VV9ls3OirLjlfueKRvGLkjLXAv3ma+z2U1fuecSYgDQb8Rjk6qRB6eyX0ViL2wiDjMRlnESNvXb9EwUzXiQa7ujhNbYLbh+mnu3ArgI4E9E5Nsi8sejctmHQwhnAWD0/1BWIw6H4/WDSSZ+HsA7AfxRCOEdAOq4CrNeRO4XkcdF5PFur37lDzgcjhuOSSb+aQCnQwiPjV7/JYZfBOdF5AgAjP5f2O3DIYQHQgh3hxDuLuRrux3icDimjCv6+CGEcyLyiojcHkJ4FsAHATw9+rsXwKdH/x+64tkk+u9WbGNQjj5WrqV930E1+lz581ETP1Q1vZQrpmecLZFoBEdRWVHOY7WoS3+6vqT21YiK4s9Zsc0D5PPPmwy/i72YMbecaAvIls3eRt+M1WovfoGySCSg1x56YfdoRUD7wtZnVn49iVBYv5Uz2ux4F4my4vbKJpuQ/fqaofp4PYDXEOxaAFOy9l6wOOZ6k9YCDGXX7ZCIS8cIwRYo6q6m14qYpuNovUo+ff3GjpVaHuF1K0vZcYn5nJm6V0nnTcrj/zsAfyYiRQAvAvi3GFoLXxSR+wCcAvALE7blcDj2GBNN/BDCEwDu3mXXB69vdxwOxzQw3ci9EMaRSSHRJllWxBLvG5y/GHe85U3qOC5TlIXlfDSxs7TiLM42IxXH5uZKV69d3DlHlWKN+a6oPmPCc0ReFtXHevkNG6lG9BVHiGUJk7SM6cwUZIlKY7VN8gqPVXdgEnjIxGbXalDQpnKHE47KW2of04BZVB+b9xyhCQCXtkiko0X1CPpG15EoNZhnoFSN/Vgo6/vJ48PRinMmctTq+Cnwubt033M2MnX3qroAELYjQoNr7jkcjhT4xHc4ZhA+8R2OGcSe1c7b4UqT6ICtq5c/E0tN9zuR8ukvaH+uXN6dDgO0n8liDfvymlK71I10W9XQS5vkp12k0snWf2723zze7s7rIV4kWtGWuGZBDxbs3JHFR9dyW+W82vfd/vHxdhodBmj/2Ya5MjjMdYcePIFDXgG9hsDU1r6ipjdVKK6h4vjcHAJsj1P3ZVPXOGhskm+dITiSI8oub8QwF2uxzyVTDnyRrodDdu1451OEWgGANUxCPx5n/fgdPj9Btj/nQhwOhyMNPvEdjhmEhAmX/6/LyUQuAngZwAEAl6Z24t3xeugD4P2w8H5oXG0/3hxCOHilg6Y68ccnFXk8hLBbQNBM9cH74f3Yq364qe9wzCB84jscM4i9mvgP7NF5Ga+HPgDeDwvvh8YN6cee+PgOh2Nv4aa+wzGDmOrEF5GPiMizIvKCiExNlVdEPiciF0TkSXpv6vLgInJcRL42kih/SkQ+vhd9EZGyiHxdRL4z6sfvjN6/RUQeG/XjCyP9hRsOEUlGeo5f2at+iMhJEfmeiDwhIo+P3tuLZ2QqUvZTm/gikgD47wB+GsCdAH5JRO6c0un/FMBHzHt7IQ/eA/CbIYQ7ANwD4NdGYzDtvrQBfCCE8HYAdwH4iIjcA+D3APzBqB+rAO67wf3YxscxlGzfxl714/0hhLuIPtuLZ2Q6UvYhhKn8AfgJAH9Prz8F4FNTPP8JAE/S62cBHBltHwHw7LT6Qn14CMCH9rIvAKoAvgXgxzEMFMnvdr9u4PmPjR7mDwD4CobiUXvRj5MADpj3pnpfACwAeAmjtbcb2Y9pmvpHAbxCr0+P3tsr7Kk8uIicAPAOAI/tRV9G5vUTGIqkPgzgBwDWQgjbKSPTuj+fAfBbiAXQ9u9RPwKAfxCRb4rI/aP3pn1fpiZlP82Jv1tq0UxSCiIyB+CvAPxGCGFjL/oQQuiHEO7C8Bf33QDu2O2wG9kHEflZABdCCN/kt6fdjxHeE0J4J4au6K+JyE9O4ZwW1yRlfzWY5sQ/DeA4vT4G4MwUz28xkTz49YaIFDCc9H8WQvjrvewLAIQQ1jCsgnQPgCUR2c5/ncb9eQ+AnxORkwA+j6G5/5k96AdCCGdG/y8A+BsMvwynfV+uScr+ajDNif8NALeNVmyLAH4RwJeneH6LL2MoCw5MKg9+jZBhed/PAngmhPD7e9UXETkoIkuj7QqAn8JwEelrAH5+Wv0IIXwqhHAshHACw+fhf4YQPjbtfohITUTmt7cBfBjAk5jyfQkhnAPwiojcPnprW8r++vfjRi+amEWKnwHwHIb+5H+a4nn/AsBZAF0Mv1Xvw9CXfATA86P/y1Poxz/H0Gz9LoAnRn8/M+2+APgxAN8e9eNJAP9l9P6tAL4O4AUAXwJQmuI9eh+Ar+xFP0bn+87o76ntZ3OPnpG7ADw+ujd/C2DfjeiHR+45HDMIj9xzOGYQPvEdjhmET3yHYwbhE9/hmEH4xHc4ZhA+8R2OGYRPfIdjBuET3+GYQfx/Lm2T0jEfT0kAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from matplotlib import pyplot as plt\n", "%matplotlib inline\n", "\n", "# a sample data point from the Olivetti faces dataset\n", "plt.imshow(X_faces[0].reshape(64,64))" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "150 samples of 4096-dimensional face image data.\n", "\n", "Johnson-Lindenstrauss lemma\n", "Target dimension: 4294\n", "Target dimension too large.\n", "\n", "Gordon's theorem\n", "Geometric complexity of face image data: 3.465\n", "Target dimension: 910\n", "Distortion for sample projection: 0.094 < 0.1\n" ] } ], "source": [ "print(\"{} samples of 4096-dimensional face image data.\".format(n))\n", "\n", "g_T_faces = geom_complexity(X_faces)\n", "m_JL_n = random_projection.johnson_lindenstrauss_min_dim(n, eps=eps)\n", "m_G_faces = gordon_min_dim(g_T_faces, eps=eps)\n", "\n", "print()\n", "print(\"Johnson-Lindenstrauss lemma\")\n", "print(\"Target dimension: {:.0f}\".format(m_JL_n))\n", "if m_JL_n < 4096:\n", " d_JL_faces = test(X_faces, int(m_JL_n))\n", " comparison = \" < \" if d_JL_faces < eps else \" > \"\n", " print(\"Distortion for sample projection: {:.3f}\".format(d_JL_faces) + comparison + str(eps))\n", "else:\n", " print(\"Target dimension too large.\")\n", "\n", "print()\n", "print(\"Gordon's theorem\")\n", "print(\"Geometric complexity of face image data: {:.3f}\".format(g_T_faces))\n", "print(\"Target dimension: {:.0f}\".format(m_G_faces))\n", "if m_G_faces < 4096:\n", " d_G_faces = test(X_faces, int(m_G_faces))\n", " comparison = \" < \" if d_G_faces < eps else \" > \"\n", " print(\"Distortion for sample projection: {:.3f}\".format(d_G_faces) + comparison + str(eps))\n", "else:\n", " print(\"Target dimension too large.\")" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD8CAYAAABXXhlaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztfWeYXMWV9qkO092Tc84Ko5xzACEkokkmgwkGAzY22Ma7Nrb3s9feXa/t9YJZGy/gNZi1yUEIg0iKSKA0ymlGMyNNzqEn9HTu+n700G+dywiNDWrh7XqfR4+qp+rWPbfCvefUSUJKSRoaGrEF09kmQENDI/rQG19DIwahN76GRgxCb3wNjRiE3vgaGjEIvfE1NGIQeuNraMQgPtXGF0JcJISoFkLUCiEe/KyI0tDQOLMQf6sBjxDCTETHiWg1ETUT0W4iulFKefSzI09DQ+NMwPIprl1ARLVSyhNEREKI54noCiI65cY3JyRIa1o6ERGJIK+T8SGUDe8iszlEoyHe4me/B10OXOPhbYXSRSBereA3S00YjpSd7nhWZxoWkXIoTqHPy+9lS8fN3UP20UgnIiJp5fe29qN/fzLqzFb+/NZG/PbkWnmnyrhalDEIOHgzU0AlxFCXjMpgEEyh9BkYRDMuFD7BqmQc6tTnMtIRN4hn8edwQkLKvc0u9BFM4O2ECb+l10BjHPo3Ke1CIU4vBZR7uXlVMFG5X5BfZx1CWSpVQcO0S9vodBARyWEzysq4kYHEBDsW2nA/v0FSanjdDra6yOP0GK78OD7Nxi8goibldzMRLfykC6xp6VR4/7eJiCjOyWnzzXRFyqEQn7yUZBeNhhlZbez3ll1TIuXUY7wPdYK65ygLwMY31RUL90bKaw/MYnXJB7HbXUW4LuU4f5ZxtxyPlA9smcjqJOaYArk+Vpf7FjZx20rs4NTcQdYu/1tYmce+k8vqzC48d2oV/t43lS82ezfaGV/CiSs7IuWe/oRIOdjKX4ShJFxob+YvIE8Rnq1wHR66e5qZtSvcgmdpvZ+Ph3vIFimnbke5byFvZ7XjRSUbOI2yCG8/uwPXedxxvF2Xcq+jfD57F+EDYxrgWyZnu9KHsuScFXz9+Sfig2K38w9W8EAK6MpX6gwfhsUVdZHyvncms7pzvrCPiIheu/VNGgs+jYw/2lvlY3KDEOJuIUSlEKIy6Bp9A2toaEQXn+aL30xERcrvQiJqNTaSUj5BRE8QEdnzi6RlhGXzZPB3REoi3vzJds471zdmRcpCYXurLNmsXcWT+DLW3JLM635VHyl7rk6LlIPb01i7bXvmR8pxZYZnWdEXKZs8+MJNWNTM2tX1ZUbKSSd5Hz3z8XVK/5B/dVb/4P1I+ZXnzgW96fxrWnVfKmiK45/rgAXj2q3wXylH+VT7klD2TOZykfcAxtXWjfd74vkdrF17Q0akbJnlZHU5Nnxdu6ejv6Cdz3vtzaDLfojLARmNaNszC/NujuNcWvJGfOV75vLxyNwAlrh/AsrBjABrZ8vHR2nVuYdY3etrlkTKRu4oZR3aNt07Hfc9yBt2mvFsYpBzJZ5JGKv4GqyJ4RJO45VZ4EZb9k1gdSfODa85b3BsW/rTfPF3E9EEIUSZECKOiG4gotc/RX8aGhpRwt/8xZdSBoQQ3yCid4jITERPSimPfGaUaWhonDF8GlafpJTriGjdZ0SLhoZGlPA36/H/FiQnFsgFs75GREQt5yawOlUlNljBTz1VFYpQ1DC2nGHWTOyH4Bp0GFQ+Sh9qXToX58g5CWV/AT89ttggc5mOg/6Q4fRVlWPzt/K623+2NlL+tw1XsLqySdBSeB7PQ39f7mHtUv4Z8mLNN7j8n5PVHym738qJlAMr+lm7VSXVkXKmqvIgolf/+7xI2Z+AcfNk8WepWFgfKR9tzGN1Ngfm0L4F8xLfaVBNfqU9Um5sS2d198zZGik/tg00mYe4ZiCpHjQOLuVrIiUJvweO4kxCGoTcYAboLSrg4/3+9DWR8p2Ny1jdxv3QJCVX4TsaWM7HW+zCyX3SCn5WMrgZ8xS0KRXTuTaHjmAcvWX8XCauPnx+0fjYQ+RpaTqtOk+b7GpoxCD0xtfQiEFEldVPsWTJxSlXERHRyfu4AYIvDSxgYj1/H+Ve3hgpd7nAYg/vy2DtEpvwLAGHwcLKhbrsDVC/+YozWbsTV4PXsrh4HyGFqw5mczFAhb0OfcgZnF3LSILaqHcrN7558JYXI+VnW6GLO360kLW7dunOSPn12umsLlSdGCn7lTE1WtaVz2iJlDtfL2J1acfB9jon4KGHivhaKVsLFWzt7QbDls347crHfA5PN5jFKWQZLTRT3oPaq3shxCxrClf3OhSDGNN7XD0755aDkfKWD6ZFyrZuvsYC8Xg2Rxcfq/5ZmOspP+Fs+uK/1ETKz76wEnRwTRzJBWD9421clB1yY73Y41Dn9nIxLlSnzK1BHfnROLb/y3+Rt75Zs/oaGhofh974GhoxCL3xNTRiEJ9Kj//XImmSj857MSyvV202OK8oHksLbuI6th0vzYyUVa8nbwmXlew9kIncBk+v/ly0HSyFzJy9h5tWxrfgXZixmlsgz83AWcPm1vGR8oDBBFM9a+gq5F5UUpHxPfn83v+8/upIOXerIqZdws8T3ngRJqT+DC4X/+r6P0XKP3zy1kjZvLCPtRt+LD9Sds/lY+XOw7JQvRzTjrFmdO3v342U/2PfBayuZzrmYvVKmJoe7Mln7VrqccZS/gyfzxN3KzcfRH9+N5d97bsg+157zwZW94cDGKukBsyt9YJu1q7/EM6LhMEZVChnD0d/yM9lrH2YG3cR5G5LMp+z2TlQW558kq998xcGImWX4jxUkMFVgr59ON/qm8jHwF0wcu8xHtnpL76GRgxCb3wNjRhEVNV59oIiWXLPA0REZOHGYpTUDHaqbQXntRyZsL4KVsN6KX1OJ2t3ffGeSPl3B89ldS8vfjxS/uJL3wZN4wZYO/9RePU5OrlWxJWPsVKt/6Sds+xxXWCVV1+4l9W9uQdiS0k5p799O9hg1ZJxuJSzwELxTotr5B5+E8+FO2Dte+WRsqquIiKKb8ez9U/hqqG4NLDYaX+BGNM5nzUjsxt9LDyHywEfHIb3WMUTUOHFP8TVYSlxuNfuNVw16Z2JeY87BDqM42FrB9ubNJez8L71ECWGSpUYCuVc9OnrgbiQdMjG6gangG0393Pp2NGBb2fFZYjDcHgjZ+cd6lQblG3+RBoVwwavyQXj6yPlg21cZPI2hztp/dWvyduoLfc0NDRGgd74GhoxiKiy+rbSQpn7o/uIiMjSy08lgzngbf9y7qOs7qH21ZHylm2wvgoZgjokV8N5Y/6XDrC699fPiJSXnQ+twcaD3IIwPgPspbs5idVZB8BB+bLA3psHudOItRRyjGUX70Nl4Y0sXu5OVJ7/8LZI+Q8bz2Pt5s6DtVi/lwevOF4DFjAxB3RIybk/Tx1EmpxdXLRqXanE+1NCeQVTuEhjTwUrat3Bn3PG1Qi9uH8tHFkCXAFCCa1KGDQLp1G1ujv2a8y71cXpzfhHRbxZywNUxLejbd/lmNukd7mTWM98PFtiLWfnh2dBVAkN8XVrSoDYkZUBK83hDTxITPz54PU7a7nFaUIT1s/QZIgV2Zv5vXoxjJQ3u53VDb4edpKqeeEhGu7UrL6GhsYo0BtfQyMGoTe+hkYMIroyfmGR/Ci8dnwHF0O8aaBD9dQjIko9ivdT4EIEdZQfck8sXxL68BdzD66S5yBHxb29O1I2V4xn7RqvhGxWeGEDq+t7qhh0KKL1UDFrRnIcZMmEbVyWFEGFxkQ+BvkXwzKwcSs6ZSGXiciajGez7+YHBcP5GLugQ1H79fBziJxKyLTn/fQDVrfmaUUVqpBotIY0lcAK0dfLLRRLx0Nt1z2EMTBvTGXtPIpzZMEWPmfd92McvQeVAKOGz1Xedqgjm1bz50wuw3pJ/COCYQyU8nYOJUDIXf/0Gqt77JdXRcp9q7l3YdKHOLTIvArR5gefKmDtcu/COUTrEA8E+4OJb0XK/+9xWFt65vCgIo49uJcvhc/FlBW1RES0+Ssvk7OqU8v4GhoaH4fe+BoaMYjosvpFRbLwm2FWP5DGrcVunIfgEi9XzWZ11kNgFeMUQ7vBUi4SqCzgzHl1rO7I+2Dp1QAbRoeMlBNggbtm8/eiGAfWtjSzN1JuH+SqrMEGsJQWN+e6khWyct5qZHUn7ihB3W6w9623cocPfz8sy5KOc9WT+VzQ5d6PGHa+VP6gIh1staXx1Gm+LENK2jBu0Eam6YoTSWUKq3NXoP/0rbAuzN7O49kd+wewvcLF2W+TGjxEibVoTFlW/hLY+eO3cToSFcec1DqsuaIfHGft0uMwt7s6S1jd8lxM2u5uXhd/L+hquA5xB4eL+fq2OPFsgRw+n1N+CmvDxoew1j01/FnmL0VqpIOvczX0cFF43bb94hHyNmh1noaGxijQG19DIwahN76GRgwiqjJ+5uRMednTlxER0abdU1mdoxUyUPI53IPL64ccG2eBDD49g2fLrft/kHsab+cylsmkBJ6shRyV2MSasdxovbO5XBzfCBpVdUr+PE5Hy17IehZDxmLrLHiFDbbyswE1pXN+IWT1zl6u/il4BqacrcsNMr5yv6ANNBpz+Jkuh6ztPMJNSFU14LkLYHq7bcs01s5UCrk4+wVuOty+CN+UuD7QZPSsy/4A9A9dwT0lc1NgAnuiBgEwUg/zZ1ZNn32pfD3blHsHFQtYaQhBE5yEZwkGTv09LHqeX9gxH52WPaqYUp83jrVrvwwHEyG/wcQ7HjJ/4hasTf8FPBDHUB/UebPHc1VzwzPhM6zqVx7+bEx2hRBPCiE6hRCHlb+lCyHeE0LUjPyf9kl9aGhofL4wFlb/j0R0keFvDxLRBinlBCLaMPJbQ0Pj7wRjYvWFEKVE9IaUctrI72oiWiGlbBNC5BHRZillxen6sRcWycJvjFjutXFu5KqvbI6U/7hnCauzJYJNMitpsnzJnHbLsJLuKZ+z+o4msGieHPDzMom3m1mOmPuHd5WzOhK4n70cbKjtbc6KD5aivHTlYVa3vwMWXaFNPGWUqxAsdsZBhT3O4WPlng2LrlAf17GVVUDsaNgPT70456nf8TYek4IGl8A6TXagfzUNGRFR3nbQ22lQfRZsBUs//C2o2zq7+ViZLOhDdnC1oppfwXw+RJOBGs5glr8MeptWc0tJNQ1V5rNglduu4Sq1JCVN+6qialbnUvJaVfdzrzvrg3ie6jvR/8wpnBU/uh1ryZ/FxZ2UA1B3xndgPHwGy85+xfHQaPlq7w6vzaNvPkyu7jOnzsuRUrYREY38n32a9hoaGp8jnPFTfSHE3UKISiFEZdDlOv0FGhoaZxx/a3jtDiFEnsLqd56qoZTyCSJ6gogoMb1IZh4IsyTumzh/+Vo9AmVkZvPT3e5OsFM2lbM1vLbsc3ES7j/O2UHvZLBy+a+BtWpfyGPWHQjAOaZkJj+tbzqKk+Xs34Etrb+Sn/4XvgeRYE8PjyOnagPSOvh1AzPBAvqSlTRci/nprmyAuGMKcK4u+DAyr4rlqCt9nj8LWTH1xx7gjjOp23FC71qKl3XAy5dL6zk40Y4r4EEU232g0d0BC7TUSi6a2PoxBlO/wcWiPSexJtRAIhUPc1XMiYchMuU8ydnoDhPoyPgWLPDaKstYO7tC18ur57K6+FqskRCPjUG3P/VepDzYgrnOsPGPnEMRbXN38IXbrxgDWu6ARiv0vzmsncWN66Zfe5TVHXwlHKUjNMYd/bd+8V8nottGyrcR0dpPaKuhofE5w1jUec8R0XYiqhBCNAsh7iSinxPRaiFEDRGtHvmtoaHxd4LTMgZSyhtPUXX+Z0yLhoZGlBDVFFr+ZEltK8OqNJshBXCCAyq7nhNcPrd3w9Ip6wDUb4MFBvIVhWLWNH7sEP8ryLH95WB0HAa1yLgF3GNORWcPruubBJnw8sW7WLsP9iAA/fxrD7K6mn9HxMTW5Qbvv2E8T8E1MLU7UsXTWNv7lDRfR3kATF8SxsqqeCHW35jH2nkmIFBm7tsGS7LbcB7g+QDXhZL5mURIsQz0NXE1mkUJhDqxGHJrWxq3Vuzogtmd59kZrC6oHL8kKumjq3+ZxekljIGth8eilwLrrN+Hc5kZ87n3ZtsUnCOtSONrZ/8eyO7OuVwN+Ic3VkXK4xdBhbdx/xTWzpqH8Ug2aNtchagbOIlnM3MnVUqpgBff9h2TWF3pxWE1dMPr/IzjVNC2+hoaMQi98TU0YhBRZfUtQ4KyPgzf0u7kTh3WIfB1gemc9Yy/EKxijx8qjuESbnUnlAAY3i7eR+elYFMTmpX4+Gnc+q/XA+urpsZMVnfFlWDp17+0IFI+8i2ushv3c1h+HezmqY76loGuUCZnGxOV1E3N1VA3pa/maaH6/BCFOufyd3fWbIyVpxZso7WPt1Mt5gauHmR13oOK2rIa7dqX8rHKrESfA9wnhVLndUXKTZugIo3vMKTyysRcpNZwNtWbhrFqb1AciQzBU+LzQX/ddTx4RUILykPPYS5aeIYrEmVQv+17i8+nypivnFLF6ra2o603qDiTpXGRw2fF+m79IhfPypHdjaw/xfxVHeUi3tAwRJXUai4uNPrDFqE+l0HfeAroL76GRgxCb3wNjRiE3vgaGjGIqMr4wZQQOS8Oy1L+YW4qa+qHbDLxaSerO7EYqiKpeOSl7+NyfP85kKvcCbxOeJT8ZOWQsQrXc5nTOVWRo/ZxeemtDqjpfOWQRz3V/FmOdinnELVc5izagHOJ+i/y964atz5tKuT63qP8rIGUQBl+Qz47y6Nou+hBBJS0mHi7Q39GUA0vt9il8etgItzzL3jOQiuXwX278JzJdQZPshD8thLa8Vw9C/i5TFEJnrPTwlWO7mLcr7QUKrbGI7wdbcWZR/Iwn081V1/3cpypFK3l62PhBcin+ErXAlaXsQfztPEIV6NlKY58LQF4XsYN8PHwz8DaTE7isfkbLwL90y2g8bvnvcHa/WrdZbjv9VztHHo5bPdr5sdGp4T+4mtoxCD0xtfQiEFEV53XZ6Lsl8NqvDgnZ/n6y8EaVX2Vp4VaWgSWdb8V7JSphgd1iN8LFWH6ha2srrEabGn6frzveiZzlizNDgvClhmcxrzNuM45AWJA5w3cM83fAroS2/i7tWcKfttbWBV5CsHadikeibl7OftqdeN3fymfwubzIAaYhtFHu5OPlf0CeDKG9vOAIF0L0La3CSKC1cnZ4/JqPHftt3ld0Ivfwg81pbmf01uYCLHOt5j3UfqP+N3wU6yJ1GN8zvqmgsa4XJ52ytcJXt/aAZGs+SoenN/5zOJIOd3Jx3vgIqj67Ef52kxqRj8Dl4LPtsfz/kMfQgSzLuXq08y5UOEdfQd6xsu/tJ+1c3Ri7aTb+HNWzwqvneAbY4uhqb/4GhoxCL3xNTRiEFFl9a3ZXsr9Ztg5Yk8tT0VktuGkc2V5Laureggn0Pf8GJlFHz10KWs343wcsU5Oamd1fz4Ma7SETrCGZkOo46AS8CHlKB+etpXKqbYJYoDJ0IdlAO9T1yxDdtVdEEe8mfyk3TygZPTthSjRPYuzbwWb8XtwNrcQI0U70rVBydg6l7OX4i2cJE+7maeT2nMI8eHiFAcpXwE/1W+8QHW44aznxdOORMpvBRRLOEPcviNdmBfazJ2zTl6D8tRs0Njq5I4+5mGMt+kQryvZBfbbFMCcWX7AHXECT0MUbLyAa2nMVWDv06q42WD9ZRhvcxD9B0P8OS1KXA7PBu5k9MA9z0XKv3gdzrA/3nIVayemGeZawRfn7iEiomfjh0/ZRoX+4mtoxCD0xtfQiEHoja+hEYOIbprs8gJZ8LN7iYjIfIJ756UdAx3BG3tZ3cA+eGapqY/sXVyOyt6nWO59j1v/tZ6AOsXejk78KVxmu2IFPPBe3cODLi6fjjOEXe/xFGAqStZBoDP3cxl/6NeQA5vquUWeKR51ITdoTDzOLQhzd6PP+kt5LPpQHsYgsRJj7FnCVY4WJRWZZSu3LjQp1l++82HF56/mKkGzV8ljYJD/P/LCJCIqvANnNu77+DPXfE8JKtrJn8WsqObU1OCpufy8wtkGuhwt/FzGXYAxTd+L8wq7QWXXPUNJB85FfEo6gbK8mAeJdSqBRCzdmCd7BQ+QmvgiaHR+kQfiLE5HnzXHcC4z6b/5Gq67AWrX+BkGOrrDdLT/5DfkrW8+Y3H1NTQ0/o6hN76GRgwiqqx+ybQk+f1X5hAR0S8qL2R1jsNgS219nKa+ZbCCMpkVkWCAs8BJeWAB3VXc88Q6CO7Hrzj6pE7tYe2GdoEVjePh/SluFRxKhpWYgY51nAW2eBXLuvH83epNg2jhaOd1aTVgv123gFUcGORikaUeLLGt1xBXX+GWzYr2J62GWyE2XYvf2Zn8QZ2VUDcVLUVKsdqTPM57fJ3CExuYy8BMiBb+DtAv7Vy0Mjk4XSpkEOOTuRn38qbxm8X1Y7w9GQZLzOOKaHUx/p56kIsEAxNOPS9mxQjPZAhpp2bqdU2AjJSxg6/NJV+tjJT/cnAmq8vNA9ve0YF1a7Hzm8lGOKtZhvhzJjaFx6Bq7ZlNoaWhofF3DL3xNTRiEHrja2jEIKJqstvdlkpP/fRyIiIKrTDEg1dyyq24Zi+re3O/Em9dUXlVVHD3tuOHEJzQYQiEYFukyPI7oR4MrOPqJblCSX+9jnti9StqRVWmFYZjEk863qcWgwWlYw7kOcs+7hU3lAd1U38HTE8fX/lH1u6+hq9EyoMTuYw8cQK8Egd8UIH1WLl8nvaBotq6mvfhK4FQ29QDmTNvPV8uzmtwNuD38bpADw4brNlQPwa6+HlF1hbQ2D2HVVHKBIxVyIqxGp7PB9WWoqjHtnBzWLoXQT9NxzAGpiCfNGlVUqD38LrE65BnoLmbnx3ZD8L7b9X0Y5Hy1laeI8CtJAlI28Xlf7cyNyVKwNGAg4+VPwF0+a/nZ1N9yeHxCbxLY8JYUmgVCSE2CSGOCSGOCCG+OfL3dCHEe0KImpH/007Xl4aGxucDY2H1A0T0HSnlZCJaRERfF0JMIaIHiWiDlHICEW0Y+a2hofF3gLHkzmsjoraR8qAQ4hgRFRDRFUS0YqTZ00S0mYi+90l9xWV5qOjrNUREFHiSBzb3J4B33lDP66xdIDN3J9rVnMfjjscXK2x6OVeFDFeCpbcp2ivPedwKrCANarTOTO7p5SuAumb6OKi5TMs5a9g6BPVed20Gq3Nsx2/XxdyaLtAKttE0BFb8q+/cwdrJLLDmKYc421jjg+WXjIc4ZcrjarRhJdz/YC9XR8oAvgdFGbAe83o46+nuwW97K6fDoohCXol2wpDW25WPe92yajOre+uX50TKq/5xW6T8bguPeydfwdx6VnExQPW2tHcogVQm8fFIKMA6cHVxdr7nBOIHqlafRERyPtZLyx2FkXLaTL4mtk6Cx2OIO6Yyz73WZejf7OFjpaboNnm5eWHapLC1a4f91OpRFX/V4Z4QopSIZhPRTiLKGXkpfPRyyD71lRoaGp8njHnjCyESiegVIvqWlHLgdO2V6+4WQlQKISq9zlP7E2toaEQPY9r4QggrhTf9M1LKV0f+3CGEyBupzyOiztGulVI+IaWcJ6WcZ0u1j9ZEQ0MjyjitjC+EEET0ByI6JqV8SKl6nYhuI6Kfj/y/9nR9uQfsdPjtcC7ropt5XPCmjciv5m/kKZeTWiDrtC2F7JvC05jRQDaEoPEZhnxzeyG3tS9AH3kpXMb3/1bJG/eNZlYn3oIMV58OJYb3EJcJVY8/YcjzlnEUMlj3XF7ps+F3yjFMzcB8zillboIKLGgzqKXilHTMhyAHplzSxtq19+H8IusVLru3rUC5sRvPWX6CM3qODNAY7OZLKaQEJXI044c3ndOregI+/+oKVuc5H5VvPLk8Uh4sNwzqdMUEewt/FmcGzk2GizD2ti5Or9iKObzqlm2sbs1ryyLl+HZOvzMZ49j9S6jYXNv4NzV1LdZ0D9f0kbscz2ntwBou3Mi9+NqWoY+Z+Xxtfrg7fO4R8I5NQz+WVkuJ6BYiOiSE+Cjs5w8ovOFfFELcSUSNRHTtmO6ooaFx1jGWU/1t9DEXjAjO/2zJ0dDQiAaiarmXnd5PX7vxTSIi+s2hczkhs8FGlvzOxuo674Pl18RvQwXmfMygWvkL2PTGCzn7bUtUrOmm4F7273GxovpesHJx2wtYXZLiBebfBBbYEs+aUfpRvCfbz+Vsac+tYN8S3uBqtCA0PuRZDhEkKY6raAIOnJVMvvkYq9vXCnFkqAzPdlP+YdbupZdXRcqdPN4IpSuGk0PFsF705BvSWLtBlyjl4khQCfppCiqpwQ3eed501Nl4bAkKKmyvyt7nTeLHSRWp+N38Qjmrq70BkyP8WAM2HuOCAsrx08Y2rk5ObMS8dy/kFqf2VqzB7hOwLnQYPpXds/AHy7BBTafkGhDKVHfN4mvTkwk6tu/gKs2khvCzdfBw/qeEttXX0IhB6I2voRGDiCqr39mXQr999RIiIorv5OyOOxtszEX/xT0N+pWUpxuWLI2Undv5eyuzTQlksY47wPRNRtm2DTHmWn7Eeb6VeWCdj22Yxur6y5T47QuU1E81PGadTwmOMXfqCVbXOYxTYGc6v65saUOk3PEyzLtMF3MNBSm/m385gVWlKCKN6iz0VP4i1i7zaqRtcvZwOryKRZ4ajMSVazi5d+E5F0+vYXWH1mDA1ZN8s4vPmar1SGo2OG6lKnkGnLiXdSpvNz/5ZKS85f7xrM5xQLmfsuSChrh6VuUA3W7holXoBuRosH/IM/UG4/FspX8BXR1zDdmaFZKvuPxDVrfhUaTvksplqdfyk3v7kzC37J3K98/gnLCoFXpNp9DS0NA4BfTG19CIQeiNr6ERg4huXP3iIpn/D98atS65DDLz4vx6Vvf2DiU4YQpUStKQhy2hCmrA4WIuB1p78Y4LJOGZcw2qIfdaBEW48743WN2jz1yG/uZD9zTo5Pq8tO0QIN1ZnEZV1itcxa0X/f+Oe/dOwbN4+XEFpSwAzcM+7hUmZ2UeAAAgAElEQVTnOaYEa1TURoGp3BNQ1EBVlHXA4LmXqZxlKPT2zuTt8scjyEX3Lh7ow6cEFRV+0JHYxL81IeXYwDrE12J8J/pouVjxNBwwqHEtynWSj3f6YfwOWpXY+Rdw3eHQSX7OocLeBZpXX72L1e35V+hCxVcxL00neEAQSz+Ed8ckfq4U2K1YgWbgmZPr+FgldChnCAt4XSAlfC7R/m//pePqa2hojA698TU0YhBRZfUT04vktAvDrH7nPM6NTFkAlcyJddz6ypsBGq1KLD01PRIRD/KQP6GL1bXWgvUq2IC/u3L5u29AvXU+t0YLKf1bWsGKB5I5C3z5oj2R8rbH5rM6q0uJuT/OwK45lJwBdkUFlsfTcPk9YHXjWrleSg3WoKrAjGmhQkqMOUc7n4uhYtTlbwN72cRTIVD5y6hrvIjfIJgDE7KUXTCLq7iJe1adfKwiUrZ9iac2bzqJOUvLR8ALZyO3yjSlwcml4Hku+vTdCQtIjwd1/iFOb7oSB29gHKtiokTeNr5fPHdCZBj+AAFBhksMATEsWCOJGTxYiMcNWkSj4mRUyp10cp/HOLYt4urC+I7wHNY+8xC523VcfQ0NjVGgN76GRgxCb3wNjRhEVGV8e36RLL3zASIiKljZxOoadsGrLHM/p6ldCWZZophFNtzI2xW/gPdY660+Vmd34LeaV09yUYnJcwllPNVxgg19mJ6C/Nl5JT8LcOyBei/A40Kw+wmDGKiqCONfgHqp/RxDvjkvnjOUzDuJT8F5QMmPUNc3k0c/77oUMnhJDo/R3rpVCRpZhXu3ncfpUANPlp1bz+qa+0G/exjnIYnbuOqzvwJ9mgzjEcrEeCcchnxb+B5Xh528CvMZmMDlZzXufd52jM2JK7gHqEk5vwmkcUJKXkV5oJSrEq2KllQ1A7Zew9XEzh1Qd154OVcJvrUO50Dp83Dd7EyeN+KtPYjgkVXE1ZG5ieGzjG13vUDOqk4t42toaHwceuNraMQgosrqO3KL5Pibw6z+YClnGy05YNH8Th6Us/hNlAeKwGr5edh7psoqXtnA6prWw9stVUlHve7hX7N2j/bOipTfbpvC6lo6FRGhV2EVTXwMHfng/0w7uEWYR1FNGlNv2frAodmUNE6lt3PPt8NtiodYFU/zparmMq+GONW8lecgCCheZbbxPJaeqxNWfbYOjLe3mEd5SN8B3taTzrlLTw7mt3Qa0no5X+TBTfoUa8BJj3AvxKbLwR7HDYLe3pncKtMyAPnJyOPm7lDUkZcqauFEgyi4G8/syeYTY+tWrP+4hMA+nabZEA1V1SERkdUKOmxb+cJ15eN+FiWWfkITp+Omb78TKT/19EW8j2lhcbPtR4+S96S23NPQ0BgFeuNraMQgohqIw+QnckQcL/g7J6kC7HF7E481Ju5DaGhnPU7Tk7K540lwL06uqxtyWZ1QWE9vmuJ08cMHWLuu5XACyvyQs2s0A6xXfPupM+IOpiphrZfzE/OSZLDV7c+UsjrneejIWw91QLebs/NCYeSy9nK2d/DL6L+uGWOVWc/ZRlceOnG18/FWp8aiGA3a93I+1+bEmPbO5/H4srdg7JqGwd7H2zgXmlMO9v7YtwzpxpSI4K4kXGfy8LUTzIMIsnoyj0H4bhqCqWTsxHJ3X8wnzdGN8XEV87EaLlDEs2wu7mS8A7F0KACxTuZwUTZrM34Xf/coqzv+BIKW5H+5LlJufJ5bsP52z4pIefJzXJRtvzQcnr7TdVoun4j0F19DIyahN76GRgxCb3wNjRhEVGX8oIOod3pYBjGP46mruiqhuknmMQZp+AjUV9aL4LFUntbL2h3MgIyVutugd1FeccFVsHpyFhnkeEWVKK7m6iXLfnhfpa6AJ5nrDX6ecOlkxLDf0syDP7Y9V4pbncst/iwnIddn74FMaF7PZV9xL+TM/jI+hYEdiNohstBH1xKDhV89nnvKf3Swupp7ENTRqmj6fIZYFY6jkOtNA3wc+5UYoKFSHBR4/Nxyz/Q6kiybuUhLwWk4wwkqadVSqwyqwwzM2Y59s1mdGI8zELMingeO8ZwGnefgWSzd/FlmLoU6tepNHnPf5FfOjpS492RQ1frjsQA/qOJrYvJtkNcbXoJr4ECFwWKzE2u65uvFrC7zQPiGJn7kc0qc9osvhLALIXYJIQ4IIY4IIX4y8vcyIcROIUSNEOIFIUTc6frS0ND4fGAsrL6XiFZKKWcS0SwiukgIsYiIfkFED0spJxBRHxHdeebI1NDQ+Cwxltx5kog+4rmsI/8kEa0koptG/v40Ef0zEf33J3cGx5Ti/+Ts2skrUHYVsiqyKimv8v8HrOKhm/JZuyXzEeTh2JHJrC6kxFsTW6D2sxmcdDLOB9s7tJHHkfOXgo8K/AksapBz+vRmFVRIoWE+xLkD4AELsrizSbddsR6rg3VXzxQutgQCYEvTLuTBK5blII7/unpYHjpe43y6Sxm6rt9yZs2+AWPlnAcLt5R9vF1/KX6H4jmP6TcpjkQDaCfH8aAiiYswBt7NPGa9Jw+Tk3xSCcDCw9kxlaM7l/PY/3HBc5Hyb9ZdHyl3GtrlvYd5ct3I52X/DsgtMpez35mH8NwWF/qYeO5J1u6QDay5UZTwP4m1lPNTsP2htSWsndkLmgfO4ePobg7vi5BhPZ8KYzrcE0KYRzLldhLRe0RUR0ROKeVHgmMzERWc6noNDY3PF8a08aWUQSnlLCIqJKIFRDR5tGajXSuEuFsIUSmEqAy6XKM10dDQiDL+KnWelNJJRJuJaBERpQohPuJtComo9RTXPCGlnCelnGdOSBitiYaGRpRxWhlfCJFFRH4ppVMI4SCiVRQ+2NtERNcQ0fNEdBsRrT3t3ewhkuPDX/3GC7gZatJJMAzOSQbm4ShUL403KEEck7icc7QbMnkwgZ8hDFRAnWVLx3XJb3I6BrYqcr1hdBLrlHTPNyGY58JM/s5rc4NeVX1HRGTxQkZ0+7msJ5WY8H6Ffn8iH49gF9R+nSH+7n7xxAL8iMO9Qvl8PIYnQHYPDPFoIUmKJ5waw35gLjdXLXoZ49Fj5jQmlsJTbaADYxyfwPsY9ChBSw3P+eUZ2yPl546ujJRV81oiIudyqEXjavmz/PR3X4qU3cuU57JzVao7C+cQw1U8mGfOXlzXYQi60nUHzluEMhddj5WydvFK3kXBjwlooFwxz96veFGW8IaJjejDsZerRf0jQ2wMLHMqjEWPn0dETwshzBTmEF6UUr4hhDhKRM8LIf6ViPYR0R/GdksNDY2zjbGc6h8kotmj/P0EheV9DQ2NvzNE1XKP3CYyHwvzJPEdnF0bUjQXuVN4vLK2DoX1coI9djXzOHL58+DO1WJQ16TtBw/UN1sJ9HEdt87zOsGW2g5ztnFwAsQF13FY023o5oEVkvah/9zreJqsAR9Y29W5PMDG248jBXjffLDicW1cJAilKkEkGjmN2YrjV5eSu2B4PA88kZSuHLS+z8dxsHT04CxLJtax3x9egJj4lgTe/3C1kspL+bv5MGejXYrnmz+DqwT/vAbsfUgRAwIz+SFxaQbUb/VdXMUrFOs/fyfY4+RKPm55L9VGym1P8Jxl/ROxdib+E/dC7PtPJTX7eoiJzgrWjIJxCv1JnIWXQhEh7egvYy/n23un4zo1rRcRkW3EGFV8VpZ7Ghoa//egN76GRgwiqqy+tIfINyF8om7v4ayWuQLWea0tnNVKSMMpvP19sMo98zlf09iG66yG9EPOHLDLGR/gBNdcbDhiVeBP4ixvVjGce1zbYD7md3PLusRW9Fm/nce6y1sAcWR3L7fMUq2x4g+BLTWe1JpqMHbJhgAbUnmV2/IQbEIYAvwNNUPzYDbEmMvdpWSpBbdN23dMYu2Kp8FqsGUft7pLmgwHqj5FFDLXcbElSzkxn/KdI6xuZxus3Qa7oQpeUcKt4j7YDEtJ26DBgScO1+WUIShKxi+5aNJwOxxnLi/exuqeq1yIdldy+s3rlbKiZFpx427W7u1350XKSSeMMd1RHP8s6Gq4lO+R8xbC+avpOwbnr6Xh9TLWU339xdfQiEHoja+hEYPQG19DIwYRXXVeUFBoMCwj+XgcBHIPKCqw9Zys/nFKOiY4MrF46kREjuO4Ti7nHlYmJfb9YDFk/EA3JyTODpVdnCET0dB2JdX2KsSsr63m8m3HQsUzzcbPEFJtEAT7f8blf/ts0JVcr1j4GbzFpmfDUnDrAS53W/swJqnxsDJzHspk7RYthyfj9uM8AkbzZSjfPBvpnl5ct4y1azoCt8RsQ9qzXh/OW+x+xdtvObe2DPlBb1sV14GJHoyHvR9j2jeBW62lH8a9Oy/gloHZmTg76urBWYPnci4/Z61CuqqD/dzfLKEWcr2ZG+6Rxa1YnE5Gud3D15VtsmLJmMVN14vXYnwa78O8W/fz9df0bQTpqLmDnzUIT/i8K2SIK3Mq6C++hkYMQm98DY0YRHRZfaJIjqOkZdw6z/IOePiBa3iWWnc3WLvMHSB5eAK3opJdYA3HpfN49i2DEBd6UtCftYGn6yIlvn9aNe+/+Rb8bn4fbHpOLWfn1bj94gJOR/VGsGuB23hs95xXFevCCSib3+fqzcIvwTwvrpuLOwkzoEbzbFRUjuN5zL3tBxFcIvEEXwaX3PRhpPzSm2DvE6bxDK0DdbDCG76ez1niOtT1zsG9LQZnnuSdmDPfebwPbzfqnr79kUj5pg/vZu1MU5QgKwa1ZVc1RBxpRV3/ND63jichrjUs5n2YlFh6xoy+5vFYLxWZmOu9SvCO8HVKaiyngYW/GGrpcZkY47pCLo70TlFEnCBXZdvbw+vAmIH5VNBffA2NGITe+BoaMQi98TU0YhBRTZNtKy+QhT+7N3zjBi6/BBIgJ0/8X25u27EAqpE5tx6MlDftnMbaqcE2q57i0cGcSjh0h6Kmc882JL5rhcwvDScg2RUIvtGjePGFWgyBLE7ifRowBB3ypSjpqXu4rJd+EdR0zV3wmAsaAnYWFcOjsG0fj/SZ2KiozmZDjk3dz/U8rkLQkdDE6RAhJW2zon0bLOHtpq4+Hinv3TeO1SWdVDzOFPL9PO4J+dIhq06bwfPBtf2xDPdGkdKOGXLbZWO8PZm8LlAK/ZtD8SA0fWAIPjoLD6qqGImIEqqVIB3FXLYueR2/W27HeJdl87Odhh7MZ+I7fBDG31EdKR9ZA/XszC/yHHsf7kad2cPnwlERVvnWPfA/5K5t1WmyNTQ0Pg698TU0YhBRVedZnSbKXhP2ZGtbzlVg9nyFvf8lVxv1N8L77YN3p0fKaYZUWycroPbqXcq9r/Jy0We3kq6LJOeKcnaCVRzO4ixfj2I2mKjE17A5OXvZvhw6leTjhiFWHtsoSniDuF9KMsajv42r8wbXQvU0+6ZqVndwI2QaaxdukFbDx8M5De98azV/zuzb6yPlln6wxIFjPGBHuwsimAgaxAVFraRkj6bABEN6ahvYY/ePuAVk/0UoJ7Sg/14u4VHZYogI6Tbe//63IPL5EiHGJffzOfPZQHDA8DnM2wExoHYCH6uWL2NcX1z4+0j5qnfuY+2SajAXtmt5LoR9m2Cx6FCmqepJLq6qq2CIZ9Ci+JfDg2zqG5t7nv7ia2jEIPTG19CIQUSV1fcnS2q9KHwKmnKAp2MaSAY735/Arels9aiLU4y7nIY0QmmvgIW3GcJJ+zcpobcv5Nep6JqrnMgn8hPc0tfADnbOBU1S8HtZBk99qGqbDocb9zAP4DEvDdaMfiUXUv12A4u9BPJCzbPcsSXjcrCRXXvxzIt+sYu1O7lzUaTsvZI7AbGMvgvAilsNCqA+F7QZcU7+DRmYhLFzNONZbDZuMeergbgwwOOSMK1HnMKaV9zIxZv965XT7ifqWZ37X8A7C8VRqzefP4xohvolay+no30+nm3i/wyxurr7UXft9nsiZUcT31oWF+7Xt5VrYoRNiTuoaIGsPKE09Z4PDYW5ie+RjmXhNRHYMjYtnf7ia2jEIPTG19CIQeiNr6ERg4iq5V7+1FR55/MriIioeoinoK77A2TVgTJWRb5syNbXLUAQw7WvL2HtvIWKLsTL32nmFMiWcVWQTS++cgdrt2YbcoQ42nkf6vnCsBIPPlDCozOkbYb8Ja7kFlz2P0Bed+Xw/l1K/AdHF+RbozWadQh1rkJ+DpFVqaRxWoWgFLlv8jMVX7JivZjJzyRMSpdqnPbBCi6fq0EizYPGAJLoM5QDOqT31Oqm/Hf5eDB16oXwOnS28SAXK2fBwu1QD1cJujx4btNO6BVNhkcZUHImGAO8ZClBRjrnGq0cUb50Fdbm9kfmczqUMydXMXehKxmPs53ed5AXwG0IgmpXzjyGxvMHELbwRLX96FHynmz+7Cz3RlJl7xNCvDHyu0wIsVMIUSOEeEEIEXe6PjQ0ND4f+GtY/W8S0THl9y+I6GEp5QQi6iOiOz9LwjQ0NM4cxsTqCyEKiehpIvo3InqAiC4joi4iypVSBoQQi4non6WUF35SP0kphXLOkrBFU9oPuUPGhCQ4wLz0wUJ+XS1Yr4EpYHHy13OWrL8c7zHXOM4KJR6Hk0pqneJYcT6n0TKoqPPyePy2ZRORZunACzAf82R8Aitextk6RyYsyzwt3FmjeLJi0fVrBNEI3G8IKqLkHcjN51aOaV/C75rvQnwKOrilZHINxk4VW4iIEpWpGTKo2FQEEpWUTh18LnypqAsqDlj2Nq7m8qajTlo5jXnvYy4S7kJMvOYtPFahJ0+Jk2gITKIGwIhvU9RmSZwbTq1DH323cZWd/XWICHO+tp/VvXsA6yC3EOKI42Gugu0vxforuJnnBVAtFhvvxdoM1XMPr09Kj5U+K7x/jtz/R3Idb/vMWP1fE9F3CQanGUTklFJ+NFrNRFQw2oUaGhqfP5x24wshvkBEnVLKPeqfR2k6KusghLhbCFEphKj0+1yjNdHQ0IgyxmK5t5SILhdCXEJEdiJKpjAHkCqEsIx89QuJqHW0i6WUTxDRE0RhVv8zoVpDQ+NT4bQbX0r5fSL6PhGREGIFEf2DlPJmIcRLRHQNET1PRLcR0drT9RVfMExzfha2h/zLe1yO35dWGikbPb3MXrwvUg5DVkps5LJY2xfwOLlZPHCj6xjUh82r0Z85hXuthYaVoBpOrqjYtwbynHUl5G5fO1cv+RRzUNMgH+IZeXg/7qnj5rbMO69T8Qjr4UEjhAWy8PDbXC3quh2/5y/DWeyejTz+/sB4JX77EB9v50KMSeYmjEHXOfzcxNKNuTAGeUxsUAKOKu5uWfs511d7o6L69HMG1HszZObOdgTNTDGkWA9ZMW4WN3+WhBa07Z6Dcl5FB2vXnolxu6qUB8B4ZcriSHlHaymrmzweZw/H98JlLvcBHkzW58Zz1mzieQySi0FXqB5/jzME5XRX4MxpYjH38GtbGz6MCfWPzQr/0xjwfI+IHhBC1FJY5v/Dp+hLQ0MjivirnHSklJuJaPNI+QQRLfik9hoaGp9PRNVyLyGzSE664ttERNTPuVymQvKtHmB17nqkPpKKJ9O8mbWsXd/3wWr1TjbEy78MrHlfPVQt1mzuqWc2gwXOSeHuUZ0bobhQ2a607VwkCFkVC6slPDBEwTNgjxsvMjBcys8Jk8BCNm3kURfcBeCrTT7eB4srosSwN6dz1aSlGjHay/+3hdUd/Uewvbnb0L8/gbOecYNKyqhlhnWUBLHA3AEvRKPFXMCB67L28LrOCyByqH3EjePrw6usj4Qmg7igaNVMU3Gdt5mrUrMqUe4416A3Ux4tM5+LkL3HlfAYWRjjkCFOYkalEhTlBh5Bpu4IrPXilEAaSXO7WTvfeog7rgV8XYV6RwLc/OIR8jY0fXaWexoaGv93oDe+hkYMIqqBOII2oqGREM3+ZH4MPDAe7yDHdn6KHZoCFkoo7Ovhdu6QET9RYe8NzI6zFiyZeooth7h11M+u+RPKxy9mdWp2pvGPg/6aO/j7M11h60KdXOTovQusoiHSNCNZZe9DBi8IcypY4JDP4PTiwr2XzEbAipaf8pRObUvwMPU3cNurtIOo8ylE2vu4ZZ3Zh3arFxxide/tmBEpl70GttSY5TWuDb/7uOKBOfQEFUefBDvXxPiUbLwD03nd0skQB2cmI8PxEycuYO086UqKq1qD85RiBerbwLMOxykKHVmINWHfzUOuq2HF66ryWd20mZBzG1/GiX93K98HX7wFDmXvt41ndXJT+H6dhmjxp4L+4mtoxCD0xtfQiEHoja+hEYOIqjrPnl8kS+55gIiI4nh8Rxoqg/xYuIHLkm1LIOsl1+HvFrchXZJDib2+jMt6ohey5MqlkEfXH5zC2w2fOlCEdQD9+3KUMwoLpyO+FkJ52jncwsrym4xIueVmrtuakAdrr/r1pZGyMTWzehiQuoL3H3wasf9Vb8Wv3Pg2a/fsI3CkVNOLEfF00unjYT23uqCKtVtTOzNS9vRwmTbhJM4azMpUDJZxVdl5C45Eytvem87q1GAk/nNxNpKWwFWwrbXwZEws5Kq+xOchJ7edq1gr9vN59mdjLmxNhqAlZUqgFafhjEJRv1kVo8ThfL6GTV48S3Ylr/PfjjHOTYQK+cjeUtbOrPRhMVhbJjWE5+zomw+Tq1ur8zQ0NEaB3vgaGjGIqLL6OVPS5fXPhFnMt+t4eqDSfwf7U/21eFZX/Do4l74KRVXGuS5ylYInLinnThIN9WAHzUpMNWM8eHcxWD5hiNunWnAlNKEPd45BzaWwZDPPPc7qdh9EVlmzm/dv78TvVdchDv6mP3HL6Oe/9atI+Zo9d7E6X52iXyoESyxbuVoxmAKWO/UAH8j+iaizDiiptgzspScTz138NhdbOuYreQcUrtootvgTleAYqXwc1XRjlIr+HVX8WbL3QJZoX8zZdNWS8YHrXouUX/zaRaxd3Q2KZtsQEIQU5yHh42OQdgR13lTUZR3g41F/BerMyYZ4eS14nkCmsv4snA6pWGkKNxdVMsvD4sKx+58iV81nF4hDQ0Pj/xD0xtfQiEHoja+hEYOIqozvyC2S424Jq/OGSrhaR5WV5KW9rK6vG95XSYchwwkD6YFlUPlYtnFzR1U+H5oL2bf0Kf7u654J2XSo2CDrKTApst6Dl69hdY/8/ouRsm8R9/Dz9kGeMw9wi+mkiQiU2d+k0G94Tkc+ApAEjhmCgORBRiwvwTmH+Fduapr+rzATPdTGTUgDdYp9qTI8QZuBEMV0+NdLnmdV39xyU6RcuA7yqLOcy6bxnYqMz492qG8W1ogapLNsRT1rd7wV3oShbp6P8JIl+yLlzU0wcw3u5+vj5ms2RsqvP3Qeq+tXrGOXreKmyVs+QHCWnKkY7zZFxUhEJG2K1+cWQ2DSREX+/wK8SI0mu0Ul8NZrrs5mdTIhPFbtP/kNees/w7j6Ghoa/3egN76GRgwiqqy+rbBIFt4fDsRh8XBu5LxLkZv4rb3cgkuN2RavqNESWjntg0XoM2g/dax7VaWUsIrHXgu8BBaq9zyeGithH6zThkrAuiXX8fdnQGFZg4Z4IKnHcd1ACb/OMxUiyPwysOInH+NRS3ovRrvLKw6yuteOwppODkAssnVy9vKe69ZFyo/s4skF4lpwHVOxGRjInHKwnsPv8Nh/PoVL9ZTCs27VFG79V9mOGPklqTxHwIEqeCjGpaKPQAuXCbImgY6OFh7PPr8IrHO8FWJQyyYem19dE950vnZUz0BrIxcl0ueBvU+OQ7veP/P+kxsgFvlSuIg3nI11EIzDIJ9/O0/vtv+7syPl5hVcbZmzIGzBuf/r/0tDx9s1q6+hofFx6I2voRGDiC6rX1Ioc3/wTSIiSjzB2R2XcspvTKVkb4Vl2SexjeuPwBow+QBnhazDeE53Fjih4RJDEDjFQYUM8exSjoBmV5GSLTeZayhs6WDFbR8msbrBcrTN2c45so6LwQ6qVlozJzSxdoN+sJsnj/JgJCqKJ0GMaTqay+pCyikzmfgamFqBmHDVu0ojZXsXp1d1RMn9kPcxdDOcZQJBPEtZBtfYuPyYp3Yn11AETkC7EEgDL15ezsWzpkoEEonr4zR600BXKuKSkCfL8Cx5igPPIK8LKP5HoXi+NrN24tmm3YsT/7Y7eXCTE9cjEEwg3iBKKCnGvrpsU6S8vpNHJmncUajQZBBlR1K/NT72EHlatJOOhobGKNAbX0MjBqE3voZGDCKqwTZJChIjsqsx5W9qMSJzOBtTWZ3qwTXuaSXYZj5X++Uoks0wF2mpZ7Eiyyted4l13DNtqAJydnwWT/c0MA+ydfa7kE07l3CRKqTIpoNzuEowrh59OCfy6zLSYZGX9m8QLA9fU8baxY9TAnY2cDWdRTnLCFWgf0exwYKwBvL0n679Lau763f3RcqBSRiP4IAhf4CServ9Cn5WEr8NajXbOVCptfQb0oG9hXbeqVx+jivDeFiO4qykdAY/J2h1Q/aNW8pTinuqIVt3L8I5gTWF5xlIT1Tm6dUMVtdzjnL24jGk4b4R6rwuL+a96ns8iOu4ApybNG/lqr6kCVj7a34O1WrcrfwsY/KyE5HyoUZubZkxJbwmWv9kOLM6Bca08YUQ9UQ0SERBIgpIKecJIdKJ6AUiKiWieiK6TkrZd6o+NDQ0Pj/4a1j986SUs6SU80Z+P0hEG6SUE4how8hvDQ2NvwN8Glb/CiJaMVJ+msI59b73SReY/ETxbeF3TfwqHiijqxpOJKnjOOPQ1wU2r20xTOF8KVylceUFsHRas4ln4zXbweZJJ/qwGOKQ3zl/W6T89OFFrC79fbDpHSvQnzCoH0UvWGKTmdMolVetZRpPxzT8AcYgux1prXKncja9ZzvkmPhBQ9xBJc1Vm5JlV5g4jSvP2x8p37jxHlZnVlRbpFhNerJ4H8lV6vIxqGcL0faKYmSffXbLUtaOpqFd+cucTXRVE08AABS2SURBVD3xZYhhpgQ85+4XZ7B2NuWy/GQec082gG1PVTKutd3Bx00dezOXNEkqcRgz9nFWv90MUaI9iHvFpXER7+ReiCPLL+KOPpWvQGR1ng+xInUtV9W6LoaFoioWEhGZHxuhv+uzzZYriehdIcQeIcTdI3/LkVK2ERGN/J99yqs1NDQ+VxjrF3+plLJVCJFNRO8JIapOe8UIRl4UdxMRWZLTTtNaQ0MjGhjTF19K2TryfycRraFweuwOIUQeEdHI/52nuPYJKeU8KeU8S3zCaE00NDSijNOa7AohEojIJKUcHCm/R0Q/JaLziahHSvlzIcSDRJQupfzuJ/VlKy+Q+f/ydSIiciTwuPe0E/KoMZWyvRc0ZuyF/O/8BW/o3AkPMWkwWlx8IeSqLceVPHJ9hsR0aaDLXs1d69T01KqnWulrXPbtnIM+vdP5IYLsQJ/5W/nYi3vw7uzdBPnOM53HkQ8q+fE+Zr7aDQHVvgvqpS3f/hVrt+h/vxMp+4u4asvchrMMNdhJsIDLrUIN0mFQc1k7IZ+nQcQnT6ZBhXkI93Z/hydbGHgf86l6OfrH8fGQQSXefAv3nvNnYc6SqkDTP9z1Imv3421XRsqmQc4IF27E/HpSDfH4lW9Z6U04RNhXVcraWXvQZ8ZBPu/OiRjI/G0YD6MXX/eNWEvx63nmxYERjW/zfz1M3ubTm+yOhdXPIaI1QoiP2j8rpXxbCLGbiF4UQtxJRI1EdO0Y+tLQ0Pgc4LQbX0p5gohmjvL3Hgp/9TU0NP7OEFXvvOSKHLnwsXAstuPHueWRGlfO7DJYRxWCxclNh7qmfS83z1O9x4Kc4yM5F9eFFCuw3B080LspAEK6ZnExQLUg9GYrsef7OL25c5DWqrWb64ZSkrk1oIqBo4rFmDIegXROY5yisvHlGCzmToDm4TLUmYYMNE5WLM76uAeh9SD4V7+iMs2dzdN1fbEQ8eye/SVPKd61EOMz6TGM/fI/72Ptfr/znEhZxHFzzty38Cxxd+Deao4EIh5j3pgnITQR422vxHOp3nhERMFE/DZabLpbwVZbMg2WmDbMTcpzaNd6IX+W5COK+JfG91zcDIg4WYm4d/0hvkfMSvCaSUtOsroce1jlu/bWN6j7WLf2ztPQ0Pg49MbX0IhB6I2voRGDiKp3Xqg7jgZ/HzZdtE/k75ykhV2Rsn8dl+GyX4L81bkAai7/FC77+pUcds+teJzV/b7z3Eh563GYSDbewGUx1cQ2kcc6JApBdLIqKh9/EpfZHD+GzJz2T1xe9G6CaajfYNaQ3IF+Zn0Z6sdDPdx0M74UsnvHNi4HqrQUFcPE03geMjsDJsGVIe4t5l8CdZnvGM4d4h/kQS4fvRFyfYohj6HVCbm76qvwBHw8dQ9r99J+nA97Mnknt/8Yue7W9yC60qVLuMnra804e+7fxJ/T8S5odivxQDMOcjFY3ACPv94B/pz2DjyLf9gY/B/9DN4MVbPdz7dWyAoZ3xiAdagF4zM8jDOhrKldvN0WGMfeW7CJ1d338h1hGgY30Figv/gaGjEIvfE1NGIQUQ62WSTzHgwH2zSmhbJ1gZ3yFHGrPmuCkjq4DqxW9h6ukmlbqrDiQ/ydllyHGy78ZmWkvG7LXNYu9Rj6yFnfwurafwMerSIdbNjO3TzuvSmbq3xU2B14tszHOa9ffz2ex6R6/BnGSipmiVOK21hd15OlkbLVhT46FvLxiJ8IFdJAN6cjvg5saUDxisvaZ7BQnIs+LW7OOie0KMFNvwB1nsfD2flQD/SupWu56Na8So3vD5FMDXpKRHThHR9GynVDPFXYgSZ4xSXswNqJG+CD2nUO1pili9NocSlqYkOQS4cSwj57D9TOPdMdrF3pjbDqO7BnHKuTFvTpaMU+MKYUL3oTolv1D/iciZGAqS0/+B15T7RodZ6GhsbHoTe+hkYMIqqn+g6Hl2ZND8cNuzN/K6v7wW/uiJRTjxucE2aC9VKdLjypvJ29HCzlcCt3Ygg24x23+z/B3idm8Xff3d/BSfJ/ll3B6m4ofj9Sfu5NWJz99Gru8PEfv70eNPVy9nigFOxm0ypDAI9h0HLxAgTK2Pm7Oaxd3yqIEmpceiIi56XQIvh6IZok1XDubzAXrGhaJWdtxaVgKd07wTq3X82deVRPqFAtZ21VRytXO9hSexZ3sLEXIxhJ1+x0Vpc+A9aF5SmIpbcjOJG1e+kgxsfczk02ExSut38y1s7kh7pZu65FeE6j2BKcgqAX/j5+JB/XjzXYvBJzW7Sea3Ou+xbEywOSs/oqex9SpsI7jY9V1TSsaVMbf07x0VLyje1brr/4GhoxCL3xNTRiEHrja2jEIKKqzrMXFsnC+8Jpso1x9TMPQN51jueeZAHFK27C0vpI+UgVtzhL24/rAhfxoA7ew0qAih4l3vxqHjgo8AqsBnvncBlcVas9dOEzkfKD+67ifTRBphUBLi+qqiFPIfesS81BUM3B4whTJozqPKXL5csPs7ovpB+IlP/1kS9FyvNuPcDard87NVJ2ZPFgIe4ByLGFf1HG9G4uF/ccgCWZMS15fAnOW5Kfg2Xa0E08wOhAF+TWy2bvZ3V/qURaaItiCWgcD+tE3CsjkT9L235Y8iVPwTnBwDEeO9+saGC9BnVy8auK5V4C/1b6kjAZak5Gs0GjqwaToWv5ODoP4nxB3RfquiciylMCt7Sey6ooviVM48k/PkTuNp07T0NDYxToja+hEYOIKqufkFkkJ18WZvX7J566nZjAY4YHGsE6Z+3F33uncI4mbqqSWuolnnK5vwzvOEc3ntnFsxlTUj3KPQZW/5KFYEVNAnWtbp4WqqEfainLs1xFNVCipO9q4WN/yXe2RMrvd42PlHMcPK7+9gOIGZhVwnMQ9A9BrRZsgnpJDeJARGTyKSJHPjcRS8zB+Js3QkTqn8ZFk4njYTXYtKmY1dkXgK0ePIoxSD/CmlHnEvC2Ji//DpWtBct98gqoLY3pwITC+weDvA+xF+ugZHV9pDxsUIM21EJsEQmG4CxdaBuK43Om0qxa4K00OBJtOYH59Lu4+tQ0BJXghBlIiW4x5EIY8EIE69nCHbeSloVF1iP3/5Fcx9s0q6+hofFx6I2voRGD0BtfQyMGEVUZ35FXJEvveICIiBZdcZDV7WqDjOit5jKzdQAii2+aElt8Jw+KEFB+2vr4cyV0QF5qWalUJHF5LnkPTCFDhpD7QeW3mI/zBLmb0ztcrngTDnPVpKquCcVzGa6gFGqe1lqoFTP28vezU3EGlIbcfCkVCCjh3QY10XAx15/GKSo8Vf1IRGQvgwztOYmgIqoqkoirFUOG4KZq7PjO1ZDVrY28oTSPXiYiim/FDYbz0Z+6HoiIPPl4NrOLj1UgEXXpB3AD8+VcpaYGfxks4/OS0II+x11Vw+rqejHGpnehgrVcwvsf2Af1odlnVPGibHXhOfumG9TJycoZiyEfRM7EsLfo4fue1jK+hobG6NAbX0MjBhFV7zwKIS31jrU81bHKHif2ck7FVaKwcs1QVw0ZWLIrl++KlDe1TGB1bdVQKWUqYd/M13MLv85xYPnyJhlSee9B0DZ/CyzOJqxqZO0KE9Dn0b4cVte7B2ojUz8fftt4iB0mL8bA7jR48ZXAaytg8Bbr7YL6SijsfXZZD2vneVtRXxmsKD2DEF3MiihhVAmqLOrgeN6JNOGbkr5VYe8NkqWq1pWF3NwtcSbEqYEWzJ91kLO5b3/hoUj5yt//I++jEWNsU8bxxtJdrN0jSyD/2ezcci9+Kuga+j7X/wYewLpN+gLWi/s9njw6TtHgWbi2mhLbMXYt5+Hvtm4u+2StVyz3ruWeku1tYTEj4DfIS6fAmL74QohUIcTLQogqIcQxIcRiIUS6EOI9IUTNyP86Fa6Gxt8JxsrqP0JEb0spJ1E4ndYxInqQiDZIKScQ0YaR3xoaGn8HGEu23GQiOkBE5VJpLISoJqIVUsq2kTTZm6WUFafqh4jIVlooc394PxER5W3kLEnXHLCRtj7OUiY1gEVrP0dhe20GFlhNE2V8pSkZVROKlYAdLn7K7DgMUSK5nvef/3XETTvYgrDWdju3aPMchbVbINHQB4zzKBjHnzPtbogMNTtLImXjaXcoB2yeyXCqv6C0PlLefgwBH/ILe1m7zv0QQUxlPGhEYQZEFZNiFVecwK0E9zwDcW1gniGTbi/YcZENem2HecAONRy40SouqKQOW1xRFynv3jqJtZPFEH1sh7imZ8JFuG5CEljxNZsWsnaqJWPyCVbFnG9svZxGXyrqUlYgzZfbx63z+k+AIbYVcl7f3YMxST2I64xBXIZzsahzt/M++iaFNTPHXn+YXN2fjZNOORF1EdFTQoh9Qoj/GUmXnSOlbCMiGvk/+5M60dDQ+PxgLBvfQkRziOi/pZSzichFfwVbL4S4WwhRKYSoDA6eOmGkhoZG9DCWjd9MRM1Syp0jv1+m8IugY4TFp5H/O0e7WEr5hJRynpRynjkpYbQmGhoaUcZp1XlSynYhRJMQokJKWU1E5xPR0ZF/txHRz0f+X3u6voRPkL0lLMNk3ssFqf63yyPl4QKuGhIBCLmWQYgvAYN8W/EE1D/uh7jMOfwsvJmcw1BXxfcavdZQHizk78UeD15cKe+inHADj23fmASV2rgXuWoo8CPI2i21XNXnHYSVnFlR58lJXJ5bUAgPrl3b+bHKvJkNkXLebIzH4a9NZe16/kHps4oHJm2Ih5xs78IY1E0xpKcuhwxqsXELSLOiZsyZirOBwR6ufhxWVLUJ2Zwj9FVhHHcM4jmNizbYh3MabwaXi/v/HRahr92KNWAx5F1Q02n3pPBzH7MSBDVo5+sld1FrpNx0GEE/Ll7O04F/8C4s9yzj+ZqwHcRzqsE2PWkGGpWznpbvGtJwvxDeC8Jg7HcqjFWPfx8RPSOEiCOiE0T0ZQpzCy8KIe4kokYiunaMfWloaJxljGnjSyn3E9G8UarOH+VvGhoan3NE1XIvrj9EJW+FHUCOm8tZnVdhf4rXcBan6Tqw7QmHwCoGy3nc8ZPXwrorcICLAeaLwcoF28DKDo3nLKq9DbyW0aLtinzErXu1F2x6QxNP25RcB56s8ULO2loGwFZb+7meLr4MakGXCfQnvcvPRirHQZ2Vt5vzdoOrcb+3X1mE+y5lzUgcRTloVKMpKkhvPsb+lqm7Wbt1v0Zugd4sgzOSDX12bYbq073MoPYbwHgHAoYjJ+Vnch1+pNZy9alzHPpwFfJnSftBfaTctBPBMKwG6zlRiTG2cI0jWZW4HwHDMVX/63g2UYJ7v1k5k7VLVlh451Ee708oNKvOTT1X8PiBogY3DxpSkT3+778jIqK7jox61PYxaFt9DY0YhN74GhoxCL3xNTRiENENxJFbJMtvDQficOdx2dSs5Cv7ypXvsrrH31kdKSc24F3lXc6DLibFQ368tmQvq/v9Xy6IlKXyuiud18za1TXBADHhCFfruLOVYBClEBIn57Szdoc/VAIr5hg8vWrQp38GFzQLMqB+a90F2TFvPlcXZjpw3d6DPA8bJSnyrxNyoEzgBxa2JJjRms18LpLWQK3YcQ6uE4a8bAkNkOuHxnO5O2cr6jovwBjYq/mZR8qyjki5p58L0H4l999dy2Dr/PJ/r2TtfMlYOw/d9XtW9/9+8pVIOfM9qJCrH+LBKk2NEOxtvacOOBI0BBxRzYz9RRjT7KwB1q6jGSa7RSU8SEfbAagBg7noI3m3QfW5BOdU/j5OyEfnRc2/eZg8zTquvoaGxijQG19DIwYRVVZfCNFFRA1ElElE3adpfqbxeaCBSNNhhKaD46+lo0RKmXW6RlHd+JGbClEppRzNICimaNB0aDrOFh2a1dfQiEHoja+hEYM4Wxv/ibN0XxWfBxqINB1GaDo4zggdZ0XG19DQOLvQrL6GRgwiqhtfCHGREKJaCFErhIhaVF4hxJNCiE4hxGHlb1EPDy6EKBJCbBoJUX5ECPHNs0GLEMIuhNglhDgwQsdPRv5eJoTYOULHCyPxF844hBDmkXiOb5wtOoQQ9UKIQ0KI/UKIypG/nY01EpVQ9lHb+EIIMxE9SkQXE9EUIrpRCDElSrf/IxFdZPjb2QgPHiCi70gpJxPRIiL6+sgYRJsWLxGtlFLOJKJZRHSREGIREf2CiB4eoaOPiO48w3R8hG9SOGT7RzhbdJwnpZylqM/OxhqJTih7KWVU/hHRYiJ6R/n9fSL6fhTvX0pEh5Xf1USUN1LOI6LqaNGi0LCWiFafTVqIKJ6I9hLRQgobilhGm68zeP/CkcW8kojeICJxluioJ6JMw9+iOi9ElExEJ2nk7O1M0hFNVr+AiJqU380jfztbOKvhwYUQpUQ0m4h2ng1aRtjr/RQOkvoeEdURkVNK+VFkkmjNz6+J6LtE9JGnUMZZokMS0btCiD1CiLtH/hbteYlaKPtobvzRPIZiUqUghEgkoleI6FtSyoHTtT8TkFIGpZSzKPzFXUBEk0drdiZpEEJ8gYg6pZR71D9Hm44RLJVSzqGwKPp1IcQ5p7vgDOBThbL/axDNjd9MREXK70Iiaj1F22hgTOHBP2sIIawU3vTPSClfPZu0EBFJKZ1EtJnCZw6pQoiPwrFFY36WEtHlQoh6Inqewuz+r88CHSSlbB35v5OI1lD4ZRjteflUoez/GkRz4+8mogkjJ7ZxRHQDEb0exfsb8TqFw4ITjTE8+KeFEEIQ0R+I6JiU8iGlKqq0CCGyhBCpI2UHEa2i8CHSJiK6Jlp0SCm/L6UslFKWUng9bJRS3hxtOoQQCUKIpI/KRHQBER2mKM+LlLKdiJqEEB/FEv8olP1nT8eZPjQxHFJcQkTHKSxP/jCK932OiNqIyE/ht+qdFJYlNxBRzcj/6VGgYxmF2daDRLR/5N8l0aaFiGYQ0b4ROg4T0Y9G/l5ORLuIqJaIXiIiWxTnaAURvXE26Bi534GRf0c+WptnaY3MIqLKkbl5jYjSzgQd2nJPQyMGoS33NDRiEHrja2jEIPTG19CIQeiNr6ERg9AbX0MjBqE3voZGDEJvfA2NGITe+BoaMYj/D0QizCICmgBNAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "d = 4096\n", "\n", "# Randomly generate a dataset\n", "X_random = np.random.randn(n,d)\n", "\n", "# a sample from the random dataset\n", "plt.imshow(X_random[0].reshape(64,64))" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "150 random samples from 4096-dimensional normal distribution.\n", "\n", "Johnson-Lindenstrauss lemma\n", "Target dimension: 4294\n", "Target dimension too large.\n", "\n", "Gordon's theorem\n", "Geometric complexity of randomly generated dataset: 3.740\n", "Target dimension: 1049\n", "Distortion for sample projection: 0.087 < 0.1\n" ] } ], "source": [ "print(\"{} random samples from {}-dimensional normal distribution.\".format(n,d))\n", "\n", "g_T_random = geom_complexity(X_random)\n", "m_G_random = gordon_min_dim(g_T_random, eps=eps)\n", "\n", "print()\n", "print(\"Johnson-Lindenstrauss lemma\")\n", "print(\"Target dimension: {:.0f}\".format(m_JL_n))\n", "if m_JL_n < d:\n", " d_JL_random = test(X_random, int(m_JL_n))\n", " comparison = \" < \" if d_JL_random < eps else \" > \"\n", " print(\"Distortion for sample projection: {:.3f}\".format(d_JL_random) + comparison + str(eps))\n", "else:\n", " print(\"Target dimension too large.\")\n", "\n", "print()\n", "print(\"Gordon's theorem\")\n", "print(\"Geometric complexity of randomly generated dataset: {:.3f}\".format(g_T_random))\n", "print(\"Target dimension: {:.0f}\".format(m_G_random))\n", "if m_G_random < d:\n", " d_G_random = test(X_random, int(m_G_random))\n", " comparison = \" < \" if d_G_random < eps else \" > \"\n", " print(\"Distortion for sample projection: {:.3f}\".format(d_G_random) + comparison + str(eps))\n", "else:\n", " print(\"Target dimension too large.\")" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Target dimension 910 for faces is smaller than the target dimension 1049 for random data. \n", "The geometric complexity of the face data (3.465) is smaller than that of randomly generated data (3.740).\n" ] } ], "source": [ "if m_G_faces < m_G_random:\n", " print(\"Target dimension {:.0f} for faces is smaller than the target dimension {:.0f} for random data. \".format(m_G_faces, m_G_random))\n", " print(\"The geometric complexity of the face data ({:.3f}) is smaller than that of randomly generated data ({:.3f}).\".format(g_T_faces, g_T_random))\n", "else:\n", " print(\"Surprisingly, the target dimension {:.0f} for faces is greater than the target dimension {:.0f} for random data.\".format(m_G_faces, m_G_random))\n", " print(\"The geometric complexity of the face data ({:.3f}) is greater than that of randomly generated data. ({:.3f})\".format(g_T_faces, g_T_random))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "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.6.5" } }, "nbformat": 4, "nbformat_minor": 2 }