{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Logistic regression\n", "\n", "To quote Daniela Witten, one of the authors of the book *An introduction to statistical learning* says\n", "\n", "> \"When we raise money it’s AI, when we hire it's machine learning, and when we do the work it's logistic regression.\"\n", "\n", "This is one of the most important methods in statistics and machine learning to learn/compute a classification method. " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "It relies on the sigmoid function\n", "$$\n", "\\sigma(x)=\\frac{1}{1+\\mathrm{exp}(-x)}.\n", "$$" ] }, { "cell_type": "code", "execution_count": 65, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 65, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAb8klEQVR4nO3de3SV9Z3v8fc3O/cQwiUh3O8XCRZbCWjFG/WGznFx2qkenU5bqTPWHu20a521pp6ZNe1Z41lnprc5M9PacjytWjunWme0FZWqtaUqoJUgCHIJbK4JAXIhBEIISfb+nj+yxRgC2YGdPPvyea2VtfNcknw2e+ezHn55nudn7o6IiKS+rKADiIhIYqjQRUTShApdRCRNqNBFRNKECl1EJE1kB/WDS0tLferUqUH9eBGRlLRhw4ZGdy/ra1tghT516lSqqqqC+vEiIinJzPafa5uGXERE0oQKXUQkTajQRUTShApdRCRNqNBFRNJEv4VuZo+ZWb2ZvX+O7WZm/2pmYTPbbGaXJz6miIj0J54j9CeApefZfiswK/ZxH/Dji48lIiID1e956O7+hplNPc8uy4Anvfs+vG+b2QgzG+fuhxKUUUSSnLvT3hmlraOL011ROiNROiNOZyRKV8TpjMYeI9EP10WidEYddyfqTjQKDkS9e507RL3HMhCN+pl1xPaNOrF9z74VeF+3B+/rjuG9V/W9z4V9r752qpw6imtn93lt0EVJxIVFE4CaHsu1sXVnFbqZ3Uf3UTyTJ09OwI8WkUTqjEQ53NJO3bFT1J84zbG2DprbOjnW1hn7vINjpzppOx2hrbOLUx0R2joinOqM9Flu0s3so8v3XzcjaQvd+ljX50vr7o8CjwJUVlbq5RcJQCTq7G1sJVzfyq4jrYQbWtnf1Mahlu4S76uYi3JDjCjMZWRRDiMKchlTnEdhbjYFuSEKc0IU5oYoyM2mMDdEfk4WOaEsskNZ5GRZ7PPYY5aRk51FTtYH64ws6/4w48yjmZHVc5kPl7O6V3xke1asMbN6Nydnlyn0XVrWa8e+9+n/64KUiEKvBSb1WJ4I1CXg+4pIArSc6uTtPU1sPHCMjQea2XKwhbaOyJntE0YUMK20iGtnlTF+RAHjR+QzrqSA8uH5Zwo8N1snxKWCRBT6SuBBM3sauAJo0fi5SHDcnZ1HWvn9jnpWV9ezYX8zkaiTEzIqxg3njgUT+djEEcwpL2Z6WRFFeYHd0kkSrN9X0syeAq4HSs2sFvgWkAPg7iuAVcBtQBhoA5YPVlgRObeGE6f51cZannv3IDsOnwBg3vjhfOW6GVw3p4yPTSghPycUcEoZTPGc5XJ3P9sdeCBhiURkQLYfOs5P1+xl5aY6OiJRPjF5BA8vm8fN88ZSPjw/6HgyhPR/LZEUFa5v5XuvVPPy1sMU5IS4a9EkvnjVVGaUDQs6mgREhS6SYlpOdfKdl3fw1DsHKMgJ8fUbZ7H8qmmUFOYEHU0CpkIXSSEvv3+Ybz7/Po2tp/nCJ6fy1U/NZPSwvKBjSZJQoYukgPbOCN96fiu/rKph7rjh/OSLlcyfOCLoWJJkVOgiSW5/00m+8m/vsu3QcR5YMoOv3zibnJDOC5ezqdBFktjGA80sf2I97vDYPZV86pLyoCNJElOhiySp13c2cP/PN1BWnMfP713ElNFFQUeSJKdCF0lCr249zAO/eJeZY4r52ZcWMqZY55NL/1ToIknmrd1NPPjURirGl/DklxZRUqDTESU++suKSBJ5/2ALf/lkFZNHFfLEPQtV5jIgKnSRJHG4pZ17Hl9PSUEOP793ESOLcoOOJClGhS6SBDq6ojzwi3dp6+ji8eULGVdSEHQkSUEaQxdJAv/wm+1s2N/MD+7+BLPLi4OOIylKR+giAVu15RCPr93H8sVTuf2y8UHHkRSmQhcJUGPraf72V1u4bGIJf3Pb3KDjSIpToYsE6JvPv8/J0xG+d8dlupxfLpreQSIBeWnzIVZtOczXbpzFLI2bSwKo0EUC0NLWyTeff5+PTSjhy9dODzqOpAmd5SISgH/9/S6OtnXw5L2LyNZQiySI3kkiQ2xPQys/W7eP/1I5iXnjS4KOI2lEhS4yxP7Xqh3k54T4bzfPCTqKpBkVusgQWhtu5LXtR/ivS2ZQVqyp4ySxVOgiQ8Td+c7LO5gwooAvLZ4WdBxJQyp0kSHy+s4G3qtt4cFPzSQ/JxR0HElDKnSRIeDu/MvvdjFhRAF/evnEoONImlKhiwyBteEmNh44xleun0Futn7tZHDonSUyyLqPzncydng+d1Tq6FwGjwpdZJBt2N/M+n3N3H/ddPKyNXYug0eFLjLIHl+7j+H52dy5cFLQUSTNqdBFBlHdsVO8vPUwdy2aTGGu7rQhg0uFLjKIfv72ftydz185JegokgHiKnQzW2pm1WYWNrOH+theYmYvmNl7ZrbVzJYnPqpIamnvjPDUOwe4qaKcSaMKg44jGaDfQjezEPAIcCtQAdxtZhW9dnsA2ObulwHXA983M01ZLhnt1xsPcqytk+W6KlSGSDxH6IuAsLvvcfcO4GlgWa99HCg2MwOGAUeBroQmFUkx/++PB7hkbDFXTBsVdBTJEPEU+gSgpsdybWxdTz8E5gJ1wBbga+4e7f2NzOw+M6sys6qGhoYLjCyS/LbVHWfLwRbuWjiJ7uMckcEXT6H39W70Xsu3AJuA8cDHgR+a2fCzvsj9UXevdPfKsrKyAYcVSRXPVNWQG8pi2cd7H/uIDJ54Cr0W6HkC7US6j8R7Wg48593CwF7gksREFEktp7si/HrTQW6aV87IIv0pSYZOPIW+HphlZtNif+i8C1jZa58DwA0AZlYOzAH2JDKoSKp4bVs9x9o6ubNSFxLJ0Or3Sgd37zKzB4FXgBDwmLtvNbP7Y9tXAA8DT5jZFrqHaL7h7o2DmFskaT1TVcO4knyunlkadBTJMHFduubuq4BVvdat6PF5HXBzYqOJpJ5DLad4Y1cDDy6ZSShLfwyVoaUrRUUS6MX3DuGO7nkugVChiyTQC5vrmD+xhKmlRUFHkQykQhdJkH2NJ9lc28Lt88cHHUUylApdJEFe3Nx9Nu+fzB8XcBLJVCp0kQR54b1DLJw6kvEjCoKOIhlKhS6SANWHT1B95AS3X6bhFgmOCl0kAV7cXEeWwa2XarhFgqNCF7lI7s5Lmw/xyRmjKSvOCzqOZDAVushF2t3Qyp7GkyydNzboKJLhVOgiF+mVrUcAuKlChS7BUqGLXKRXtx7mskkjGFuSH3QUyXAqdJGLcKjlFO/VtnBzRXnQUURU6CIX47Vt3cMtt8xToUvwVOgiF+HVbUeYXlrEjLJhQUcRUaGLXKiWtk7e2t3EzfPGat5QSQoqdJEL9Ied9XRFnZs0fi5JQoUucoFW76hndFEun5g0IugoIoAKXeSCRKLO6zsbuG52GVmamUiShApd5AJsqjlGc1snSy4ZE3QUkTNU6CIXYPWOerIMrp1VFnQUkTNU6CIXYHV1PQumjKSkMCfoKCJnqNBFBujI8Xa21h3XcIskHRW6yAD9oboegCVzVOiSXFToIgO0ekcD40ryuWRscdBRRD5ChS4yAF2RKGvDjVw/p0xXh0rSUaGLDMDmgy2cON3F4pmlQUcROYsKXWQA1u5qBOCqGSp0ST4qdJEBWLu7kYpxwxlVlBt0FJGzqNBF4nSqI8K7+49x9SwdnUtyUqGLxGn9vqN0RKJcNWN00FFE+hRXoZvZUjOrNrOwmT10jn2uN7NNZrbVzF5PbEyR4K3d3UhOyFg0bVTQUUT6lN3fDmYWAh4BbgJqgfVmttLdt/XYZwTwI2Cpux8wM11xIWlnbbiRT0weSWFuv782IoGI5wh9ERB29z3u3gE8DSzrtc+fAc+5+wEAd69PbEyRYDWf7GBr3XGu1umKksTiKfQJQE2P5drYup5mAyPN7A9mtsHMvtDXNzKz+8ysysyqGhoaLiyxSADe2tOEOyyeqfFzSV7xFHpfl8N5r+VsYAHwJ8AtwN+Z2eyzvsj9UXevdPfKsjLddlRSx9pwI0W5IeZP1OxEkrziGQysBSb1WJ4I1PWxT6O7nwROmtkbwGXAzoSkFAnY2nAjV04fTU5IJ4ZJ8orn3bkemGVm08wsF7gLWNlrn+eBa8ws28wKgSuA7YmNKhKM2uY29jW1cZXGzyXJ9XuE7u5dZvYg8AoQAh5z961mdn9s+wp3325mLwObgSjwE3d/fzCDiwyVdeEmAP1BVJJeXOdfufsqYFWvdSt6LX8X+G7iookkh7W7Gykdlsfs8mFBRxE5Lw0IipyHu7M23MTimaN1u1xJeip0kfPYeaSVxtbTLNbdFSUFqNBFzmNtuPt2uYt1Qy5JASp0kfNYG25k6uhCJowoCDqKSL9U6CLn0BmJ8se9RzU7kaQMFbrIOWyuPUarppuTFKJCFzmHteEmzOCT03X/FkkNKnSRc1gTbmTe+OGM1HRzkiJU6CJ9aOvoYuOBZp2uKClFhS7Sh/X7mumMuMbPJaWo0EX6sDbcSG4oi4VTNd2cpA4Vukgf1oYbuXzKCApyQ0FHEYmbCl2kl6Ox6eY0fi6pRoUu0stbu7tvl6vL/SXVqNBFelkTbqQ4L5v5E0qCjiIyICp0kV7W7W7kiumjyNZ0c5Ji9I4V6aHmaBv7m9p0uqKkJBW6SA/rdsdul6tClxSkQhfpYW24ibLiPGaN0XRzknpU6CIx7s663Y0snqHp5iQ1qdBFYqqPnKCxtUPDLZKyVOgiMWt2afxcUpsKXSRm3e4mppcWMV7TzUmKUqGLEJtubk8TV83UZBaSulToIsB7Ncc42RHR/VskpanQRei+3N8MPjlDR+iSulToIsC6cBOXji9hRKGmm5PUpUKXjHfydBcba5p1doukPBW6ZLx39h2NTTen4RZJbSp0yXhrdjWSm63p5iT1qdAl4725q4FFU0eRn6Pp5iS1xVXoZrbUzKrNLGxmD51nv4VmFjGzzyYuosjgOXK8nZ1HWrlGsxNJGui30M0sBDwC3ApUAHebWcU59vs28EqiQ4oMljdjl/tfrUKXNBDPEfoiIOzue9y9A3gaWNbHfl8FngXqE5hPZFCt2dVA6bBc5o4dHnQUkYsWT6FPAGp6LNfG1p1hZhOATwMrzveNzOw+M6sys6qGhoaBZhVJqGjUWRNuYvHMUrKydLtcSX3xFHpf73TvtfzPwDfcPXK+b+Tuj7p7pbtXlpWVxZtRZFDsOHyCxtbTXDNL70VJD9lx7FMLTOqxPBGo67VPJfB0bFKAUuA2M+ty918nJKXIIFgT7v5f4tW6oEjSRDyFvh6YZWbTgIPAXcCf9dzB3ad98LmZPQG8qDKXZPfmrkZmlw9jbEl+0FFEEqLfIRd37wIepPvsle3AM+6+1czuN7P7BzugyGBo74zwzt6jXD1Twy2SPuI5QsfdVwGreq3r8w+g7n7PxccSGVzr9x3ldFeUa2ZruEXSh64UlYy0ZlcjuaEsrpimy/0lfajQJSO9sauRBVNGUpgb139SRVKCCl0yTsOJ02w/dFxXh0raUaFLxnl9Z/fpitfN1h9EJb2o0CXjrK6uZ0xxHvPG63J/SS8qdMkoXZEob+xsYMmcMcQuhBNJGyp0ySgb9jdzor2LJZdouEXSjwpdMsrq6gZyQqb5QyUtqdAlo6zeUc/CqaMozs8JOopIwqnQJWMcPHaK6iMnWDJnTNBRRAaFCl0yxh+qu+deWXKJCl3SkwpdMsbvt9czaVQBM8qKgo4iMihU6JIRTp7u4s1wIzfOLdfpipK2VOiSEd7Y2UBHV5Rb5o0NOorIoFGhS0Z4ddsRRhbmUDllZNBRRAaNCl3SXmckyu+2H+GGueVkh/SWl/Sld7ekvT/uOcrx9i5urigPOorIoFKhS9p7ddthCnJCXKu7K0qaU6FLWotGnVe3HuHa2aXk54SCjiMyqFToktY2H2zh8PF2bqrQ2S2S/lToktZefK+OnJBxk8bPJQOo0CVtRaPOi5sPcd3sMZQU6GZckv5U6JK2qvY3c/h4O7dfNi7oKCJDQoUuaeuF9+rIz8nixrkabpHMoEKXtNQVibJqyyFumFtOUV520HFEhoQKXdLSW3uaaDrZwe3zxwcdRWTIqNAlLb3wXh3D8rK5fo4uJpLMoUKXtHOqI8KqLYe5Zd5YXUwkGUWFLmnnN+8fovV0F3dWTgw6isiQUqFL2nmmqoapowtZNG1U0FFEhlRchW5mS82s2szCZvZQH9s/Z2abYx/rzOyyxEcV6d/+ppO8vecod1RO0sxEknH6LXQzCwGPALcCFcDdZlbRa7e9wHXuPh94GHg00UFF4vEfG2rJMvjM5ROCjiIy5OI5Ql8EhN19j7t3AE8Dy3ru4O7r3L05tvg2oMFLGXKRqPPshlqunV3GuJKCoOOIDLl4Cn0CUNNjuTa27lzuBX7T1wYzu8/MqsysqqGhIf6UInF4c1cDdS3t3Fk5KegoIoGIp9D7Goj0Pnc0W0J3oX+jr+3u/qi7V7p7ZVmZzg+WxPrZun2UDsvjhrljgo4iEoh4Cr0W6HnIMxGo672Tmc0HfgIsc/emxMQTic+ehlZWVzfw51dOJi9b555LZoqn0NcDs8xsmpnlAncBK3vuYGaTgeeAz7v7zsTHFDm/J9/aT07I+NwVU4KOIhKYfu9a5O5dZvYg8AoQAh5z961mdn9s+wrgm8Bo4EexU8W63L1y8GKLfOh4eyf/XlXD7fPHU1acF3QckcDEdRs6d18FrOq1bkWPz/8C+IvERhOJz39U1XKyI8LyxdOCjiISKF0pKimtKxLliXX7WDBlJB+bWBJ0HJFAqdAlpT2/qY4DR9u479rpQUcRCZwKXVJWJOr8cHWYueOGc7MmgRZRoUvqenFzHXsbT/K1G2bqvi0iqNAlRUWizg9+H2ZOeTE3V4wNOo5IUlChS0p6acshwvWtfPWGmWRl6ehcBFTokoLaOyN8+zc7mDtuOLdeOi7oOCJJQ4UuKeena/Zy8Ngp/u4/zSWko3ORM1ToklLqT7Tzo9Vhbq4o56oZpUHHEUkqKnRJKd9/ZScdkSh/c9vcoKOIJB0VuqSMP+5p4pdVNSxfPI2ppUVBxxFJOip0SQmnOiL89bObmTyqkK/fOCvoOCJJKa6bc4kE7buvVLO/qY2n/vJKCnP1thXpi47QJemt33eUx9ft5fNXTuGTM0YHHUckaanQJak1tZ7mr57ayKSRhTx06yVBxxFJavq/qyStSNT5q6c3cvRkB89+5SqK8vR2FTkf/YZI0vqn31azNtzEdz47n0sn6F7nIv3RkIskpX+vquGR1bu5a+Ek7qyc1P8XiIgKXZLPb7cd4aHntnDNrFL+ftmlQccRSRkqdEkqb+1u4oFfvMulE0pY8ecLyM3WW1QkXvptkaTx2rYj3PP4O0weVcjj9yzUH0FFBkiFLknh2Q21fPnfNjBnbDG/vO9KRhXlBh1JJOXoEEgC1RWJ8r1Xd7Li9d0snjma//P5SobpyFzkgug3RwJTf6Kdr/5iI3/ce5TPXTGZb95eQV52KOhYIilLhS5Dzt159t2D/M+XttHeGeGf7ryMz1w+MehYIilPhS5DateRE/z9i9t4c1cjC6aM5Nt/+jFmjikOOpZIWlChy5CobW7jn1/bxXPv1lKUm83Dy+bxuSumaIJnkQRSocug2lRzjJ+u2cuqLYcIZRn3Xj2Nr1w/U2exiAwCFbok3JHj7azcVMez79ay4/AJivOy+dLiqSxfPI3xIwqCjieStlToctEiUWdTzTFW76hndXU9W+uOA/DxSSN4eNk8Pn35RJ2KKDIE9FsmAxKJOrXNbWw/dIJNNcfYVNPMltoWTnZECGUZCyaP5K+XzuGWeWOZUTYs6LgiGSWuQjezpcC/ACHgJ+7+j722W2z7bUAbcI+7v5vgrDIE3J3mtk4Ot7Rz5Hj3R01zG3saTrK7oZV9TW10dEUByAkZFeOG89kFE1k4bRTXzCyjpDAn4Gcgkrn6LXQzCwGPADcBtcB6M1vp7tt67HYrMCv2cQXw49ijJIi70xV1ItHYY8TpikY/XD7zGKW9M0p7Z4T2ziinOiO0d0bOPLZ3RjjVEeV4eyctpz78OB57bGrtoCMS/cjPDmUZU0YVMr2siCVzxjC9rIhZ5cVUjBtOfo4uBBJJFvEcoS8Cwu6+B8DMngaWAT0LfRnwpLs78LaZjTCzce5+KNGBX9/ZwMMvbqP7R4F/sMH56DKctY+f2cc/utzziwbydb2291x79j7n+J7neR6RM0UdJdor48UqzA1RUpBDSUEOwwtymDSqkEsLchhdlEv58HzGluRTPjyf8uF5lA/PJyek2/6IJLt4Cn0CUNNjuZazj7772mcC8JFCN7P7gPsAJk+ePNCsAAzLy2ZOeexCFPvIA90jPx8ud6/rZ58z2+0c+59je69vYD1+aL9f2yvL2Vm7P8kOGaEsIzur52PWh8s9tmfZB8tZ5GVnUZATIj8nFHvMIv+D5dwQ+dlZZKugRdJOPIXe15UfvY8X49kHd38UeBSgsrLygo45F0wZyYIpIy/kS0VE0lo8h2m1QM85wCYCdRewj4iIDKJ4Cn09MMvMpplZLnAXsLLXPiuBL1i3K4GWwRg/FxGRc+t3yMXdu8zsQeAVuk9bfMzdt5rZ/bHtK4BVdJ+yGKb7tMXlgxdZRET6Etd56O6+iu7S7rluRY/PHXggsdFERGQgdKqDiEiaUKGLiKQJFbqISJpQoYuIpAnz3te9D9UPNmsA9gfyw+NXCjQGHSJAmfz8M/m5g55/Mj//Ke5e1teGwAo9FZhZlbtXBp0jKJn8/DP5uYOef6o+fw25iIikCRW6iEiaUKGf36NBBwhYJj//TH7uoOefks9fY+giImlCR+giImlChS4ikiZU6L2Y2R1mttXMomZW2WvbfzezsJlVm9ktQWUcKmb2P8zsoJltin3cFnSmoWBmS2OvcdjMHgo6z1Azs31mtiX2mlcFnWcwmdljZlZvZu/3WDfKzH5rZrtijykzo44K/WzvA58B3ui50swq6L4X/DxgKfCj2ATa6e5/u/vHYx+r+t89tfWYFP1WoAK4O/baZ5olsdc85c7FHqAn6P597ukh4HfuPgv4XWw5JajQe3H37e5e3cemZcDT7n7a3ffSfe/3RUObTobAmUnR3b0D+GBSdElD7v4GcLTX6mXAz2Kf/wz4z0Ma6iKo0ON3romw092DZrY59l/TlPmv50XI1Ne5JwdeNbMNsYndM035BzOuxR7HBJwnbnFNcJFuzOw1YGwfm/7W3Z8/15f1sS7lz/k8378F8GPgYbqf58PA94EvDV26QKTl6zxAi929zszGAL81sx2xI1lJchlZ6O5+4wV8WVpOhB3vv4WZ/V/gxUGOkwzS8nUeCHeviz3Wm9mv6B6GyqRCP2Jm49z9kJmNA+qDDhQvDbnEbyVwl5nlmdk0YBbwTsCZBlXszfyBT9P9B+N0F8+k6GnLzIrMrPiDz4GbyYzXvaeVwBdjn38RONf/2pNORh6hn4+ZfRr4AVAGvGRmm9z9ltjE2M8A24Au4AF3jwSZdQh8x8w+TveQwz7gy8HGGXznmhQ94FhDqRz4lZlBdz/8wt1fDjbS4DGzp4DrgVIzqwW+Bfwj8IyZ3QscAO4ILuHA6NJ/EZE0oSEXEZE0oUIXEUkTKnQRkTShQhcRSRMqdBGRNKFCFxFJEyp0EZE08f8Bpzm3iv05OM8AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pylab as plt\n", "import numpy as np\n", "import sklearn.datasets\n", "x = np.arange(-12, 12, 0.1)\n", "f = 1 / (1 + np.exp(-x))\n", "plt.plot(x, f)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# The Math: Defining a Likelihood Function\n", "\n", "The logistic regression problem is based on a statistical viewpoint of machine learning. The basis is the assumption for a binary classification to put things into a category $0$ and $1$. We want to understand the binary outcome given explanatory variables. This can be achieved using a Bernoulli distribution\n", "$$\n", "0, 0, 1, 0, 1,0,1,1,0\n", "$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "and a basic statistics course tells us that the probability density function of a Bernoulli distributed variable is\n", "$$\n", "P(Y=y\\mid X=x,p)=p^y(1-p)^{1-y}\n", "$$\n", "with $y\\in\\left\\lbrace0,1\\right\\rbrace$ the response, where $1$ has probability $p$ and $0$ probability $1-p$, and an observed input $x$. We now want to find the $p$ that generates the above numbers with $p^{y_1}(1-p)^{1-y_1}$ for the first number and $p^{y_2}(1-p)^{1-y_2}$ for the second and so on. " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "If we assume independence of the variable the joint probability distribution is the product of such terms. It can be shown that an important function in this context is the sigmoid function\n", "$$\n", "\\sigma(z)=\\frac{1}{1+e^{-z}}\n", "$$\n", "used to model the probability density. Assuming a linear relationship\n", "$$\n", "y=\\mathbf{X}^{T}\\beta+\\varepsilon\n", "$$\n", "with $\\mathbf{X}=[x_1,x_2,\\ldots,x_n]$ with $x_i\\in\\mathbb{R}^d$. " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "The components in the error vector $\\varepsilon$ \n", "are computed from\n", "$$\n", "\\varepsilon_i=\n", "\\left\\lbrace\n", "\\begin{array}{cc}\n", "1-p_i&y_i=1\\textrm{ with probability } p_i\\\\\n", "-p_i&y_i=0\\textrm{ with probability } 1-p_i\\\\\n", "\\end{array}\n", "\\right.\n", "$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "We can write this also as\n", "$$\n", "\\sigma(\\beta^Tx)=\\frac{1}{1+e^{-\\beta^Tx}}\n", "$$\n", "where this represents the unknown probability $p$. Using the sigmoid or logistic function as probability $p$ we obtain now\n", "$$\n", "P(Y=y\\mid X=x,p)=\\sigma(\\beta^Tx)^y(1-\\sigma(\\beta^Tx))^{1-y}.\n", "$$ \n", "We now consider the likelihood function\n", "$$\n", "L(\\beta)=\n", "\\prod_{i=1}^{n}P(Y=y_i\\mid X=x_i)\n", "$$\n", "for independent variables or equivalently\n", "$$\n", "L(\\beta)=\n", "\\prod_{i=1}^{n}\\sigma(\\beta^Tx_i)^{y_i}(1-\\sigma(\\beta^Tx_i))^{1-y_{i}}.\n", "$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "We then obtain the log-likelihood function as\n", "$$\n", "L_{\\log}(\\beta)=\n", "\\sum_{i=1}^{n}{y_i}\\log(\\sigma(\\beta^Tx_i))+\\left(1-y_{i}\\right)\\log(1-\\sigma(\\beta^Tx_i)).\n", "$$\n", "and we are interested in maximization of this quantity. We consider the gradient of $L_{\\log}(\\beta)$ via" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "\\begin{align}\n", "\\frac{\\partial L_{\\log}(\\beta)}{\\partial \\beta_j}&=\\frac{\\partial }{\\partial \\beta_j}\\left(\t\\sum_{i=1}^{n}{y_i}\\log(\\sigma(\\beta^Tx_i))+\\left(1-y_{i}\\right)\\log(1-\\sigma(\\beta^Tx_i))\\right)\\\\\n", "&=\\sum_{i=1}^{n}{y_i}\\frac{\\partial }{\\partial \\beta_j}\\log(\\sigma(\\beta^Tx_i))+\\left(1-y_{i}\\right)\\frac{\\partial }{\\partial \\beta_j}\\log(1-\\sigma(\\beta^Tx_i))\\\\\n", "&=\\sum_{i=1}^{n}{y_i}\\frac{x_i^{(j)}\\sigma'(\\beta^Tx_i)}{\\sigma(\\beta^Tx_i)}+\\left(1-y_{i}\\right)\\frac{(-x_i^{(j)}\\sigma'(\\beta^Tx_i)) }{(1-\\sigma(\\beta^Tx_i))}\\\\\n", "&=\\sum_{i=1}^{n}{y_i}\\frac{x_i^{(j)}\\sigma(\\beta^Tx_i)(1-\\sigma(\\beta^Tx_i))}{\\sigma(\\beta^Tx_i)}+\\left(1-y_{i}\\right)\\frac{(-x_i^{(j)}\\sigma(\\beta^Tx_i)(1-\\sigma(\\beta^Tx_i)) }{(1-\\sigma(\\beta^Tx_i))}\\\\\n", "&=\\sum_{i=1}^{n}{y_i}x_i^{(j)}(1-\\sigma(\\beta^Tx_i))+\\left(1-y_{i}\\right)(-x_i^{(j)}\\sigma(\\beta^Tx_i)\\\\\n", "&=\\sum_{i=1}^{n}{y_i}x_i^{(j)}(1-\\sigma(\\beta^Tx_i))+y_{i}x_i^{(j)}\\sigma(\\beta^Tx_i)-x_i^{(j)}\\sigma(\\beta^Tx_i)\\\\\n", "&=\\sum_{i=1}^{n}{y_i}x_i^{(j)}-x_i^{(j)}\\sigma(\\beta^Tx_i)\\\\\n", "&=\\sum_{i=1}^{n}x_i^{(j)}\\left({y_i}-\\sigma(\\beta^Tx_i)\\right)\n", "\\end{align}" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "and based on this the gradient is given as\n", "$$\n", "\\nabla L_{\\log}(\\beta)=\\mathbf{X}\\left({y}-\\mathbf{p}\\right)\n", "$$\n", "with \n", "$\\mathbf{p}_i=\\sigma(\\beta^Tx_i)$." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "Now, we also need the Hessian matrix, which can be computed entry-wise via the following formula\n", "\n", "\\begin{align}\n", "\\frac{\\partial^2 L_{\\log}(\\beta)}{\\partial \\beta_k\\partial \\beta_j}\n", "&=\\frac{\\partial}{\\partial \\beta_k}\\sum_{i=1}^{n}{y_i}x_i^{(j)}-x_i^{(j)}\\sigma(\\beta^Tx_i)\\\\\n", "&=\\sum_{i=1}^{n}-x_i^{(j)}x_i^{(k)}\\sigma(\\beta^Tx_i)(1-\\sigma(\\beta^Tx_i))\\\\\n", "&=\\sum_{i=1}^{n}-x_i^{(j)}x_i^{(k)}\\mathbf{p}_i(1-\\mathbf{p}_i):=(-\\mathbf{X} D\\mathbf{X}^T)_{kj}\\\\\n", "\\end{align}\n", "with\n", "$$\n", "D:=\n", "\\begin{bmatrix}\n", "\\sigma(\\beta^Tx_1)(1-\\sigma(\\beta^Tx_1))&&&\\\\\n", "&\\sigma(\\beta^Tx_2)(1-\\sigma(\\beta^Tx_2))&&\\\\\n", "&&\\ddots&\\\\\n", "&&&\\sigma(\\beta^Tx_n)(1-\\sigma(\\beta^Tx_n))\\\\\n", "\\end{bmatrix}.\n", "$$\t" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "This leads us to the following Newton iteration\n", "$$\n", "\\beta=\\bar{\\beta}+(\\mathbf{X}D\\mathbf{X}^{T})^{-1}\\mathbf{X}\\left({y}-\\mathbf{p}\\right)\n", "$$\n", "where $\\bar{\\beta}$ indicates the iterate from the previous Newton iteration. " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "If we decide to include a regularization term in the log-likelihood\n", "$$\n", "L_{\\log}(\\beta)=\n", "\\sum_{i=1}^{n}{y_i}\\log(\\sigma(\\beta^Tx_i))+\\left(1-y_{i}\\right)\\log(1-\\sigma(\\beta^Tx_i))-\\frac{\\lambda}{2}\\Vert{\\beta}\\Vert^2\n", "$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "\\begin{align}\n", "\\frac{\\partial L_{\\log}(\\beta)}{\\partial \\beta_j}&=\\frac{\\partial }{\\partial \\beta_j}\\left(\t\\sum_{i=1}^{n}{y_i}\\log(\\sigma(\\beta^Tx_i))+\\left(1-y_{i}\\right)\\log(1-\\sigma(\\beta^Tx_i))\\right)-\\frac{\\partial }{\\partial \\beta_j} \\frac{\\lambda}{2}\\Vert{\\beta}\\Vert^2\\\\\n", "&=\\sum_{i=1}^{n}x_i^{(j)}\\left({y_i}-\\sigma(\\beta^Tx_i)\\right)-\\frac{\\partial }{\\partial \\beta_j} \\frac{\\lambda}{2}\\sum_{i=1}^{n}\\beta_i^2\\\\\n", "&=\\sum_{i=1}^{n}x_i^{(j)}\\left({y_i}-\\sigma(\\beta^Tx_i)\\right)-\\lambda\\beta_j\n", "\\end{align}\n", "and as a result\n", "$$\n", "\\nabla L_{\\log}(\\beta)=\\mathbf{X}\\left({y}-\\mathbf{p}\\right)-\\lambda \\beta\n", "$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "Now, we also need the Hessian matrix, which can be computed entry-wise via the following formula\n", "\n", "\\begin{align}\n", "\\frac{\\partial^2 L_{\\log}(\\beta)}{\\partial \\beta_k\\partial \\beta_j}\n", "&=\\frac{\\partial}{\\partial \\beta_k}\\left(\\sum_{i=1}^{n}{y_i}x_i^{(j)}-x_i^{(j)}\\sigma(\\beta^Tx_i)-\\lambda \\beta_j\\right)\\\\\n", "&=\\sum_{i=1}^{n}-x_i^{(j)}x_i^{(k)}\\sigma(\\beta^Tx_i)(1-\\sigma(\\beta^Tx_i))-\\lambda \\delta_{kj}\\\\\n", "&=\\sum_{i=1}^{n}-x_i^{(j)}x_i^{(k)}\\mathbf{p}_i(1-\\mathbf{p}_i):=(-\\mathbf{X} D\\mathbf{X}^T)_{kj}-\\lambda \\delta_{kj}\\\\\n", "\\end{align}\n", "with\n", "$$\n", "D:=\n", "\\begin{bmatrix}\n", "\\sigma(\\beta^Tx_1)(1-\\sigma(\\beta^Tx_1))&&&\\\\\n", "&\\sigma(\\beta^Tx_2)(1-\\sigma(\\beta^Tx_2))&&\\\\\n", "&&\\ddots&\\\\\n", "&&&\\sigma(\\beta^Tx_n)(1-\\sigma(\\beta^Tx_n))\\\\\n", "\\end{bmatrix}.\n", "$$\t" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "we obtain the following Newton scheme\n", "$$\n", "\\beta=\\bar{\\beta}+(\\mathbf{X}D\\mathbf{X}^T+\\lambda I)^{-1}\\left(\\mathbf{X}\\left({y}-\\mathbf{p}\\right)-\\lambda \\bar{\\beta}\\right)\n", "$$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "Let's try this on the famous IRIS dataset [see here](https://en.wikipedia.org/wiki/Iris_flower_data_set)\n" ] }, { "cell_type": "code", "execution_count": 91, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(150,)\n" ] } ], "source": [ "iris = sklearn.datasets.load_iris()\n", "x_train = iris.data[:, :2]\n", "y_train = (iris.target != 0) * 1\n", "x_train = x_train.T\n", "y_train = y_train.T\n", "print(y_train.shape)" ] }, { "cell_type": "code", "execution_count": 92, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAdoAAAFiCAYAAABRfRm3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nOzdd1hU19YH4N+ZxjD0oqAYBFQUBUWKogJWsKKxRCzX9hmjUROTmJirURNjSa435d7E2FI0aqLX2A32aBSN2MAolkjAigUpUmeYtr4/0AkjEBlgHMT1Ps88Cds9+6zDDCxmn33WFogIjDHGGDMPkaUDYIwxxuoyTrSMMcaYGXGiZYwxxsyIEy1jjDFmRpxoGWOMMTPiRMsYY4yZkcQcg7q6upKXl5c5hmaMMcZqpTNnzmQSUb3H282SaL28vHD69GlzDM0YY4zVSoIgXC+vnaeOGWOMMTPiRMsYY4yZESdaxhhjzIw40TLGGGNmxImWMcYYMyNOtIwxxpgZcaJljDHGzIgTLWOMMWZGnGgZY4wxM+JEyxhjjJkRJ1rGGGPMjDjRMsYYY2bEiZYxxhgzI060jDHGmBlxomWMMcbMiBMtY4wxZkacaBljjDEz4kTLGGOMmREnWsYYY8yMONEyxhhjZsSJljHGGDMjTrSMMcaYGXGiZYwxxsyIEy1jjDFmRpxoGWOMMTPiRMsYY4yZESdaxhhjzIw40TLGGGNmxImWMcYYMyNOtIwxxpgZcaJljDHGzIgTLWOMMWZGnGgZY4wxM+JEyxhjjJkRJ1rGGGPMjDjRMsYYY2bEiZYxxhgzI060jDHGmBlxomWMMcbMiBMtY4wxZkacaBljjDEz4kTLGGOMmREnWsYq6cSJE4jp0QMvuLqic0gIdu7caemQGGPPAE60jFVCQkIC+nXvjvDTZ/CTWIqRaVcxZcQIrF2zxtKhMcZqOYGIanzQkJAQOn36dI2Py5il9OnSBd3OnsNwGxtD2xl1MaaJRUi9fRsiEf/NytjzThCEM0QU8ng7/3ZgrBLOnD2LLnK5UVuQVIasBw+Qm5troagYY88CTrSMVYJ3o0a4qFEbtd3Q6SCTyWBnZ2ehqBhjzwJOtIxVwvT338f7GjUuPEy2N7VaTC9WYcprr0EikVg4OsZYbca/IRirhJdeegk5WVkY9957KC7Ih14kwtTXXsPcDz+0dGiMsVqOF0MxZgKdTofMzEw4OTlBJpNZOhzGWC1S0WIo/kTLmAnEYjHc3NwsHQZj7BnC12gZY4wxM+JEyxhjjJkRJ1rGGGPMjDjRMsYYY2bEiZYxxhgzI060jDHGmBlxomWMMcbMiBMtY4wxZkacaBljjDEz4kTLGGOMmREnWlanHTp0CKNfGooBPXpg6dKlUCqVlg6JMfac4UTL6qxPPv4YYwcMgO/+/eiVmITNs2ejR6dOUKlUlg6NMfYc4UTL6qTMzEwsmDcPmxQ2GGdji/7WCnwvV0B67Rp++OEHS4fHGHuOcKJlddJvv/2GYFs7NBT/tUGVIAgYSMD+bdssGBlj7HnDiZbVSU5OTrin1eLx/ZbvEsGpXn0LRcUYex5xomV1UqdOnaCxt8MPyiJDsk3VarBaq8H4ya9aODrG2POEEy2rk0QiEXbs24c1To7oVlSAIcUqvFiQj/mffYaQkBBLh8cYe45IntyFsWdT8+bNkZyailOnTiEvLw9hYWGwtbW1dFiMsecMJ1pWpwmCgHbt2lk6DMbYc4ynjtlTlZOTg/v371s6DMYYe2o40bKn4saNG+gVGYnGDRqg6QsvoGNgIM6dO2fpsBhjzOw40TKz02g06BEejsDzyUhydsU5Z1cMvn4D0Z07Iycnx9LhMcaYWXGiZWa3e/duOBUW4jWFDeSCAIkgYJjCBh1FYqxdu9bS4THGmFlxomVmd/36dfhBKNPeQqPB9bQ0C0TEGGNPDydaZnbBwcE4otVAW6pKExHhsESC0LAwC0bGGGPmx4mWmV2HDh3gFxyMV5SFSFQX44JGjbeVhVC61cegQYMsHR5jjJkVJ1pmdoIgYPOuXYh4803MtLXBVJkUnuPG4eDx45DJZJYOjzHGzEp4vOh6TQgJCaHTp0/X+LiMMcZYbSUIwhkiKlPjlT/RMsYYY2bEJRgZq2HZ2dnYsGED7ty5g8jISHTv3h0iEf9Ny9jzihMtYzUoISEBMdE9ESGVwFOjxRtffonGbdti2969fD2asecU/5nNWA3R6/UY/dJL+EgixZdyBd6xs8cuhS0Kz57FihUrLB0eY8xCONEyVkMuXrwITV4eesnlhjaJIGC8SIxN339vwcgYY5bEiZaxGiISiaArZxW/7uG/McaeT/zTz1gN8fPzg52rK7YplYa2YiKs1OsQ+3//Z8HIGGOWxIuhGKshgiBg3ebN6NO9O34uVsJTq8MB0iOkc2e8/PLLlg6PMWYhnGgZq0Ft27bFnzdvYvPmzbh79y7+LzIS7du3hyCU3VSBMfZ84ETLWA2zsbHB6NGjLR0GY6yW4Gu0jDHGmBlxomXPpPfeew9uCgWsRSI0sLXFxx9/bOmQGGOsXDx1zJ45b06bhjVLluDfDk4IsndEgroY7856D1qtFrNnz7Z0eIwxZoR372HPHEepDCsdHNDJ6q/CELuVSsxSFSGjqMiCkTHGnme8ew+rE4qKipCn1aCjzMqoPcLKCg9K3b/KGGO1BSda9kyRy+VQiES4oNUYtf+uUcNWKrVQVIwxVjFOtOyZIhKJ0CsmBpNzsnFZU5Jsz6nVePNBDl4aO9aywTHGWDl4MRR75mzcsgWD+/dH/927odXrIRWJMCg2FitWrrR0aIwxVgYnWvbMEYlE2Przz9Bqtbh79y7c3d0hkfBbmTFWO/FvJ/bMkkgkaNSokaXDYIyxv8XXaFm16XQ6zJkzBwMGDMBPP/1k6XDM6s8//8TmzZuRlJRk6VAYs5gbN25gy5YtSEhIQHVvES0uLsbevXuxY8cO5OXl1VCEtQvfR8uq5ddff0X/7t2hAOArkeK0Rg1He3tcy8iATCazdHg1RqPRYPzIkdgdF4dQGxtcUKng07IltuzeDScnJ0uHx9hTQUR47bXXsHbtWvj6+uLu3btwdXVFXFwcGjZsaPJ48fHxGDx4MFxdXSGTyXD16lUsXboUI0eONEP05lfRfbScaFm1uMisMEQmw2x7B4gEAfl6PYZk3od1m9Y4VYfeAx8vWoTdi/+Nb61Lyj7qiDBbWQTq0R1r6/ineMYeWbVqFRYuXIipU6dCoVCAiPDzzz+jqKgIv/zyi0ljFRYWwtPTE6NGjULLli0BAOnp6fjiiy9w6tQpNG3a1BynYFZcsILVuJMnT6JQo8Y7dvYQPdwGzk4kwix7B/xZx6ZWVy9fjnekMliLSn5kxIKAd+XW2LpjB1QqlYWjY+zpWLlyJXr16gWFQgGgZA/mXr164dSpU7hz545JY+3cuRNeXl6GJAsAHh4eaNeuHdasWVOjcVsaJ1pWZenp6ZAJAuSP7bXqLBJBb4aZEkvKLyyEs8j4x8X24XkXFxdbIiTGnrr8/HzY2toatUkkElhbW6OgoMCksfLy8mBjY1OmXaFQ1LlrtZxoWZX1798fBOBgsfEnuh+KCmHn4mKZoMykZ69e+N9j57lTqURA8+ZwcHCwUFSMPV19+vRBQkKCUdvly5ehUCjQpEkTk8aKiorCuXPnkJ+fb2jTaDRISkpC3759ayTe2oJv72FVJhaLMfKVVzBp5UqMVtjATyrFLqUSR9XF2L17l6XDq1HzPv4Y4QcO4LayCOFEuCAI2KrT4Oevv7Z0aIw9Ne+++y46dOiAb775Bq1atUJGRgaOHz+O9evXQyQy7XObt7c3Xn/9dXz66afo1KkTZDIZTpw4gXbt2qFHjx5mOgPL4MVQrNrWrVuH96ZPhyovD25NmuB/P/0EPz8/S4dV47Kzs/H1ypU4e/w4fFq0wMQpU+Dp6WnpsBh7qvLz87Fq1SocPXoUL7zwAiZOnAhfX98qj3fo0CGsXbsWKpUKgwcPxsCBA01O2rUFrzpmjDHGzIhXHTPGGGMWwImWVdvp06cxbvhwdG/fHrNmzMDdu3erPFZubi4WLViAqA4dMGzAABw8eLDcfpcvX8bk8ePRvX17THv1VaSlpVX5mIwxZk48dcyqZevWrZg0ahQmSmTwFYvxC+mxXyzCb2fOmHz9Mjc3Fx2DgtA8JwcvCiLc1uuxVKvFjAXzMXXaNEO/48ePIyYqGmOkUgSJxDih12G9TosD8fFo06ZNTZ8iY4xVCl+jZTVOr9ejScOG+ESrR5iVlaH9o8IC6IcMxtJvvjFpvI8XLcKJTz7BEuu/7q27ptWiX0Eebty9Czs7OwBARFAQhl69jkEPb5oHgNWFBTgeHISfTaxOwxhjNYWv0bIad/PmTRQXFKD9YzWNB8hkOLRvn8njHYqLwwDB+C3pJZGgmbU1EhMTAQBarRbHf/8dMdbWRv1etFbg8LFjJh+TMcbMjRMtqzJ7e3sUarUofGxWJF2ng0sVCla4ursjXaczatMR4XZxMVxdXQGU3Ltrb22Nu4/1S9dp4cKFIxhjtRAnWlZlTk5O6NOrF+arlCh+mGzv6XT4t06LSdOnmzzeK6+/jmU6LdK0GgAlSfYLZRE8mzZFq1atAJTUVn15wgS8ry5GoV4PAMjV6zFfo8YrU6fW0JkxxljN4UTLqmX5998jNzgIYbk5GFisQre8BxgyeXKVtrnq3LkzZv/rXxhQWID+KiU65OfiuLcXfvr5Z6N+8z/+GO7RUQjLzcEgtQodH2QjYNAgzJg5s4bOijHGag4vhmI14urVq7h16xb8/f2rvT9rYWEhkpKS4OLi8rcVptLT05GamormzZvDzc2tWsdkjLHq4lXHjDHGmBnxquNajIhw9epV3L9//6ke8+bNm0/cQ1Kj0SA1NRW5ublPKTLG6r779+/j2rVrMMcHHVb7cKK1sAMHDqCFZ2N0DAhA0xc80bdbd9y7d8+sx0xMTESwnx+CWrRAqyZNEBkSgj///LNMv69XrMAL9euja9u2aOzujldGj+ZNzhmrhjt37iA6OhpNmjRBSEgImjVrhkOHDlk6LGZmPHVsQSkpKejQti3+Y2WNLlZWUAH4r7IIpzxfwPGzZyE8tqF6TcjOzkYLb2/MEUnworU1dAC+VyqxWm6FP65fg+zhPbFxcXGYPGwYvrFWoKVUhhy9Hu+qlGjUPwYrvv++xuNirK4jIrRt2xaNGjVC7969IZFIkJycjHXr1uHMmTPw8fGxdIismnjquBZa+dVXGC6VoatcDkEQYC0ImGGtwP1r12CuP1TWrVuHSIkEgxQKiAQBUkHAywoFPDRq/Fxqde9/Fi3CP6UytJSWJF4nkQj/klvjfxs38jQyY1Vw/PhxZGdnIyYmBlKpFIIgICAgAGFhYVixYoWlw2NmxInWgm6lpaHZY20iQUBTmQy3bt0yzzGvX0czjbZMe1O93uiY6TdvoZlEatTHSSSCvVSKzMxMs8TGWF1269YtNGzYsMxMVf369XHjxg0LRcWeBk60FtSuSxcceOyHLk+vx6nCQgQHB5vlmO07dsQBsdhoEYaaCIe1WrRr1+6vfp06Yr9GbfTcixo1NBIJb3bOWBWEhITg8uXLZdY5XLx4EZ06dbJQVOxp4Gu0FpSbm4vQgAB0yM9HrESGbL0en+s06Dg0Fl+uNM9UkkajQWRoKOpfv4HxEgmKibBUp4VD+/bYunu34a/tK1euIDw0FONEYvSQSpGi1eJfGjVmf/IJJkycaJbYGKvrJkyYgCNHjqBXr16wsbHBb7/9htu3byMxMRG2traWDo9VE1+jrYUcHBxw7MwZ1Bs9GjNsFfjKowEmLl6M/y5fZrZjSqVS7IuPR9upUzDXwR4fu7qg76xZ2Lhjh9GUlq+vL+JPnkRGr56YJpdjq18LfLV+PSdZxqphxYoVePPNN3Hs2DFs2bIFISEh+O233zjJ1nH8iZYxxhirAfyJljHGGLMATrSsQkVFRXjrrbfQPigIAwcOxB9//FFuv5SUFAwaNAjt2rbFW2+9haKionL7Xbx4EXNnz8bMd9/FyZMnzRm6wbVr17Bg/nzMmD4dv/zyC1fiYYw9dZxoWblu374NT1dX7F+yBN3/TIV+71609fPDunXrjPr98MMPaN28ObR79qBHahoOLFmCF1xcytye9OV//oMuoaHIXPIVlMtXYHD37nhn2jSznsNPP/2E4FatcPWzz0Fff4NXBw3CyMGDoX+4vR5jjD0NfI2WlatreDgUZxKxwsnZsEhqa1ER3i8qxH1lEUQiEfR6PeopFHjfWoHBChsAJdVvJuVkIz+wDQ4fPw6g5P7BAF9f7LFzQCOJBADwQK9Hr4I8bP7lF6PbimpKYWEhXnBzw3qFDfwfFt1QEWFQUSHmfvsNBg8eXOPHZIw93/gaLTPJuZMnMdHW1mgl8gBra6g1apw4cQJASc3kYrUaL1orDH0EQcBEW1ucP3PG0BYXF4ce1gpDkgUAR5EIg0VibN20ySzxHz58GC3l1oYkCwByQcAIAFt+/NEsx2SMsfJwomXlEolEKH5stkMHQA/A2toaACCXy6F72F6amkoqXD0ilUqhLqdss1oQIJXJyv5DDZBKpSimslPEasBQz5kxxp6GJyZaQRCsBEEYIQjCLEEQ5j56PI3gmOV07NEDn+XnQ10q2a4qKICdQoHAwEAAgL+/PxxtbPBtQb6hj5oIn+bnoX23boa2/v3740hRES6UqjSVrtVik0aNYSNGmCX+zp0746ZejyOlqvDk6PVYrddj+LhxZjkmY4yVR/LkLtgOIBfAGQDF5g2H1RY/bNiAAB8ftL93Fz3kclzUaPCnXocd+/YZ9Vu/cydievTATpUK/lIpDqhUkDs7Yd9PPxn6uLq64us1axA7Zgw6W1tDRsB+ZRHmLVyIli1bmiV+mUyGjdu3Y1C/fggRABci7FMpMX7SJERFRZnlmIwxVp4nLoYSBCGZiPxNGZQXQ9UNer0e33//PeLi4tCkSRPMmTOn3Ao2RUVFmD9/PlJSUtC7d2+MGzcOIlHZyZKsrCxs374darUa/fr1Q6NGjcx+Dvn5+di+fTtyc3MRHR2NZs0e38aBMcZqRkWLoSqTaFcC+JKIzlf2YJxoGWOMPW8qSrQVTh0LgnAeAD3sM04QhDSUTB0LAIiIWpsrWMYYY6yu+LvFUP0AxADoDaApgOiHXz9qZzUgJSUFL/XrBweFAo1cXTH7n/9EcXHVLoWrVCr06t4ddhIJZIIAT2cXbKrG7TNLly5FAzs7yAQBjlIphg0dWqbYAxFh2dKl8GvcGLZyObq2b4/4+PgyY+Xl5eH1SZNQ38EBzra2GDd8OO7cuVPl2N6YNg2uVlaQCQLqWVtjzpw5VR6rNjt58iS6desGGxsbeHt74/PPP69ywY3CwkK89dZbqF+/Puzt7TF8+PBy90G9fPky+vfvDzs7OzRs2BBz586FWq0uZ0TGWGVUZup4LRGNelJbaTx1XDkZGRkI9PPDGJ0ew+XWyNLr8S9NMezDw7Fx506Txwts3hySa9ewyMERDcUSbCsqxKL8POw6eBBdunQxaaw1a9bg1bFjscjBEdFya6RoNZj+IAc+4eHYe+iQod9H8+fjx08+xQKZDL6SksVQH2qKsevQIYSGhgIoScYRISHwSLuKaTIrWAkCvlMpsd/eDmcvX4ZCoagojHK9OnEiNn/zLf7j6IQgmQwJ6mK8+SAHU2fNwvz5800aqzZLTk5GREQEYmJiEBQUhHv37mHz5s0YNGgQPv74Y5PGIiL06NEDRUVF6NOnD+RyOY4cOYKzZ8/iwoULsLe3B1BSESwwMBCdO3dGhw4dkJubi+3bt6Nly5b4ke8/ZuxvVecabSIRBZX6WgzgPBFVuFyUE23lLPjwQ1z6z3+xuFTBh2IidMjNwZHERDRv3rzSY507dw7tAwNxxq0BHEotRFqQ+wCnfbxx6tw5k2Jr4uaGsapi/J+tnaHthlaLbvfvIT0zE87OzlCpVPCoVw87FbZoXKoYxfeFBUjs2AGbdu0CABw6dAiTBw7EPoWt0f21o1VKjPr0E4wdO9ak2BykUqxycEJ7KytD236VEu8UFeK+UmnSWLXZqFGjUFRUhJ49exracnNzMX/+fNy6dcuQHCvj5MmTePHFFzF37lyjhWrffvstRo0ahalTpwIAZs+ejYSEBAwdOtTQR61WY/bs2UhKSoK3t3cNnBljdZPJlaEEQZgpCEI+gNaCIOQ9fOQDyEDJLT+smpJPnULYY3/nWAkCghQ2uHDhgkljHT58GD4SiVGSBYBwKznuXL1qcmzZ2dnoaCU3avOUSGArCEhKSgIA3LlzBwpBMEqyANBBZoXzpRJ7cnIy2oslRkkWAMI0Wpx/OFZlFRUVIV+rRbvHik50kFnhQal7ZuuCc+fOwdfX16jNwcEBLi4uuGria5qcnIymTZuWWQ3u4+OD33//3fD12bNn0aRJE6M+MpkMTZo0Mfk9yRgrUWGiJaKPiMgOwL+JyP7hw46IXIho5lOMsc7yDQhA0mMVk7REOKdSlvkF+yQdOnTAVa0WRY9dvzutLkY9Dw+TY3NwcEDiY9fl7uh0yCdCQEAAAMDNzQ0FOh1u67RG/RI1avi2aGH42tfXF0l6fZmdc85KJWjeqpVJcSkUCtiIxfhdozE+ploN+zpW8alFixZlEmpBQQGysrLg6elp0li+vr64du1amdfg1q1b8PPzM3zt5+eH69evG/XRarW4du2aye9JxliJv/tEGyQIQhCAnx79f+nHU4yxzpo4ZQp+1mmxrqgQxUS4q9NhurIIgaGh8Pc36dZlhISEwNPDA6/kZOOmVgsNEbYWFWFFYQEWfvaZybHN+PBDzM97gAMqJfRESNFoMD47E6EhIahfvz6AkqQ3efJkTFWpkKLRQE+EX1RKLNao8c7cv4qH9ejRA2J3N8xTFiFHr0OhXo8lRYU4LxFj+PDhJsfWf8gQvJqThbNqNYgIJ4uL8caDHPxj4kSTx6rNZsyYgX379iExMRF6vR737t3DqlWrMGrUKDg5OZk0VqdOndCgQQNs2LABeXl5UKvV+OWXX5CSkmI0dT916lScOnUKx44dg1arRU5ODtasWYOOHTtyomWsqoio3AeAQw8fxwFoAJxGSXUoDYCjFT2PiBAcHEyscs6ePUs9OnQksUhEtlZWNHn8y5Sfn1+lsXJzcyksMJBkgkACQPUVClq+fHmVY5s3bx45y2QkAGQtCBTdpQsVFxcb9dHpdLRw3jxyc3QkkSBQm6bNKC4ursxYGRkZNOqll0gulZJULKYBUdGUmppapbh0Oh2NiI0lW7GYBIDsxRKaOGFClcaq7Q4cOECBgYEkEonIycmJ3nvvPdJoNFUaKzs7m8aMGUNyuZzEYjFFR0fT5cuXy/Q7ffo0RUREkFgsJltbW5o8eTIVFhZW91QYq/MAnKZycmJlFkNtALCQHhasEATBH8DbRDS2oufwYijTabVaiMVio91yqkqv10OtVkMulz+5cyUUFRVBLpeXW+3pESKCVquFVCp9YmxEBLFYXGOxmbpq+VlU0++PyrwGNXlMxp4HJhesKKUFlaoKRUTJgiAE1mh0DBJJZV6KyhGJRDWWZAFUKpEJgvDEJAvgb5N1VTwPSRao+ffH0z4mY8+zyvwkXRIE4RsA61BSKeofAC6ZNSpWa1y5cgWnTp2Ch4cHIiMjy/0lrdfrER8fj5s3byI0NNSk25JY7ZSYmIhvv/0W9evXx/Tp08utcW0pu3btws6dO9GyZUtMmTKlxv94Y6ymVWbqWA7gVQCRD5uOAFhGRBXeS8FTx88+nU6HV0aPxs5t29BJYYMrWg0k9eph16FD8Ci1ivn27dvo07UrNBkZ8JVIcayoEH3798c369bV2PQwe7qioqJw+PBhNG/eHLm5ucjIyMD333+P2NhYi8alUqng7++P27dvw9fXF+np6SguLsbRo0fRujVXhGWWV+Wp44cJ9fOHD/acWL58OS7E7cIxe0coRCKQVIbP72di/PDh2HPkiKHf+OHD0TXjPqZb20AQBCglUozavQdLly7Fa6+9ZsEzYFXx6aef4uTJk1i4cCEcHBwAAKdOncLYsWMxcOBAyCx4C1VsbCyICIsXL4ZMJgMRYcuWLejVqxdu375tsbgYe5K/u71n48P/nhcE4dzjj6cXIrOEtcuWYZpYAsXDaTlBEDDZWoGEU6eQkZEBALh//z5+O3ESU60VhgUz1iIR3pBIsGbZMovFzqpuxYoV6N27tyHJAkBoaCjs7e2xcuVKC0YG/Prrr0bJXhAExMTE4P79+yYX8GDsafq7T7TTHv6339MIhNUuhYWFsH/s2pcVALlYjKKiIgAlU3lWYhGsHnuuvSAy9GHPluLi4nIXmCkUCjx48MACEf1Fp9PB2traqE0ikUAsFiM3N9dCUTH2ZH9XGerR1irdAciI6Hrpx9MJj1lK38GD8YNWbVRJaH+xCi716qFx48YAgEaNGqG+mxv2Plb68AetBn0HDXqq8bKaER0djUOHDkGn0xnabt++jfT0dLz88ssWjKykutXhw4eN2pKSkiCTyfgaLavVKrMY6kMA4QAao6RgRTyAeCI6W9FzeDHUsy87OxuR7dqhYXYOuuv1SBGJsE2jxpa4OERGRhr6HT16FC/27o0XpTI00+txUCTCLSdHHDl5Ei4uLhY8A1YVRUVF8Pb2hkwmQ3h4OHJzc3Ho0CGMGjXK4lPHly5dQkhICLy9vREYGIibN28iISEBX3zxBSbWsapg7NlU5d17Sg1gDWACgLcBeBBRhUtKOdHWDYWFhVi3bh1OHj6MRj4+GP/KK+XW2L158ya+XbkSN1NTERoZiVGjRsHGxsYCEbOaoFKpMGvWLOzevRu2traYNWsWBg4caOmwAAB3797Fm2++icTERDRs2BCLFy82bMfImKVVZ5u82QA6AbAFkATgKEo+0Va4azcnWsYYY8+b6lSGGgRACyAOwGEACX93Dy1jjDHG/vLEkipUsul7dwAnAUQBOC8IwlFzB1abnTlzBqOGvISOrVtj8viX8eeff1Z5rLy8PCycPx8Rbduib5cu2LRpU5OgzFgAACAASURBVJmtzABg9+7daN2sGdwVCvg1boyNGzdW5xQq5ebNm3jrtdfQqU0bxPbvj6NHn+uXvUZpNBosXboUkZGR6Ny5M5YvXw7NY1v/1XY7duyAn58fXFxc0Lx5c2zdurXcfjt37kTfvn0RFhaG999/H9nZ2WaP7ciRIxg0aBDat2+Pt956C+np6WX6aLVafP311+jSpQsiIiLw5ZdfQv3Y1pBAyW1ss2bNQvv27RETE4O9e/dWOS4iwsaNGxEVFYUOHTrgo48+Qn5+fpXHY8+Gykwd+wOIANAZQAiAmyiZOp5b0XPq8tTx3r17MWrwYEyWyhAgliBep8WPeh0OHjtm8tZ2hYWF6BQUDK+MDMSKxMjW67FEr8WgiROxcPFiQ7/169dj/MiRmGRji05WcpxRF+O/BflY+NlneOONN2r6FAEA165dQ8fgYPTXE3qIJUjRavGlToNPV6zA8BEjzHLM5wURISYmBteuXUPnzp0BlNwj2qxZM2zduvWZKOK/atUqTJo0Cb169UKzZs2QmpqK3bt34/PPP8err75q6PfRRx/hq6++QnR0NBwcHHD69GncvXsXp06dgqOjo1liW7NmDaZPn46ePXvCzc0NycnJ+P3333Hq1Cm88MILAEpeg9jYWJw/fx5du3aFSCTCkSNH4O7ujj179hjKOmZmZiI4OBg+Pj5o27YtsrKysG/fPsycOROvv/66ybG9/fbb2LRpE6KiomBjY4OEhAQUFxfjt99+K3PrEnv2VOca7aMp46MAThHRE//srquJlojQplkzvJ3zAN3lf/1QfFNYgPOdOmJTXJxJ4y1duhTbZs/Gd/K/Cj5k63SIyM3BpbQ0NGjQAADg6eSMKUT4h81f9WbjlEV4T1mE+0plDZxZWZPGjYPV1m2YUeqYSWo1XoUeV+/c4YLz1XDo0CGMGTMGM2fONJSp1Gq1WLRoETZs2IDw8HALR/hk9evXR58+fdCxY0dD28mTJ7F161ZkZWUBAHJyctC4cWO89957Rvvnrl69GgMGDMC7775b43FpNBp4eHjglccW7m3duhXe3t5YunSpIdYBAwZg9uzZhs0wdDodFi9ejOXLl6Nnz54AgDlz5iA+Ph4jSv1xmZGRgU8++QTp6ekmLfq7efMmWrVqhXnz5hmeR0RYvnw5pkyZYvHbp1j1VZRoKzN13JeIFhPRb5VJsnVZXl4e0m7cRFcr451x+sqtcaQK06rxe/agL8HoE4yzWIx2trZISEgwtN19kIO+1sZFBHrJrZGtUqGgoMDk41bGkYMH0e+xcnttZTKguBjXr/Nt1NURHx8Pf39/o1rQEokE/v7+iI+Pt2BklZeVlYWgoCCjtqCgIGRnZ0Ov1wMoucTSuHHjMpvUt2nTBgcPHjRLXGlpaZBIJGVWx7dt2xa//vqr4etHr0HpHafEYjH8/f2N7tU9dOgQAgONNyurX78+6tevj3PnTCuQd/z4cfj5+RklZ0EQEBAQYLbvB6sdeNsLEygUCkgkYtx7+IvkketaLdxdXU0ez+2FF/B4yiIi3NBo4ebmZmiTi8W4odUa9UvX6SAVBLNtE+dWvz6ua3VGbXl6PfI1Gjg7O5vlmM8LNzc35OTklGnPyckxet1rM5lMhszMTKO2zMxMyGQyw7Srm5sb7t+/b0i8pfs9mq2paS4uLsjPz4fqsSIqmZmZRt9bNze3cq8V5+TkGMXm7u6O+/fvG/XRarXIysoy+bVyc3NDZmZmmTUYWVlZZvt+sNqBE60JpFIpxo4ZgznFKuQ//OVxT6fDAq0ak956y+TxXn71VazVqJGoLgYAaImwVFkEWT1XdOjQwdCva8+eeDc3B1kPq/Xk6vV450EOQoKDzbZF2OQZM7BYp8HNhwleqdfjA5USMX37lfmEwkwTGxuLK1euIDExEUQEIsKZM2eQlpaGIUOGWDq8SomIiMC6desMMyqFhYVYu3at0fs2ICAAnp6eiIuLM1SaunnzJg4dOoTJkyebJS5XV1dER0dj8+bNhoVN9+/fx65du4yuqQ4cOBC3bt3CiRMnDK/B77//jkuXLmH48OGGflOnTsWBAwcMmxZoNBrs2LEDbdq0gY+Pj0mxRUREQCKR4MCBA4Y/PlJTU3H8+HFMmDChuqfOarNHb7KafAQHB1NdpVQqaeywYeRobU2tnV3I0dqaZr39Nun1+iqNt2nTJnJ3cqLmjo7kbmNDHQMD6dq1a0Z9NBoNhbVpQ1aCQD4SCckFgQKaNKHCwsKaOKUK/WvhQnJUKCjAyZmcrK1pSN++lJeXZ9ZjPi+OHz9OPj4+1LBhQ3J3d6emTZvSyZMnLR1WpRUXF1NgYCBJpVJyc3MjqVRK/v7+pFQqjfrdvn2bIiIiyNnZmby8vMjV1ZXWrFlj1tgePHhAMTEx5ODgQD4+PuTg4ECLFy8u0y8pKYmaN29O7u7u1LBhQ/Ly8qL4+Pgy/VasWEFOTk7k7e1Njo6O1L17d8rIyKhSbGlpaRQSEkKurq7k6elJbm5utG3btiqNxWofAKepnJxY4WIoQRB2omSj94oSdP+K/q2uLoYq7d69e7h58yaaNm1a7dWTGo0G58+fh729PZo2bVphv6tXr+LIkSNo164d/Pz8qnXMyiooKMDly5fRsGFDNGzY8Kkc83lBREhOToYgCGjVqtUzsdr4campqTh69CjCwsLQvHnzCvtdvXoV2dnZ8Pf3h5XV49tQmEd6ejru3LlT5rpoaUSEixcvQqfTwd/fv8IZIpVKhQsXLsDV1dVQ67s6UlJSUFBQgICAAF5YWIeYvOpYEITOfzcgER2u6N+eh0TLGGOMlWZyZai/S6SsZul0Ovzxxx+wt7dHo0aNqj3enTt3kJ2dDV9fX6NVlaXp9Xrs378fVlZW6NKlS7WPyZ5PeXl5uH79Ojw9PY32sDWn9PR0/PTTT+jYsSPatWtXrbF0Oh3Wr18PrVaLUaNGGa0EZ6zGlDefXPoBoBmATQAuAkh79Pi759Tla7Q17eeffybP+vXJy96eXKytqVtYGN26datKY2VmZlJMjx7kKJdTEwcHauDkRD+sW1em39dff00OEik5iURkL4jIxcqKdu7cWd1TYc8RnU5HM2bMIDs7O2rcuDHZ2dnRm2++SVqt1qzH9ff3J6lUSk5OTiSVSsnZ2bnKPy9ffPEFWVtbk7W1NdnY2JBcLqfZs2fXcMTseQJTr9E+8rDc4vsAPgcQA2AcSqac36/oOTx1XDmXLl1CZGgolskV6GBlBTURliiL8GsDd5y+cMHka3bRERHwTL6Af1orYC0IOKdW4/+URdh6YD/CwsIMxwxq1QpfODqjt1wOPYAfiwoxPz8Pd3JyYG9vb4YzZXXN4sWL8d1332H8+PFwcHBAXl4evvvuO4wYMQKzZ882yzF79+6NkydP4q233kK9evWgVCrx3Xff4fbt22VuwXmSGzduoGnTphg2bBg6duwIQRCQmJiIVatW4fjx4wgODjbLObC6rcoFKwBYE9EvKEmu14noAwDdajrA59HXS5dipFSGDg8Xh8gEAW9aK5B3+zZOnDhh0lgpKSk4l5SE2Q+TLAC0lskwSSLBss8/N/SbOXMmesmt0cfaGoIgQCwIGGVjCz+JBB988EGNnRur27788ksMGTLEMF1sb2+PoUOH4ssvvzTbMePj4zFs2DDUq1cPAGBtbY3Ro0fjwYMHSE1NNWms8ePHw9vbG+Hh4RCJRBAEAcHBwfD39+dbbViNq0yiVQmCIAKQIgjCVEEQBgKob+a4ngt3rl+H92NtgiDARyrDnTsV7kJYrrt378JTLof0sU/BPmIx0ktVcrp76xaal7PKsblEyhWfWKVlZGSUKdjwqEDFk2bJqkqn05U5pp2dHaRSKS5dumTSWOnp6eWuovfw8DD50zFjT1KZRPsGAAWA1wEEAxgFYIw5g3pedIyKwp7HEmOOXo+ThQUmb2bdunVrpBQV4bbOuILUHr0e4VFRhq+7RUdjm1IJfalfhsVE2KtSYsCAAVU4C/Y8Cg0NRVJSklFbUlISQkJCzHabklwuR2JiolFbamoq9Ho9unfvbtJYMTExSEpKgrZUxTW9Xo9Tp04hMjKyRuJl7JEnXqM1dBQEewBERE/c04mv0VZOfn4+wgID4Z+dg1ixBNl6Pf6r06LnmNH45IsvTB7v4wUL8O2//403xRI0EouxTavFIbkVTp47Z5huU6lUaOzqilZaHSbZ2qKYCJ/n5yPXxRmpJn6KZs+vY8eOoV+/foiKikLTpk2RlpaGffv2YcuWLWZbxf71119jypQpiIqKQuvWrZGeno7NmzcjPDzc5K3rdDodnJ2d4erqir59+0IsFmPv3r24fv06MjMzeScdViUVXaOtzKrjEADnAVx7+PgdQPDfPYdXHVdeVlYWzf7nPym4RQvq2q4drV27tspVpoiINm/eTNGdwimoeXN6+4036M6dO2X63L9/n3r37En15XJyUyho6Esvlanow9iTJCUl0fDhw8nf359iY2Pp9OnTZj/m999/T66urqRQKMje3p6mTJlS5bEePHhAgYGBZGNjQzY2NtSiRYsqr2BmjKh6q47PAZhCRPEPvw4HsJSIWlf0HP5Eyxhj7HlTnVXH+Y+SLAAQ0VEAT5w+ZowxxljlEu1JQRBWCILQRRCEzoIgLAXwqyAIQYIgBD3x2c8QjUaD//3vf5gyYQI+/OAD3Lhxo9x+mZmZ+PSTTzB5/MtYuXIlCgsLn3KkFTty5Ai6RkairZ8f3nnnnTLbhQEllwuOHDmCN6ZOxTtvvokzZ86UO5ZKpcKaNWsw5eUJ+GjRIpNXQptTTk4O/vvf/2LSpElYunQp8vLyqjyWXq/HRx99hMDAQHTo0AFbt24tt19hYSG+/vprTJo0CZ988kmZbeJMtWLFCgQHByMkJATffPNNuX3UajV+/PFHTJo0CQsWLMCtW7fK7Xfp0iX069cPAQEBGDduXLVjq0lr1qyBl5cX3NzcMHjwYCiVyjJ9iAgHDx7E66+/jnfeeafMQqtHsrOz8fLLLyMgIAB9+vRBcnJyuf1u376NRYsWYdKkSVi3bh2Ki4tr9JzKk5WVhc8++wyTJk3C8uXLK9wr+uLFi5g5cyamTJmC3bt3l9lG0FKICPHx8Zg2bRqmT5+OU6dOlduvuLgY69atw6RJk7Bo0SLDzkbsb5Q3n0zG12gP/c3jYHnPeRav0RYWFlJ4cDC1c3KiD+wd6P+cnMnZxob27dtn1O/8+fPk5uhIQ52d6UN7R4p2dqZmL7xQ7rXQp23OnDlkLQg01saG3rd3oBZSKbnb2VNubq6hj16vp9deeYW87ezpXXsHetPBkRrY2NDHCxYYjfXgwQNq26IFRTo70zx7Bxrl5Ewutrb022+/Pe3TKuPKlSvk5uZGHTt2pNjYWGrXrh298MILdOPGDZPH0mg05OPjQy4uLjRw4EDq2bMnyeVyGjFihFG/e/fuUZMmTSgoKIhiY2MpPDycXF1d6dy5c1U6h7CwMLKzs6OYmBiKiYkhW1tbioiIMOpTUFBAoaGh1LJlSxo6dCh169aNnJyc6ODBg0b9NmzYQDKZjNq1a0exsbHUokULsrGxoYsXL1YptprUv39/kslk1K1bNxoyZAg1aNCAbGxsKCsry9BHr9fThAkTyMPDgwYOHEj9+vUjZ2dn+vTTT43G+vPPP8nGxoaaNWtGsbGxFBYWRjKZjFavXm3ULz4+npycnKhr1640dOhQ8vf3p8DAQKOfg5p28eJFql+/PoWHh1NsbCyFhISQl5cXpaenG/X75ptvyNHRkXr37k2DBg0iLy8vGjJkCOl0OrPFVllvvPEGNWjQgAYMGEAxMTHk6upK8+fPN+qTm5tLgYGB5O/vT0OHDqWuXbuSk5NTubsePY9Q1Wu0VfEsXqP99JNPsG/hInxjrYDo4e0Jx4pVeFsk4OqdO4YaqF3bt0fPKykYpfhrN5B5hQUQBr6IFatWWSR2oGSXnfr29ljv4ooQWUkBDB0RhmZlwmvgi9iwYQMAICEhAbHR0dhrYwe7hzuV3NPp0CPvARIvXoSXlxcA4L1330XKyq/xmbXCcLvGz8oifOXshN9TUiy600yvXr1gZ2eHqFK3Le3cuRO2trZYv369SWPNmzcPX331FebOnQuZTAagZP/SDz74AMnJyYYdaSZOnIi0tDS89NJLhuceOXIEV69exbFjx0w65qZNmzBmzBjMnz8ftra2AEpWoM+ZMwcbN25Ev379AACLFi3Ctm3bMH78eMP3Ozk5GTt27EBaWpphpxlnZ2fExMSgU6dOhmOsW7cOeXl5OHv2rEmx1aS7d+/C09MTb7/9tuF9pdPpsHjxYnh4eODo0aMAgMOHD2PEiBF49913IZfLAZR8cv3oo49w8eJFeHh4ACi5pUgqlWLMmDGG78eJEyewadMm5ObmAij54ODr64uoqCgEBgYa2lavXo3o6GjMmzfPLOfatWtXuLu7o2vXroa2rVu3ws3NDatXrwYAPHjwwPD9cHd3B1Ayi/bZZ5/h008/xYsvvmiW2Crj9OnT6Nu3L2bOnAmFQgEAyM3NxaJFi3Dq1CnDrmJz587F/v37MXbsWMNrcPbsWezfvx9Xrlx5JnegqklVvkYrCIKbIAjfCoKw++HXLQVBGG+OIC1px4YN+IdIZEiyANDJSg4rtRrnz58HUJLMEhITEWutMHruGCs5dm7f/lTjfdyPP/4Id7HYkGQBQCwIGG9ji9/27DG07di2DQMFkSHJAoCbWIxoawV27dplaNv50yaMkkiMfnD6yK1x786dCqfUnwadTodffvmlzL2OkZGRiIuLM3m8jRs3okuXLoYkCwD16tWDr68vlixZYmjbsWNHmWN27NgRiYmJJk9br1y5EmFhYYYkC5QUXggNDcWyZcsMbVu2bEF4eLjRa9CqVSsUFxfj8uXLAEq2n8vPzzeU2HykW7duuHLliklx1bQFCxbA1dXVkGQBQCwWo3v37kZTvtu3b0doaKghyQIlfzwEBARg9+7dhrZLly6hW7duRt+P0NBQqFQqnDt3DgCQlpaGBw8eoE2bNoY+giAgIiKiwksC1aVSqXDs2DGEh4cbtUdGRuLnn382fH3w4EE0bdrUkGQBQCqVon379tiyZYtZYqusHTt2IDg42JBkAcDBwQGBgYFGP1fbtm1DRESE0WvQpk0bPHjwAGlpaU815mdJZa7RrgawF8CjMipXUFLEok6Ry+UoeuzTvZ4ISp3O8AtAIpFAJBJB9Vi/ItJDLns6e2xWxM7ODsq/pvsNikgPcalKUHJraxSV81dnkQCjX3RWcisoHxtLC0Ct1z+1/UTLIwgCJBIJ1Gq1UXtxcXGV4pLJZOVex1apVEaJ0MrKqsx1Po1GAwAV7pBUEWtr6wqPWfr+TblcXuaYRAS1Wm14raysrEBERoUXgJLvR0V7qz4tdnZ2UKvVZd6TxcXFRr+o5XJ5mdfzUb/S70mxWFymn1arBREZ9pu1srKCRqMpc91TpVIZjVWTxGJxubE9/p78u/MsneAsoaLYSr/XHvV7/D2p1+uh0WjM9v2tCyrzk+hKRBsB6AGAiLQAdGaNygL+MXEilul1KCj1A/qjUol6Hh6G6UO5XI5+PXvivyql4ZeHlgifq9UYOW6sJcI2iI2NRZEgYHNRkaEtT6/Hf/LzMWTcOEPbsOHDsUWjxtVSv5h/V6sRrzSuDPWPV17BF1qt0R8VK5VKBAUGGv1F/rSJRCIMHToUcXFxhtdAr9dj165dGDlypMnjvfHGGzh48CBycnIMbX/88Qdu3ryJ6dOnG9pGjx6NPXv2QKcreesTEfbs2YPevXubXNxg1qxZOH36tNEikvT0dCQlJWHWrFmGtrFjx+LAgQNGv9gOHz4Mb29v+Pj4AAAaNmwId3d3o++HTqfD9u3b0b59e5Piqmnz5s1DQUGB0WK7wsJC7Nq1C9HR0Ya2kSNHIiEhwWgBV2pqKlJSUhATE2NoCw8Px/bt2w1/4ADAnj174OzsjCZNmgAAGjVqhJYtW+LgwYOGPmq1GgcOHMCYMeYpaCeVSjFgwADs2rXL6DXYvXu30Xuye/fuuHfvHi5evGhoy83NxdGjRzFq1CizxFZZw4YNw+nTp3H37l1D240bN5CcnIxBgwYZ2saMGYP9+/cbJeWDBw+iZcuWhil+Vo7yLtyWfgD4FYALgMSHX4cBOPx3z3kWF0PpdDqaNG4c1VMo6CUXF2rv7EKN3d3p0qVLRv3u3btHQX5+1MrRkYa5uJKnnR317datVhR82LRpEylEImojlVE/a2uyFQQK8fcvs9Di6xUryMHamvq5uFJPF1dyUiho27ZtRn00Gg2NGDSI3G1sKNbFldo6OVHzxo3p2rVrT/OUypWdnU2hoaHk6elJnTt3poYNG1Lnzp0pPz+/SuP17duXpFIpBQYGUrNmzUgqlZZZiKNUKqlPnz7k7u5OkZGR5O3tTa1bt6a7d+9W6ZivvvoqSaVS8vf3p1atWpFUKqXXX3/dqI9Op6OxY8eSs7MzRUREUPPmzalx48aUkpJi1C8pKYlsbW3Jzc2N2rVrR7a2tuTh4WHWxT+VtWjRIpJKpeTl5UVBQUFkZWVFjRo1KrOd3pIlS8jOzo7CwsIoKCiIHBwcaNeuXUZ9CgsLydPTk2xsbCg0NJTc3d3JxsaGTp48adQvLS2NfHx8yNfXlyIiIsjFxYVGjhxp1i387t+/T23btiUvLy+KjIykBg0aUI8ePaiwsNCo35EjR8jFxYUCAwOpY8eOZGdnRwsXLjRbXKb47rvvyM7Ojtq1a0chISFkb29PmzdvNuqj1Wpp5MiR5OLiQhEREeTr60s+Pj6UlpZmoahrF1SjYEUQgC8B+ANIBlAPwBAiOlfRc57FxVCP/PHHH4iPj4e7uzt69uxZ7rSgXq/Hr7/+irS0NLRp08bkusTmVFBQYFhyP2rUqAprwN6/fx+7d++GRCJB3759K9y0+/z58zhx4gQaNWqEqKioWrMxNj28FeGPP/5AQEAA2rdvX62FGImJiViyZAkcHBwwc+ZM1K9f/r4Zp0+fxtmzZ+Hj44MuXbpUa3o2JSUFn376KQRBwIwZM+Dt/fgWEyUuXryI3377DQ0aNEDPnj0hKWdTCK1Wi88++wwXLlxAdHR0lT7dm8v9+/fxyiuvICMjA1OmTMGIESPK7Xfv3j3s2bMHMpkM/fr1g52dXbn9/ve//yEuLg5+fn6YPn260fX1R7RaLfbv34/09HSEhYXB39+/Rs+pPESEw4cPIyUlxfB7obz3ZFFREeLi4lBQUICoqCg0atTI7LFVVmZmJnbt2gWxWIy+ffvC0dGx3H7JyclISEiAh4cHoqKiyn1PPo8qWgxVqVXHgiBIADQHIAD4g4g0f9f/WU60jDHGWFWYvOpYEIRQQRDcAcN12WAACwF8KgiCs9kiZYwxxuqQv5v3WgFADQCCIEQC+BjAGgC5AFaaP7TaiYiwYtkyNPPwgFgkQlDz5tixY4elw2LV8GjRjUwmg4ODA6ZMmYL8/LJVRnfu3Ik2bdpALBbDy8sLS5cuNdveq4/Ex8ejQ4cOkEgkaNCgARYuXGhYkPUIEeGrr75C48aNIRaL0bZtW6NbtR7Jy8vDpEmTYG9vDysrKwwYMKDcDdPPnz+PqKgoSKVSODs74+233y53lfT69evRokULiMVi+Pr6Yu3atVU+z2vXrmHQoEGQy+Wws7PDhAkT8ODBgyqP96y7fPky+vTpA6lUCkdHR7z++uu1qgIdM02FU8eCIPxORG0e/v9XAO4T0QcPvz5LRIEVDVqXp46//M9/sPT99/EvmRVaS2U4UqzCP4tV+Pann9C7d29Lh8dMlJWVhVatWqFTp06IjIyEUqnEzp07IZPJcPDgQcM1tn379mHkyJEYNmwYWrZsievXr2PDhg2YNm0a3nzzTbPEdvbsWXTr1g2DBg1C27ZtkZGRgZ9++gm9evXCZ599Zuj373//G0uXLsWwYcPg6emJCxcuYP369di4caPhGj0RGe4D7tevH+RyOQ4fPoyEhARcvHjRcC3u1q1bCAwMRHR0NDp06IC8vDxs3boVjRs3xubNmw3H3LBhA9544w2MGDECzZo1Q2pqKn788Ud8/PHHGD16tEnnmZ+fDz8/P4SEhKBLly5Qq9WIi4uDSqVCQkLCc1cEISMjAwEBAejcuTPCw8NRWFiIHTt2wMHBAXtK3RPPah+Tr9EKgpAMIJCItIIgXAbwChEdefRvRFTh6oK6mmj1ej0audbD9xIJWkr/WoDxs7IIP3h74UgFNYNZ7fXJJ59g+/btRrdX6PV6zJs3Dzt37kRISMnPTKdOndCqVSvD10DJLTnLli3DnVKVw2rSiBEjoNFojCpg5eXl4cMPP8SNGzfg6OgIrVYLd3d3vP7660a3XZ04cQJpaWk4dOgQAOD48eN46aWXMHv2bKMFXKtXr0ZsbCymTZsGAJg5cyYSExMxePBgQx+NRoM5c+bg5MmThgpBLVu2RM+ePeHn52fod+XKFWzbtg1//vmnSee5bNkyrF69GuPH/1UHh4jw0UcfYfXq1Wbb37a2WrBgAQ4cOGC0aEyn0+H999/HgQMH0Lp1hRunMQurSmWo9QAOC4KwHYASwKNt8pqiZPr4uZOXl4eCwgKjJAsA7WRWuJSSYqGoWHUkJyejcePGRm0ikQhNmjTBpUuXDG1//PGHIck84uHhAaVSaSj/V9MuXLhQ5pj29vZwcXHB9evXAZSUKnyUbEtr2rSpoXoUUFJVycfHp8wq6caNGxtVaTp//rxRJSeg5D5Rb29vo/H+/PNPw72rpY+Zmppq8nT6hQsXyrwGgiDAx8fH6DV4XiQnJ5dZgS4Wi8u8J9mzo8JES0QLAUxHSWWocPrrp0cE4DXzh1b72Nvbw87WFhc1xhVUTqqL0dLX10JRserw9/c3JK1H9Ho9UlNTjT6tURlhnwAAIABJREFUtWjRoswntfT0dCgUigpvjaquVq1alTlmXl4esrKyDMnQ2dkZUqnUqNAAUJIIW7RoYfjaz88PaWlpZSomXb9+3ejWl9atW+Pq1atGfTQaDa5evWr0/Xg0Xfz4MZs2bWryVG95rwERIS0tDS1btjRprLogICCgzGug0+mQmpr6XH4/6oK/vQmQiBKIaCsRFZZqu0JEieYPrfYRiUSYOXcuXlOpcEZdDDURflEpMU9djJkLFlg6PFYF48aNQ0pKCvbu3QulUomsrCysXbvWcM3wkblz52LLli04d+4cNBoN0tLS8P3332PmzJlmu7f43XffxS+//IITJ05ArVbj1q1b+O677/Dyyy8bkrtEIsGMGTOwatUqXL16FRqNBr///ju2bduGOXPmGMYKCwuDt7c31q1bh+zsbBQVFWHPnj24du2aUcWkyZMnIzExEYcOHYJKpUJGRgZWrVqFbt26GX2CnTt3LtavX49Lly5Bq9XiypUr+PHHHzF37lyTz3PkyJFIT09HXFwcCgsLkZOTgx9//BFubm5l6ks/DyZMmIBLly7hwIEDUCqVyMzMxJo1axAUFISAgABLh8eqorwqFtV9PIuVoSpLr9fTimXLqJmHB0lEImrr60s7duywdFisGlJSUigmJoZkMhnZ29vTlClTyq0ytXPnTmrdujWJxWLy8vKiZcuWkV6vN2ts8fHx1KFDB5JIJNSgQQNauHBhmQpHer2elixZQo0bNyaxWEyBgYFlqioRlWxxNmnSJLK3tycrKysaMGAApaamlul37tw5io6OJun/t3fnUVFdWdvAnyMyCBQWKDgwqDgAQQREQcQBMUCLdAQVZ2OrJCZmmQF9k9j9xry97MTuxAzGNit+bUdFo3bUJGoUxKGjIDIqzjgrKIgYI6jMxfn+AKutLkxE61oMz28tV8Lh1Lm7wOWue++5e5uaSjs7O7lgwQJZUVGhN2/jxo3S3d1dmpiYyD59+sh169Y98fu8cuWKHDt2rDQ3N5cqlUq+/PLL8s6dO0+8XnOXm5urrVimVqvl66+/rldlipoesE0eERGRcp64TR4RERE9ORaoJELdxqakpCRYWlpi9OjROi3yHrZ27Vrs27cPXl5eeOuttxqs8VpdXY3ExEQUFRUhKChIZxORUmpqarB06VKcPn0aYWFhmDZtWoPz7t69q31GNTw8HF26dGlw3uXLl7F//36o1WpEREQ02KGotrYWK1euxOHDh+Hv74+5c+c2WPu5oqICCQkJuH37NkaMGKHtPKSkmpoaJCUloaCg4JnVOjYGKSUyMjJw/Phx9OrVC8OHD38m7REf1EDv2rUrwsLCWOv4tzR0Pflp/7Tke7TU8nz00UfSxsZGDh48WPr6+kpbW1u5b98+nTm//PKL7Ny5s7SxsZEBAQHSwcFBqlQqeeLECZ15Z86ckc7OztLDw0MOHTpU2traytmzZ+t1UDKkh7v3BAQESJVK1WD3nqSkJKlWq2X//v1lYGCgVKlU8tNPP9Vbb+HChdLGxkYGBQXJfv36SXt7e5menq4zJz8/X3bo0EHa2trKgIAAaWdnJ21tbfW6O2VnZ8tOnTpJLy8vOWTIENm+fXs5f/58Re9tX7x4Ubq6uko3Nzdt957Jkycr2r3HGO7fvy9Hjhwpu3TpIocPHy67d+8ufXx8ZHFxsWLHrKmpkZMnT9Z273Fzc5Ourq4N3utvjcB7tET6MjMzERERgQULFsDW1hZA3TOza9asQX5+vrYh9/Dhw3H79m3MnTsXJiYmkFJi69atyM3NRV5eHoC6D60+Pj7w8vLS7patqKjA8uXL8b//+7+K9UN1cnKCl5cXoqKiIISARqPB8uXL0aVLF+zbtw9AXVcnZ2dnxMbGonfv3gDqnsFdunQp9uzZA19fXwB1FbBmzZqFuLg47Vl9Tk4Otm3bhqtXr2rPXLy9vWFhYYGZM2eiTZs2kFIiPj4ev/zyi7bfam1tLVxdXREWFqbdwX3//n189tlnWLFihU6vWUMaPHiwtqsMUNeP9ssvv8Rrr72G1157TZFjGsPbb7+NAwcO6PwOtmzZgg4dOmDTpk2KHHPFihVYsWIF5s6dq+2a9KBLUmpqqiLHbE54j5aoAevWrUNQUJA2yQKAm5sbnJyckJSUpB3LysrCmDFjtI/yCCEwevRoFBQUaBuWnzt3DoWFhRgyZIj2dRYWFhg5ciTWrFmjSPwFBQW4ceMGIiIitM+vmpiYYMyYMUhPT9fOS0hIQI8ePbRJFqh7BjcwMBDr16/Xjq1ZswbDhg3TuXTu4+ODdu3a4dChQ9qxM2fOICoqSnuZUgiBqKgonDt3TvusbkZGBoQQ8PPz077OysoKwcHBWLt2rYF/EnWuXbuG06dPIyQkRDtmZmaG0NBQxY5pLOvXr8eoUaN0fgcRERHYtm0bqqt/tcHaE1uzZg1CQ0N1WhOGhITg9OnTuHbtmiLHbAmYaKlVq6ioaLCfqZmZmU4h/draWpibm+vMedCr+MG8yspKmJmZ6d0jMzc3b7AovyFUVlZCCKF3j8zc3FynOMXjvs+Kigq99wnUfWD4rZ9HQ8c0NzfXK2BhZmaG8vLyx3yHjVNZWQlTU9Nn+jswlqqqKr3fgZmZGTQajV7jCUOprKzUO2abNm1gamqKyspKRY7ZEjDRUqsWFRWFjIwMVFX9p9pXcXExzp49q1NjuFevXtq6wQ+kpqZCrVZrG3d7enpCCKFTJq+2thYpKSkYO3asIvH36NEDKpUKaWlpOuP79+9Hn4eqlYWHh+P06dPas2+g7h/NzMxMREVFacfGjRuHw4cP6/xDnZ+fj+vXr2Po0KHaMScnJ72fx7///W84Ojpqk9ygQYNQXFyMK1euaOdoNBocPnxYp5ayIbm6ukKtVuPYsWPaMSklkpOTER0drcgxjSUyMhIHDx7UGUtJScGQIUNgYWGhyDGjo6ORnJysU2bz2LFjsLW1fSab3Jqthm7cPu0fboai5kKj0cgpU6ZIZ2dnGR0dLUeNGiVtbW3lV199pTPv5MmT0tLSUrq7u8sJEybIAQMGSDMzM7llyxadeXv27JFqtVqOHDlSxsTESHd3dzlo0CBFiw1s2rRJmpmZSX9/fzlhwgTZp08faWVlJXNzc3XmLV++XNrZ2cmIiAgZHR0tnZyc5IsvvqizMam6ulqOHj1a9ujRQ44bN06Gh4dLtVotN23apLNWSkqKtLCwkP369ZMTJ06UPj4+0tzcXG8T2datW6VarZZhYWFy/PjxsmfPnjIsLExWVVUp9vNITk6Wtra2Mjg4WMbExEhPT0/p6+srS0tLFTumMRQUFMju3btLPz8/OXHiRBkUFCQ7deokT58+rdgxS0tLpa+vr/T09JQxMTFyxIgR0tbWViYnJyt2zOYE3AxF1DApJfbu3Yvt27fDysoK06dPh6enp968W7duYf78+Th69Ci6d++Ojz/+GG5ubnrz8vPzsWbNGhQWFmL48OEYO3as9jKzUs6cOYO3334bV69ehZ+fHz7++GN07NhRb96JEyfwzTffoKysDGPGjEFISIjepV2NRoOEhAQkJiZCrVZjxowZOvd2H7h27Rrmz5+P06dPw93dHZ988glcXFz05l26dAlr1qzBzz//jLCwMERGRipWtvKBgoICrF27Fvn5+QgKCsL48eMbvCTe3N27dw/ffPMNjh49Cjc3N8yYMQN2dnaKHrOyshJbtmzBoUOH4OzsjD/84Q+PfEystWl0m7ynwURLREStDXcdExERGQETLT0TNTU1WLFiBQYNGgQ/Pz8sWbIEZWVlxg4LQF3ruffffx++vr4ICgrCqlWr9NrJNcaDx2OsrKxgY2PzVI/25OXl4bXXXkO/fv0QHh6OXbt2PfFaj+vWrVt455134O3tjeHDh2Pjxo0N9pjNycnB5MmT4eXlhZiYGGRnZyseG1FzxERLz8SUKVPw5Zdfwt/fH8HBwfjuu+/w/PPPo6amxqhxlZeXY+jQodizZw9GjhwJHx8ffPTRR3jllVeeaD1nZ2ccP34ckyZNwrx58+Dv7485c+YgPj6+0Wvl5+fD398fly5dwujRo+Ho6IjZs2fjyy+/fKLYHkdJSQkGDRqE9PR0hIeHw8PDA3/84x+xcOFCnXmHDh1CSEgINBoNfv/738PExAShoaH46aefFIuNqLniPVpSXHZ2NiIiIrBo0SLtpqDa2lp8/vnnWLJkiVEfu1i9ejWWLVuGV199VbspqKKiAu+//z4yMzPRq1evx17r7Nmz8PLywqJFi+Dg4KAd37hxI44cOYKSkpJGxTZv3jxcuHBB5+dz48YNLFu2DNevX1fkEY6lS5di69atmDlzpnbs7t27+POf/4zz58+jU6dOAIAhQ4bAzc0NAQEB2nlZWVnIyclBRkaGweMiag54j5aM5vDhw/D09NTZedumTRs899xzSElJMWJkwMGDB7XPvz5gYWEBDw8PHD58uFFrLVu2DO3bt9dJsgDg5+f3RAUEkpOT0a9fP52xzp07Q6VS4dy5c41e73EcOHBAr7m4SqVCz5498fCH58zMTPj4+OjM8/X1RVZWVoOXmYlaMyZaUlyXLl1QXFysN37r1i107drVCBH9h6Ojo15sUkrcvHmz0bENGjQIpaWlOsUvAKCoqOiJko+joyOKiop0xiorK3H79m3tmaWhOTk54ebNmzpjtbW1uHHjhs7Pw8HBQS+2oqIi2Nvb6z0uRNTaMdGS4iIjI/HLL7/gwIEDqK2thZQSR44cwalTpzB9+nSjxhYbG4vMzEycOnUKUkpoNBrs2bMHJiYmGDFiRKPWevHFF9G2bVts3LhRW47u2rVr+P7773Vq7z6uN998E4mJibh+/TqAukvamzdvRmhoqGKJdu7cuTh48CAuXLgAoK7l344dO+Ds7KxzBjtv3jxs2bJFezm8tLQU3377LebNm6dIXETNWkNVLJ72DytD0X87c+aM7N+/v7Szs5MODg6yd+/e8tChQ8YOS0pZV82pW7dusnPnzlKtVsvAwEB5+fLlJ1orNTVVWlpaSjMzM6lWq6Wpqans3bv3E8e2cuVK2aFDB+nk5CRVKpUcO3asXvs7Q/vhhx9k165dZdeuXWX79u3lyJEjZWFhoc4cjUYj3377balSqaSLi4tUqVTyrbfeanGt6IgaA6wMRU3BlStXUFVVhd69ezepS4y1tbU4d+4cLC0tG6xu1Fjbtm1Deno6Xn/9dXTu3Pmp1qqsrMSFCxdgb2+vd/9XKRqNBufOnYONjQ0cHR0fOe/u3bu4evUqXFxcYGNj80xiI2qqWBmKiIhIQdx1TK3W+fPnkZ2drbdJ6UncuXMHGRkZuHHjxq/Ou3r1KjIzM59pUY7CwkJkZGQ0+jEiosa6d+8eMjMz2YP2MTHRUouVl5eHgIAABAYGYvz48XBycsLmzZufaC0pJf74xz/CxcUFU6dORZ8+fTBt2jS9Hqe3b99GeHg4vL29MWnSJDg6OmLFihWGeDuPVF5ejsmTJ8PNzQ1Tp06Fs7MzFi1axMdsSBEff/wxHB0dMXnyZPTt2xcvvPACP9z9hra/PYWo+ZFSIiIiAm5ubpg5cybatGmDK1euYM6cOXB3d9d7VvS3rFy5Et9++y3ee+89tG/fHhUVFVi/fj3i4uJ0KjVNmTIFtbW1+OCDD2BqaoqioiIsXrwYvXv3RlhYmKHfJgDgjTfewKVLl/CXv/wFFhYWuHPnDr766iu4uLggNjZWkWNS6/T9999j2bJleOedd9CxY0dUVVVh8+bNmDVrFrZu3Wrs8JosntFSi5SWlobS0lKEhYVpG5F3794dw4YNw8qVKxu93vLlyzFmzBi0b98eQF1Ri5iYGKxbt057STo/Px9paWmIiorSFufo1KkTwsPDsXz5cgO9M13l5eXYuHEjYmJitJWi1Go1XnjhBfz9739X5JjUen3xxReIiIjQtmA0MzPDuHHjkJSUhFu3bhk5uqaLiZZapJs3bzZYPMHOzg4FBQWNXq+4uBj29vY6YzY2NtBoNLh//752jq2trV7vWXt7+9+8p/uk7t27B6CuetN/H/O/C08QPa2ioiK9PscWFhZQqVT4+eefjRRV08dESy1SQEAAzp8/j9LSUu2YlBLHjx9/ouIRQ4cO1etOc+rUKbi4uECtVgMAPDw8UFJSgsLCQp15R48eRXBwcOPfxGPo2LEjOnXqhNzcXJ3x7OxsDB06VJFjUusVHByMo0eP6oxdvXoVGo0Grq6uRoqq6eM9WmqROnfujNdffx3Lly/H888/D5VKhfT0dFRWVuoUzH9cixcvxrBhw1BeXg53d3fk5eVh7969+Oabb7Rnze3atcOSJUuwaNEihIeHw97eHjk5OTh37hw2bNhg6LcIABBC4PPPP8eMGTMQGhoKJycnnDlzBhkZGUavI00tz8KFC+Hv7w+NRgMvLy/cuHEDSUlJ+Oyzz/Su5NB/8IyWWqzFixfjs88+Q0FBAdLT0xEZGYmUlBRYWVk1eq3nnnsO6enpcHZ2RnJyMtq2bYukpCSMGjVKZ96cOXOwYcMGlJaWIiUlBQMHDkRWVpZiJROBuhKXCQkJEEIgOTkZPXr0QEZGBtzc3BQ7JrVOzs7OyMrKgpeXF1JSUlBRUYGtW7di2rRpxg6tSWPBCiIiIgNgwQoiIiIj4D3aVkhKiaSkJGzZsgWmpqaYMmUKhgwZ8sTr5eXlYdWqVbh8+TKCgoIwffr0J7o8q4Rjx45h9erVKCkpQUREBKKjo9G2re5f+9raWuzYsQPbtm2DlZUVZsyYgQED9D6UoqSkBGvXrkVmZiZ69+6N2NhYo7f5U0JZWRkWLlyI3bt3w9raGgsXLsS4ceOMHRaAul3W8fHxOHz4MFxdXREbGwtnZ2djh0X0q3hG28pIKREbG4uXXnoJJSUlKCoqwvjx47Fo0aInWi85ORk+Pj5IS0sDAPzjH/9A//79m8RW/1WrViEkJASXLl1CRUUF/vSnPyEyMhLV1dXaObW1tRg/fjzi4uJQVlaG/Px8jBo1Cp9//rnOWgUFBejXrx82bNgAIQR++ukn9OvXD0eOHHnWb0tR9+7dQ48ePfDdd9+hf//+cHBwwLRp05pE4Yvi4mL4+vpi9erVAIDU1FR4e3sjNTXVyJER/Treo21lUlJSMGHCBLz77rvaAgelpaX44IMPkJWVhZ49ez72WlJKeHh4ICQkRKdX6YYNGzBgwAAsXbrU4PE/rpKSEjg7O2PBggXa7jkajQbLli3De++9h6lTpwIAtm/fjjfffBPz58/X7pq8ffs2PvzwQ1y4cEHbLWf27Nm4fv06xo4dqz1GamoqcnNztR8yWoKXXnoJe/fuxTvvvKMt9FFQUIAPP/wQeXl5T92J6Gm88cYbOHHiBCZNmqQdy87ORmpqKo4fP96kukFR68R7tAQA+PHHH+Hn56dNskBd4QVvb28kJiY2aq1r167h5s2b8Pb21hkPCgrCjh07DBLvk3qw+/bhxGBiYoKAgABs27ZNO7Z9+3YMHDhQ59EEOzs7uLu7Y+/evdqxnTt36j2XGhAQgGPHjuk8q9vcJSUlYcSIEdokCwBdu3aFo6MjVq1aZcTI6v7uBgUF6Yz5+voiLy9PsYIgRIbARNvKWFtb6xXCB+pK+TX2vmq7du1QXV2NmpoanfGysjJYWlo+VZxPy9LSEuXl5Xrj5eXlsLa21n5tbW39yHkP/zwsLS31OvFUVlZCCNGinh80NzdvsONQWVmZtjCHsTT0O62pqUFNTY3OB0eipoaJtpWZMmUKMjMzdcoQXrx4EWfPnkVUVFSj1urYsSMCAwORmJio7RRTVVWF3bt3Y9asWQaNu7GGDRuG+/fv69xDvXPnDg4ePKhTsGLGjBk4dOiQTp3WkydP4saNGwgPD9eOzZw5Ezt37tTe35VSYufOnYiMjES7du2ewTt6NubMmYOEhASdbiyZmZkoLS3Fyy+/bMTIgFmzZiExMVFbW1pKiYSEBAwfPhy2trZGjY3oV0kpDf7Hz89PUtMVHx8vVSqV7N+/v/T29pZqtVomJiY+0VrXr1+XXl5eslu3bjIwMFDa2trKF198UdbU1Bg46sbLysqSnTt3lh4eHnLgwIFSpVLJv/71r3rzli9fLq2treWAAQOkp6en7Nixo0xJSdGZU1lZKceNGyc7dOggBw8eLJ2cnOTAgQPlrVu3ntXbeWZCQ0Olqamp7Nu3r3RycpLm5uZy06ZNxg5LVldXy6lTp0o7OzsZGBgoXVxcpLe3tywsLDR2aERSSikBZMkGciI3Q7VSd+7cwZ49e2BqaoqwsLCnutQrpURycjLy8/MxcOBA9OnTx4CRPp2qqirs3bsXJSUlGDFixCM389y6dQv79u2DpaUlQkNDH3kp8tSpU8jJyYGrqysGDRrUYjfg5OTk4J///Cfs7e2xYMECo98KeNjZs2eRlZWFbt26ISgoqMX+Dqj5edRmKCZaIiIiA+CuY9JTU1MDjUZj7DAUV1NT0+AGHyKiZ4GJthW6fPkyIiMjYWlpCSsrK0yaNKlF9i69ceMGvLy8YGFhAWtra9jb22Pjxo3GDouIWhleOm5l7t27B3d3dwQEBCA4OBg1NTVITExEQUEBcnJyYGJiYuwQDaZLly5wdHTExIkTYWlpibS0NHz77bc4ePAgAgICjB0eEbUwvHRMAIBNmzahS5cuCA8Ph7m5OaysrDB27FhUVlYiKSnJ2OEZzLZt21BaWorZs2ejffv2MDU1xdChQzF48GDMnz/f2OERUSvCRNvK5Obm6hVhF0LAxcUF586dM1JUhpeWlgYXFxe9M/SePXvi2rVrRoqKiFojJtpWxtvbG5cuXdIZq62txcWLF9G3b18jRWV4ISEhuHLlik4DAQA4c+YMXF1djRQVEbVGvEfbypSXl8PLywu9evXCiBEjtPdoq6qqkJaW1qKeSezRowfatWuHCRMmwNraGocOHcKuXbuQnZ3doj5UEFHT8Kh7tOxH28q0a9cOKSkpePfdd7FkyRKYmppi4sSJ+PDDD1tUkgWAEydOICIiAn/7299QU1ODrl274scff2SSJaJnime0REREBsBdx0REREbAS8cKunPnDn744QeUlZXhd7/7XbPbhCOlRHp6OjIzM+Ho6IjIyEiYmZkZO6xGKS0t1T7qExoa2qTqMDdlubm52Lt3L9RqNcaMGQOVSmXskIiaLV46Vsju3bsxedw4DLZoB5WU2FNRjnlxcXh/8WJjh/ZYqqqqMG7cOGRnZ8PDwwM3b95EaWkp9u/fj549exo7vMdy4MABREdHo2fPnrCyssLx48cxe/ZsfPTRRy3ufrShSCkRFxeH+Ph49OvXD3fv3sXly5exfft2vabrRKSLTQWeobKyMrh07oxVZhYYaG4OALil0eD39+9hY2JCs/gH65NPPkF8fDzmzJmDtm3rLnzs27cP165dw6FDh4wc3W+rqqqCk5MTpk6dCg8PDwDA/fv38cknn+Drr79GWFiYkSNsmnbt2oVXXnkFcXFx2o49J0+exObNm5GXl9eimtwTGRrv0T5DSUlJ6Gtmrk2yANDRxARTTUywce1aI0b2+OLj4xEaGqpNsgAQHByMkydPorCw0IiRPZ4DBw6gQ4cO2iQLAFZWVhgyZAjWr19vxMiatnXr1mHYsGE6bfH69u0LlUrVLD5gETVFTLQKqK6uRkOf+83rv9ccVFdX61VVEkLAxMSkWbyH6upqnQ8JD7Rt2xZVVVVGiKh5qKqqarDetampKX9uRE+IiVYBoaGhyCq7j/MPJaSy2lpsgkT0xIlGjOzxxcTE4ODBg3j41kJWVhacnJz0Sjg2RcOHD0d+fj7y8vK0Y9XV1UhLS0NMTIwRI2vaJkyYgNTUVJ0PU1euXEFhYSGGDRtmxMiImi/eo1XIuvh4vPnqq4g2M4dNbS1+kLUYGRWF/7d2bbPYiHP37l2MHDkSd+/ehbu7O4qLi5Gbm4vdu3fDz8/P2OE9li1btiA2NhZ+fn6wtrbGkSNHEBQUhA0bNqBNG37GbIhGo8HEiRORkZEBX19f3Lt3D9nZ2Vi9ejWio6ONHR5Rk8bNUEZw8eJFbPjmG9y/exejX3gBQ4YMaRZJ9oHq6mps27YNaWlpcHZ2xrRp09ChQwdjh9UoV69exfr161FSUoJRo0YhODi4Wf0OjEFKif3792P37t1Qq9WYNm0aXFxcjB0WUZPHREtERKQg7jomIiIyAiZaoiZs//79cHd3h0qlgr29Pd5+++0nXqumpgZffPEFfH194eHhgYULF+LOnTsGjJaIGsISjERNVEpKCkaNGoXnn38eEyZMwM2bN7Fq1SqcP38e33//faPXmz59Oo4fP46wsDBYWFggJSUFP/74IzIzM2FhYaHAOyAigGe0RE3W3LlzERwcjDFjxsDJyQn9+/dHXFwcdu7cidu3bzdqrVOnTmHPnj149dVX4eHhgR49emDatGkwMTHBpk2bFHoHRAQw0RI1Wfn5+fDx8dEZc3BwgI2NDfbv39+otTIyMvDcc8/pNIUQQsDDwwOpqakGiZeIGsZES9REWVlZ6ZW7rKioQGlpqU5pycfh5OTUYOnMoqIidOvW7aniJKJfx0RL1ETFxcXhu+++w5UrVwDUNUWIj4+Hk5MTPD09G7VWSEgI2rRpg4SEBFRXV6O2thZHjhzBsWPHMHPmTAWiJ6IHuBmKqImKi4vD2bNn8emnn8LU1BSVlZXo3r070tLSGr2WiYkJ9uzZgxkzZmDhwoUwNTWFg4MDduzYga5duyoQPRE9wIIVRE1cRUUF0tLS4OrqapAKTTdv3kRFRQWcnZ1ZJYvIgB5VsIJntERNnIWFBYKDgw22noODg8HWIqLfxnu0RERECmKiJSIiUhATLRERkYKYaImIiBTEREtERKQgJloiIiIFMdESEREpiImWiIhIQUy0RERECmKiJSIiUhBLMNIjlZQxMV3FAAAFkElEQVSU4Ouvv0ZqaipcXFzwyiuvoHfv3sYOi4ioWeEZLTWouLgY/fv3x7/+9S+oVCrk5ubC398fe/fuNXZoRETNCs9oqUFLlixBt27dMGnSJO1Yz549MWfOHFy4cIFdX4iIHhPPaKlBu3btQmBgoM6Yp6cnSktLtY3IiYjotzHRUoNUKhXu3r2rM1ZdXY3KykpYW1sbKSoiouaHiZYa9NJLLyExMRHl5eUAgNraWiQkJCAoKAj29vZGjo6IqPngPVpqUGxsLHJycrBo0SL06dMHhYWF6Ny5M3bu3Gns0IiImhUhpTT4ogMGDJBZWVkGX5eevby8PGRlZcHZ2RkDBgzgJigiokcQQmRLKQf89zjPaOlXubi4wMXFxdhhEBE1W7xHS0REpCAmWiIiIgUx0RIRESmIiZaIiEhBTLREREQKYqIlIiJSEBMtERGRgphoiYiIFMRES0REpCAmWiIiIgUx0RIRESmIiZaIiEhBTLREREQKYqIlIiJSEBMtERGRgphoiYiIFMRES0REpCAmWiIiIgUx0RIRESmIiZaIiEhBTLREREQKYqIlIiJSEBMtERGRgphoiYiIFMRES0REpCAmWiIiIgUx0RIRESmIiZaIiEhBTLREREQKYqIlIiJSEBMtERGRgphoiYiIFMRES0REpCAmWiIiIgUx0RIRESmIiZaIiEhBTLREREQKYqIlIiJSEBMtERGRgphoiYiIFMRES0REpCAmWiIiIgUx0RIRESlISCkNv6gQxQCuGnxhIiKipqublNL+vwcVSbRERERUh5eOiYiIFMRES0REpCAmWiIFCCH+JIQ4JYQ4LoTIEUIEGHj9YCHEj487boDjRQkhnnvo65+EEAMMfRyilqitsQMgammEEIEAIgH0l1JWCiE6AjAzclhPKwrAjwBOGzsQouaGZ7REhtcFwC0pZSUASClvSSkLAEAI4SeEOCCEyBZC7BZCdKkf/0kI8bkQIlUIcVII4V8/7l8/drT+v26PG4QQwkoI8bUQIrP+9WPqx/8ghPhOCJEohDgvhPjoodfMFkKcq4/nH0KIvwshBgN4AcDH9WfnPeunxwghMurnDzXED46oJWKiJTK8JADO9QnoSyHEcAAQQpgCWA5gvJTSD8DXAD546HVWUsrBAObWfw8AcgEMk1L6AlgE4MNGxPEnAPullAMBjEBdorSq/54PgIkAvABMFEI4CyG6AngPwCAAoQDcAUBKmQpgO4D/kVL6SCkv1q/RVkrpD+BNAO83Ii6iVoWXjokMTEp5TwjhB2Ao6hLcv4QQ7wLIAtAXwB4hBACYACh86KUb619/UAhhI4RQA1ABWCuE6A1AAjBtRChhAF4QQiyo/9oCgEv9/++TUpYAgBDiNIBuADoCOCClvF0/vhlAn19Z/7v6/2YD6N6IuIhaFSZaIgVIKTUAfgLwkxDiBIAZqEtIp6SUgY96WQNfLwbwbylltBCie/2aj0sAGCelPKszWLcxq/KhIQ3q/i0QjVgbD63x4PVE1ABeOiYyMCGEW/0Z6AM+qKuUdhaAff1mKQghTIUQng/Nm1g/PgRASf0ZZ3sA1+u//4dGhrIbwDxRf/oshPD9jfkZAIYLIWyFEG0BjHvoe3dRd3ZNRI3EREtkeNaou9x7WghxHMBzAP5PSlkFYDyAvwkhjgHIATD4odf9IoRIBfAVgNn1Yx8BWCKEOIS6S82NsRh1l5qPCyFO1n/9SFLK66i7B5wOYC/qdhiX1H97E4D/qd9U1fMRSxBRA1iCkagJEEL8BGCBlDLLyHFY199jbgvgewBfSym/N2ZMRM0dz2iJ6GH/J4TIAXASwGUAPxg5HqJmj2e0RERECuIZLRERkYKYaImIiBTEREtERKQgJloiIiIFMdESEREpiImWiIhIQf8ft2Ppmz+KDYIAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "x_min, x_max = x_train.T[:, 0].min() - .5, x_train.T[:, 0].max() + .5\n", "y_min, y_max = x_train.T[:, 1].min() - .5, x_train.T[:, 1].max() + .5\n", "\n", "plt.figure(2, figsize=(8, 6))\n", "plt.clf()\n", "\n", "# Plot the training points\n", "plt.scatter(x_train[0,:], x_train[1,:], c=y_train, cmap=plt.cm.Set1,\n", " edgecolor='k')\n", "plt.xlabel('Sepal length')\n", "plt.ylabel('Sepal width')\n", "\n", "plt.xlim(x_min, x_max)\n", "plt.ylim(y_min, y_max)\n", "plt.xticks(())\n", "plt.yticks(())\n", "\n", "intercept = np.ones((1,x_train.shape[1]))\n", "x_train_i = np.concatenate((x_train,intercept), axis=0)\n", "# x_train_i = x_train" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Let's define the sigmoid function and some training data" ] }, { "cell_type": "code", "execution_count": 93, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [0]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]\n", " [1]]\n" ] } ], "source": [ "y_train = y_train[:,np.newaxis]\n", "print(y_train)" ] }, { "cell_type": "code", "execution_count": 96, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(3, 150)\n" ] } ], "source": [ "print(x_train_i.shape)" ] }, { "cell_type": "code", "execution_count": 98, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from numpy.linalg import inv\n", "from numpy import linalg as LA\n", "\n", "n=y_train.shape[0]\n", "\n", "def sigmoid(z):\n", " return 1/(1+np.exp(-z))" ] }, { "cell_type": "code", "execution_count": 103, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 103, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXUAAAD7CAYAAACVMATUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAZC0lEQVR4nO3dfYwdV3nH8e+Tu0HEvCiVvGoojndbKUIkiECySh2CkJU1FUms8A9/BJkikKqt14BCAUWlkSKBFFWVKsRbk3QLqkAxIAqBojSpKG9tUEXQ2nFCXLdSANuYhGZBImlq1NbO0z/mbryevS8ze889c+bM7yONfO/MeOaZcycn4zPPOcfcHRERycMFTQcgIiLhqFIXEcmIKnURkYyoUhcRyYgqdRGRjKhSFxHJSOVK3cx6ZvaImd0/YNtuM3vGzI70lzvChikiIlXM1Nj3VuAY8PIh2x9y972ThyQiIltVqVI3sx3ATcCdwAdCnHj79u0+Pz8f4lAiIp1x6NChX7r77LDtVZ/UPw7cBrxsxD7XmtmjwJPAh9z96KgDzs/Ps7q6WvH0IiICYGYnRm0f26ZuZnuBp9390IjdDgNz7n4l8Cng60OOtWRmq2a2ura2Nu7UIiJSU5UXpdcBN5vZceBLwPVmdu/GHdz9WXd/rv/5AeBCM9tePpC7r7j7grsvzM4O/deDiIhs0dhK3d0/7O473H0euAX4jru/Y+M+ZnaJmVn/8zX94/5qCvGKiMgIdbJfzmNm+wHc/R7gbcCymZ0BfgPc4hr+UUQkOmuq7l1YWHC9KBURqcfMDrn7wrDt6lEqSTh4EObn4YILij8PHmw6IpF22nLzi0goBw/C0hKcPl18P3Gi+A6wb19zcYm0kZ7UpXG3336uQl93+nSxXkTqUaUujTt5st56ERlOlbo0bufOeutFZDhV6tK4O++EbdvOX7dtW7FeROpRpS6N27cPVlZgbg7Mij9XVvSSVGQrlP0iSdi3T5W4SAh6UhcRyYgqdRGRjKhSFxHJiCp1EZGMqFIXEcmIKnURkYyoUhcRyYgqdRGRjKhSFxHJiCp1mZgmuBBJh4YJkIlogguRtOhJXSaiCS5E0qJKXSaiCS5E0qJKXSaiCS5E0qJKXSaiCS5E0qJKXSaiCS5E0qLsF5mYJrgQSYee1DOnHHKRbtGTesaUQy7SPXpSz5hyyEW6R5V6xpRDLtI9qtQzphxyke5RpZ4x5ZCLdI8q9Ywph1ykeypnv5hZD1gFfu7ue0vbDPgEcCNwGniXux8OGahsjXLIRbqlzpP6rcCxIdtuAC7rL0vA3RPGJXIe5duLVFOpUjezHcBNwGeG7PJW4PNe+AFwsZm9IlCM0nHr+fYnToD7uXx7Vewim1V9Uv84cBvw/JDtrwR+tuH7qf46kYkp316kurGVupntBZ5290Ojdhuwzgcca8nMVs1sdW1trUaY0mXKtxeprsqT+nXAzWZ2HPgScL2Z3Vva5xRw6YbvO4Anywdy9xV3X3D3hdnZ2S2GLF2jfHuR6sZW6u7+YXff4e7zwC3Ad9z9HaXdvgG80wq7gGfc/anw4UoXKd9epLot56mb2X4z29//+gDwE+AJ4G+AAwFiEwGUby9Sh7lvavqOYmFhwVdXVxs5t4hIW5nZIXdfGLZdPUplpAMHYGameEKemSm+i0i6NJ66DHXgANy9oRvZ2bPnvt91VzMxichoelKXoVZW6q0XkeapUpehzp6tt15EmqdKXYbq9eqtF5HmqVKXodbnM626XkSapxelMtT6y9CVlaLJpdcrKnS9JBVJlyp1Gemuu1SJi7SJml9ERDKiSr3F9uwpOgWtL3v2NB3R1mkSDEleiJs0xo3u7o0sV199tcvWLS66F1NGnL8sLjYdWX333uu+bdv517FtW7FeJAkhbtJANzqw6iPqVo390lI2aAT7voZ+0i2bny9mMyqbm4Pjx2NHIzJAiJs00I2usV8keZoEQ5IX4iaNdKOrUpfGaRIMSV6ImzTSja5KvaUWF+utT5kmwZDkhbhJI93oqtRb6lvf2lyBLy4W69tGk2BI8kLcpJFudL0oFRFpEb0ozViMtFnlj4u0i4YJaKmDB4txWE6fLr6fOHFuoK2q/5obd4wQ5xCRuNT80lIx0maVPy6SHjW/ZCpG2qzyx0XaR5V6S8VIm1X+uEj7qFJvqRhps8ofF2kfVeotFSNtVvnjIu2jF6UiIi2iF6VbECM3u8o5lCMunaAbPaxR4/JOc0l1PPUYY3tXOYfGGJdO0I1eGxpPvZ4YudlVzqEccekE3ei1jWt+UaVecsEFgyeZMIPnn493jhhxiDRON3ptalOvKUZudpVzKEdcOkE3enCq1Eti5GZXOYdyxKUTdKOHN6rBfZpLqi9K3Yt3NHNz7mbFn9N4Z1PlHDHiEGmcbvRa0ItSEZF8TNymbmYvNrMfmtmjZnbUzD4yYJ/dZvaMmR3pL3dMGrjAgQMwM1O8M5qZKb7X2Q7p5NyLSCSjHuP7T/EGvLT/+ULgYWBXaZ/dwP3jjuUtaX5JwfLy+am768vycrXt7unk3ItIOIRsfjGzbcD3gWV3f3jD+t3Ah9x9b9VjqflltJkZOHt28/peD86cGb8d0sm5F5FwgqQ0mlnPzI4ATwP/tLFC3+DafhPNg2Z2xZDjLJnZqpmtrq2tVbqArhpUYW9cP247xBkPXWOui6SlUqXu7mfd/XXADuAaM3tNaZfDwJy7Xwl8Cvj6kOOsuPuCuy/Mzs5OEnf2er3R68dth3Ry7kUknlp56u7+a+B7wFtK65919+f6nx8ALjSz7aGC7KL1uUCHrR+3HdLJuReRiEY1uPfb22eBi/ufLwIeAvaW9rmEc0MOXAOcXP8+bNGL0vGWl917veLlY693/kvQKtvd08m5F5EwmPRFqZm9Fvgc0KN4sv+yu3/UzPb3/6dwj5m9F1gGzgC/AT7g7v866rh6USoiUt+4F6Uz4w7g7o8Brx+w/p4Nnz8NfHqrQYqISBga+2WAEJ1pqnQMmvQYMSbaCHEdyQjxw2p2E0ndqLaZaS6ptqmH6ExTpWPQpMeIMdFGiOtIRogfVrObSALQ2C/1hOhMU6Vj0KTHiDHRRojrSEaIH1azm0gCNElGTSHG7Dcbvq1qcY87RoyJNkJcRzJC/LCa3UQSoEkyagrRmaZKx6BJjxFjoo0Q15GMED+sZjeRFlClXhKiM02VjkGTHiPGRBshriMZIX5YzW4ibTCqwX2aS6ovSt3DdKap0jFo0mPEmGgjxHUkI8QPq9lNpGHoRamISD7Upt5i49KdlQ6dqBSS+1OIQZox6jF+mkvKzS8pGJfurHToRKWQ3J9CDDI1qPmlncalOysdOlEpJPenEINMjZpfWmrc5BOanCJRVWYv6UIM0hhV6okal+6sdOhEpZDcn0IM0hhV6okal+6sdOhEpZDcn0IM0pxRDe7TXPSidLxx6c5Kh05UCsn9KcQgU4FelIqI5KNzL0pjDJkdKwVYeeg1taXAQgyUH0KIjhCxxqiX6kY9xk9zmUbzS4whs2OlACsPvaa2FFiIgfJDCNERItYY9XIexjS/ZFWpz80N/u9lbi7cMdabKctLr5fetXRKWwps3A0U6zrGnadKHDH+g5NNxlXqWbWpxxgyO9YY4xqWu6a2FFiIgfJDGHeeWGPHt+V3S0in2tRjDJkdKwVYeeg1taXAQgyUH0KIjhCxxqiXWrKq1GMMmR0rBVh56DW1pcBCDJQfQoiOELHGqJd6RrXNTHOZVp56jCGzY6UAKw+9prYUWIiB8kMI0REi1hj18gK61KYuIpK7TrWph5JTrrskKEZe9p49xc21vuzZ08wxJL5Rj/HTXFIdJiCnXHdJUIy87MXFwTfY4mLcY8hUoOaXekKMUz7uGBruusNiDIQfIu82Vu6u1Kbml5pCjFM+7hga7rrDNBC+TJkq9ZKcct0lQcrLlilTpV6SU667JChGXvbiYr310zqGNGNUg/s0l1RflLrnlesuCYqRl11+0bmVF5whjiHBoRelIiL5mPhFqZm92Mx+aGaPmtlRM/vIgH3MzD5pZk+Y2WNmdtWkgYuISH1V2tT/B7je3a8EXge8xcx2lfa5AbisvywBdweNsi/WmP0hhJgHIYlrCRFElZ5WMc5T5Rwp9Aqr0uknRO+2GDdYW270tsRZxai2mfICbAMOA79fWv/XwNs3fP8P4BWjjlW3TT3WmP0hhJgHIYlrCRFElZ5WMc5T5Rwp9Aqr0uknRO+2GDdYW270tsTZR4hJMoAecAR4DviLAdvvB9644fu3gYVRx6xbqccasz+EEPMgJHEtIYKoMqtIjPNUOUesGVBGGXT+9WVdiJlcYtxgbbnR2xJn37hKvdaLUjO7GPga8D53f3zD+n8A/tzdv9///m3gNnc/VPr7SxTNM+zcufPqE4N61g0Ra8z+EELMg5DEtYQIokrPxBjnqXKOFHpRhiivWGU+Tltu9LbE+cIpA/YodfdfA98D3lLadAq4dMP3HcCTA/7+irsvuPvC7OxsnVNHG7M/hBDzICRxLSGCqNLTKsZ5qpyjLb3CQvRui3GDteVGb0ucVY16jO8/xc8CF/c/XwQ8BOwt7XMT8CBgwC7gh+OOqzb1FlyL2tQ3xzptalNXm/oYTNqmDrwWeAR4DHgcuKO/fj+wv//ZgL8Cfgz8iDHt6b6FSt093pj9IYSYByGJawkRRJWeVjHOU+UcKfQKq9LpJ0Tvthg3WFtu9LbE6eMrdXU+EhFpkc6N0tiWVNJOSSUHOEQcsY4R4lpy0aVrDWHUY/w0l2mM/ZJQs5esS6W9MkQcsY4R4lpy0aVrrYgQeerTWKZRqSeUSirrUskBDhFHrGOEuJZcdOlaKxpXqWfVpp5QKqmsSyUHOEQcsY4R4lpy0aVrrahTbeptSiXtjFRygEPEEesY43TpRu/StQaSVaUeY/4BqanKjxLjhwsRR6xjhLiWXHTpWkMZ1TYzzWVak2QkkkoqG6WSAxwijljHCHEtuejStVZAl9rURURy16k2dUlUiHHMY+UqxxgIP5VrzSn/O5V+DikY9Rg/zSXlOUoloBBjrsTKVY4xaE8q15pT/ncq/RwiQc0v0qj5eRg0xPLcHBw/XnyemYGzZzfv0+vBmTPVjhFCiDjG7ZPKtcY6TwwxriWh8hrX/KJKXaYrxDjmsXKVYwyEn8q15pT/nUo/h0jUpi7NCjGOeaxc5RgD4adyrTnlf6fSzyERqtRluqrkGS8tDf676+tj5SqHiGPcPqlca07536n0c0jFqAb3aS56UdohIcYxj5WrHGMg/FSuNaf871T6OUSAXpSKiORDbepdl0JubYgYrriieCm1vlxxRTNxhDhPCr+J5GvUY/w0FzW/RJBCbm2IGC6/fHBu9+WXx40jxHlS+E2k1VDzS4elkFsbIoZxaYCx4ghxnhR+E2k15al3WQq5tSFiCFGpp5L/ncJvIq2mNvUuSyG3NoUYYsYRYzx1kRFUqecshdzaEDFcfnm99dOKI8R5UvhNJG+jGtynuehFaSQp5NaGiKH8srTOS9KQcYQ4Twq/ibQWelEqIpIPtanL9IXIu04lt1s55DJMW+6NUY/x01zU/JKJEHnXqeR2K4dchkno3kDNLzJVIfKuU8ntVg65DJPQvaE8dZmuEHnXqeR2K4dchkno3lCbukxXiLzrVHK7lUMuw7To3lClLpMJkXedSm63cshlmDbdG6Ma3Ke56EVpRkLkXaeS260cchkmkXsDvSgVEcnHxG3qZnapmX3XzI6Z2VEzu3XAPrvN7BkzO9Jf7pg0cBERqa9Km/oZ4IPu/mpgF/AeMxs06MZD7v66/vLRoFHmKEaHnVhCdBxK5VpCOHAAZmaKzIiZmeJ7bDmVp9Qzqm1m0AL8PfDm0rrdwP11jtPpNvUYHXZiCdFxKJVrCWF5+fzrWF/K85BOU07lKZsQsk3dzOaBfwFe4+7Pbli/G/gqcAp4EviQux8ddaxOt6nH6LATS4iOQ6lcSwgzM3D27Ob1vR6cORMnhpzKUzYJ1vnIzF4K/DNwp7vfV9r2cuB5d3/OzG4EPuHulw04xhKwBLBz586rTwy68bogRoedWEJ0HErlWkIIMaHHpHIqT9kkSOcjM7uQ4kn8YLlCB3D3Z939uf7nB4ALzWz7gP1W3H3B3RdmZ2crX0R2YnTYiSVEx6FUriWEXq/e+mnIqTyltirZLwZ8Fjjm7h8bss8l/f0ws2v6x/1VyECzEqPDTiwhOg6lci0hLC3VWz8NOZWn1Deqwb3fNPNGwIHHgCP95UZgP7C/v897gaPAo8APgDeMO26nX5S6x+mwE0uIjkOpXEsIy8vuvV7xgrLXi/uSdF1O5SnnQZ2PRETyoQG9UpVTHnEKedkiAsBM0wF00sGDRRvr6dPF9xMnzrW57tvXXFxbceAA3H33ue9nz577ftddzcQk0mFqfmlCTnnEKeRli3SIml9SdPJkvfUpG1Shj1ovIlOlSr0JOeURp5CXLSIvUKXehJzyiFPIyxaRF6hSb8K+fbCyUrShmxV/rqy07yUpFC9Dl5fPPZn3esV3vSQVaYRelIqItIhelJa0Kj28LcG2Jc5YVB7SpFHdTae5NDFMQKuGmW5LsG2JMxaVh0wZGibgnFalh7cl2LbEGYvKQ6Ys2HjqoTVRqbdqmOm2BNuWOGNReciUqU19g1alh7cl2LbEGYvKQxrWqUq9VenhbQm2LXHGovKQpo1qcJ/m0tR46q0aZrotwbYlzlhUHjJF6EWpiEg+1KYusi7EuO/KQZfEaTx16YYQ477nNA6+ZEvNL9INIcZ9Vw66JEDNLyIQZtz3nMbBl2ypUpduCDHuu3LQpQVUqUs3hBj3XTno0gKq1KUbQoz7ntM4+JItvSgVEWkRvSgVEekQVeoiIhlRpS4ikhFV6iIiGVGlLiKSEVXqIiIZUaUuIpIRVeoiIhkZW6mb2aVm9l0zO2ZmR83s1gH7mJl90syeMLPHzOyq6YQrIiKjVHlSPwN80N1fDewC3mNml5f2uQG4rL8sAXcjk9OEDCJS09hK3d2fcvfD/c//BRwDXlna7a3A5/tT6P0AuNjMXhE82i5Zn5DhxAlwPzchgyp2ERmhVpu6mc0DrwceLm16JfCzDd9Psbnilzpuv/3cDDvrTp8u1ouIDFG5UjezlwJfBd7v7s+WNw/4K5tGCjOzJTNbNbPVtbW1epF2jSZkEJEtqFSpm9mFFBX6QXe/b8Aup4BLN3zfATxZ3sndV9x9wd0XZmdntxJvd2hCBhHZgirZLwZ8Fjjm7h8bsts3gHf2s2B2Ac+4+1MB4+weTcggIlswU2Gf64A/BH5kZkf66/4M2Ang7vcADwA3Ak8Ap4F3hw+1Y9YnXrj99qLJZefOokLXhAwiMoImyRARaRFNkiEi0iGq1EVEMqJKXUQkI6rURUQyokpdRCQjjWW/mNkacKKRkxe2A79s8Px1tCVWxRlWW+KE9sSaQ5xz7j6092ZjlXrTzGx1VFpQStoSq+IMqy1xQnti7UKcan4REcmIKnURkYx0uVJfaTqAGtoSq+IMqy1xQntizT7Ozrapi4jkqMtP6iIi2elEpW5mPTN7xMzuH7Btt5k9Y2ZH+ssdDcV43Mx+1I9h00hnKU3uXSHWVMr0YjP7ipn9e3/i9GtL25Mo0wpxplKer9oQwxEze9bM3l/ap/EyrRhnKmX6J2Z21MweN7MvmtmLS9vrl6e7Z78AHwC+ANw/YNvuQesbiPE4sH3E9huBBylmmdoFPJxwrKmU6eeAP+p/fhFwcYplWiHOJMqzFFMP+AVFznRyZVohzsbLlGLKz58CF/W/fxl416Tlmf2TupntAG4CPtN0LBPS5N41mNnLgTdRTPCCu/+vu/+6tFvjZVoxzhQtAj9293IHwsbLtGRYnKmYAS4ysxlgG5tnjKtdntlX6sDHgduA50fsc62ZPWpmD5rZFZHiKnPgm2Z2yMyWBmxPaXLvcbFC82X6e8Aa8Lf9prfPmNlLSvukUKZV4oTmy7PsFuCLA9anUKYbDYsTGi5Td/858JfASeApihnjvlnarXZ5Zl2pm9le4Gl3PzRit8MU/zS7EvgU8PUowW12nbtfBdwAvMfM3lTaXmly70jGxZpCmc4AVwF3u/vrgf8G/rS0TwplWiXOFMrzBWb2IuBm4O8GbR6wrpH7dEycjZepmf0WxZP47wK/A7zEzN5R3m3AXx1ZnllX6hRT8d1sZseBLwHXm9m9G3dw92fd/bn+5weAC81se+xA3f3J/p9PA18DrintUmly7xjGxZpImZ4CTrn7w/3vX6GoPMv7NF2mY+NMpDw3ugE47O7/OWBbCmW6bmiciZTpHuCn7r7m7v8H3Ae8obRP7fLMulJ39w+7+w53n6f4Z9h33P28/xOa2SVmZv3P11CUya9ixmlmLzGzl61/Bv4AeLy0WxKTe1eJNYUydfdfAD8zs1f1Vy0C/1barfEyrRJnCuVZ8naGN2k0XqYbDI0zkTI9Cewys239WBaBY6V9apdnlYmns2Nm++GFSbPfBiyb2RngN8At3n/tHNFvA1/r32MzwBfc/R9LcaYyuXeVWFMoU4D3AQf7/wz/CfDuRMt0XJyplCdmtg14M/DHG9YlV6YV4my8TN39YTP7CkVT0BngEWBl0vJUj1IRkYxk3fwiItI1qtRFRDKiSl1EJCOq1EVEMqJKXUQkI6rURUQyokpdRCQjqtRFRDLy/0dMHhtgBD9BAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "pos = np.flatnonzero(y_train == 1)\n", "neg = np.flatnonzero(y_train == 0)\n", "\n", "D = np.zeros((x_train_i.shape[1],x_train_i.shape[1]))\n", "\n", "plt.plot(x_train_i[0, pos], x_train_i[1, pos], 'ro')\n", "plt.plot(x_train_i[0, neg], x_train_i[1, neg], 'bo') " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "Now the main loop!" ] }, { "cell_type": "code", "execution_count": 104, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "29.99478823090076\n", "7.231943804354122\n", "6.3873372606361905\n", "5.795082269692554\n", "5.363623020737249\n", "5.04027181566039\n", "4.792588266652513\n", "4.599516692098444\n", "4.446788647042965\n", "4.324374516665336\n", "4.224994422875631\n", "4.143209952122904\n", "4.074850144530657\n", "4.016637841958841\n", "3.9659404074961837\n", "3.92060006404582\n", "3.878816649876752\n", "3.8390657999775204\n", "3.8000416763968246\n", "3.7606171133996624\n", "3.7198163661107704\n", "3.676797098632501\n", "3.6308391436938368\n", "3.581338113814782\n", "3.5278022773319115\n", "3.4698513260401\n", "3.4072158225320224\n", "3.3397362743804377\n", "3.267360970884204\n", "3.1901419514510843\n", "3.108228750583675\n", "3.0218598645924066\n", "2.9313521798094753\n", "2.8370888548339765\n", "2.739506326409844\n", "2.639081185392073\n", "2.5363176374819325\n", "2.4317361329463267\n", "2.325863547419039\n", "2.2192250605540367\n", "2.112337653335366\n", "2.005704965349225\n", "1.8998131463106525\n", "1.7951273109378239\n", "1.6920882567148827\n", "1.5911092102440705\n", "1.4925725027269765\n", "1.396826209720987\n", "1.3041809005169949\n", "1.214906711991938\n", "[[ 7.63301243]\n", " [ -9.16127903]\n", " [-12.61143099]]\n" ] } ], "source": [ "lmda = 1e-2\n", "theta = np.random.rand(x_train_i.shape[0],1)\n", "I = np.identity(theta.shape[0])\n", "# iterator 50 steps\n", "for x in range(0, 50):\n", " h = sigmoid(x_train_i.T.dot(theta))\n", " error = y_train-h\n", " tmp = (-1)*y_train*np.log(h) - (1-y_train)*np.log((1-h)) # add the regularization in this objective function\n", " J = 1/n*np.sum(tmp) # evaluate the log-likelihood function but its not used further\n", " np.fill_diagonal(D, h*(1-h))\n", " H = 1/n*((x_train_i).dot((D.dot(x_train_i.T)))+lmda*I)\n", " dJ = 1/n*(x_train_i.dot(error)-lmda*theta)\n", " step = np.linalg.solve(H,dJ)\n", " print(LA.norm(step))\n", " if LA.norm(step)<1e-2:\n", " break\n", " #update theta\n", " theta = theta + .1*step\n", "\n", " \n", "print(theta)" ] }, { "cell_type": "code", "execution_count": 106, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[-15.6679642 29.36680913]\n", " [-37.65503387 7.37973945]]\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWgAAAD8CAYAAABaZT40AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO2df7RVZZ3/X597QcrfXRE0geDbkGmaiXdBLZfZBKWpyeB3/C7AMWuaQWYxs3SVq1Fp1be+2bdZM+PYjDRKNmWGMlaSaDalNhq1Er7AqGRg4IhKoai3IKJA7v18/zj7XM7dd+9znr3PPmf/eD6vte7i7Oc853ke8d73ffPs570/oqoYhmEYxaMn7wUYhmEY0ZhAG4ZhFBQTaMMwjIJiAm0YhlFQTKANwzAKigm0YRhGQTGBNgzDSICIvE5E1onIEyLylIh8JmjvE5EHRWRr8Ocb2p7LzkEbhmG4IyICHKGqe0VkLPBj4CrgEmBAVb8gItcCb1DVv21nLnPQhmEYCdAae4PLscGXAnOB24P224E/aXeuMe0OkJajju3T4984Ka/pU9E7ZiDvJXSNo/f05r2EriLj/sCGX0yOff+st7zQxdV0mN/uynsFidmwk1dU9fh2xphx4pG658CgU99nBv7wFPCHhqblqrq8fiEivcAG4I+AZaq6VkQmqupOAFXdKSIT2lkv5CjQx79xEjeseCCv6RNz5BvuynsJXeUDDx6T9xK6Rs/0nwMwbf6neP6lvlHvT5k4wLpbP9vtZXWMnke+mPcSEiOf4bl2x9hzYJAbz5vq1HfuXVv+oKr9ce+r6iDwDhE5FlglIqe1u74obIvDGIVP4tzIDX9xP4ePOzCi7fBxB7jhL+7PaUXZU0ZxLjKq+hvgEeB84CUROREg+LPtf6qYQDvgm3v2ibp7Blg4ZyO3XrOSKRMHEFGmTBzg1mtWsnDOxhxXaBQNETk+cM6IyOuBOcAWYDVwRdDtCuDedufKbYvDKCY+uedGca6zcM7GygqyuefMOBG4PdiH7gHuVtX7ReSnwN0i8lHgeeDSdicygW6BT+7ZJ3H2DRPn7FDVJ4EzI9pfBWZnOZdtcTTBJ3H2jSj3bBhFwwTaAMw9Vxlzz+XFBDoGc8/VxdyzURZMoCPwTZx9cs++ibO553JjAu05Pomzb5g4lx8T6BC+uWef8M09G+XHBNpjzD1XF3PP1cAEugFzz9XF3LNRRkygA3wTZ5/cs2/ibO65OphAe4hP4tzInQ/NYNr8TzHmvf/EtPmf4s6HZuS9pMwxca4WFvXGP/fsE3X3fOdDM7jyH+azb/9hADz/Uh9X/sN8gMo+e8MoP+agPcNX97z0touGxbnOvv2HsfS2i3JaUfaYe64e3gu0uefq0rj3/MKu6Pqdce2GUQS8FmjfxNkn9xy+MTh5wq8j+8W1lw1zz9XEa4H2CZ/EOYoqV0sxca4u3t4k9M09+0Tcg/ihthf9wq43MHnCr7nhL+63G4RGofFWoH3Cd/dcp4rVUsw9VxsvtzjMPVcXn0IpJs7VxzuB9k2cfXLPPomz4Qe2xVFhqirOd619M59cNZMXBo5kct9ePjdvHQtmPZP3srqKuWc/8MpB++aeq8hda9/M4jvezfMDR6EIzw8cxeI73s3K516X99KMDrHkh8vzXkJueCXQPlFV9/zJVTPZd2DsiLZ9B8ZWKhHYCnPP/uCNQJt7rgYvDBwZ3e5JItA3cS6iexaRySLynyKyWUSeEpGrgvY+EXlQRLYGf7b9TemFQPsmzlV1zwCT+/ZGt1ckEWgcoojiHHAQ+LiqngK8E1giIqcC1wIPq+p04OHgui28EGifqLI4A3xu3joOP+y1EW1VSQS2wjf3XFRUdaeqbgxe/xbYDJwEzAVuD7rdDvxJu3NV/hSHb+656tRPa3xy1Uxe+PWRlgisKAVwz+NFZH3D9XJVHbUoEZkKnAmsBSaq6k6oibiITGh3EZUXaJ+ounuus2DWMyyY9YxX557NPbfPngMT+f6Ojzn2XvSKqvY36yEiRwLfBq5W1T0i0vYawzhvcYhIr4j8l4iM+rek1PhnEdkmIk+KSCFKVZh7ri4mztWlAO65JSIylpo4r1DVe4Lml0TkxOD9E4Fd7c6TZA/6Kmp7LVF8AJgefC0C/rXNdbWNb+Lsi3uG7ouzD6WyikJJxFmArwCbVfXGhrdWA1cEr68A7m13LieBFpFJwIXAbTFd5gJf1xqPAcfWf5MYnccnce429VJZz7/Uh6oMl8rqlkj75p5LwtnA5cB7ReTx4OsC4AvA+0RkK/C+4LotXPegbwI+ARwV8/5JwAsN1zuCtp3pl5Ye39yzT3TbPTcrlWU3JrOlDO4ZQFV/DMRtOM/Ocq6WDlpELgJ2qeqGZt0i2jRirEUisl5E1v/21wMJlmnEYe65s+RZKsvcs+GyxXE2cLGIbAdWUrP13wj12QFMbrieBPwqPJCqLlfVflXtP+oNfSmX3Bxzz9UljxuDeZXK8k2cy+Keu01LgVbV61R1kqpOBeYDP1TVPwt1Ww18KDjN8U5gd/08YDfxTZx9cs95ndqocqmsomDiHE/qc9AishhAVW8BHgAuALYB+4CPZLI6IxafxDlP8iiV5Zt7NuJJJNCq+gjwSPD6loZ2BZZkubCk+OaefSLvM89VLJVVFMw9N8eexVFCzD1XF3PPRiOViHqbe64u7bjnJTddwpfvO5vBoR56e4b4yw/+hGVX39P6gznhmzibe25N6QXaN3H2yT23K8633HsO9ROgg0O9wTWFFmlfMHF2w7Y4SoRP4twuX77vbEYfz5egvXj45p4NN0ot0L65Z59o98bg4FD0t3Zcu9E9zD27Y9+tJcHcczJ6e4YSteeJuWcjjtIKtLnn6pLFsbq//OBPGP20AQ3ai4Nv4mzuORmlvEnomzj75J6zOvNcvxFYplMcVcfEOTmlFGif8Emcs2bZ1fcUWpB9c89Gckq3xeGbe/aJvBODRucw95wOc9AFxtxzjfd//Eoe3vjW4evZM7bwg3+8tSNz3fnQjK48d8PcM/TvWsvF21fRt3+AgXF9rJ46j/UTZkX2WQtn5bTMXCmVgzb3XF3i3PMhcZbhr4c3vpX3f/zKzNfQreopvolzlHvu37WWhVvv4Lj9Awhw3P4BFm69g/5dayP7+EppBNo3cfbJPTfb2jgkzo3ICEedFc2qpxjpiNvauHj7KsYNjXyM67ihA1y8fVXTPr5RGoH2CZ/EuUh0o3qKb+45jr4YV9zYHtfHJ0oh0L65Z58o0o3BvKqnVJVmNwYHxkVXVGpsj+vjE6UQaJ8w9zyS2TO2EBU4qbVnS6erp5h7PsTqqfPY3zNyO2l/z2GsnjqvaR/fKLxAm3uuLi7u+Qf/eGuDSNe+OnWKY+Gcjdx6zUqmTBxARJkycYBbr1mZySkO38S51bG69RNmcef0y3l1XB8KvDqujzunXz7iFEdjH18p9DE738TZJ/ecZGujU0fqorDqKe3jeuZ5/YRZo47VxfXZsGbRhizWlhUi8m/ARcAuVT0taOsD/h2YCmwH/peqtrU/VngH7Qs+ibNv+OaePeFrwPmhtmuBh1V1OvBwcN0WhXXQvrlnn0h6Y9AlPBLVB7pb7NXwJzGoqj8Skamh5rnAe4LXt1Or3/q37cxTWIH2CXPP8dTDI/XzyfXwCByquB3V56N/txBVeG1wTOznuoG559IyXkTWN1wvV9VWv30mqupOAFXdKSIT2l1EIQXaJ/fsmzgndc/NwiN1oY3qc+Dg6G/t8Oc6jW/iXHT3POWo51j23kVOfb+0hldUtb/DS2pJ4fagfRJn30hz5tklPJIkSJJl6MQ4RNHFuUu8JCInAgR/7mp3wMIJtE/45p7T4BIeSRIk6VboxDf3bACwGrgieH0FcG+7AxZKoM09V5e0iUGX8EhUn8PGHGRs78GmnzOywUf3LCJ3AT8FThaRHSLyUeALwPtEZCvwvuC6LQq5B+0D5p7daNxnjjuNEden1ec6hbnn6qOqC2Lemp3lPIURaJ/cs2/i3O7zNlzCI3F9un2szjdx9tE9d5NCbHH4JM6+UaSHIRnZYuLceQrjoH2hqu75rrVv5pOrZvLCwJFM7tvL5+atY8GsZ1p+LhwwueCdP+OBx04rZbikku55E7VM3G7gGGr/gD891xV5Re4Cbe65/Ny19s0svuPd7DswFoDnB45i8R3vRk74JQunx38uKmByy73nUH9Af17hEiNgE3Af8FpwvTu4Bpa8bO65GxRii8MXquqeP7lq5rA419l3YGzLSiRRAZNw9ZSyVDSppHt+mEPiXOe1oN3oCrkKtE/uuariDPDCwJHR7S1CIa6hkaKHSyopzlBzzBFoTLuRPbkJdO8YK2dTFSb37Y1ubxEKcQ2NWEWTnIjxFFbppHvYFkcXqLJ7BvjcvHUcftjIfwu7hEKiAibh6ilFD5dU1j1D7YbgyJ2rUVVPjM5iAm20zYJZz3DL5T9iSt9vE1Uiiapgsnjumo5UNDFScDrwQYaddFTVE6Oz5H6Ko+pU3T3XWTDrGRbMeibxuefogMk92S2sg1TaPdc5vfZlZ57zwRx0B/FFnOv4FErxQpyN3GnpoEXkdcCPgHFB/2+p6qdDfd5D7clNzwZN96jqZ7NdqlFkshLnJTddwpfvO5vBoR56e4b4yw/+BGBE23ve8Qu2/nJi4jCLS2UWYzTmnvPDZYtjP/BeVd0rImOBH4vI91T1sVC/Napa/AOrXcI395wFS266ZERQZXCoN7hmRNvDG99K0jCLS2UWV8w9G92i5RaH1qifoxobfGmTjxiekZV7/vJ9ZxMOqtSuo9oO4RJmaVaZxYjH3HO+OO1Bi0iviDxOrULAg6q6NqLbu0TkCRH5noi8LWacRSKyXkTW7371d20su9iYe07H4FD6WyJpQzFJQzA+uWcT5/xx+olQ1UFVfQcwCZgpIqeFumwE3qSqZwD/AnwnZpzlqtqvqv3HHHdEO+suLL6Jc5Y3Bnt7hlJ/Nm0oJkkIxidxNopBIsuiqr+hVkr8/FD7nvo2iKo+AIwVkfFZLdIoJlmf2qjdEAzvnmlM2yHShmKKHoLJE3PPxaClQIvI8SJybPD69cAcYEuozwkiIsHrmcG4r2a/3GLjm3vOmmVX38PiuWvo7RkElN6eQRbPXTOqbfaMLYnDLFGhmCQhGHPPRh64nOI4EbhdRHqpCe/dqnq/iCwGUNVbgD8F/kpEDgK/B+arqt1IrDCdOvO87Op7WHb16KBKVFtSXCqzGOaei0RLgVbVJ4EzI9pvaXh9M3BztksrF+aeq4tP7tnEuVhY1DsDfBPnKPfsEgKJCqGEnXFUn7NP256q6orLfC3/W13F2aXySLjPdGBri88YXmMCbSQiTpxbhUCahVDqohnX57b7z+bgYO/w2C5VV1zmy4wmlUeGBTeqz/qGMaI+02XMPRcPexZHm/jmnqNwCYHEhVBq7c371MW5sa3ZXK7ztcLZPbtUHonqE8aqlZQGETlfRJ4WkW0icm2n5jGBNpyJuzHoEgKJC6E0tmcZVHGZLzPiKozsdujjOlaHMffsTnBgYhnwAeBUYIGInNqJuUyg28Dccw2XEEhcCKWxPcugist8zUh0YzDu2+AYhz6uY3UQE+fEzAS2qep/q+oBYCUwtxMTmUCnxDdxbnasziUEEhdCqT+trlmfMb2Do9qazeU6XxyJT21EVB5hbNDerE+Y8GeMonIS8ELD9Y6gLXPsJqHRklZnnus355qd4qjfmGt2qiKuT5pTHC7zZUb9pl6zUxxRfQpwisMr93zUBIbes8Ct72e+OF5EGm/jLlfV+l9W+OYGdOgBcpJXnmT62yfpjauvymXudjH3XF18OvMM5RHoL61ZtEFV+9sZo//kibruVjeB7v3jL8bOJyLvAv63qp4XXF8HoKr/t531RWFbHEZTfBJn3yiLOBeQ/wdMF5FpInIYMB9Y3YmJbIsjIb65524TFXiB5tsnWdF19/xdYAO1fxwLcBZwYXfGMnFOj6oeFJG/Br4P9AL/pqpPdWIuE+gE+CbO3XbPUYGXj/7dQlThtcExw21pK6E0Ixdxbtzh1IbrpCKd5ViGE8FTOx/o9Dy2xWFEksfWRlTg5cDBMcPiXKcSlVA2JGzPcCxzz+XBBNoR39xzHiSpbpK0EkozcrkxGHdvPs09+yzHMgqFCbQxirxuDCapbpKkbyGJOqjVrD2jscw9lwsTaAfMPXeHqMDLYWMOMrb34Ii2LCuh5Has7qyE7RmMZeJcPuwmYQt8E+c8j9XFBV6i2rK4QZjrmef6zbssTnFkOZZRKEygjWGKcOY5rupJJSuhXEh2ItpiLHPP5cS2OJrgm3v2Cd8Sg0Y5MQdtAM3ds0u1lKg+P/nZ1MTPwnCZK3dcqqekJRw4mQoMpJirYY2vjuujf+paAC7evoq+/QMMjOtj0xtO5/Rfbxq+Xj11HusnzBoxTP+utSM+E9XH6Bwm0DGYe67hUi0lqs9HvnAZBwd7SFLRxGWuLGjLPbtUT0lLVODk2YZr17lCazxu/wCX/+JrqMJYBofbzn3x0eGDHsftH2Dh1jsAhgW4f9daFm69g3FDB2L7GJ3Ftjgi8E2cm7lnl2opUX1qVVCSVTRxmatd2t7acKmekhaXkIrLXBFrHKODw+JcJ/x/Z9zQAS7evmr4+uLtq4bFOa6P0VlMoD2n1Y1Bl2opSUIjzSqauMyVOy7VU9LiGixpNVcba+nbPxD5Oq6P0VlMoEP45p5b4VItJUlopFlFE5e52iGTG4Mu1VPS4hpSaTVXG2sZGNcX+Tquj9FZTKA9xuVYnUu1lKg+tSooySqauMyVOy7VU9LiElJxmWs27O8JbTlJL68xsvhu+P/O/p7DWD113vD16qnzRo0T7mN0FhPoBsw9j2bhnI3ces1KpkwcQESZMnGAW69ZOeKmXVSfr167gsVz19DbUxPq3p5BFs9d0/QUh8tcacnsWN3pwAc55FKPCa6zOMVxIdDPISctwLTkcy15eTl3Tr+cV8f1odROcdzxlg/zjZM/PKLt0RPOHXF95/TLR9z8Wz9h1qhxwn2MzmIVVQJ8E+cihFK6hW9nnqsYSilSRZVuYg7aQ3wSZ9+oojj7jJ2Dxj/33Am6GWZJQuncc9oQzNeBZ+FmFgGw+ei3suyMj43osuSJGzllz5bh66g+RrEwB+0ZnXDP9YDJ8y/1oSrDAZM7H5rRtM+HP38Zt9x7DoNDtTPT9TDLkpsuyXyNpaAeMKkfk6sHUza1+FwgzlDbthbglD1bWPLEjcNd6uIsTfoYxcN7gfbJPXdqayNtmGVQk4dZklA695w2BPPs6Ka6ANepi3OzPkbx8FqgfRLnTtLNMIsrpRNn6GwIxiglXgu0T3TyxmA3wyyVJmUIxipbVRdvBdrcc3akDbP0SvIwiwuldM+QOgSz+ei3Rvwt1tqT9DGKh7cC7ROdPlaXNszyteuTh1kqTYoQzJIfLmfZGR8bFuD6V/iEhksfIxtE5FIReUpEhkSkP/TedSKyTUSeFpHzWo7lY1DFJ/fs25nn0rrnFPh05rlMQRUROQUYAm4FrlHV9UH7qcBdwEzgjcBDwFtUdTBuLO8ctE/i7Bs+ibNRXFR1s6o+HfHWXGClqu5X1WeBbdTEOpaWQRUReR3wI2Bc0P9bqvrpUB8BvghcAOwDPqyqBSuD4R9p3XM4UHLBO3/GA4+d1naVkyU3XTIqlHL2ads7MhfgFvpIU8Ek/JmzgCmhufqA7bRfxLXhjDNQezbHmfDqd/u4ef+i4Son03Zv45wX19DDEEP08PTRb2Hi/l2JK6GEK6i4VF1xIaoyC+BbtZaTgMcarncEbbG03OIIxPcIVd0rImOBHwNXqepjDX0uAP6GmkDPAr6oqk3/pvPY4vDNPacR6HBVkxp1lalx+LgDiR9itOSmS4KKKo2ncZUeGWJIe0e0pZlrlHsOVz6B2g23xj3dcAWTKFw/I7Q+TtFPMpEOi3MMBxF60dDf7Mi/6f09h7V80FG4gkracVzGfY1eRGqFBFqNfd6kG5l715a2tzjOmvomXbt0qVPfsYuufA54paFpuaoO7ymJyEPACREfXaqq9wZ9HmHkFscy4Keq+o3g+ivAA6r67bh1tNzi0Bp76+sOvsLfinOBrwd9HwOOFZETW41tdI607jkqUBIOk6SpclILn4yOSowU52zmAtxCH2kqmMR9xuVWjst8jTiIM8CYkDhD62opUURVUEkzjsu4YxkcIc5px+4gr6hqf8PXiA1/VZ2jqqdFfN3bZMwdwOSG60nAr5otwmkPWkR6ReRxYBfwoKquDXU5CXghtJBR1l1EFonIehFZv/vV37lMnRk+ued2bgy6BkqSVjlpJ3zSaq7IvWeX0EeaCibt3FPP+cByq0oorpVSklZUSdI/3Pe8SZWKoq8G5ovIOBGZBkwH1jX7gNNPjaoOquo7qCn+TBE5LdQlqhbEqG9HVV1e/410zHFHuEydCT6Jc7u4BkqSVjlpJ3zSbK7YG4MuoY80FUxcPxNFO5/NgFaVUFwrpSStqJKkf2PfsoqziMwTkR3Au4Dvisj3AVT1KeBu4OfAfwBLmp3ggISnOFT1N8AjwPmhtxJbd6MztHusLipQEv5dm6bKSS18Mjoq0SPh78/25wLcQh9pKpjEfcZFfF3ma2Ta6CYFwn9jB5HIEEojLpVQoiqopBnHZdzX6OWgjNzeqkq1FlVdpaqTVHWcqk5U1fMa3rtBVd+sqier6vdajdVSoEXkeBE5Nnj9emAOEH7CymrgQ1LjncBuVd2Z6L+qQ5h7TkZUoGTx3DVtVzlZdvU9kaGU269fkXqupsfqXEIfaSqYRH2mH5gX+ty0iD5JT3F8iBEiXQ+XfP3kj46slnLyn/PoCecySE8g4D1sPvqtiSuhRFVQaVV1xYWocb9x8oe54y0fjh27rO45a1xOcbwduB3opSbod6vqZ0VkMYCq3hKc9LiZmrPeB3ykfucyjm6d4vBJoC2UUl18CqVEiXMOpzgKUVGl5TloVX0SODOi/ZaG1wosyXZp7WPiXF1MnA0fqGxFFZ/EudO4VEvpJpmKczjMMh3YSvKKJi7jkKJPxIOSokIf3Qx4dHJ+29oYSWUF2ic66Z7DwZV6tRQgV5HOhHCYZTcjQyj1iibQXKRdxvkOtb3owSZ97qW20Tx0qG3/dw6jf/raYQEMhz6O2z/Awq13AHRFpPOe3zcq+SwOc8/Z4VItpZtk6p6jwixhXCqauIwzxOjjF2EGOSTOAeHwRlToo5sBj07Ob+55NJUUaJ/o9N6zS7WU0uJaqaRVvw5XPGkMb8SFPpKGR7JYS5bzmzhHUzmB9sk9d+PGoEu1lG6R+Y1B12+VVv06/C3XGN6IC30kDY9ksZY85veNSgm0T+LcLVyqpXSDjpzaiAqzhHGoaOI0Tg+1g6rNqB9kbSAc3ogKfXQz4NGJ+c09x2M3CUtKt47V1W8EFukUR2bUb/y1e4rDdRzHPq9+ty/2hET9dV6nOPKe3zcqI9DmnjvHwjkbq3OsLszpuB2jy2ochz6fevkLTd9fP2FWroKY5fzmnptTqS0OX/AtlOITPoVSTJxbUwkH7ZN79k2cO54YDAdMoiqhhKuluGx7uFRdCW1xfPWEj3Lp7hUdqYzSyW2ItHOdtH03pzzxMq/fd5DfHz6GzWcczy+n+vOz7ELpBdoncfaNrohzOGASfvbzemAjI8IjLcMr4aorceOEgiof2v1VehgafsZSL0OcsmfL8LVrKKSbYZK0c/35wU9zxroXGTNYexbQ4fsOcsa6FwFMpBuwLY4S4Zt77jguARMYFR5pGV6Jq5zS4pHYvQ3iXCeryiidCrOkneuUJ14eFuc6YwaVU554OfM1lplSC7S55+rSlYchtRMwafbZDldOSVsZpRNhljRznTfpRl6/72Dke3HtvlJqgfYJc88doJ3f780+2+HKKWkro3QiTJJ0rvqNwd8fHr27GtfuK6UVaJ/cs2/i3LVHiboETGD0T0mr8Epc5ZQWP231B+43klVllE6FWdLOtfmM4znYO/I32cFeYfMZx2e+xjJTyl9XPomzb3T1Oc9RAZMsTnHUK6c4nOJ4ddPIUMq03dvaPsXRzTBJkrkaj9XVbwTaKY7mlFKgfcI399x1MgyYjOBCoktcNYyz5IfLYebIt9dPmMU3p1+WcLLRdDPMknauX049xgS5BaXb4jD3XF18qpLiGxZKSYc56AJj7jmGcLjEtepJ2rGWAa80XI9ndIG3hGta8sPlkQEPqN5zLkyc01MqgfbJPfsmzs7uOSpc4lL1JO1YYXEmuF7GIZFOuKa6OIcDHpf/4muowtjgyf5WrcQozRaHT+LsG4m2NqLCJS5VT9KOFRbnqPYUa4oKeIzRwWFxrtPNaimdwEf3LCJ/LyJbRORJEVklIsc2vHediGwTkadF5LxWY5VGoH3CN/eciLiASJrQSVZjJRin/jCkJKGRblVLMTLjQeA0VX078AvgOgARORWYD7wNOB/4kog0fUp4KQTa3HN1SXxjMO5bIc23SFZjpRgnSWikrNVKfHTPAKr6A1WtRyIfAyYFr+cCK1V1v6o+C2xj1DmekZRqD9oHzD23YDYj93vBrepJ2rHGE73NMT75mhofJbp66rwRe9AAB6V3xB40dLdaSpYUUZz3HD3I997n/M+j8SLS+Dir5aqa5lmwfw78e/D6JGqCXWdH0BZL4QXaJ/fsmzinOlYXFS5Je4rDZawltD7F4TBO+DnPcQGPqDa7QZgLr6hqf9ybIvIQcELEW0tV9d6gz1LgILCi/rGI/k2f3FJogfZJnH2jrTPPWVVBcR0rfKQu7Tgh4gIeZRfkIrrnrFHVOc3eF5ErgIuA2apaF+EdwOSGbpOAXzUbpxR70D7gm3v2CZ+qpBggIucDfwtcrKr7Gt5aDcwXkXEiMo1ayYZ1zcYqrIM291xdcksMuoRJovrg8Lkmc3123LXebFUkcc8VrqhyMzAOeFBEAB5T1cWq+pSI3A38nNrWxxJVHWwyTnEF2ifMPXcBlzBJVJ/vUNs5HGzyuRZz+RI4SSrOVa2ooqp/1OS9G4AbXMcq5BaHT+7ZN3HOzT27hEmi+k89rG8AAAyGSURBVAwBYY/TKhgTMU7ZAydZYxVV3CicQPskzr6R68OQXMIkSQIqzfrGvFflwEnSG4NWUcWNwgm0T/jmnnPFJUySxBs06xvzXlkDJ53AKqq4USiBNvdcXXJ/lGhU9ZRwmCSqTw8QDuO2CsbMpmsVTYpAmmN1VlHFDft1lRM+uefcxRncQilxfVp9LsSSl5fTP330o0SreIMw7Zlnq6jiRmEE2if37JM4FwqXMElcH8cQSv3MczcrmpQVq6jSmkJscfgkzr5RCPdsdAQfEoN501KgRWSyiPyniGwWkadE5KqIPu8Rkd0i8njw9anOLLf8mHuuLpYYNLLGZYvjIPBxVd0oIkcBG0TkQVUNK80aVb0o6QLMPVcXJ/ecNt2X5lkcLinB6cDWDOZy4NKtK0ZU8F5zwjmpC8ZGlc/q5BaLuefu0FKgVXUnsDN4/VsR2UztEXlmBRPik3t2Fuc06b40Ja6ixrmX2rPEhhraGh8wmWCupO750q0rOPfFR4cfb9bLEOe++ChAYpGOKp/VyeSiiXP3SLQHLSJTgTOBtRFvv0tEnhCR74nI21zG88k9+yTOzqRN96UpcRU1ziCHxDkOh7nSbG2c8+KaUc+elKA9KVHlsyy5WA2cT3GIyJHAt4GrVXVP6O2NwJtUda+IXEDtCQbTI8ZYBCwCmNJnh/arivONwXbSfVmVper0Z2PoifnNENfejLiEYieSi+aeu4uTgxaRsdTEeYWq3hN+X1X3qOre4PUDwFgRGR/Rb7mq9qtq//ijjmxz6eXB3HMM7aT7sipL1eZn094YHIr50Ytrb0ZcQtGSi+XH5RSHAF8BNqtq5K9PETkh6IeIzAzGfTXLhRrlINGxurTpvjQlrqLG6aX1T0DaclotWHPCOaNKaWjQnpTVU+d1Jblo7rn7uGxxnA1cDmwSkceDtuuBKQCqegvwp8BfichB4PfA/IYqAl7jk3tOfOa5nXRf0pMVrinBBKc42jlWV78RmMUpjrjyWVneIDRxzgeXUxw/JrqWVmOfm6k9pNpowCdxTk076b6s5koxdhZnnr85/bLUx+rCWHKxmhQiSWiUH0sMVhdzz/lRmGdxVA1zz45kFUJx4bvABmqbvQKcBVyYfBhLDBrdwgTaaJvU7jmrEIoL32VkCEUbrlOItC+Ye84X2+LoAD6557a2NrIKobiwIWF7DD65ZxPn/DGBzhifxLltsgqhuBB3pijBWSOfxNkoBibQRmravjGYVQjFhbhzSE3PJ/mLuediYAKdIeaeE5JVCMWFsxK2hzD3bLgiIv9HRJ4MHr38AxF5Y8N714nINhF5WkTOazWWCbSRikyO1Z0OfJBDjvmY4LoTpzguBPo55JgluLYbhKMw99w2f6+qb1fVdwD3A58CEJFTgfnA24DzgS+JSLji5QjsFEdG+OSeMz3znFUIxYULsWN1LTBxbp/Qw+SO4NCdjrnASlXdDzwrItuAmcBP48Yygc4An8TZN3wS56ozeLCPvb9e4Nj7E+NFpPFg5nJVdf5mEJEbgA9Ru+X9x0HzScBjDd12BG2xmEAbiai75xWbFrD04c/z/O4pTDnmeW6YfT2XnX5X+xN0MrjSzVCMA92ugtIKc88jeEVV++PeFJGHgBMi3lqqqveq6lJgqYhcB/w18Gmib0k3PUdkAt0mPrrnFZsWsOi+L7PvtSMAeG73VBbd92WA9kS6k8GVFGN30j13uwqKkS2qOsex653UYlKfpuaYJze8Nwn4VbMP201Cw5m6e1768OeHxbnOvteOYOnDn29vgk4GV7oZinGgaFVQzD1nh4g0Fiu5GNgSvF4NzBeRcSIyjdqzE9c1G8scdBv45J4bbww+v3tKZJ+4dmc6GVxJOHan9567WQWlFSbOmfMFETmZWkG154DFAKr6lIjcTa2e60FgiaoONhvIBDolPolzmCnHPM9zu6dGtrfFMUQLZhbBlQRjd+PG4MC4Po6LEGOrglJ+VPV/NnnvBuAG17Fsi8NoSfhY3Q2zr+fwsb8b0Xb42N9xw+zr25uok8GVboZiHOhWFZRWmHsuNuagU+Cze4ZDNwIzP8WRVfWUNsbu1rG6blRBMcqPCbTRlLhQymWn35XNsbownQyudDMU40DeVVDMPRcf2+JIiE/u2bcqKT6FUkycy4E56AT4JM6+0a44X7p1RcsCsEULphjFxwTaiMQ399wOl25dwbkvPjocE+tliHNffBQ4VL27SMEUc8/lwbY4HDH3XF3adc/nvLhmVIZXgvY6RQumGOXABNoYhbnnZPQw1LK9KMEUc8/lwgTaAZ/cs2/inMWNwaGYH6PG9rgASjeDKSbO5cMEugU+ibNvZHVqY80J54x6JJkG7XWKEkwxyoXdJDSG8c09Z0X9RmCzUxx5B1PMPZcTE+gmmHuuLlmfef7m9MtGHasLk3cwxSgftsVhAOaeq4y55/JiAh2DT+7ZN3G2xKBRFkygI/BJnH3DJ3E2yo8JtOf45p59wtxz+TGBDmHuubqYezbKhgm0x5h7ri7mnquBCXQDPrln38TZJ/ds4lwdTKADfBJn3/BJnI1qYQLtIb65Z58w91wtTKAx91xlzD0bZaalQIvIZBH5TxHZLCJPichVEX1ERP5ZRLaJyJMiMqMzyzXapfDueRNwE/CZ4M9N+S6nTJh7LhYico2IqIiMb2i7LtDJp0XkvFZjuDyL4yDwcVXdKCJHARtE5EFVbbSdHwCmB1+zgH8N/iw8PrnnUojzfcBrwfXu4BpSFXv1yT2bOBcLEZkMvA94vqHtVGA+8DbgjcBDIvIWVR2MG6elg1bVnaq6MXj9W2AzcFKo21zg61rjMeBYETkx4X9T1/FJnEvBwxwS5zqvBe0J8UmcjULyT8AnYMSTaOcCK1V1v6o+C2wDZjYbJNHT7ERkKnAmsDb01knACw3XO4K2naHPLwIWBZf7xy668mdJ5i8I44FX8l5EQkqx5rPgrMg3dsOGz7Ah2WiLWnfpDLn8XX+pvY+X4fvjTe0O8OzmTd9fOGPy+NY9AXidiKxvuF6uqk6/9UXkYuCXqvqEyIhiaCcBjzVc13UyFmeBFpEjgW8DV6vqnvDbER8JP8Oc4D9weTDeelXtd52/KJRx3bbm7lHGdZdxzWlQ1fOzGktEHgJOiHhrKXA98P6oj0Utq9k8TgItImOpifMKVb0nossOYHLD9STgVy5jG4ZhlA1VnRPVLiKnA9OAunueBGwUkZmk0EmXUxwCfAXYrKpxdyJWAx8KTnO8E9itqjtj+hqGYVQSVd2kqhNUdaqqTqUmyjNU9UVqOjlfRMaJyDRqhyrWNRvPxUGfDVwObBKRx4O264EpwYJuAR4ALqC26b0P+IjDuGW9i1PGdduau0cZ113GNZcOVX1KRO4Gfk7tdNySZic4AES16RaIYRiGkROWJDQMwygoJtCGYRgFpesC7RIdLxoi8joRWSciTwRr/kzea3JFRHpF5L9E5P681+KKiGwXkU0i8njoLGphEZFjReRbIrIl+N5+V95raoWInBz8Hde/9ojI1XmvyzhE1/egg4ThiY3RceBPQtHxQhGcZDlCVfcGRw5/DFwVpCYLjYh8DOgHjlbVi/Jejwsish3oV9WihyeGEZHbgTWqepuIHAYcrqq/yXtdrohIL/BLYJaqPpf3eowaXXfQjtHxQhFE2PcGl2ODr8LfXRWRScCFwG15r6XKiMjRwLupHUdFVQ+USZwDZgPPmDgXi1z3oJtExwtHsFXwOLALeFBVC79mas+D+wQwlPdCEqLAD0RkQ/B4gKLzP4CXga8G20m3icgReS8qIfOBu/JehDGS3AS6RXS8cKjqoKq+g1r6Z6aInJb3mpohIhcBu1Q14TMsCsHZqjqD2lMSl4jIu/NeUAvGADOAf1XVM4HfAdfmuyR3gi2Zi4Fv5r0WYyS5CLRDdLywBP90fQTILNffIc4GLg72c1cC7xWRb+S7JDdU9VfBn7uAVbR44lcB2AHsaPhX1beoCXZZ+ACwUVVfynshxkjyOMXhEh0vFCJyvIgcG7x+PTAH2JLvqpqjqtep6qQgbjof+KGq/lnOy2qJiBwR3Dwm2CZ4P1Dopx4GMd4XROTkoGk2tbRYWViAbW8UkkSPG82IyOi4qj6Qw1pcORG4PbjT3QPcraqlObZWMiYCq4IHzYwB7lTV/8h3SU78DbAi2C74b9wed5A7InI4tQfLX5n3WozRWNTbMAyjoFiS0DAMo6CYQBuGYRQUE2jDMIyCYgJtGIZRUEygDcMwCooJtGEYRkExgTYMwygo/x+4ZgRgkLypkQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "#plot boundary to separate data points\n", "u = np.linspace(np.ndarray.min(x_train[0:]), np.ndarray.max(x_train[0:]), num=2)\n", "v = np.linspace(np.ndarray.min(x_train[1:]), np.ndarray.max(x_train[1:]), num=2)\n", "z = np.zeros((len(u), len(v)))\n", "for i in range(0, len(u)):\n", " for j in range(0, len(v)):\n", " z[i, j] = theta[1]*v[i] + theta[0]*u[j]+theta[2]\n", " \n", "print(z)\n", "plt.contourf(u, v, z)\n", "plt.colorbar();\n", "\n", "# plot_x = [np.ndarray.min(x_train[1:]), np.ndarray.max(x_train[1:])]\n", "# plot_y = np.subtract(np.multiply(-(theta[2][0]/theta[1][0]), plot_x), theta[0][0]/theta[1][0])\n", "# plt.plot(plot_x, plot_y, 'b-')\n", "\n", "plt.plot(x_train[0, pos], x_train[1, pos], 'ro')\n", "plt.plot(x_train[0, neg], x_train[1, neg], 'bo') \n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Lets try the gradient scheme" ] }, { "cell_type": "code", "execution_count": 107, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.0833254467443032\n", "0.3967390410267632\n", "0.253880339708813\n", "0.25048477189513413\n", "0.24875367042035454\n", "0.24703230127762024\n", "0.24531647830518077\n", "0.24360661720059154\n", "0.2419031435733602\n", "0.2402064686980325\n", "0.23851698862620221\n", "0.2368350840079264\n", "0.23516111998892622\n", "0.23349544615354997\n", "0.23183839651017277\n", "0.2301902895168571\n", "0.22855142814505544\n", "0.22692209997909155\n", "0.22530257734910977\n", "0.22369311749517018\n", "0.22209396276014948\n", "0.22050534080913142\n", "0.21892746487298323\n", "0.2173605340138604\n", "0.21580473341043416\n", "0.2142602346606768\n", "0.2127271961001355\n", "0.21120576313367662\n", "0.20969606857877435\n", "0.20819823301850546\n", "0.20671236516249464\n", "0.205238562214157\n", "0.2037769102426727\n", "0.20232748455822058\n", "0.2008903500891137\n", "0.19946556175954258\n", "0.19805316486676203\n", "0.19665319545662385\n", "0.19526568069645572\n", "0.19389063924438601\n", "0.19252808161428037\n", "0.19117801053555344\n", "0.18984042130718579\n", "0.18851530214536463\n", "0.18720263452422098\n", "0.18590239350921883\n", "0.184614548082806\n", "0.18333906146199952\n", "0.1820758914076327\n", "0.18082499052504308\n", "0.17958630655603255\n", "0.1783597826619659\n", "0.17714535769793052\n", "0.17594296647789862\n", "0.17475254003088403\n", "0.17357400584811292\n", "0.17240728812124492\n", "0.17125230797172394\n", "0.17010898367135255\n", "0.16897723085419317\n", "0.1678569627199479\n", "0.16674809022894221\n", "0.1656505222888921\n", "0.16456416593362194\n", "0.16348892649390875\n", "0.16242470776066062\n", "0.1613714121406106\n", "0.16032894080474347\n", "0.15929719382965482\n", "0.1582760703320574\n", "0.15726546859664645\n", "0.15626528619753827\n", "0.1552754201134938\n", "0.1542957668371465\n", "0.15332622247843722\n", "0.1523666828624722\n", "0.1514170436220143\n", "0.15047720028480005\n", "0.1495470483558989\n", "0.14862648339529821\n", "0.14771540109091502\n", "0.14681369732722288\n", "0.14592126824967191\n", "0.14503801032508984\n", "0.14416382039823333\n", "0.14329859574465875\n", "0.142442234120085\n", "0.14159463380639192\n", "0.14075569365443052\n", "0.13992531312377457\n", "0.13910339231956764\n", "0.13828983202660788\n", "0.13748453374079525\n", "0.13668739969807733\n", "0.13589833290102055\n", "0.135117237143117\n", "0.13434401703095225\n", "0.1335785780043396\n", "0.1328208263545274\n", "0.1320706692405822\n", "0.13132801470404334\n", "0.13059277168194633\n", "0.12986485001829992\n", "0.12914416047410615\n", "0.1284306147360019\n", "0.12772412542360304\n", "0.12702460609562358\n", "0.12633197125484505\n", "0.125646136351995\n", "0.12496701778861388\n", "0.12429453291895595\n", "0.12362860005099698\n", "0.12296913844659543\n", "0.12231606832086513\n", "0.12166931084080908\n", "0.12102878812326151\n", "0.12039442323218734\n", "0.11976614017537782\n", "0.11914386390058782\n", "0.11852752029115265\n", "0.11791703616112013\n", "0.11731233924993759\n", "0.11671335821672114\n", "0.11612002263414523\n", "0.11553226298197716\n", "0.1149500106402896\n", "0.11437319788237457\n", "0.11380175786738604\n", "0.11323562463273583\n", "0.11267473308626291\n", "0.1121190189982005\n", "0.11156841899296063\n", "0.11102287054075319\n", "0.11048231194905993\n", "0.10994668235398021\n", "0.10941592171146204\n", "0.1088899707884366\n", "0.1083687711538692\n", "0.10785226516973838\n", "0.10734039598196037\n" ] } ], "source": [ "lmda = 1e-2\n", "theta = np.random.rand(x_train_i.shape[0],1)\n", "#iterator 500 steps\n", "for x in range(0, 140):\n", " h = sigmoid(x_train_i.T.dot(theta))\n", " error = y_train-h\n", " tmp = (-1)*y_train*np.log(h) - (1-y_train)*np.log((1-h))\n", " J = 1/n*np.sum(tmp)\n", " np.fill_diagonal(D, h*(1-h))\n", " step = 1/n*(x_train_i.dot(error)-lmda*theta)\n", " print(LA.norm(step))\n", " if LA.norm(step)<1e-2:\n", " break\n", " #update theta\n", " theta = theta + .1*step" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Speed up with Newton?" ] }, { "cell_type": "code", "execution_count": 111, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.007267524438509921\n" ] } ], "source": [ "lmda = 1e-2\n", "#iterator 500 steps\n", "for x in range(0, 100):\n", " h = sigmoid(x_train_i.T.dot(theta))\n", " # print(h)\n", " error = y_train-h\n", " tmp = (-1)*y_train*np.log(h) - (1-y_train)*np.log((1-h))\n", " J = 1/n*np.sum(tmp)\n", " np.fill_diagonal(D, h*(1-h))\n", " H = 1/n*((x_train_i).dot((D.dot(x_train_i.T)))+lmda*I)\n", " dJ = 1/n*(x_train_i.dot(error)-lmda*theta)\n", " step = np.linalg.solve(H,dJ)\n", " print(LA.norm(step))\n", " if LA.norm(step)<1e-2:\n", " break\n", " #update theta\n", " theta = theta + .3*step" ] }, { "cell_type": "code", "execution_count": 112, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(3, 1)\n", "[[-16.95298391 31.84029078]\n", " [-41.18466395 7.60861074]]\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWgAAAD8CAYAAABaZT40AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO2df5hV1XnvP+8MP4xKtBOEoEDgWkzqr5sABXt5TG0g0SiKeps8oDE2sUXy0FZt8lgVa6qV3CRN1KS1Ipq0xChcbqoBiV4jpDQmzwUEajQGE6EqEImIE0G0gsy894+9z3Bmz97nrL3PPmf/OO/neeaZs9dZe601A/Od76z1vvsVVcUwDMPIHx1ZL8AwDMMIxwTaMAwjp5hAG4Zh5BQTaMMwjJxiAm0YhpFTTKANwzByigm0YRhGAkSkU0T+Q0RW+dddIvK4iDzvf/6dRucwgTYMw0jGVcCWquvrgDWqOgFY4183hAm0YRhGTERkNHAecG9V8yxgif96CXBho/MManSApAw7tkuPO370gPbOQd0ZrKY1vPzGa1kvoWl8QEdkvYSG2fzScZHvTTppRwtX0iLe2J31CpzZtIs9qhr9D+TAxFFH676DPU59t3W//SzwdlXTYlVdXHV9B3AtMKyqbaSq7gJQ1V0i0vAPRWYCfdzxo1l4/yP92o7+naUZrab53Pzv3+EYGt6SyiXrDv5F1ktIhROvm8P27mED2seO7GbD3bdksKLm0rH2G1kvwRm5mZcaHWPfwR5uO3ucU99ZS597W1Unh65FZCawW1U3ichZja6rFrbFYTREWcQZ4NaLNnDkkHf6tR059CAL/3RVRitqHkUS5xwyDbhARF4ElgEfEZHvAq+IyCgA/3PDf6LkRqDL7p6N/DNn6jYWXfZjxna9gaCMHdnN3V9YxiUzNme9NCNHqOr1qjpaVccBs4EfqeqngJXA5X63y4EVjc6V2RZHNSbOxaRM7rnCnKnbmDN1Gx0TfpH1UpqGueem8WVguYhcAWwHPtHogLkQaMPIEybOhiuquhZY679+DZie5viZb3GYey4mZXTPhpE3MhfosmLiXEzMPRt5IlOBLrN7NgzDaJTMBLrMCSnmnouJuWcjb9gWh2Fg4mzkExPolDH3bBhGWphAp4iJczEx92zkFRNowygpJs7FxwQ6Jcw9F5Myu2ej+FgmoVGTMonz0vUncuNDU9jRfTRjuvaz8HMPccmErFfVHMw9lwMT6BQos3suC0vXn8i8+z7MWwcHA7C9exhXfm02gD0MycgttsXRIGUW5zK55xsfmtInzhXeOjCEBffOzGhFzcPcc3kwgTbagh3dR4e37y5XEQUT53JhAt0A5p6Lw5iu/eHtI37b4pUYhjsm0MYAyibO0B7VUsw9lw8T6ISU2T2Xkb5qKSO7EbFqKUXiyUn5CoUUkSNEZIOI/ExEnhWRm/32LhF5XESe9z83vH9mURwJKLM4l9E9V7j0Uw9z6aceznoZTcHcc0s5AHxEVfeLyGDgJyLyKHAxsEZVvywi1wHXAX/dyETmoI22oMwJKWUW57y5ZwD1qBxqDPY/FJgFLPHblwAXNjqXCXRMzD0bhiEinSLyFF7l7sdVdT0wUlV3AfifRzQ6j21xxMDEOR8EMwJvvWgDc6Zui+xv7rmYpO2ejx0/ggseuMat89Irh4vIxqqWxaq6uHKhqj3AB0XkWOAhETk11cX6mEAbhSIsI3DefR8GqCnSRrHIwdbGHlWdXK+Tqr4uImuBc4BXRGSUqu4SkVF47rohbIvDEXPP+SA0I/DgYG58aEpof3PPRtqIyHG+c0ZE3gXMAJ4DVgKX+90uB1Y0Opc5aKNQRGYEhrSbOBeTHLjneowClohIJ57JXa6qq0Tk/wHLReQKYDvwiUYnMoF2wNxzfhjTtZ/t3cNC2w2jFajq08CHQtpfA6anOZdtcdTBxDlfhGYEDnmHWy/a0K/N3HMxKYB7binmoI1CUTkIjBPFUSZMnNsLE+gamHvOJ3OmbmvbsDqjvXDe4vADs/9DRAY8XUY8vikiW0XkaRGZmO4yjTQpsjjXo8zibO65/YizB30VsCXivY8DE/yPucBdDa4rc8rsno14PLB6IuNn38Sgj9zO+Nk38cBq8x9Ga3ASaBEZDZwH3BvRZRbwHT9HfR1wrB+oXUjKLM7mnuPxwOqJXPm12Wx/pQtVYfsrXVz5tdktF2lzz+2Jq4O+A7gW6I14/wRgR9X1Tr/NMArNgntn8taBIf3aWl0qy8S5fakr0CIyE9itqptqdQtp05Cx5orIRhHZuPe1N2Mss3WYey4mzdp7jiqJVbZSWUY+cXHQ04ALRORFYBnwERH5bqDPTmBM1fVo4OXgQKq6WFUnq+rkY95zVMIlG0kwcU5GVEmsVpXKMvfc3tQVaFW9XlVHq+o4YDbwI1X9VKDbSuDTfjTHGcDeymP3ikSZ3bORjIV/uoojhx7s11a2UllZYOLsRuI4aBGZB6Cqi4BHgHOBrcBbwGdSWV0LKbM4m3tOTqUk1oJ7Z7Jj9+8wZsRvWfinq1pSKqvM7tlwI5ZAq+paYK3/elFVuwLz01yYYdSjVTHPl8zY3PLahWUWZ3PP7tizODD3bBhGPmn7VG8T52Li4p7n33Ex9zw8jZ7eDjo7evmz83/KnVc/2ILVNYa5Z6NC2wu0UU7m33Exi1acSSUCtKe307+mECJdRkyc49PWWxzmnouJi3u+5+FpDAzPF789v5TZPRvxaWuBNoqH68FgT2/4f+2o9jxQZnEuk3sWkTEi8m8iskVEnhWRq/z2LhF5XESe9z83nM2U3/+tTcbcc7np7Ah/KkFUu2HE4BDweVX9PeAMYL6InAxcB6xR1QnAGv+6IdpSoE2ci0mcsLo/O/+nDHzagPrt+cPcc3FQ1V2qutl//QbeUz5PwHto3BK/2xLgwkbnskNCo5RUDgKLEMVh4twa3tw7jCcfPtO1+3AR2Vh1vVhVFwc7icg4vPqE64GRlQxqVd0lIiMaW3EbCrS552KSJCnlzqsfzKUgG4Vgj6pOrtVBRI4G/hW4WlX3iYQ9M64x2nKLo4yYOBcTc8/FREQG44nz/apacQGvVJ6D73/e3eg8beWgy+yey8zHPn8lazZ/oO96+sTn+OHX727KXA+snpjJczfaka61qxh93+0M2fMbDg5/Lzsvu4bus2aG9pkEkzJa5gDEs8rfArao6m1Vb60ELge+7H9e0ehcbeOgyyzOZXbP5yw60xdn6ftYs/kDfOzzV6Y+V6urp7Sze+5au4pxd97E0Fd3IaoMfXUX4+68ia61q0L75IxpwGV4j15+yv84F0+YPyoizwMf9a8boq0ctFE8DotzNdLPUadFreopabvodhZngNH33U7ngbf7tXUeeJvR993e56LD+uQBVf0J4UVKAKanOVdbOGhzz8Wk1XvPVj2ldQzZ85u67VF92om2EOiyYuKcLq2qntLu7hng4PD31m2P6tNOlF6gy+ye24HpE58jLOHEa08Xq57SGHGiNnZedg09Q4/o19Yz9Ah2XnZNzT7tRqkFuszi3C7u+Ydfv7tKpL2PZkVxXDJjM3d/YRljR3Yjoowd2c3dX1iW6v5zmd1zHLrPmsmL82/hwHGjUBEOHDeKF+ff0i+Ko7pPu2KHhEauCNvaaFZIXRjNrJ5SZnFOEvPcfdbMAWF1UX02XXDypqRrKzKlddDmng3DKDqldNAmzsUk6mDQJXkkrA9kU+w1DHPPRhJKKdBGeagkj1TikyvJI3C44nZYnyu+cgmq8E7PoMj7jMYxcW4updviMPdcTKLcc63kkVp9Dh4a1CfOUfe1ijK7Z6O5lE6gjeJRK+bZJXkkTiJJq5NOyizO5p6bT6kE2txz+XBJHomTSJJ20olhNJPSCLSJczGplzHokjwS1mfIoEMM7jxU875mY+7ZaBQ7JDRyTeVAr1Y0RlSfevc1ExNnIw1KIdDmnouJ6/M2XJJHovpYxIZRZEqzxVFGTJyLiblnIy0K76DL7J6LxNL1J3LjQ1PY0X00Y7r2c+tFG5gzdVvd+4IJJuee8XMeWXdqLpJLDOAZYA2wFzgGuq5YVTc9ux0QkW8DM4Hdqnqq39YF/G9gHPAi8ElVbehUutAOusziXCT3vHT9icy778Ns7x6GImzvHsa8+z7M0vUnhvavuOewCiaLVpzZsoomzaBU7vkZ4GE8ccb7HKx60sb8C3BOoO06YI2qTsD7tXZdo5MUWqCNfHDjQ1N46+Dgfm1vHRzMjQ9NqXlfWIJJsFBFVsklSSiVOIMnMe/0b6pUPWl3VPXHQHegeRawxH+9BLiw0XkKu8Vh7jk/7Og+2rm9eu/ZNWnEKppkxN7w5qJWOjnqzWf5/U0nu3YfLiIbq64Xq+riOveMVNVdAKq6S0RGJFlnNYUV6LJSNHEGGNO1n+3dw0LbqwkeDI4Z8Vu2v9JVf/wCJJeUzj0DHEOoSLdJpZM9qjo560UUcoujzO65iNx60QaOHNL/b+Ejh7zDrRdtqHlfWIJJsHqKVTTJkOnUrXpi9OMVERkF4H/e3eiAhRPoMotzEd0zwJyp21h02Y8Z2/UGgjK26w0WXfbjflEcYWF1YRVM5s16oqkVTZpBKd0zwGnUrXpi9GMlcLn/+nJgRaMD2haHkQpzpm6LDKurFfMcnmDyYIoray6lFWcOxzybIA9ERJYCZ+HtVe8Evgh8GVguIlcA24FPNDpPoQTa3LNhGHlAVedEvDU9zXnqCrSIHAH8GBjq9/+eqn4x0OcsPDv/gt/0oKrekuZCTZyLSdyMwfl3XMw9D0+jp7eDzo5e/uz8nwL0azvrg7/i+V+PjJ3M4lKZJQ7t4J6NbHFx0AeAj6jqfhEZDPxERB5V1XWBfk+oqv0tZCRm/h0Xs2jFmVRioXt6O/1r+rWt2fyBvmvXSikulVkMDxPn/FD3kFA9KvFSg/0PrXFL6ph7LiZx3fM9D08jmKjiXYe1HcYlmcWlMkscyuyejfzgFMUhIp0i8hRe2Mjjqro+pNsfiMjPRORRETklYpy5IrJRRDbufe3NBpZt5J0kD0Pq6U0eVFQvmcWlMosrZRZnc8/5wuknQlV7VPWDwGhgioicGuiyGXifqv534B+A70eMs1hVJ6vq5GPec5TTAs09tw+dHb2J762XzOJSmcUw8kYsy6KqrwNrCTwkRFX3VbZBVPURYLCIDG90cSbOxSTpo0S9A8Hg7plGtB3GJZnFpTKLC+aejVZSV6BF5DgROdZ//S5gBvBcoM97RUT811P8cV9Lf7lGmbnz6geZN+sJOjt6AKWzo4d5s54Y0DZ94nOxk1nCkmLiJsGYOButxiWKYxSwREQ68YR3uaquEpF5AKq6CPhj4HMicgj4L2C2qjZ0kGjuuZg0+iD+O69+kDuvHpioEtYWF5fKLIaRJ+oKtKo+DXwopH1R1et/BP4x3aWVExPnYmLu2ciCXGYSltk9l5Wl60/kxr/545pJIGFJKEFnHNZn2qkvJqq64jJfagQqjzAdOK1OnwnA83XuMdqa3Al0mcW5rO556foTmffdM2smgdRKQqmIZlSfe1dN41BPZ9/Y1X2iEk5c5nOlrnuuVB6pPNBvr38NhwU3rE/104bD7mkB5p7zTeGeZmfkjxsfmlI3CSQqCcVrr92nIs7VbbXmcp3PBaetjZDKI7zjt9fqEyR4T5Mxcc4/uRJoc8/FZMdvIyqqVCWBRCWhVLenmajiMl9qRFQe6dce1cd1LKMtyY1AmzgXk44Jv3BKAolKQqluTzNRxWW+ejgfDB7j0B7Vx3WslDH3XAxyI9BGcXFJAolKQqk8ra5Wn0GdPQPaas3lOl9qTMd7Qk01g+n/4MmwPkGC9zQJE+fikAuBNvdcTCphdS5JIFFJKNUHdlF9/vm6+2NXXXGZr+bXFies7jTgfA6732P869Pq9Jlc5x4jt4jIOSLySxHZKiLXNW2eBvNJEjPh9NF628qrABPoImIxz8WkqO55ygUnb2q0iOvk40U3znXrKzcTOZ+ftPcr4KPATuBJYI6qpv7NzdxBmzgbhlEwpgBbVfU/VfUgsAyY1YyJMo2DNnEuJs10z2FVT4BUK6HUoqXu+QfAJrytcgEmAec1b6yiuufUGDaC3rOiKlUFuPkbw0WkOlJ9saou9l+fAOyoem8nMDWVNQbIXaKK0b6EVT254iuXoArv9Azqa2tWJZSWi3P1j79WXccVaYex2l6c47OnxpZKMMAemlTEJLMtjpffKO/D7sw9JyOs6snBQ4P6xLlCI5VQcsOmmO2tGstwYScwpup6NPByMybKfA+6bJg4JydOdZMklVBq0fKDwSi/lcSH1RnL3HPqPAlMEJHxIjIEmA2sbMZEJtBGbohT3aTwlVDC/kiu1d6qsYy6qOoh4M+Bx4AteI9gfrYZc5lAp4i558YIS3gZMugQgzsP9WtLUgmlFpmE1U2K2Z5wLHPPzUFVH1HVk1T1RFVd2Kx57JDQyA2VQ79WRnFkFvNcOQhMI4ojYqwnbzJxLjom0Clh7jkdoqqelLISynkkD6tr5lhGbrAtjhQwcS4mljFo5B1z0EZdwpJHgo42rM9Pfz4udkUTl7kyx6V6SlKCCSfjgO54cz056Rd0rV3F6PtuZ8ie33Bw+HvZedk1AP3aXp/8hxy78d/79ek+q3/4Ytg4wT5G8zCBbpCyu+ew5JFgokhYn898+VIO9XQQp6KJy1ypfW1J3bNL9ZSkhCWcvFB17ThX19pVjLvzJjoPvA3A0Fd3Mf6bC0CVjp5DfW0jHl3WF+gx9NVdjLvzJoA+AQ4bJ9jHaC62xWGEUtnaCEseCSaKhPXxqqDEq2jiMlcaNLS14VI9JSkuiSV15npy0i8Yfd/tfaJaoePQO33iXCH4r9N54G1G33d733XYOME+RnMxgW6AMrvnClEJIdXtcZJGalU0cZkrc1yqpyTFNUmlzlxD9vwm8RKq740ap5HxjXiYQCekzOJcfTDoUi0lTtJIrYomLnM1SsMHgy7VU5LimlgSMVflYPDg8PcmXkL1vVHjNDK+EQ8TaKMmLtVSwvp4VVDiVTRxmStzXKqnJMUlSSViruqojZ2XXUPP0CP6vd87aDC9nf2PnIL/Oj1Dj+g7TIwaJ9inFfz++U+0dL48YYeECWgX9wzRySPVh3ZRfeJGcbjM1dDXlkZYXeVwrhlRHGEJJ+OIHcVROcBrNIojahw7IGwdmVVUOWr8u/T3/vZ3M5m7Ucoq0BbzXEzKHPNccc8iVzZeUeX9I3XD3W7Pg+78o280PF8a2BZHTMoqzoaRN9p5a6OCbXHEoMzi3Kh7bmUyS1wK5Z5jJsFU3PNJN36Wdz+9rq993+ln8Ktbv92vr0sfI1+YgzYappJgsv2VLlSlL8HkgdUTa/b5ky9dyqIVZ9LT68VMV5JZ5t9xcWprK5w4P8zhMLpKYsoz4d2D4izQ9/Hup9dx0o2f7evr0idPmHv2MIF2xNxzNEmTWXo0fjJLqUmYBFMR3moqAhynj5E/TKAdMHGuTSuTWeJQKPcMsZJg2uFgsKiIyCdE5FkR6RWRyYH3rheRrSLySxE5u95YJtBGw7QymaXUNDMJpiAUXZx9fg5cDPy4ulFETsYrj3UKcA7wTyLSWWsgE+g6mHuuT9Jklk6Jn8ziSuHcMzgnwQTd877Tzwj5LnrtcfoY6aCqW1T1lyFvzQKWqeoBVX0B2ApMqTWWCbTRMJfM2MzdX1jG2JHdiChjR3Zz9xeWDUhmCfb5lxvuZ96sJ+js8IS6s6OHebOeaDiKo5DiDF60xvkcdszH+NdVURxhWxu/uvXbfQJc+QhGaLj0yZqcuefhIrKx6mNuCmOeAOyout7pt0ViYXY1MPfsTlQllHp9LpmxOfWwukJzGomyEl2ENk9inAV64Ah6nz/ZtfueWokqIrIaCHsoyQJVXRF1W9iyai3CBDoCE+diUlj37IAdDOYHVZ2R4LadwJiq69HAy7VuqCvQInIE3mb3UL//91T1i4E+AnwDOBd4C/gTVc1ZGQzDlWBCybln/JxH1p3a8PMx5t9x8YCklGmnvtiUuQC3pI8kFUyC90wCxgbm6gJepPGCsN8BXvCGmczJ7Dv9DPbMuHjA8zGO3rKZEY8th95e6Ohg36lTOGLXS7GfoRGsoOJSdcWFWhVehu7ZBWO7YOGFcOnU2GMXiJXAAyJyG3A8MAHYUOuGus/i8MX3KFXdLyKDgZ8AV6nquqo+5wJ/gSfQU4FvqGrN73Sen8XRzu45WNXEo6IyHkcOPThgj7ke8++42K+oUv1XntIhvfRqZ7+2pHP1c8/ByifgHbhV7+kGK5iE4XqPUP95zpOJJ9K+OFdT2T+uPjzqFUFUA9/Z/t/pnqFH8OL8W2qKa7CCStJxXMbt7RwEInQcqvoHOnIILP5UqEin8SyOSePep+sXLHDqO3hu8vlE5CLgH4DjgNeBp1T1bP+9BcBngUPA1ar6aK2x6h4Sqsf+yrr9j+B/xVnAd/y+64BjRWRUjK/JaAEuWxthCSXBrbMkVU685JOBqRL9xTn5XAO2NlySPpJUMIm6x+WZYy7zVfPCwCZh4A9tR0CcK/2qcamEElZBJck4LuN29BzqL84Abx2EBd+PNXYeUdWHVHW0qg5V1ZEVcfbfW6iqJ6rq++uJMzhGcYhIp4g8BewGHlfV9YEuTqeTIjK3cip66I0el6lbTpndswuuCSVxq5w0knySqKKKS9JHkgomjTz8MZsHR/ZRrxKKa6WUuBVVYvXf3h1r7LLj9FOjqj2q+kG8Te0pInJqoIvT6aSqLlbVyao6edCwmvHZmVBmcXY9GHRNKIlb5aSR5JN6c4UeDLokfSSpYOJ6TxiN3JsC9SqhuFZKiVtRJVb/sV2xxi47sWyNqr4OrMXLgqkm9umkkU/CEkqCv2uTVDnxkk8Gpkp0SPAvqcbnAtySPpJUMIm6x0V8XearYm9EcknwV12vSGi/alwqoYRVUEkyjsu4DO6EIQGTduQQ76DQ6KOuQIvIcSJyrP/6XcAM4LlAt5XAp8XjDGCvqu5KfbVNxNyzR1hCybxZT9RMQnHhzqsfDE1KWXLD/Q3NFRlW55D0wXl4B3cVcRVgfIJ7JgMXBe4bH9InZhRHVHLJC3/1VQ4cNwoV4cBxo3jhmq+w++Oz0Y4Or19HB/tOP6NfH5eDve6zZvLi/Fv63bf747Njj+MyLv98OXz7cnhfl/f9eV9X5AFhO+MSxXE6sAToxBP05ap6i4jMA1DVRX6kxz/iOeu3gM+oas3z8bxFcZRVoC3muZhYzHN/ihTFkSZ146BV9WngQyHti6peKzA/3aW1jrKKs2HkjaIlpGRN22cSllmc03LPLtVSWk1q7jmYzDIBeJ74BWFdxsGtz4EnRjF5zymRiSFhSR+tLOSa9fztRNsLtFGbYOJKpVoKkJlIpyrO1ckse+mfhFKpaAK1RdplnO/j7bX21OizAnoZxNAe7/hm6Ku7GHfnTcDhCtvBpI+wPs2kkfnNPcenrZ9mZ+65Pi7VUgpLWDJLEIeKJk7j9HJYnKPo8RI4qgkmhoQlfSRJHklK1vO3G20r0CbObrhUS2klqR4MRiWzxO3nOk5CqhM9opI+4iaPpLGWOPObe05G2wq04YZLtZTC4lqppF6/Jlc8qU70iEr6iJs8ksZaXOc3cU5OWwq0uWd3XKqltIrUw+rCklmChFQ0STROB16gag16Bw32HiJURTAxJCzpI0nySFKynr/dsENCoyaVg8CsoziaEvNcOfhrNIrDdZyIPrpPBjyCMypCovI6qyiKuPObe26MthNoc8/xcamWUlgSVjBJPE5In42Tnu137ZLxl2VYW9bztxNttcVh4lxMLGOwmJh7bpy2c9CG0Y9ggklYJZRgtRSXbQ+Hqiuv/I/ZnP5P01tSGaWZ2yBhc534dYeT0/vXe89/3t5dqooqIvL3eE9xOQhsw3v0xev+e9cDV+AFXf6lqj5Wc6x6z+JoFq1+Foe552LSVPccVnUldBH0f4RcsMpKkKiqK4FxgtVKejs6kd6eplRGSVIJxYWouTq/dUltsb1/Pcz9rveQ/grlqajyMeBHqnpIRL4CoKp/LSInA0uBKXglr1YDJ6lqZIR8W21xlBET5wZwSTCBgc/3rJe8ElU5JTBO8AmlHQFxDuuTtDJKs5JJouaqWxllwff7izOUqaLKD1W1knG0Du/xy+BVnlqmqgdU9QVgK55YR9IWAl1m92w0QCMJJrXubfIfpUkrozQjmSVyzHqVUaLez09FleGV6k/+x9yE43wWqJS2cqo8VU3p96DLLM7mnhvkGJKLdK0tVpcCsg3gUhll6KsDH8fejGSWqLnqVkYZ2wUvhYhxEyuqPCe7OWPIP7h231Nri0NEVgNh39AFqrrC77MArzjs/ZXbQvrX/J/SFg7aMEJxSTCBgT8l9ZJXoiqnBMYJ/mT2dnQ2rTJKs5JJQquluFRGWXih1y/ufTlBVWeo6qkhHxVxvhyYCVyqhw/6YleeKrVAm3suJi0LqwuruhJWCeVCaldZCRJVdeVCalYreeHq/9W0yijNOCAMzhWrMsqlU71+JayoIiLnAH8NXKCqb1W9tRKYLSJDRWQ8XqrShppjlTmKo6wCbeJcTCzmOTlpRHHE0ZxNf/LzRqI4tgJDgdf8pnWqOs9/bwHevvQh4GpVfTR8FI/S7kGXVZwNI29YQkp/VDXyt4CqLgQWuo5VSoEuszibe2Zgcolr1ZOkY90J7Km6Hs7AAm91xglzz2EJHpDdczaM/FFKgTaKRyxxDlYvcal6knSsoDjjX9/JYZGuM06UOAcrk4z/5gJQ7Xtof6urpSTB3HNzKd0hobnnkhOWXOJS9STpWEFxDmtPsKawBI+OQ+/UrahitBelEmgT52IS62AwKm45STxzWmPVGCfqYDBO0kirqqXExdxz8ymVQBttQFSCSJKqJmmNlWCcOEkjraqWEgcT59ZQGoE291xMYofVhSWXuFQ9STrW8Ih7q9sjxtl2xVcjpw5L8HCpqGK0F3ZIaGRGopjnsOolSaM4XMaaT/0ojpBxtl3x1ZoHe1GVScLa8nZAaO65dZRCoM09txlpVUFxHSsYUucwTvek+qIaVZkkb4JsZCb80IQAAA2FSURBVEfhtzhMnIuJZQwWE3PPraUUDtownHBJSgnrg8N9PhVxbmVFk1aRujiXtKJKmhRaoM09F5NM3LNLUkpYn+/jPcynp8Z9AcKSUPKecNJyghVVXur2rsFEuorCb3GUERPnJuCSTBLWp5fD4hx1n0/FPbeyokmrSN09l7iiSpoUVqDL7J6NJuCSlBInQaVG31ZWNCks+a+okgsKKdBlFmdzz03CJZkkToJKoG/1wWBUYkkeE05caMrBYFTllCZWVCkihRRow4iNS1JKWJ8OoLP2fcGojVZWNGk2TYvaKHhFlVZROIE291xMMg+rC6ueEqyMEtbnQrxazDEqqrSyoklhKXdFlb8TkadF5CkR+aGIHF/13vUislVEfikiZ9cbq9BRHGXCxLkFuCSlRPWpE1YXJCoJpUg0Peb50qmlEOQQ/l5V/wZARP4SuAmYJyInA7OBU4DjgdUicpKqBo+h+yiUgy6zezaMPGEJKclR1X1Vl0dxuPbvLGCZqh5Q1ReArcCUWmPVFWgRGSMi/yYiW0TkWRG5KqTPWSKy17f0T4nITe5fjhtlFmdzz8WkzBmDBsNFZGPVx9w4N4vIQhHZAVyK56ABTgB2VHXb6bdF4rLFcQj4vKpuFpFhwCYReVxVg/87n1DVYv9NZ6SKkzgnze5L8iwOlyzBCcDz9edKQ5zH3nULIx5bDr290NHB7rM/yfbPJfM2aWYultU9Hz/sPXzxDz/t1PcCrt1Tq2isiKwGwsJyFqjqClVdACwQkeuBPwe+yOE679XUrNpdV6BVdRewy3/9hohswVP9ltkHc88lJWl2X5ISV2HjrMD78eitattYdU8j5bTqMPauWxjx6LLDP7G9vYx4dBlAbJG2zMXWo6ozHLs+APwAT6B3AmOq3hsNvFzr5lh70CIyDvgQsD7k7T8QkZ+JyKMickqccWth4lxMnNxz0uy+JCWuwsbp4bA4RxEyVxruecRjywfYKfHb45Jm5mJZ3XMrEZEJVZcXAM/5r1cCs0VkqIiMx/t7bUOtsZyjOETkaOBfgasDm+AAm4H3qep+ETkX7wkGE0LGmAvMBRjynmDAqdF2NJLdl1ZZqmbfG0VvxG+GqPYapJW5aOKcGl8Wkffj/fp/CZgHoKrPishyvN2HQ8D8WhEc4OigRWQwnjjfr6oPBt9X1X2qut9//QgwWEQG1KJQ1cWqOllVJw8aFoz+H4i552LifDDYSHZfWmWpYt6b2sFgR8SPXlR7DcqWuVh0VPV/quqpqnq6qp6vqr+uem+hqp6oqu9X1UfrjeUSxSHAt4AtqnpbRJ/3+v0QkSn+uK+5fkFGeYgVtZE0uy9JiauwcTqp/xNQNVeaURu7z/7kgNMh9dvjkkbmornnfOKyxTENuAx4RkSe8ttuAMYCqOoi4I+Bz4nIIeC/gNmqWvN0sh7mntsAl5JTaZW4ihon2OYYxdEolYPANKI4ospn2QFh8ZEGdTQxR41/l/7e3/5u6HsmzsXEYp6LSRHcs8iVm2qFvbkw4fTRetvKAWkcoVww/tqG50uDQmUSGkYWmDgbWZG7Z3GYey4mid1zWkkoLvwA2IS32SvAJOC8Js1lGCmQO4EuKybOIaSVhOLCD+ifhKJV1zVE2tyzkSW52uIos3s2QkgrCcWFTTHbDSMH5EagyyzO5p4jSCsJxYWos/AaZ+Tmno2syY1AG21IWkkoLoQ9pqZGu4mzkQdyIdDmnotJw2F1aSWhuDApZrth5IDMDwlNnItJKjHPaSWhuFA5CHSI4jD3bOSFzAXaaHNcylClxXm0dVidiXPxyHSLw9xzMbGMQcNoDeagjVhUxPn+Z+awYM2X2L53LGOP2c7C6Tdw6WlLG5+gmYkrdcZutTinWQWlHuaei0lmAv0BHZHV1E2nzO4ZPHGe+/A9vPXOUQC8tHcccx++B6AxkW5m4kork2IcsCoohgu5iOIoE2UW54p7XrDmS33iXOGtd45iwZovNTZBMxNX6ozdavecZhWUeph7zgYR+YKIaPWz8UXkehHZKiK/FJGz641hWxxGbLbvHRur3ZlmJq60MinGgbSqoNTDxDkbRGQM8FFge1XbycBs4BTgeGC1iJxUq6qKOegUaQf3DDD2mO2hfaLanWlm4kqNsbM4GLQqKKXnduBa+ueqzgKWqeoBVX0B2ApMqTWICXRKtIs4AyycfgNHDn6zX9uRg99k4fQbGpuomYkrEWNvu+KrKQwenzSqoNTD3HM2iMgFwK9V9WeBt04AdlRd7/TbIrEtDiM2lYPA1KM4mpm4EjF2VgdyVgWl9fQc6mL/b+c49r52uIhUP/9wsaourlyIyGog7M+dBXgVpz4W8l7YgwVqVkwxgU6BdnLPFS49bWk6YXVBmpm4Ehg765jn7rNmWlhdftlTq6KKqs4IaxeR04DxwM/8Mq2jgc1+rdadwJiq7qOBl2stwrY4jLYka3FuJibO2aGqz6jqCFUdp6rj8ER5oqr+BlgJzBaRoSIyHq8C5oZa45mDbpB2dM/GQMbedUvdArCtTEwx8oeqPisiy4FfAIeA+bUiOMAEuiFMnItJ2u557F23MOLRZYc3GHt7GfHoMuBw9e5WJaaYe84Xvouuvl4ILHS937Y4DKNBRjy2fMDpj/jtFVqZmGKUBxPohJh7LiZN2Xvu7a3b3orEFHPP5cME2uiHiXMCOiJ+jKram52YUmZxXnlJ+/6VYQKdgDK7ZyM+u8/+5IBgVvXbK7QiMcUoH3ZIGJMyi7O552RUDgJrRXE0MzHF3HN5MYE2jBTY/rmbBoTVBWlGYoqJc7mxLY4YmHsuJmVOSjHKjQm0YeJcUMw9lx8TaEfK7J4NI0+YOB/GBNqBMouzuediUmb3bBzGBNooJSbOxcTcc39MoOtg7tkwjKwwga6BiXMxMfdcTMw9D8QE2jCMzDFxDscEOgJzz8XE3LNRJkygjdJg4lxMzD1HU1egRWSMiPybiGwRkWdF5KqQPiIi3xSRrSLytIhMbM5yW4O55wx5BrgDuNn//Ey2yzGMuIjI34rIr0XkKf/j3Kr3rvd18pcicna9sVyexXEI+LyqbhaRYcAmEXlcVatV7ON49bUmAFOBu/zPhcPEOUOeAR4G3vGv9/rXULeQrLnnYlJi93y7qn6tukFETgZmA6cAxwOrReSkWmWv6jpoVd2lqpv9128AW4ATAt1mAd9Rj3XAsSIyKtaXYxhrOCzOFd7x2w2j+MwClqnqAVV9AdgKTKl1Q6yn2YnIOOBDwPrAWycAO6qud/ptuwL3zwXm+pcHBs+98udx5s8Jw4E9WS8iJoVY8ySYFPrGXth0M5tq331yE1aUiEJ8rwMUYc3va3SAF7Y889glE8cMd+x+hIhsrLperKqLY0z35yLyaWAj3g7Eb/E0cV1Vn4pORuIs0CJyNPCvwNWqui/4dsgtwWeY43+Bi/3xNqrqZNf580IR121rbh1FXHcR15wEVT0nrbFEZDUQVg5nAd4W79/haeDfAV8HPoujTlbjJNAiMhhPnO9X1QdDuuwExlRdjwZedhnbMAyjaKjqDJd+InIPsMq/jK2TLlEcAnwL2KKqt0V0Wwl82o/mOAPYq6q7IvoahmGUlsD520VAZSt3JTBbRIaKyHi8oIoNtcZycdDTgMuAZ0TkKb/tBmAsgKouAh4BzsXb9H4L+IzDuHH2c/JEEddta24dRVx3EdecZ74qIh/E2754EbgSQFWfFZHlwC/wouPm14rgABDVmlsghmEYRkZYJqFhGEZOMYE2DMPIKS0XaJfU8bwhIkeIyAYR+Zm/5puzXpMrItIpIv8hIqvq984HIvKiiDzjp8lurH9H9ojIsSLyPRF5zv+//QdZr6keIvL+qnTkp0Rkn4hcnfW6jMO0fA/aP+EcVZ06DlwYSB3PFX4ky1Gqut8POfwJcJWfNZlrROSvgMnAu1V1ZtbrcUFEXgQmq2rekyf6EJElwBOqeq+IDAGOVNXXs16XKyLSCfwamKqqL2W9HsOj5Q7aMXU8V/gp7Pv9y8H+R+5PV0VkNHAecG/WaykzIvJu4MN44aio6sEiibPPdGCbiXO+yHQPukbqeO7wtwqeAnYDj6tq7teM9zy4a4HerBcSEwV+KCKb/McD5J3/BrwK/LO/nXSviByV9aJiMhtYmvUijP5kJtB1Usdzh6r2qOoH8bJ/pojIqVmvqRYiMhPYrap1nmGRS6ap6kS8pyTOF5EPZ72gOgwCJgJ3qeqHgDeB67Jdkjv+lswFwP/Jei1GfzIRaIfU8dzi/+m6Fkgtr79JTAMu8PdzlwEfEZHvZrskN1T1Zf/zbuAh6jzxKwfsBHZW/VX1PTzBLgofBzar6itZL8ToTxZRHC6p47lCRI4TkWP91+8CZgDPZbuq2qjq9ao6WlXH4f35+iNV/VTGy6qLiBzlHx7jbxN8jMOpsrlEVX8D7BCR9/tN0/GyxYrCHGx7I5fEetxoSoSmjqvqIxmsxZVRwBL/pLsDWK6qhQlbKxgjgYe83+MMAh5Q1f+b7ZKc+Avgfn+74D9xe9xB5ojIkcBH8dORjXxhqd6GYRg5xTIJDcMwcooJtGEYRk4xgTYMw8gpJtCGYRg5xQTaMAwjp5hAG4Zh5BQTaMMwjJzy/wE9UwCvbaMTgAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "#plot boundary to separate data points\n", "print(theta.shape)\n", "u = np.linspace(np.ndarray.min(x_train[0:]), np.ndarray.max(x_train[0:]), num=2)\n", "v = np.linspace(np.ndarray.min(x_train[1:]), np.ndarray.max(x_train[1:]), num=2)\n", "z = np.zeros((len(u), len(v)))\n", "for i in range(0, len(u)):\n", " for j in range(0, len(v)):\n", " z[i, j] = theta[1]*v[i] + theta[0]*u[j]+theta[2]\n", " \n", "print(z)\n", "plt.contourf(u, v, z)\n", "plt.colorbar();\n", "\n", "# plot_x = [np.ndarray.min(x_train[1:]), np.ndarray.max(x_train[1:])]\n", "# plot_y = np.subtract(np.multiply(-(theta[2][0]/theta[1][0]), plot_x), theta[0][0]/theta[1][0])\n", "# plt.plot(plot_x, plot_y, 'b-')\n", "\n", "plt.plot(x_train[0, pos], x_train[1, pos], 'ro')\n", "plt.plot(x_train[0, neg], x_train[1, neg], 'bo') \n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 134, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "text/plain": [ "LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,\n", " intercept_scaling=1, l1_ratio=None, max_iter=100,\n", " multi_class='auto', n_jobs=None, penalty='l2',\n", " random_state=None, solver='lbfgs', tol=0.0001, verbose=0,\n", " warm_start=False)" ] }, "execution_count": 134, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import numpy as np\n", "import pylab as pl\n", "from sklearn.linear_model import LogisticRegression\n", "\n", "logreg = LogisticRegression()\n", "\n", "# fit the model with data\n", "logreg.fit(x_train_i.T,np.ravel(y_train))" ] }, { "cell_type": "code", "execution_count": 135, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 3.38822960e+00 -3.16424935e+00 2.82754169e-04]]\n", "[-8.32412399]\n" ] } ], "source": [ "print(logreg.coef_)\n", "print(logreg.intercept_)\n", "theta=logreg.coef_\n", "theta=theta.T" ] }, { "cell_type": "code", "execution_count": 136, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(3, 1)\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWgAAAD8CAYAAABaZT40AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO2de5Qc1XWvvz0jaZCw0ANJBmmwRRI/riPr2niMZZNrEyGbRwCFPFgyQZg8Ljdcbi4icAWGC47BZiGCDdgmxlqEy9ssYiwLsAjmERljDLaEeRqS2EaQ0cjoCYhI6DGz7x9dPeqpqeo+VV3VXXV6f2v10vTpU+ccPWbrN7v2r7aoKoZhGEbx6Gr3AQzDMIxoLEAbhmEUFAvQhmEYBcUCtGEYRkGxAG0YhlFQLEAbhmEUFAvQhmEYCRGRG0Vko4g8Hxr/GxH5VxF5QUSubHYfC9CGYRjJuQk4pnZARH4fWAjMVdXfBa5qdhML0IZhGAlR1UeBraHhM4ErVHVXMGdjs/uMaXaBtEycPFWnz+xt1/ap6R4T/jvxg4HtW9p9hMz57fWDieY/t2du7GcfHPtss8fJjHGHTmn3ETLnP4em1/38pV++sFlV609qwGEHv0Pf3O32b+JXW99+AXi7Zmi5qi5vcNl7gf8mIl8Orj1PVX+W6rABbQvQ02f28uXbV7Vr+1S8Y8q3232EXPjiD29hEv5903/vwtcTzZ/32vdYPzhaNMzq7mfV9I9ndaymmHn7H7f7CLmwdvv/qPv54Sd+4JVm93hz9yBfPXq209yF337pbVXtS7jFGGAKMA/4KHCXiPyWNvE8DUtxGF6SNDgDnD/xCsbLjhFj42UH50+8IqtjGRE0Cs4loh/4rlb4KTAETGtmQQvQHc4Xf3hLu49QGE6asJJlk5Yyq7sfYYhZ3f0sm7SUkyasbPfRAH/Vs0d8D5gPICLvBcYBm5tZsG0pjrLha3rDR9Ko5yonTVhZmIDcCZRVPYvIt4EjgWki0g98AbgRuDEovdsNfK6Z9AZYgO5oTD2XB1PPxUJVPxvz0alZ7mMpDgdMPZeHZtSzYRQNC9Adiqnn8uCrei5reqOVWIBugKnn8mDq2fANC9AdiKlno92YenbDAnQdTD2XB1/Vs6/pDcMNC9Adhqlno92YenbHAnQMpp6NdmPq2bAA3UH4qp59TW/4iKnnZFiAjsDUs9FuTD0bYAHaKDlJ1POKHQuZ99pPeNfAK8x77Ses2LEwx5MZYUw9J8es3iF8Vc++pjdcWbFjIee/cSU7dQIA6wd7Of+NSkeioj17w9SzUcUUtFFakqjnZdsvGA7OVXbqBJZtvyDrYxkRmHpOhwXoDqDT1TPAwODMROPtwtSzUYsF6Bp8TW/4SNLKjZndA4nGjeww9ZweC9CeY+q5Qhm6pZh6NsLYTcIAU8/lIU3dc/VG4LLtFzAwOJOZ3QOcP/GKwt0gNIxaLEB7jKnnkRS5W4qv6tnSG81hKQ5MPZcJcw0anYQFaE8x9Wy0G5/Vs4jcKCIbg/6D1bG/F5GXRORZEVkhIpOb3afjA7Sp5/JQVc++OQJ9TW94zk3AMaGxB4E5qjoX+Dfg881u0vEB2kd8Vs9VR+D6wV6UrmFHYNmDtG/4rJ4BVPVRYGto7Aequjd4+wTQ2+w+HR2gTT2XD98cgaaeveUvgPubXcSqODzDV/VcTW+UxRHYyRRVPb+5+5080P+3jrPPmCYia2oGlqvqcpcrReQiYC9we9IzhunYAG3quZzM7B5g/eDonxzL6Ag09VxoNqtqX9KLRORzwPHAUaqqzR6io1McRjmoLa0rgyOwkymqem4FInIMcD5woqruaDTfhY4M0L6qZ1/TG7WcNGElyyYtZVZ3P8IQs7r7WTZpaWENKHGYei43IvJt4CfA+0SkX0T+EvgGMBF4UESeFpHrm92nY1McRjmIMqYU2RHYyXSSelbVz0YM/2PW+zgraBHpFpGfi8h9EZ+JiHxNRH4ZFGkflu0xjUZ0gnr2BVPPhitJUhxnAy/GfHYs8J7gdQbwzSbPlRu+pjd8pCi2bt+MMXnQSeq5lTgFaBHpBf4AuCFmykLgFq3wBDBZRA7O6IxGA0w950fWxhhTz0YSXBX0NcBSYCjm81nAf9S87w/GCoWp5/JQFPXsmzEmD0w950fDAC0ixwMbVXVtvWkRY6NqAEXkDBFZIyJrtm/bGnGJkRRTz/mSpTHG1LORFBcFfQRwooisA+4E5ovIbaE5/cAhNe97gVHOAVVdrqp9qto3ccrUlEdOh6nn8lAU9QzWKstoLw0DtKp+XlV7VXU2sAh4RFVPDU27BzgtqOaYB7yhqhuyP65Ri6nn/MnKGOOrerb0Rr6kroMWkb8GUNXrgVXAccAvgR3An2dyuoww9VweiqSewVplGe0lUYBW1dXA6uDr62vGFTgry4MZ9TH13DrMGBONqef88d7qbeq5PBRNPWeFr+kNI3+8D9A+YurZjQu3XcbsgZc5ZOBVZg+8zIXbLmv3kbzB1HNr8PpZHKaeO5cLt13GrTs/R7UCdJAxwXu4fMrFLTuHqWejGUxBlwxf1XPW6Y07dp7K6PJ8CcaNZjD13Dq8DdCmnjubQboTjeeBqWejWbwN0EZ5yOPmYDeDicYNN0w9txYvA7Sv6tnX9EYenDL+NkY/bUCD8fwx9Wxkgdc3CY3ik1dpXfVG4B07T2WQbroZ5JTxt7X0BqFvmHpuPRagS4Kp5+RcPuXitgRkU89GVniX4vA1veEjvhpTfMTUc3vwLkD7SKer50WbbuOQgVeHX4s25ZdHbrZ7iqlnd6auvo+5f3kUfQt/l7l/eRRTV4/qpjc85yPwkTYcMRYROUdEXhCR50Xk2yKyXx77eBWgTT2XB1f1vGjTbfx4zyep1DRXXj/e88lcgnTW3VOMeKauvo/Z111Cz6YNiCo9mzYw+7pLRgTp2jlFQkRmAf8b6FPVOUA3lSd9Zo5XAdpHOl097wvOtUgwni3Ndk/xVT3nkd7ovfVqune9PWKse9fb9N56dd05BWIMMF5ExgATiHj+fVabeIGp5/JQ1Nxzlt1TjPqM2/ybhuNxc9Jy4EHv4PT/83tOc//hR0wTkTU1Q8tVdTmAqq4XkauAV4GdwA9U9QeZHjbAFHSB6XT13Gqse8po8ro5uHvaQQ3H4+a0iM3V7k/Ba3n1AxGZQqVR9qHATGB/EcnlGQJeBGhTz+UhqXo+YuyjRBlOKuPZ0kz3FF/TG3nRv/gcBntG3lcb7NmP/sXn1J1TEBYAL6vqJlXdA3wX+EQeG3kRoH3E1HOFO6efWhOkK68jxj7KndOzFywnTVjJsklLmdXdjzDErO5+lk1a2rEP68+ztG7rkcez7qxL2TX9YFSEXdMPZt1Zl7L1yOMj5xSMV4F5IjJBRAQ4Cngxj41Kn4M29ew/eQTjONJ0TzH1nI6tRx4/IiDXm7P2xA+sbdGxGqKqT4rId4CngL3Az4Hl9a9KR+kDtI/4qp6LenPQGI0ZU+qjql8AvpD3PqVOcZh67gxczCNRc5o1nbhg6tnIE1PQRktIq56r5pFqfXLVPAL7Om5HzTn39a8AsIee2OuMaEw9F4fSKmhf1bOv6Y20uJhHoubsoWc4OMdd1yymno28KW2ANspDM7lnF/NIEiOJmU7qY+q5WFiALhCmnkfjYh5JYiTJynRi6tloBaUM0L6mN3yk2coNF/NI1Jyx7GIsu+peZ4zE1HPxsJuEBcHUczTVG3rLtl/AwOBMZnYPcP7EK0bc6Iub0+i6tJh6NlpF6QK0qefykFXds4t5JG6OVWwYZaaUKQ7fMPVcHnxVz5beKCalUtCmnovLlsePY+Dus9m95SDGHfgbrhq8HBzU64odC0ekIeaPe4hHdi/IPC1hpGPq6vvovfVqxm3+DbunHUT/4nMa2rON7ChVgPYRH9TzlseP45Wb/g7dPR6A3Vtmcr40NoVEGUxu3fk5qg/oN3NJa4hTz9WOJtWH5le7ngAWpFtEaVIcpp6Ly8DdZw8H5youppAog0m4e0rW5pJm8DW9EYdL1xMjX0oToH3EB/UMsHtL9IPVG5lCXE0jZi7Jj3q5Z5euJ0a+lCJAm3ouNuMOjP6GbWQKcTWNFKGjSaepZ3DremLkSykCtI/4op4BZv7xtci4nSPGXEwhUQaTcPcUM5fkR6PKDZeuJ0a+FP4moann4nPgJ1YBlVz0ni0HOVdfRBlMiljF0YnqGfbdCLQqjvZR+ABtlIMDP7GKH913R6WFZgKiDSYXZ3YuIxrXumeXridGfhQ6xeGrevYpveE7naqejWLQUEGLyH7Ao0BPMP87QbuX2jlHAiuBl4Oh76rqpdke1SgyWdm6L9x2GXfsPJVBuulmkFPG3wYwYmze2B+zbui3E6dBwqaYIqRP2oG5BrNBRLqBNcB6Vc3lxwyXFMcuYL6qviUiY4HHROR+VX0iNO9HeR3SJ0w9x3PhtstGGFUGGRO8Z8TYj/d8kqRmFpfOLGFMPRsNOJtKN+8D8tqgYYpDK7wVvB0bvLTOJZnga3rDR7JSz3fsPJWwUaXyPmpsH2lNMUUywbQKU8/ZICK9wB8AN+S5j9NNwkDKrwV+B7hOVZ+MmPZxEXkGGADOU9UXItY5AzgDYNpBs1IfuqyYeq7PIN2pr01riokbN/XsH/t3beIjE7/lOn2aiKypeb9cVZfXvL8GWApMzOp8UTjdJFTVQVX9ENALHC4ic0JTngLerar/Ffg68L2YdZarap+q9k2cMjV2P1PP5SEr9QzQzWDqa9OaYopggmkVpp4Tsbkaq4LXcHAWkeOBjaq6Nu9DJKriUNXXgdXAMaHxN6tpEFVdBYwVkWlZHdIHTD03pnJDMJw905ixfaQ1xcRdZ+rZaMARwIkisg64E5gvIrflsVHDAC0i00VkcvD1eGAB8FJozkEiIsHXhwfrbklzIFPP5SFL9Qxw+ZSLWTz+ZrrZCyjd7GXx+JtHjR0x9lFmdfcjDDGru59lk5Y6mWKWTVqa+DrDCKOqn1fVXlWdDSwCHlHVU/PYyyUHfTBwc5CH7gLuUtX7ROSvg8NeD/wJcKaI7AV2AotUNfcbiWXB1LM7l0+5mMunjDaqRI0lxaUzi6/q2dIb5aRhgFbVZ4EPR4xfX/P1N4BvNHsYU8/lIWv1bBhlRFVXU0n75oJZvXOmU9SziwkkyoQSVsZRcz7asybV8zpc9suMh15CbngcNm6HGRPRv/oELHh//TnzZiNPrKt/TQaYei4vhQnQpp7LQ1g9u5hA6plQqkEzbk4lyI4ZXtul64rLfmFSpzceegm56mFk197K+9e2w1UPV25lVgNu1JyVz+2r6I66xuh4Cv0sjrLTKerZxQQSZ0KpjNefMzhKRzQ2qrjslxVyw+P7Am91bNfeilquNye8TuiaLDD1XG4KEaBNPZcbFxNInAmldjxLo4rLfrU0dXNw4/bG43FzXNcyOpJCBGgf8VU9R90cdDGBxJlQasezNKq47JcZM2LMZLXjcXNc10qBqefy0/YAbeq5/LiYQOJMKNWn1dWbU6mBHjlWby/X/ao0W1qnf/UJtGdkGkZ7xlRu+tWbE14ndI1htD1AG+UhrrTOxQQSZ0KpvWEXN+fqyUtGrL14/M0NDScu+2XGgvej5x2FvnMiKlR+Pe+okTf7ouYs/GD9a5rA1LMfSLv8JL/1gbn65dtXeamgOym9UXbMmFIODj/xA2tVta+ZNfre90796bc+6zS3+/evbXq/LGirgvYxOPuKj8HZV3wLzp1MYeqgfcFX9dwqogwvQC6dUNqunq95BLn3eRhS6BL0hDmwZH7Ta/V1fZ2NR5/Mq2deku15jZbTtgDdPWZru7Y2EtIq9RxleDn39a8AsIee4TGXDiqF55pHkFqjypDCyucqNw6TBulRaw0x4/47ASxIlxy7SZghpp6bI8rwsoee4eBcJYtOKO1Wz3Lv89G9Y+59PrO1ZjxwV8rTGUXBArRRl1bmnht1RUk7t5AMxdycjxtPtdZQ8rUKyE1//1i7j9A2LEBnhKnn5knS3aSZTijtVs8AdIU1b4PxVGvZt3fZsb9BI5ZWV25EGV7Gsoux7Box5tJBpejoCXOie8ecEO4ml36tjUefnPJ0xaGT1TNYgM4EU8/ZEGV4+crkc/nK5HP964SyZH7FqNIllcDcJejCD6ar4lgyn43HLkK7uoK1uth47CK7QegBVmZnRNKuuue4ridZBeRCpDeqLJmPpi2rC/HqmZd4F5A7XT2DKeimMfVstBszpviLKWijIS7dUqLm/GxXX+KOJi57pSUz9ezSPSUtYfPKh2Yh699IvNfU1ffRe+vVjNv8G3ZPO4j+xecAjBh7ve9TTF7zwxFzth55fMN1wnPyoOjqWUSOAa4FuoEbVDWXmyIWoJvAV/Vcm95w6ZYSNeec168Jnr3s3tHEZa+249I9JS1R5pWn+ut2XYlSz1NX38fs6y6he9fbAPRs2sChX7sIVOka3Ds8NuP+O4fX7tm0gdnXVVIk1QActU54TicSNNC+Dvg00A/8TETuUdVfZL2XpTiMurh0S4maU+mCkqyjicteaclKPbt0T0m9dozhJOlevbdePRxUq3Tt3TMcnOPW7t71Nr23Xl13nfCcPCi6egYOB36pqr9W1d3AncDCPDayAG2MIHxz0KVbShLTSL2uKS57tR2X7ilpcTWpBHvF5Z7Hbf5N6iPUXhu3TjPrl4hpIrKm5nVGzWezgP+oed8fjGWOpThS4mt6I8zM7gHWD/ZGjjeaE0W9jiYue6Uh08qNGRMrqYao8WbpErcg3WCv3dMOomfThlRH2D3toIbr1M7JmjzV8+6XtzHwZ3e7Tt9c53GjUc6gXJ7bbAraGCaqtM6lW0rUnOqD8kcS3dEkyV7txqV7Suq1YwwnUXvVq9zoX3wOgz37jRgbGjOWoe76HV0Ge/YbvpkYt054TofSDxxS874XaE5FxGAKOgWdop5h3825epUVcXOSVnG47JWUzOueF7y/EtjyqOJYMr+ydqiKg6gqjjoZleoNvGarOOLW6eQbhAE/A94jIocC64FFwCl5bNS2jirvmdurX73n7Lbs3Sw+BmhfH8hfKGNKRvha91wvvfEPPzqj6Q4nc8eN01XT3+k095CB/rr7ichxwDVUyuxuVNUvN3O2OExBJ8TH4OwrPgZnoxio6ipgVd77WIA2MlHPrTSzeE1KE8x7/+9fcMCzTwy/f3PuPP7tSzcmntMuSlBa1xbsJmECTD1HUzWYrB/sRekaNpis2LGw7pwlr1/NrTs/N1wzXTWzXLjtsqbPVEr1XDXBvLYdUSq/XvUwPPTS8JSo9EY18AoMvw549gne+3//ItEco3hYgO5wslDPac0sQ4wlqZnFZ9KaYKqBd8R1wXiSOe3C1HM8FqAdMfUcTyvNLF7TwATj681BIx4L0B1MVpUbcUaSsJnFlXpmFqfzlDG9AfEGlCxMMAXF1HN9LEA7YOq5PmnNLF3sIamZxWfqmWDqqec3586LNLi8OXdeojlG8bAAbTRNVCeUcNeTqDnXTD6HxeNvHnYddrOXxeNvbqqKo7TqGSommPOOQt85ERUqv553VMMqjn/70o3DAbj6CldouMxpNaaeG2NGlQb4qp7NmFIefM09JwnQRTOqtApT0IY3+BicfcXUsxsNjSoish/wKNATzP+Oqn4hNEeodBc4DtgBnK6qT2V/XCMLGqnnsKFk/riHeGT3gqafj3HhtstGmVI+2rMml70AN9NHmg4m4WtOmANzZo7ca9Yk5On1I+ek6D+4/ewn6Xv22uH3b86dx+YFfzTq+RjvePEpZjxwFwwNQVcXb845nP02vJL4GRrhDiouXVdcCK/bN/U4AE5ct4Kpu7aytWcq98w+iTUzPpZ4bZ9pmOIIgu/+qvqWiIwFHgPOVtUnauYcB/wNlQD9MeBaVa37J12GFEcnpjfCXU0qKLX1yuNlR+LO2hduuyzoqFJbjat0sTeoh25ur1HqOdz5hOCGW21ON9zBZNTu7tcgIDpyLDwncdfuc+9GajuqsC9/XPuj75AIolr39zHYsx/rzrq0bnANd1BJu47LunvoRgTG6L6KnV1d47jjPYsjg7SlOGLQCm8Fb8cGr3BUXwjcEsx9ApgsIgdne1QjCxqp5yhDSdhMkqbLScV8MtoqMTI4Z7MXuJk+0nQwibtGdPTYqDn3Pu94+uCaUHCurhP+pu0KBeeo/V06oUR1UEmzjsu6YxkcEZwBeoZ2c+K6FYnW9h2nHLSIdIvI08BG4EFVfTI0xanDgIicUe1Q8MaW/0x75pbgq3puhKuhJGmXk2bMJ432isw9u3Q+SdjBJNE1UTRzbQY06oTi2iklaUeVJPOn7tqaaG3fcQrQqjqoqh+i8mDqw0VkTmiKU4cBVV2uqn2q2jfpwP2Tn9ZoCpfKDVdDSdIuJ82YT1J1VHExfXRF/bPN6JooElybR+VGo04orp1SknZUSTJ/a8/URGv7TqIqDlV9HVgNHBP6qGUdBlpBp6pniDaUhP+vTdPlpGI+GW2VqJhV0u8VV7nh0vkkSQeTRteojB4bNeeEsK6pT5y5ZCg0NiTS8Pfh0gklqoNKmnVc1t1DN3tl5E9Vu7rGcc/skxKt7TsNA7SITBeRycHX44EFwEuhafcAp0mFecAbqpquKZqRC651z1GGksXjb65rQnHh8ikXR5pSrpl8TuZ7AW6mjyXz0YUfRLsqAU67BD2sN/k1Cz+IXnj0yOsO6x01x/UGYVU9x5lLXv7bK9k1/WBUhF3TD+blc5ax8dhFaFdXsF8Xb86dN2KOy429rUcez7qzLh1x3cZjFyVeJ27dLT1TUWBLz1Rue9/p3Pre00eMxd0g7GRcqjjmAjdT6RzQBdylqpeKyF8DqOr1QaXHN6go6x3An6vqmnrrFrWKw1f17KMxxde6ZzOmjKZTqzga1kGr6rPAhyPGr6/5WoGzsj2akRU+BmejXJgxJR3WUaUGX9Vzs7h0S2k1mannsJll3mzkiXXJG8K6rAMN57y88yjm3npUXWNI2PTR6kau7d6/LIjIecDfA9NVdXOaNSxAe06z6jlsXKl2SwHaHqSbJmxmeW071JpQXtsOVz1cuVFWL0i7rHPFD0AE2TsUP2fZgxyqD9M1WFmnZ9MGZl93CbCvw3bY9BE1J0/S7N+J6llEDgE+DbzazDr2LI4AU8/RuHRLKSuRZpbwHIeOJk7rDOq+4Bw3Z+/QcHCuEjaGRJk+0phH0tLu/UvE1cBSIsqNk2AK2qiLS7eUVpNZeiPOzJJ0nus6Kak1esSZPpKaR7I4i8t4idXzNBGpLXRYrqrLXS4UkROB9ar6TKV+Ij0WoPFXPWdxc3Bm9wDrB3sjx0vPjImV1ILLvCzWSUmt0WP3tIPo2TS6gjWpeaSZs7Rz/2b41axu/vDvJrtNPr1/c70qDhF5CIj6TV8EXAh8Js0Zw1iKw6iLS7eUVpJlaV2kmSU8J2RUSb1Ot6BjuurOGRozlqHukeuEjSFRpo805pG0JNm/xOq5Iaq6QFXnhF/Ar4FDgWdEZB0V095TIpLqfzBT0J6SVWld9UZg0ao4MmHB+ytBMlRZQdIqDtd1Yue8NVwNAdStkKh+3a4qinbvX3RU9TlgRvV9EKT7rIojJb6mN7LkpAkrCxGQczGmLHg/GgrAqe7qOK4TNSdsTHFx/LUzILrs77N6biUdH6B9xIwp5cFX16BRQVVnN3N9RwdoU8/lITdbd9hgEtUJJdwtxSXt4dJ1Zd5s5j5+V0s6o+SZhgjvdcfU48DhmRp9G5+0jioN6OgA7SOmnhMQZTB5bfu++uQhhZXPwX3PI4O6b04j80q460rcOiufoyeY0rNpA4de83lkaLDmuiEOePaJ4feuppRWmlmi9jply60AdYNt38YnOeXfb6VnaDcAB+7ayin/3vi6TqNjqzhMPZeHvNSzk8EE9gXV6lgD80ps15XwOqE5XbXBOWZO2s4oeZlJovZy6Yxy4roVw8E5yXWdRscGaB8x9ZyQZgwm9a7NuXNK2s4oeZhZ4tZs1Bkl7nPrqDKSjgzQpp7LQ66PFG1kQEl7bTNdVxxI2xklDzNJ3JqNOqPEfW4dVUbSkQHaR0w9J8fJYELFZDJirIF5JbbrSnid0Jyhru7cOqPkZWaJ2sulM8o9s09iV9e4xNd1Gh13k9DUc3nI/YH8UQaTWZMgooqDJFUcS+ZX1o2o4qisUzGmvN73KSav+WHmVRytNJOM2GvTBudqjOrnVsVRn4YdVfKiXR1VfAzQvqpn65hSHvI2pmTRUWX/Q8frf/m733Gau/b05wvRUaWjUhw+BmejXPgYnI386LgUh290pHoOm0tcu56kXev0W5BXtg2/1XdPgZtOa/pMUWYSKN9zLszWnR8dE6BNPXtClLnEpetJ2rWC4Dzi9t4r29DTb9kXpB3PVKueowweh37tIlCt21HF6Cw6KsVhlIN66jnSXOLQ9STtWqOCM4HppEZRpzlTlMGja++ehh1Vioap53zpiADtq3r2Nb1RlziDSBrTSVZrOawTzj0nMY20qluKUTw6IkAb5aFh5UacQSSN6SSrtVKsk8Q0UtRuJaae88cCdEnpSPVMjLnEoetJ2rX03VOiTSfvnuK8TlTlRpTBw6WjitFZeH+T0Nf0ho841T1HmUvSVnG4rHXTaZUbgvWqOFKcKc5MEjVWxBuEpp5bg/cB2kc6VT0PE9G9JNe1bjqtcZeVmHXq1T3HdSYpYkA22oPXKQ5Tz+XBV9egYTSDKeiS0fHquRlczCRRcyC5Meahl9iz/Fn6Nn+t0KmKNGSV3vC1o4qIfAi4HtgP2Av8T1X9aZq1vA3Qpp7LQ0vUs4uZJGrOFT8AEWTvUPx1MXv17DLDSRyed1S5Eviiqt4vIscF749Ms5DXKQ7fMPWcHidTStScQd0XnGOuc9mr6IYTV7JSz553VFHggODrScBA2oW8VNCmnstDy3LPLqaUJAaVenNjPjPDyT5K0FFlmoisqXm/XFWXO167BHhARK6iIoKT14AGeBmgfcTUc5PMmFhJT0SNN5oTt17CvYpqOHEly9K6rT1TOTAiGOfZUWXmxAP5wqdOazwROJGlm+s9blREHgKi/hV9RlcAAAyWSURBVEIvAo4CzlHVu0XkZOAfgQUpjuxfisPUsxGFkyklak63oGO66l4X5td/dnHLOpqUlbJ3VFHVBao6J+K1Evgc8N1g6j8Bh6fdxxR0CfBVPbe0tM7FTBI3p9F1IVrZ0aRVZG1M8byjygDwKWA1MB/497QLeRWgTT0bdXExpcTMcTXGVI0pcSYUYx9rZnzMl4Ac5r8D14rIGOBt4Iy0C3kVoI3yYMaU8mC27mSo6mPAR7JYq2EOWkQOEZF/EZEXReQFERnVSFBEjhSRN0Tk6eB1SRaHS4Kv6tnX9IaPWDsrI2tcFPRe4FxVfUpEJgJrReRBVf1FaN6PVNV+pjMaMkI9p3X3pXkWh4tLcN5s5Il12bTTasC7vnnpiA7eG48+mVfPTKdtotpnZZFiMfXcXhoGaFXdAGwIvt4uIi8Cs4BwgDYyxnv1nNbdl6bFVdQ6yx4EVWRQ942tfG5fB5UEeyVVz+/65qXMuP/OfXsNDTHj/jsBEgfpqPZZ5lz0g0RldiIyG/gw8GTExx8XkWdE5H4R+d0MzuaMr+kNH6lVz6ndfSlaXEWus3doX3CujoWvS9lOqxEzHrgrspXWjAfuSrxWVPusLJyLpp7bj/NNQhF5B3A3sERV3wx9/BTwblV9K/Cefw94T8QaZxDc0Zw+c3LqQ3cC3qtnaM7dl1VbqgyuTZV7HhpKNl6HOIeiORfLj5OCFpGxVILz7ar63fDnqvqmqr4VfL0KGCsi0yLmLVfVPlXtm3Tg/k0evYKp5/IwqnLDpVVU3m2p8r42jq6Yb7248TrEORTL7lw03Ko4hIpV8UVV/WrMnIOCeYjI4cG6W7I8aCfREeqZJtx9KVpcRa4zpgvtHploGNXeqsFeaSs3Nh59cmQrrY1Hn5x4raj2Wc06Fy29UQxcUhxHAIuB50Tk6WDsQuBdAKp6PfAnwJkishfYCSxS1YZNKJrF1HN5iKx7bsbdl7SywtUlOG82tKCKo3ojMIsqDh+di0YFlyqOxxh97yQ85xvAN7I6VCfTKep5mCbcfVntFR7LXVkEvHrmJanL6sJk6Vw09VwcSvuwJFPP5cFX16AZU4y8Mat3geg49QzZmVBcuOYR5N7nYUihS9AT5sCS+fnsVVJMPReLUgZoU8+ekJUJxYVrHkFqTShDCiufq+yVIkibejZaQWlTHL7hq3qul97IyoTigtz7fKQxRO59PvO9yoqp5+JRugBt6tkjsjKhuDAUc+svbrwOpp6NVlG6AG2Uh4Y3B7MyobjQFVOIFDfeYZh6LialCtC+qmdf0xuNyMqE4rTXCXMijSF6wpxE65h6NlpJKW8SGsXHqbQuKxOKC0vmV/ayKo5RmHouLhag20ynqudhsjKhuLBkPtpEQDb1bLSa0qQ4fE1v+IivxhQfMfWcPSLyp0H3qSER6asZ/7SIrBWR54JfG6oFU9BtpMzqecWOhSzbfgEDgzOZ2T3A+ROv4KQJK5tfOE/jShNr56Ge8+qCYrSd54E/Ar4VGt8MnKCqAyIyB3iASvOTWEoRoE09F4sVOxZy/htXslMnALB+sJfz37gSgLNWNPFPKk/jSitNMQ4UoQuKqed8UNUXAYIHfNaO/7zm7QvAfiLSo6q74tYqTYrDN8qsnpdtv2A4OFfZqRNYtv2CptbN07jSzNp5qOe8uqAYpeGPgZ/XC85QAgVt6rl4DAzOjB4fqvvTWmPyNK600hTjgHVBaT2De6fy1rbPOs5eOk1E1tQMLFfV5dU3IvIQENUR4SJVrZvrC1oCLgM+0+gUhQ/QPlJm9Qwws3uA9YO9o8YPmbGtuYVnTKykHqLGmyXl2nlVbuyedhA9mzZEjrcCS280ZLOq9sV9qKoL0iwqIr3ACuA0Vf1Vo/mFTnGYei4m50+8gvGyY8TYhJ7dfPmv7mtq3TyNK600xbiQRxcUo9iIyGTg+8DnVfXHLteYgm4xZVfPwHC1Rm0VxxXnPcIpC55qbuE8jSutNMU40M4uKKae80VETgK+DkwHvi8iT6vq0cD/An4HuFhELg6mf0ZVN8atVdgAbeq52Jw0YeVwoM607jlP40rCtfM2pmTZBcUoDqq6gkoaIzz+JeBLSdYqdIrDN3xQz0a5MfVcLgqpoE09G+0mqXp+1zcvbdgA1owpRlIKGaB9xFf1bLbuIDjff2dNt5YhZtx/J7Cve7cZU4w0FC7FYerZaDdJ1fOMB+6K7NYy44G7ht+bMcVIQ+ECtFEeTD0HDA01HG+3McXUczkpVID2VT37mt7wkVSVG10x30Y143EGlFYZU4xyUqgAbZQHU8/72Hj0yZHdWjYeffLw+3YaU0w9lxe7SZgzpp7LQ9q65+qNwHpVHO00phjlpTAB2tf0ho+Yeh7Nq2deMqqsLkw7jCmmnsuNpThyxNRzebB2VkYRKUSANvVcHkw9G0brKESA9hFTz+XBV/Vs6Y3y0/YAbeq5PJh6NozW0vYA7SOmno12Y+rZD9oaoE09lwdf1bOv6Q3DD0xBZ4ypZ6PdmHr2h7YF6IHtW9q1tWEApp6N4mMKOkN8Vc++pjd8xNSzX1iANjoSU8/l4ejer7b7CG2jYYAWkUNE5F9E5EUReUFEzo6YIyLyNRH5pYg8KyKH5XNco9W0XD0/9BKy6EZk/rXIohvhoZdau3+JMfVcDETkT4NYOSQifaHP5orIT4LPnxOR/eLWAbdncewFzlXVp0RkIrBWRB5U1V/UzDkWeE/w+hjwzeDXjsHX9EZLeegl5KqHkV17K+9f2w5XPVx5UlyGjWRNPZeHkqrn54E/Ar5VOygiY4DbgMWq+oyIHAjsqbdQQwWtqhtU9ang6+3Ai8Cs0LSFwC1a4Qlgsogc7Pq7MYpJq9Wz3PD4vuBcHdu1F7nh8Zaeo4yYei4Oqvqiqv5rxEefAZ5V1WeCeVtUdbDeWqIafpJtnckis4FHgTmq+mbN+H3AFar6WPD+YeB8VV0Tuv4M4Izg7Rwq/9OUjWnA5nYfIiGlOPNH4CNxn62Fta08SxOU4s86RBnO/G5Vnd7MAiLyz1R+ry7sB9T2KFuuqssT7rcaOK8aB0VkCZV/4zOA6cCdqnplvTWcHzcqIu8A7gaW1Abn6scRl4yK/MFvcHmw3hpV7Rt1VcEp47ntzK2jjOcu45nToKrHZLWWiDwERLXDuUhVV8ZcNgb4PeCjwA7gYRFZq6oPx+3jFKBFZCyV4Hy7qn43Yko/cEjN+15gwGVtwzCMsqGqC1Jc1g/8UFU3A4jIKuAwIDZAu1RxCPCPwIuqGpexvwc4LajmmAe8oaobkp7eMAzDYx4A5orIhOCG4aeAX9S7wEVBHwEsBp4TkaeDsQuBdwGo6vXAKuA44JdUpPufO6ybKJ9TIMp4bjtz6yjjuct45sIiIicBX6eSZ/6+iDytqker6jYR+SrwMyop4FWq+v26ayW5SWgYhmG0DnMSGoZhFBQL0IZhGAWl5QHaxTpeNERkPxH5qYg8E5z5i+0+kysi0i0iPw9q1UuBiKwLbLBPi8iaxle0HxGZLCLfEZGXgn/bH2/3mRohIu8L/oyrrzeDWl2jILQ8Bx04DA+utY4DfxiyjheKoJJlf1V9Kyg5fAw4O3BNFhoR+VugDzhAVY9v93lcEJF1QF+1HKkMiMjNwI9U9QYRGQdMUNXS+P9FpBtYD3xMVV9p93mMCi1X0I7W8UIRWNjfCt6ODV6Fv7sqIr3AHwA3tPssPiMiBwCfpFKOiqruLlNwDjgK+JUF52LR1hx0YB3/MPBkO8/hQpAqeBrYCDyoqoU/M3ANsBQYavdBEqLAD0RkbfB4gKLzW8Am4P8F6aQbRGT/dh8qIYuAb7f7EMZI2hagG1jHC4eqDqrqh6i4JA8XkTntPlM9ROR4YKOqluUZFrUcoaqHUXlK4lki8sl2H6gBY6g4wr6pqh8G/hO4oL1HcidIyZwI/FO7z2KMpC0B2sE6XliCH11XA5n5+nPiCODEIJ97JzBfRG5r75HcUNWB4NeNwArg8PaeqCH9QH/NT1XfoRKwy8KxwFOq+lq7D2KMpB1VHC7W8UIhItNFZHLw9XhgAVDoJ8mr6udVtVdVZ1P58fURVT21zcdqiIjsH9w8JkgTfIaCP/VQVX8D/IeIvC8YOooGFt6C8VksvVFInJ9mlyGR1nFVXdWGs7hyMHBzcKe7C7hLVUtTtlYy3gmsqPw/zhjgDlX95/YeyYm/AW4P0gW/xu1xB21HRCYAnwasi0EBMau3YRhGQTEnoWEYRkGxAG0YhlFQLEAbhmEUFAvQhmEYBcUCtGEYRkGxAG0YhlFQLEAbhmEUlP8PExcM0DDkDbcAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "#plot boundary to separate data points\n", "print(theta.shape)\n", "u = np.linspace(np.ndarray.min(x_train[0:]), np.ndarray.max(x_train[0:]), num=2)\n", "v = np.linspace(np.ndarray.min(x_train[1:]), np.ndarray.max(x_train[1:]), num=2)\n", "z = np.zeros((len(u), len(v)))\n", "for i in range(0, len(u)):\n", " for j in range(0, len(v)):\n", " z[i, j] = theta[1]*v[i] + theta[0]*u[j]+logreg.intercept_\n", " \n", "plt.contourf(u, v, z)\n", "plt.colorbar();\n", "\n", "# plot_x = [np.ndarray.min(x_train[1:]), np.ndarray.max(x_train[1:])]\n", "# plot_y = np.subtract(np.multiply(-(theta[2][0]/theta[1][0]), plot_x), theta[0][0]/theta[1][0])\n", "# plt.plot(plot_x, plot_y, 'b-')\n", "\n", "plt.plot(x_train[0, pos], x_train[1, pos], 'ro')\n", "plt.plot(x_train[0, neg], x_train[1, neg], 'bo') \n", "plt.show()" ] } ], "metadata": { "celltoolbar": "Slideshow", "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.6" } }, "nbformat": 4, "nbformat_minor": 2 }