{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Gaussian mixture model\n", "\n", "Gaussian mixture model (GMM) is a probabilistic model created by averaging multiple Gaussian density functions.\n", "It is not uncommon to think of these models as a clustering technique because when a model is fitted, it can be used to backtrack which individual density each samples is created from.\n", "However, in `chaospy`, which first and foremost deals with forward problems, sees GMM as a very flexible class of distributions.\n", "\n", "On the most basic level constructing GMM in `chaospy` can be done from a sequence of means and covariances:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "execution": { "iopub.execute_input": "2021-05-18T10:57:26.262571Z", "iopub.status.busy": "2021-05-18T10:57:26.262113Z", "iopub.status.idle": "2021-05-18T10:57:26.271457Z", "shell.execute_reply": "2021-05-18T10:57:26.271735Z" } }, "outputs": [ { "data": { "text/plain": [ "GaussianMixture()" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import chaospy\n", "\n", "means = ([0, 1], [1, 1], [1, 0])\n", "covariances = ([[1.0, -0.9], [-0.9, 1.0]],\n", " [[1.0, 0.9], [ 0.9, 1.0]],\n", " [[0.1, 0.0], [ 0.0, 0.1]])\n", "distribution = chaospy.GaussianMixture(means, covariances)\n", "distribution" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "execution": { "iopub.execute_input": "2021-05-18T10:57:26.274489Z", "iopub.status.busy": "2021-05-18T10:57:26.274179Z", "iopub.status.idle": "2021-05-18T10:57:26.363215Z", "shell.execute_reply": "2021-05-18T10:57:26.362921Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA5cAAAFyCAYAAABhgy79AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAuJAAALiQE3ycutAAAoXklEQVR4nO3df6jl+X3f99euRknW6Y7RVNGy0UpRGkkfrUEYy9SspXa9LSgqcWzFFIvFlROIa5I4rmXiUpdQKbJjr21RaEIKCqSBEIEhieJGChgZ0yQoRFnb6drUbsJHiCj+Qap1lC3ZxbJrSZ7+Mfd4z9y5984593x/fD6f7+MBQpo7Z+Yc7d7z+byf9/M95zxw+/btAAAAwCkeXPsBAAAA0D9xCQAAwMnEJQAAACcTlwAAAJxMXAIAAHAycQkAAMDJbkzxl5RSHkzy8SQ/m+SVSb4qyZ+utf6Hc7d7T5Ink3z+7L7fX2v1WSgAAACdm/Lk8tla6w/WWt+f5EtJ/tT+b5ZSXpvkmSTvq7V+MMmbknzLhPcPAADASiaJy1rr79RafyhJSimvSPKHkvzLczd7V5Lnaq1fPvv1J5O8e4r7BwAAYF2TXBa7U0p5Z5I/n+SfJ/kn5377kSQv7v36xbOv3eXXX/qiy2QBAAAa8JqHX/nAobed9A19aq0/XWv9E0n+XZIPn/vt55Pc3Pv1zbOvAQAA0LlJ4rKU8lWllD++96XPJHlLKeVGKeV1Z1/7qSRvO7tsNrnzxj4fm+L+AQAAWNcDt2+ffhVqKeWPJPlQkl84+9LXJPnhJF+Z5Jla6xNnt3tPkqeSvHB2u3veLdZlsQAAAG045rLYSeJySuISAACgDau95hIAAIBtEpcAAACcTFwCAABwMnEJAADAycQlAAAAJxOXAAAAnExcAgAAcDJxCQAAwMnEJQAAACcTlwAAAJxMXAIAAHAycQkAAMDJxCUAAAAnE5cAAACcTFwCAABwMnEJAADAycQlAAAAJxOXAAAAnExcAgAAcDJxCQAAwMnEJQAAACcTlwAAAJxMXAIAAHAycQkAAMDJxCUAAAAnE5cAAACcTFwCAABwMnEJAADAycQlAAAAJxOXAAAAnExcAgAAcDJxCQAAwMnEJQAAACcTlwAAAJxMXAIAAHCyG1P8JaWUNyT50SQ/n+RVSW4m+Z5a65fO3e7ZJL+196Wna62fm+IxAAAAsJ5J4jLJrSQfrbV+NElKKZ9I8m1J/va5232i1vrBie4TAACARkwSl7XW55I8t/elB5O8dMFN31pK+f4kvyfJv6m1fmSK+wcAAGBdU51c/q5Sytcn+UKSj13w2x+qtf5MKeWBJH+3lPL7aq1/Y+rHwHx+7YXfnPXvf+zWQ7P+/QAAHMf8x6EeuH379mR/WSnl65J8V5I/W2v9rfvc9s8l+S9rrd+6//Vff+mL0z0gTjb3YnIdFiAAgOmY97jKax5+5QOH3nayk8tSyjuTfFOS70jyQCnlW5P870kerbX+ainlLUneUWv9m2d/5E1JPj3V/TONFheX8656jBYiAIB79TDj7dt/vOa7fkz1brFfm+QfJPm5JP9H7rzm8rkk/z7JM0meSPJikm8spbw2ye/NnXeV/cAU9880elt0LnLZ/weLEgCwBSPMc+ft/j+Z59o36WWxU3BZ7PJGXIQOYYECAHpmhmMJx1wWKy43bquL0mUsVgBAi8xsdzOzLUdcchCL1P1ZuACApZnRDmNOW4a45L4sWtdnIQMApmQuuz5z2fxWebdY+mEBO835f34WNQDgGGYxRiUu4URiEwC4iphkK8QlTExsAsC2icnl/NoLv2nWaoi4hJmJTQAYm5iEO8QlLExsAkDfxCRcTFzCyvY3KKEJAO0Rk3AYcQkNcaoJAG0QlHA8cQkNc6oJAMsQk3A6cQmdcKoJANMSlDAtcQmdcqoJAMcTlGMxA7VFXG7QY7cesrAORmgCwMXMPLAccQmDcfksAFsnKGEd4hIG51QTgC0QlNtjrmmPuNwol8Zuk9AEYCRmGWiLuISNEpoA9EZMsmN2aZO43DCnl+wITQBaZVaBfohL4C7eEAiAtQlKrmI2aZe43Dinl9yPU00AlmAegf49cPv27bUfw11+/aUvtvWANsKCzrGEJgCnMn9wLPPH8l7z8CsfOPS2Ti6Ba3GiCcB1CEquy7zRPnFJEpfHchqhCcBVzBiwDS6L5S4Wf6YiMgG2zUzBlMwV6znmslhxyT1sBkzNhgCwDWYI5mCOWJfXXAJNcdkswLgEJXMyN/RFXHIPr79kTkITYAxmBeA8l8VyKZsGSxKaAO0zG7Aks0EbvOaSydhEWIPNBKAdZgHWYBZoh7hkUjYV1mJjAViHvZ812f/bIi6ZnE2GtdloAOZlr6cF9vv2iEtmYdOhFTYegOnY32mF/b1N4pLZ2IBojY0I4Hj2c1pjP2/X4p9zWUp5Q5IfTfLzSV6V5GaS76m1func7d6T5Mkknz+77/fXWsUkcG27AcmmBHA1QUmr7OHjmOpzLm8l+Wit9aNJUkr5RJJvS/K3dzcopbw2yTNJSq31y6WUv5PkW5L8xESPgQXsnvw2KFrj8zMBLmbPpmX27LFMEpe11ueSPLf3pQeTvHTuZu9K8lyt9ctnv/5kkndHXHbpsVsP2axoltAEts4eTQ/s0eOZ6uTyd5VSvj7JF5J87NxvPZLkxb1fv3j2NTolMOmBy2aBrbAn0xP78pgmjctSytcl+TNJnq61/s65334+ydfs/frm2dfomMCkF04zgVHZh+mNfXhcD071F5VS3pnkvUm+I8mXSinfWkq5UUp53dlNfirJ20oprzj79ZO593STDlkg6M2vvfCbv/sfgB5Zx+iVuXFsk3wUSSnla3PnNZQ/d/alB3PnNZgfT/JMrfWJs9u9J8lTSV44u9097xbro0j6ZYOjZzY7oHX2WXpnr+2Tz7lkNTY+RmDzA1pib2UE9tZ+iUtWZRNkFDZCYC32UkZiP+2buGR1NkVGY2MElmD/ZDT2z/6JS5pgg5zOv/r8i/e/0YEef/XNyf6urbJRAlOyX57OPtkee+U4xCXNsGFez5Sb5DFsqMexcQKnsEcex97YD/vjWMQlzbGBHmatjfMqNtXD2EiBQ9gP789e2Df74XjEJU2yoV6txc30MjbZy9lUgYvYAy9n/xuHPXBM4pJm2Vwv1tPGehkb7r1ssrBt9rx79b7f2esuZ88bl7ikaTbbe/W+2V7EBvwyGy5si33uZfa3bbDPjU1c0jwb78tG3HgvYjO2+cLI7Gt3bGFPs5+9zL62DcfE5Y05HwhcZrcY2Yy34/zAscXNef/73YYMY9j6PraFmORi9jEuIi5Z1WO3Htr8xrxV+wPJlkPT5gz92fK+JSZJ7F1cTlyyOoHJlk81nWZCP7a6VwlK9tmruIrXXNKMrW7aiY37KlsKzR0bN7Rjq3uTfelyW9yXEnvTlnnNJV3yOkwussVTTaeZsL6t7UVikqvYiziUuKQ5LpPlKlt7rabXZsKytrT/CEoOYf/hGC6LpVlb2uATm/ypthCaiU0e5rCl/cZec5qt7DWJ/YaXuSyWIbhMlmNs5UTTJbMwna3sL4KSY9lfuC4nl3TBAMB1jRyaO4YAOM4W9hT7yfS2sJ8k9hTudczJpbikG4YBTjX6YGAggMuNvofYP+Y1+v6R2EO4nLhkaAYEpjD6oGBIgDvsGUzBnsGWiUuGN/KwYFBY3shDg4GBLRp5j0jsE0sbeY9I7BPcn7hkE0YeHgwO6xl5iDBAMDr7AlOzJ4C4ZGNGHSYMEusbdagwUDAa+wBzGHUPSOwDHEdcsjkGC+Y24pBhuKB3I6791v12jLjuJ9Z+jicu2aQRh4zEoNGaEYcNgwY9GXGtt863x1oPLxOXbJrBg6UYPmA51naWYm2Hu4lLNs8QwtJGG0YMIrRitPXcWt42azncS1zCGUMJSzOYwOms3azB+g0XE5ewx5DCWgwqcBzrNWsZab22VjM1cQkXGGloMbD0x+ACl7M+sybrM1xNXMIlDDCszRADd1iPaYE1Ge5PXMJ9jDLUGGj6Zqhhi0ZZfxNrcO9GWYOtv8xNXMIBRhlwDDf9G2XASQw5XM6aSyusuXAccQlHGGHgMeyMY5Shx8BDMsb6umOdHYM1Fo4nLuFIowxAhp9xjDIAJYagLbKm0qJR1lVrKktbPC5LKY8l+YEkb6+1Pn7JbZ5N8lt7X3q61vq587cTl6xphIHIMDQeAxE9GGH93LGOjmeEddQaylqOicsbE93nk0l+Isk7rrjNJ2qtH5zo/mAWj916qPsB6fFX3zQYDWb377P34Wj/uWVIGkfva+aOdXNcva+diTWTfkx2WWwp5akkf73W+pZLfv/vJ/nZJL8nyb+ptX7kots5uaQFIwxLBqWxjTAsJQamXo2wRu5YK8c1wjppjaQFa5xcHuJDtdafKaU8kOTvllJ+X631byx4/3Cw3WLe8wDlBHNsTjNZQ89r4j5r4/h6XxsTayJ9Wiwua60/c/bft0sp/yjJH00iLmla75G521wNUuPa/3fb+zAlNNvU6/p3EWvhNvS+Flr/6NlscVlKuZHk0Vrrr5ZS3pLkHbXWv3n2229K8um57hum1vtrMZ1ibsMop5nJy0FjyFpHz+vdeda+bel5/bPeMYJJ4rKU8s4kTyf5j0sp/1OSDyf56iTPJHkiyYtJvrGU8tokvzfJq5J8YIr7hqWMcIppyNoGp5lcR69r22Wsd9vS+1pnfWMUPucSrqH3IczQtT29D17nGcSm0/t6ts/atk09r2/WMnqw+OdcTklc0pOehzJD2Hb1PIhdxHB2nJ7XrctYz7ar5/XM2kUvxCUsqOdBzUC2bT0PZVcxsN2r53XqKtaw7ep5/bJG0RtxCSvoeXgzoNHzoHaVrQ5xPa9H92O9ouf1aqtrEn0Tl7CiXoc6AxtJ30PbIUYd7Hpdd45hjaLn9WnUtYdtEJewsp4HPQMcOz0PcsfocejreY05hvWInV7Xox7XFzhPXEIjeh0ADXTs63WoO9XaQ2Gv68cUrEHs63UNWnsNgamIS2hIzwOiAY/zeh3y5nSdAbLndWEu1hvO63W9EZWMRlxCg3odJg18XKTXoY/2WGO4SK9rjLBkROISGiYyGU2vQyDrsZ5wmV7XE1HJyMQlNE5gMqJeh0KWYf3gfnpcQ0QlWyAuoRMik1H1OCQyPWsFh+h1vRCWbIW4hM70GJmGRg7V6+DI9VkfOFSP64OoZGvEJXSox8BMDJEcp8dBksNYCzhGj2uBqGSrxCV0rMfINFRyHT0Ol9zNc5/r6PG5LyzZMnEJAxCZbEmPw+YWeY5zih6f56ISxCUMQ2CyVT0OoaPynOZUPT6fRSW8TFzCYEQmW9bjYNozz12m1OPzV1jC3cQlDKq3yDSkMoceh9WWeZ4yhx6fp6ISLiYuYWC9BWZieGVePQ6xa/J8ZG69PSdFJVxNXMIGiMz+fepXXrrw629//cMLP5Ix9TbgzsFz7voue37u81y9W4/POWEJ9ycuYUN6i0zD7h2HDK6XMdBeX4/D7yE8r67vlOfiztafkz0+r0QlHE5cwgaJzH5MMcyet/XhdiotD8lbfs5MYY7n3b6tPgdbfs5cRFTC8cQlbFRvgZlsc2Cee8jd2eqwu4ZjBuwtfs+vYann2c7Wnm+iErZDXMLG9RaZWxu2lx56d7Y2/LItaz2vkm09t3qLykRYwqmOicsbcz4QYB27jbSXyNwNK1uLzKWdH763NBAzljVDcst6C0tRCcsTlzAwkclVLhrQBSctEpPrEpXAoVwWCxvRS2DujB6YPQ3LgpMl9fTcSMZ+fohKIPGaS+AKIrMdvQ3R+0YeqFlOz8+BnVGfC8IS2BGXwH2JzPWNMFjvG3XIZhqjfb8nY37Pi0rgPHEJHKynyBSY/Rlx+OZqo39P74z2vS0qgcuIS+AoPQVmMl5kbmUY3zfaYL5FW/y+3Rnt+7ensBSVsDxxCVyLyFzPlgf180Yb3Hvm+/Jio3yPikrgEOISOInIXIdB/nCjDPdr8v12PSN87/UUlYmwhLWJS2ASPUWmwOQyI8TAoXz/zGuE76WewlJUQhsWj8tSymNJfiDJ22utj19ym/ckeTLJ55PcSPL+Wus9dy4uoT0ic1kCoV1zxIV/333oPSxFJXBdx8TljYnu88kkP5HkHRf9ZinltUmeSVJqrV8upfydJN9y9meAxj1266FuAnM3QPUcmW9//cOCo1H+vWxTz2EpKoElPTjFX1Jr/fEkv3HFTd6V5Lla65fPfv3JJO+e4r6BZTx266GuNv7HX32zq6HqvJ6HWRhJr8/F3tbAnvYX4HKTxOUBHkmyf4zw4tnXgM70GJm96nWohVH0+hzsad3rbU8BrrZUXD6fZH+lu3n2NaBTPQ0Evf0Ef1+vwy30rsfnXk9rXU97CHC42eKylHKjlPK6s1/+VJK3lVJecfbrJ5N8bK77BpbT04DQ0+C1r8chF1hWL2tbT3sGcLyp3i32nUmeTvLNSf5qkg8n+eokz9Ranzi7zXuSPJXkhbM/5t1iYUC9vPFP0t+b/ngzGVhGTz/Q6SkqgT75nEtgVQJzPgIT5tVLWPYSlYmwhN6JS6AJInMeAhPmISynJSphDOISaIrInJ7AhGn1EJaiEljDMXG51LvFAhvW0xs49DK89TAIA9PpYW3qaa0H5uHkElhcLyeZPZxiOsGE07X8w5oeojJxWgkjc1ks0AWROQ2Bebpf/OwL979RZ976h2+t/RC6ICxPIyphfMfE5Y05HwjAVR679VAXgfn4q282H5jcMWIkXtex/yzEaDtEJdArJ5dAE3qIzKTdU8wtnV4KyGWNHJ0tnloKS6A1LosFutVDZArM5QjJ9owUm63FZethKSphm8Ql0D2ReT0jBKag7EuvsSksDycqYdvEJTCM1iNTYE5HVPatt8hsJS5bjspEWALiEhhM64GZtBWZvcWlqBxLD5EpLO9PVAI74hIYUuuRKTCPIyrH1nJkthCXwhLoxTFx+eCcDwRgSo/deqjpoefxV99sZmBsYXi+irAcn3/Hl2tlnbhIy2ss0D4nl0CXnGIepsUTTNGxPS2dYq79g5dWw1JUApdxcgkMr4dTTO4lLLfJv/c7Wl0XWl5Lgb6IS6BrLUdmC5fJrn1Ks09gbJt//21qdf0E+iQugSG0PCAJTLhjy4G59jpwkZbXTaBPXnMJDKfV12Ou+TrMtV97OXpUfObTz0/6973xzY9M+ve1Zq3XYK71g5bWwlJUAsc45jWXN+Z8IABreOzWQ00G5uOvvrlaYL799Q+vHpg9mzoeT72/0eMTgD6JS2BIu5/MtxaZawYmh1s6Jo910ePrKTh/8bMvNPUOsnNyaglsictigeG1Fpg7a0TmGqeXPVwS23pMXkcPsbl0YK5xWWxLcSksges45rJYcQlsRouRuYXAbDkuR4zKi7QcmiMHprAERuA1lwAXaPG1mC6TXcdWonJn//9vy6HJPIQlsBQfRQJsSoufi7n06caWP5rkM59+fnNheV5r/wxaPtkG4DjiEtikrQfmFrUUVC1oKTIF5nxaW+uAsYlLYLO2PHRt7fSylYhqUSuRuVRg+kgegPmIS2DTWgpMp5fTayWceuCfFQCnEpfA5m01MLd2eslh1gxMl8cC9E1cAmS7gTkyp3DXt+Yp5hKBuZVLY1ta14BtEJcAZ1oaxJYKzCVOL5f+HEOmIzBPs/bHDLX20UvA+MQlwJ6WApPrc2o5nbVOMV0iC9AfcQlwTiuBOdLpJf0bMTC3cHoJsCRxCXCBrQUmHEJgAnAVcQlwiVYCcwlOLznUGpfJ9h6Ya55eet0lsCRxCXCFFgJzhNNLb+ozntECc24ujwW24MZUf1Ep5RuSvDfJZ5M8muT7aq2/fe42zyb5rb0vPV1r/dxUjwEAWM5nPv183vjmRxa7v1/87Auz/aDiU7/y0rAn+L/2wm828YMyYHyTnFyWUr4iyUdyJyifSfKFJN99wU0/UWt9au8/whJoXgtD2RKnl6MO1sxr6ctk5zzBHPnyWIAlTHVZ7BNJPldr3a2an0zy7gtu99ZSyveXUt5fSvn2ie4bgA4seWnskqdp3CEwD7NWYHrtJbCEqeLykST7q+WLZ18770O11h9L8kNJvrmU8p0T3T/ArJxewv0teYr5i599YbbI/NSvvDRrZApMYFRTxeXzSfannptnX7tLrfVnzv77dpJ/lOSPTnT/AHTAG/tsg1PM+xOYwIimistnkzxaStkF5pNJPlZKuVFKeV2SlFLeUkr5jr0/86Ykn57o/gFm5/SyLy6NXddIp5hzEZjAaCaJy1rrF5J8e5L/pZTygSS/P8n/mjuR+ffObvZikm8spXyglPLDSV6V5EemuH8A+uG1l9sywinmiIEJMIcHbt++vfZjuMuvv/TFth4QwJ4WfuK/xDA695uaLPmZhUt/HiOXWzL25/ohxpwn+0t/pm0LV2MA7XvNw6984NDbikuAIwnMaQjM7VoqMgXm/QlM4H6OicupXnMJAM1yeWxblno95lyvxZzz3WSXvky2hR+WAeMQlwBcaO439ln6nWMFZntGiMw5CEygV+ISoENLXzo3F4FJsmxkTk1gArzMay4BjtTKELbU8DnSay93vAazbUv8EGCOH2zMddq/5A+TvAYTOM9rLgHoxtKnl4kTzNYtcZI5x6WyI5xitvLDM6BPTi4BjtTK8LXkwDn36WWyzglm4hSzB3P/MKCXU0wnmMAanFwCbMAor7vcWeMEM3GK2YO5TzJ7OcV0ggm0zsklwJFaGrpGO71M1jvBTJxi9qK3k8yeTzGdYAJOLgHo1lonmIlTzF7sTjLn+mHA1CeZPZ9itvTDNKB9Ti4BjtDaoLX0xxUsdXqZrHuCmTjF7M1cPxho/RTTCSYwt2NOLsUlwBG2HpfJtgIzEZm96SEyBSbQE3EJMJPW4jIZ+/RyR2RyrDkiU2AKTNgicQkwgxbDMhn/9HKnhcBMRGZvthSZAhOYg7gEmIG4vNsagZmITK6n5cgUmEDLxCXADMTlvbYemInI7M3UkdliYCbLRKbAhG0QlwATazUsk3XjMlkvMBORubT/9zOfPvrPvOqNb57hkZxuysgUmMDIxCXAxMTl1dYMzERkzuE6IXmoVoKzxVNMgQm0RlwCTKjlsNwRmHeIzGnMGZYXWTs2Rz/FFJjAKcQlwER6CMukjbhM2gjMRGSeaum4PG+t2JwqMlsLTG/yA5xCXAJMRFwer5XATNqKzKSf0Fw7LvctHZqtnWL2FJjiEsYkLgEm0EtYJm3FZdJWYCZtRWYPgdlSXO5bMjRbOsUUmMCaxCXABMTl6VqLzKSd0Gw5MluNy52lIlNgXo/AhLGIS4AT9RSWOwLzOCLzcq3H5c4SkdlKYHqDH2At4hLgBD2GZdJuXO60GplJG6HZUmT2Epc7c0emwLwegQljEJcA19RrWCbtx+WOyLxaK5EpMO82WmCKS+BQ4hLgmnqOy6SfwEzajsxk/dBcOzJ7i8udOSNzisD0+kugN+IS4Bp6D8ukr7jcaT0yk/VCU2Bez+iB6fJYYEniEuBII4Rl0mdc7ms9NLcYmQLzXiMFprgE7kdcAhxhlLDc6T0wk/YjM1knNNeKTIF5r1MD0+WxQC+OicsH53wgAK0bLSxH8fbXP/y7/2nVW//wrUkC4RhvfPMjk72xzDGW+lzJqbUcxWu/phdgDk4ugU0bNS5HOL28TMunmksGwxqnmC3H2mVaPr1MXB4LtM/JJcABRg3LZJlBcS0tn2oueZq5xilmjyeYPQYxQK+cXAKbNHJY7hv5BPO8Fk80lzrJXPoUs7dga/n00msvgdat8oY+pZRvSPLeJJ9N8miS76u1/va527wnyZNJPp/kRpL311rvegDiEpjbVsJyZ0uBua+12FwiNJeMzJ4Cs+W4TFwaC7Rt8ctiSylfkeQjuROUzyT5QpLvPneb1yZ5Jsn7aq0fTPKmJN8yxf0DHGprYZmMfYnsVVq7fHaJS2aXvFT2VW98c5eXybKuLa7BsCVTvebyiSSfq7Xufjz+ySTvPnebdyV5rtb65StuAzCbLQ81Ww3Mnf3QXDs2l4rMpQjMcWz1KgdgOjcm+nseSbK/Ir149rVjbwMwiy2H5c4uMA2Qd19GuNbls7vAnOty2V1gLnGp7H5gtna57Bbi91O/8tLqPzQBSKY7uXw+yf6PxW+efe3Y2wAws8dffXPzJ5n71j7RnPskc413lG0l6Fp5HNzND/tgXFPF5bNJHi2l7KaVJ5N8rJRyo5TyurOv/VSSt5VSXrF/m4nuH+BSBpmLicx7rXn57JyRudbHlqwZd8LyelzZAJxiyneLfSrJtyf55SSvSfIXkvxnSZ6ptT5xdpv3JHkqye4aIO8WC8xKWB7OUHm1pS+fnfPdZZf+6JKdpS6ZXTIsW3i32GS6d4xNlnuNtneOhT6s8lEkUxGXwFSE5WnE5uWWDM0RI3Nnqthc65RyqtNgcQm07Ji4nOoNfQCaIixPd37AFJsvW/INgeZ8458l3/TnIveLwovi0+WuAO1ycgkMR1guR3C+bKnTzJFPMnsy5WtYt3pymTi9hB44uQRgEVcNoVsLz6VOM0c+yexFa2EJ0Aonl8BQnFr2ZfQAXeI000nmsqZ+193WTi0TJ5fA3byhD7BJwnI8I8Xn3KEpMufXYlgm4hKYl7gENkdYbk+v4dlzZCbbDc05PidUXN4hMKFt4hLYFGHJeb2EZ8+huZXInCMqk2lfaykugTmJS2BTxCWHaDk4e47MndFic66oTNoOy0RcAncTl8BmCEtO0VpwjhCZSd+hOWdUJtO/O2zvp5aJuITWiUtgE4QlU2spNucMzaUic6f12Jw7KHdaD8tEXAL3EpfAJohL5tZCbI5ymnnemsG5VEzu6yEsE3EJ3EtcAsMTlqxhzdgcNTIvMlV4rhGRF5k6LBNxCSxHXALDE5esTWhyP3NEZTJWWCbiElonLoGhCUtas1Zozh2ZidC8jrmiMpkvLBNxCVxMXALDEpa0Tmhu15xRmYwZlom4hNaJS2BY4pKerBGaS0Tmjti8Y+6oTOYNy0RcApcTl8CwxCW9Gj00k23F5hJBmcwflcm6YZmIS2iduASGJCwZxRZCMxkrNpeKyX1bCMtEXELrxCUwJHHJaLYSmft6Cc41YnLfEmGZiEvg/sQlMBxhyei2GJoXWTo+147I85aKyqSNsEzEJbROXALDEZdsycjvOMvFlozKpJ2wTMQltO6YuLwx5wMBAI63P/gvGZr7gSM0l7F0VCZthSUwFnEJAA1rITQTsTmlNYJyR1gCcxKXANCJtUIzcap5qjWDckdYAnMTlwAsaq3XV432ut1dKKzx+syLQklw3quFoExEJbAcb+gDdGG0MBhZr2/OMcL32FpvBHSZrQVnKzG5r/Ww7HW9gC3xbrHAcEYY/EezhaGw9++71mJzZ5TobDEmd1qPyp0trCPQO3EJDKf3Ib9nhr+X9fx92GpoXqSV+Gw5Hi/TS1TuWF+gfeISGFLPg30PDHnX0+P3ZU+hyeGEJTAHn3MJwJUMddPZ/2fZS2iejxCx2bfeohIYl5NLoCu9DO8tEZLr6fH7VWj2YYSgtDZBH1wWCwytx4F9CQa19vX4vSs22zFCUO5Yr6Af4hIYXo9D+pQMZv3r+XtYcC5rpKjcsYZBP8QlsBk9D+iHMIBtwwjfx4JzOiPG5D7rGvRFXAKb0/twbthiX+/fzzuC8zCjx+Q+ax30Z/G4LKW8L8kfTPJAkl+ttf61C27zhiSfSPK5sy99rtb69PnbiUvgVK0O5oYqrqPV7+dTbTE8txSRl7EOQn8WjctSyn+a5K/UWt9x9utnk3xXrfW5c7d7Q5Knaq1/66q/T1wCc5h7QDcwsYRRQ/MivcangLycdRL6tPTnXH5Tkn+29+t/muTdSZ676LallD+Q5CuT/GSt9VMT3D/AfRlqGMH57+ORY/M6kTZHkIrFaViDYRsOistSyk/mzmWv5308ySNJfnXvay8mee0Ft/13Sf5SrfWXSim/P8kvlFL+RK31/z7yMQMA2VZsHkIItklYwnYcFJe11j922e+VUn4wyf5qfjPJ8xf8Hb+R5Jd2/7uU8lySp5KISwCYwP4Qv/XQpA3CErblwQn+jn+Y5B17v/7Pk3wsSUopX1lKuXX2v/9kKeWte7d7U5JPT3D/AMA5j9166K7/wJJ838E2Tflusa/PnXeL/ezu3WJLKR9IcrPW+t+XUv6LJH8uyc8neV2Sf1tr/aHzf5c39AGA+TnZZC6iEsbicy4BgKOITU4lKmFMS79bLADQOW8OxHWJSmBHXAIA97goGAQn+0QlcJ64BAAOcllMiM5tEZXAZcQlAHASp5zjE5TAIcQlADC5q2JEeLZPTALXIS4BgEUJz/aISWAK4hIAaMahkSNCr0dEAnMSlwBAd46JpK2EqHAE1iYuAYChzRVdU0WrKARGIS4BAK5BFALc7cG1HwAAAAD9E5cAAACcTFwCAABwMnEJAADAycQlAAAAJxOXAAAAnExcAgAAcDJxCQAAwMnEJQAAACcTlwAAAJxMXAIAAHAycQkAAMDJxCUAAAAnE5cAAACcTFwCAABwMnEJAADAycQlAAAAJxOXAAAAnExcAgAAcDJxCQAAwMnEJQAAACcTlwAAAJxMXAIAAHAycQkAAMDJbkzxl5RSnkjyY0n+Va31z15ymweT/OUkX0zyB5L841rrR6e4fwAAANZ18sllKeVGkrck+Sf3uel/neQNtdYPJnlfkh8rpTx66v0DAACwvpPjstb6pVrr30py+z43/aYk/2z3Z5I8l+S/OvX+AQAAWN9Bl8WWUn4yyR+84Lc+Xmv9wIH39UiSF/d+/eLZ1+7ymodf+cCBfx8AAACNOCgua61/bIL7ej7Jzb1f3zz7GgAAAJ2b9d1iSylfWUq5dfbLf5jkHWdfv5HkbUk+Mef9AwAAsIxJ4rKU8t8meTLJ20op37v3W+9L8hfP/vffT/LLpZS/nOSvJfn+Wuv/M8X9AwAAsK4Hbt++3/vwLK+U8q4k/02S/yt33on2X9Ra//q6jwru75CP5YEWlFK+Icl7k3w2yaNJvq/W+tvrPiq4WinlsSQ/kOTttdbH1348cIhSyhuS/GiSn0/yqtx5adj3nL3BJTTr7KMkP57kZ5O8MslXJfnTtdb/cNmfmfWy2BO8LskHa63/c5Lvyp2PLXn9yo8JrnTEx/LAqkopX5HkI7kTlM8k+UKS7173UcFBnkzyE0m8+R89uZXko7XWH6u1/o9J/pMk37byY4JDPVtr/cFa6/uTfCnJn7rqxk3GZa31f6u1/uu9L305yf+31uOBQxzxsTywtieSfK7WunsH708mefeKjwcOUmv98SS/sfbjgGPUWp+rtX5070sPJnlprccDh6q1/k6t9YeSpJTyiiR/KMm/vOrPHPRusXM44uNNvjfJj9RavbMsq5voY3lgbQd9NBQA0yqlfH3uXC3ysbUfCxyqlPLOJH8+yT/Pfa7QWy0uD/l4k1LKdyZ5Ra31RxZ4SHBfE30sD6zNR0MBLKyU8nVJ/kySp2utv7P244FD1Vp/OslPl1L+YpIPJ/nOy27b5GWxSVJK+R+SpNb6I6WUryqlfPXajwlgEM8mebSUsgvMJ+On6ACzOTv5eW+S70jypVLKt678kOC+zhrsj+996TO58/4il2r13WL/uyR/KckvnX3pVpIP1Fr/wWoPCg5w9rE835bkP0ry47XWv7LuI4KLlVKeSvLtSX45yWuS/AXvFkvrzgb0p5N8c5K/muTDtdZ/v+6jgquVUr42d17b/nNnX3owyXO11u9d7UHBAUopfyTJh5L8wtmXvibJD9da/8/L/kyTcQkAAEBfmr0sFgAAgH6ISwAAAE4mLgEAADiZuAQAAOBk4hIAAICTiUsAAABOJi4BAAA4mbgEAADgZP8/4R8/n4MU2f4AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import numpy\n", "from matplotlib import pyplot\n", "pyplot.rc(\"figure\", figsize=[15, 6], dpi=75)\n", "\n", "xloc, yloc = numpy.mgrid[-2:3:100j, -1:3:100j]\n", "density = distribution.pdf([xloc, yloc])\n", "pyplot.contourf(xloc, yloc, density)\n", "\n", "pyplot.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Model fitting\n", "\n", "`chaospy` supports Gaussian mixture model representation, but does not provide an automatic method for constructing them from data.\n", "However, this is something for example `scikit-learn` supports.\n", "It is possible to use `scikit-learn` to fit a model, and use the generated parameters in the `chaospy` implementation.\n", "For example, let us consider the [Iris example from scikit-learn's documentation](https://scikit-learn.org/stable/auto_examples/mixture/plot_gmm_covariances.html) (\"full\" implementation in 2-dimensional representation):" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "execution": { "iopub.execute_input": "2021-05-18T10:57:26.365938Z", "iopub.status.busy": "2021-05-18T10:57:26.365618Z", "iopub.status.idle": "2021-05-18T10:57:26.465942Z", "shell.execute_reply": "2021-05-18T10:57:26.465593Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[5.006 3.428 ]\n", " [6.5464 2.9495]\n", " [5.9171 2.778 ]]\n", "[[[0.1218 0.0972]\n", " [0.0972 0.1408]]\n", "\n", " [[0.3874 0.0922]\n", " [0.0922 0.1104]]\n", "\n", " [[0.2755 0.0966]\n", " [0.0966 0.0926]]]\n" ] } ], "source": [ "from sklearn import datasets, mixture\n", "\n", "model = mixture.GaussianMixture(3, random_state=1234)\n", "model.fit(datasets.load_iris().data)\n", "\n", "means = model.means_[:, :2]\n", "covariances = model.covariances_[:, :2, :2]\n", "print(means.round(4))\n", "print(covariances.round(4))" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "execution": { "iopub.execute_input": "2021-05-18T10:57:26.469079Z", "iopub.status.busy": "2021-05-18T10:57:26.468588Z", "iopub.status.idle": "2021-05-18T10:57:26.613160Z", "shell.execute_reply": "2021-05-18T10:57:26.612792Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA5MAAAFyCAYAAABoaI6HAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAuJAAALiQE3ycutAAAkBUlEQVR4nO3dfax1WUHf8d8zL7xEZuxM6UzNgDAUWEAztmhiRkjGqQklhCCByIQaU6nUVDDxJf4BoYIYZYokTVtNYy0xTSQxgRIUowkTY2NIi6DJaEsirum0EJhamI6jzjRCMwy3f9xz5T73nnvP2edtr7X355OQmXvPvs/dPOfsvdZ31j77XDk6OgoAAAAMcc3YOwAAAEB/xCQAAACDiUkAAAAGE5MAAAAMJiYBAAAYTEwCAAAw2HXrblhKuTbJfUkeqrW+acnjn0zylVPfemOt9Ytb7yEAAADNWTsmk7wryZ8kecYFj3+s1vrurfcIAACA5q0Vk6WU1yd5KMlnk9x9wWZ3lFLeluQpST5Xa/3ATvYQAACA5qyMyVJKSfKdtdYfLaW86ZJN31dr/VQp5UqSD5VSnlZrff/ZjR5+/ImjzXcXAACAXbnlhuuvbPqzV46OLm+7Uso/z/Fq45eTfFuS5yX5j0n+Ta31yxf8zFuSfFet9Q1nHxOTAAAAbdgmJleuTNZa33Py74uVybtrre8tpVxXSnl2rfULpZQXJXl5rfWXF5u+IMkDm+4UAAAAbRtyN9dXJXlNkttLKT+U41i8N8mdSR5L8upSym1JnprkphzfsAcAAIAJWnmZ6665zBUAAKAN21zmes0udwQAAIB5EJMAAAAMJiYBAAAYTEwCAAAwmJgEAABgMDEJAADAYGISAACAwcQkAAAAg4lJAAAABhOTAAAADCYmAQAAGExMAgAAMJiYBAAAYDAxCQAAwGBiEgAAgMHEJAAAAIOJSQAAAAYTkwAAAAwmJgEAABhMTAIAADCYmAQAAGAwMQkAAMBgYhIAAIDBxCQAAACDiUkAAAAGE5MAAAAMJiYBAAAYTEwCAAAwmJgEAABgMDEJAADAYGISAACAwcQkAAAAg4lJAAAABhOTAAAADCYmAQAAGExMAgAAMJiYBAAAYDAxCQAAwGBiEgAAgMHEJAAAAIOJSQAAAAYTkwAAAAwmJgEAABhMTAIAADDYdWPvAEAPHnr0y3v5c59189P38ucCAOybmARmb1+huM3vFpkAQOvEJDAbY0bjUGf3VVwCAK0Rk8Dk9BSN6zr5/yQqAYBWiEmge1OMx4uISgCgFWvHZCnl2iT3JXmo1vqmJY/fk+SuJI8s/tx31lqPdrSfAEnmFY6XEZUAwNiGrEy+K8mfJHnG2QdKKbcluTdJqbU+WUr5YJLXJfnITvYSmC3xeLmHHv2yoAQARrFWTJZSXp/koSSfTXL3kk1emeT+WuuTi68/nuS1EZPAQOJxOKuUAMAYrlm1QSmlJPnOWuv7L9ns1iSPnfr6scX3AFZ66NEv//X/2Jy/PwDgkFbGZJLvSfIXpZS3J3l1kjtKKW8vpZz+T+BfSnLjqa9vXHwP4JzT8SiAdsvfJwBwKCtjstb6nlrrT9Va35vkt5J8evHvT5RSnr3Y7L4k37q4SU9yfCOej+5lj4EuiUcAgGm5cnS03g1XSymvSvJPk9ye5N8neSDJvbXWOxeP35Pj91M+uviRpXdzffjxJ9zhFWZCOI7H+ycBgHXccsP1Vzb92bVjclfEJEybgGyDmAQA1rFNTA75aBCAc8Rjm3xkCACwb2IS2IiIBACYNzEJrE1AAgBwQkwClxKQAAAsIyaBcwQkAACrrPycSWA+fA7ktHguAYB9sjIJMyc4AADYhJiEGZpiQH7mkcc2/tkXP/PGHe4JAMA8iEmYkZ4jcptY3ObP7jk0fc4kALBPYhImrseA3Gc4DnV6X3oOSwCAXROTMFG9RGRL4bjKyb6KSgAAMQmT0npA9hSOl/nMI48JSgBg9sQkTEDLETmVgDzLKiUAMHdiEjrWYkRONR4v0uoqpZvvAAD7JiahM60F5NziEQCAY2ISOtFSRArIq7W6OgkAsE9iEhrXSkQKyMu1FJQucQUADkFMQqNaiEgBCQDARcQkNGbsiBSQfbMqCQAcipiEBowdkImIBABgGDEJIxo7IgXktFiVBAAOSUzCCMaMSAEJAMAuiEk4IBHJvliVBAAOTUzCAYwVkQISAIB9EZOwRyJyuU98/vGNfu5l33zDjvdkGqxKAgBjEJOwB2NEZGsBuWkwrvtnCstjQhIAGIuYhB2aa0TuIxyH/M4WwvLFz7xx7F0AADgoMQk7MLeIHCMeL/OJzz/eRFAemlVJAGBMYhK2dOiQHCMiW4vHZeYWlEISABibmIQNHTIiBeR6xgpKl7gCAHMkJmGgKUdkjwF51hxWKK1KAgAtEJOwpqlG5BQCckyHXpUUkgBAK8QkrDDFiBSQfRKSAEBLxCRcQESyyiFXJYUkANAaMQlLHCokDxGRArJ/QhIAaJGYhFOmEpECcv/cwRUAmDsxCRGRDOPyVgAAMcnMiUhaJiQBgJaJSWbrECEpIg9v358xeahVSSEJALROTDI7va9GCsjxCEkAgK8Tk8yGiJy+fa5KCkkAgKuJSWah50taReT4hCQAwHlikkkTkfOxr1VJIQkAsJyYZJJ6vqRVRA6375vu7JuQBAB6JCaZnF5XI0Vkew6xKikkAYBeiUkmo9eITNoMyU9/9tGrvr7j9ptH2pPL9Xx5q5AEAHp25ejo6KC/8OHHnzjsL2QWeg3JsSPybDCuo6WoFJIAANu55Ybrr2z6s2KSronI4TYJyGXGjkohCQCwvW1i0mWudGvfITmliNxVQE6dkAQAWN/KmCylXJPkN5L8fpLrk7wkyQ/UWv/yzHafTPKVU996Y631izvcV0hiNXJdUw7IfaxKCkmgV4e6g3mvnHthf9ZdmfxkrfVnk6SU8sEk35/k589s87Fa67t3uG9wjtXI1aYckUmfIWkiA1xECO7fNn/Hzt9wuZUxWWv9WpKTkLw2yXOSvH/JpneUUt6W5ClJPldr/cAud5R563E1UkTu3q5D0moksCuicJrWfV6d65mrtd8zWUp5RZIfTvJ7SX53ySbvq7V+qpRyJcmHSilPq7Uui04YxGrk5eYQkYmQBA5LHDLEOq8XYwJTNPhurqWUdyS5vdb6g5ds85Yk31VrfcPZx9zNlXVZjbzcmBF5yDu59nhZa2LSAK0SibTImMGY9no311LKS5I8r9b6m4tvPZjkVaWU65J8U631C6WUFyV5ea31lxfbvCDJA5vuFFiNvNhcViITIQkMIxTp1UWvXeMJrVvnMtf/l+SflFJeuvj6pUl+LMldSe5NcmeSx5K8upRyW5KnJrkpybt2vrfMQm8hObeIPNSqpJAETohE5mrZa984Q0sGX+a6LZe5chGXtV5MSG7PHVuhXWIRtmMMYht7vcwVDsFq5HKtRGTSb0hajYQ2CEbYn7PHl3GJQxGTjE5ILicktyck4bAEI7RBXHIoYpLR9BaRyWFCsqWITA4Tkj1e1poYnJkfsQh9Epfsi5hkFL2FpNXI/bEaCe0RjTBtp49xYxrbEJMcnJBcTkhuT0jC+gQjkAhLtiMmORh3a71YKyHZa0QmQhIuIhqBdQlLhhKTHERvq5GJkNyXHlcjE4Mq7RONwC4JS9YhJtk7IXmxFkKy14hMrEYyX8Kxf/sYuw7hUP8Bj7acnHOMiZx15ejo6KC/8OHHnzjsL2RUvYXkoSIyGT8kDxWRiZCETYnGdvQaf60Ro9NgfJyWW264/sqmP2tlkr3w/sh2icjVDJIcmmg8HFE4rnX//kVn26xUcsLKJDsnJFcbY1XykBGZCEm4iHDcPYE4P2KzLcbOvlmZpBlCsk1WI9djMGTXhON2RCIXuei1ITLHYaVyvsQkO9Pb+yOTaYfkoVcik35D0uDHLgjHYYQi+7DsdSUwD+ehR79sTJ0ZMclOCMk2jBGQSb8RmQhJNiMcVxOLtEJgHpZVynnxnkm20uNlrSfGjsldvm9yShGZCEnaIhyXE4tMjcDcPeNsH7xnklEIyfGMFY8nRCRTJRyvJhiZk9Ovd2G5Gy57nT4rk2xESO7HstXKscPxrJ4vaU2EJFcTj4IRVhGW2zP2tm2blUkxyWBCcp6sRtK7uYejaITtCcvNGYfbJSY5GCE5P71HZGIAm6s5xqNg3L1Wx459nZtZn7AcznjcJu+Z5CDmODGbMxFJT+Z2fhKNl2s1AHdpm/+PQnQ3To5DUbk+76GcHiuTrOVQEzWrkuPb5yRDSLIrc4pH4egcPgbBOZyoXJ8xui0uc2WvhOQ8iEhaNod4nFs0Ojf3TWxeTFSuZqxui5hkb3oPycSEZZWpRGRicJqSqcfj1MPReXeeBObVROXljNnt8J5J9mIKIcnFRCQtmXI8TvEcJxZZ5uzrYu5x+ZlHHhOUTJ6VSZaaSkia8Jy378FdSLKOqcbjVMLRuZN9mWtgisrzjN/tsDLJTk11kjdnhxi8RSSXmeJ5pfdwFIyM4fTrbk5haZXyPHd2nQYxyWh6n4j1QEQylqnFY8/nK9FIq05em3OJSkHJFIlJrjKVy1vnbooRmQjJlk0pHns8PwlGejan1UpBydSISf7alCaDc3SoAVhEkkzrfNFTPIrG5NOffXTsXdjKHbffPPYuNG0OYSkomRI34CHJYSeGh5y4zWHiNeWITIRkS6YQkMKxHb1H4SEIz+kGZeKmPIkxvhVuwMNWphqSU3bIwVVEzpd4PIwpRaNA3K2L/j7nFJlze18l9EZMztwUJotzceiBVETOzxTOB63HY8/hKBTbMcfInGJUutyVKRCTM3boiWPrk7wWjTFoish56T0gWz6v9BaOYrF/Z5/DKcblJz7/+KSCEnrnPZMzNpeY7G1CN9YgOeZ/HRWShyMe96OX84xgnK+pheVUgnLOK5PG/nZs855JMTlTY0wox5wEtjzRG3tAFJHT13NAthiPLZ9PEsHIalMJy7HHz10Qk7RATDLI3ELyRCsTwBYGv7EHLwPIfonH3Wrl3HGWaGQXeg/LFsbUbYw9Ho/JXKAd7uYKazg94BxqctjaIDf2oGXg2J9eA7K1eGwxHEUj+/Tpzz7afVDSH/OB6bAyOTNjTThbmzCetovJY2vReNbYEZkYOPahx4Bs7VzQWjzOLRwffOBLo/3u57/w1tF+d6t6jcrWx+DLtDA+j8GcoC1WJllLjxPPQ+h5ELpMKwOUAWN3ej2GWwrIVuJxStE4ZhBuY+h+zyE+T16XvUYlfTAvmBYxCRMjIqelx4BsJR6F43Z6jcR9Wfb3MdXAdOkrsC4xORM9TkhZXysBmYjIXejxeG0hIFuIx57CUSxub8qBKSjZB3OE6RGTMzD2xLSFSeZUicjpGPs4HaqF41o8Xk4sjuP03/tUwpL9aGkMPwTzhGkSk9CZ1gYfg8PmBORwYwZkq+EoGtvVe1hanWRXzBWmS0xOXAuT1Rc/88YmJqE9ay0gEwPDplo4JocY+9gVj18nGvt28vz1GJXsXovjOmxCTELDWhxsRORwPQXkXOOxpXAUjdMmKpkb84Zp8zmTE9bSBHbsCWpPWgzIxGAwVEvH3ypjH59jBGQr8Sgc562XoGz5UtceP96r1XF+H8wd+uBzJmmeS10v1/LAYiBYn4Bcz1zjUThy1oMPfKmboGQ3Wh7vd838YR7EJIyk9QHFILAeAbmeQwfk2PEoHJkKq5K70/q4v0vmEPOxMiZLKdck+Y0kv5/k+iQvSfIDtda/PLPdPUnuSvLI4s99Z63VJa0jaXGCO/fVyV4GEQPAai0eXxcRkPs3pXD88wcfGHsXrnLT81849i7sndXJeehlDrAtc4j5WXdl8pO11p9NklLKB5N8f5KfP3mwlHJbknuTlFrrk4ttXpfkIzveXzo3t6DsafAwAFxOQK4mHtvVWiSua9l+zyEwWa2nVcme5gLbMI+Yp5UxWWv9WpKTkLw2yXOSvP/MZq9Mcn+t9cnF1x9P8tqISZaYclD2OGA4+V9MQK42h4DsJR57DcahTv//FJb71eolrkKyPeYS87X2eyZLKa9I8sNJfi/J7555+NYkp2cyjy2+B0tNJSh7HSSc9C/XS0TOISDF49XmEozrOvv3IS53R0hur9c5wlDmFPO2dkzWWn87yW+XUt6R5BeT/OCph7+U5KWnvr5x8T1G0MtE+PRJtoewnMKg4IR/sV6OGwG5ey3Go2jczJ8/+ICgnDAh2RZzCpL1bsDzkiTPq7X+5uJbDyZ5VSnluiTfVGv9QpL7kryjlHLt4lLXu5J8cF87zfScnHRbicqpDQJO+MsJyMtNNSBbi0fhuFuCcnstrkoKybaYV3DiytHR5TdcLaX8nSTvS/JHi2+9NMl7knxjkntrrXcutrsnyd1JTmYES+/m+vDjT7jD6571MkFexz4n0VM/2TvRX6yXY2SMiDxUQM41HoXj4fQSlK3dyVVIbm7q84oT5hfTc8sN11/Z9GdXxuSuicn962WizH44yS/Xy3EhIHejhYAUjuPqISaF5OV6ichkHiFpfjFd28Tk2u+ZBNrlBL+cgLzcISJyTgEpHumZkNzMHCIyMc/gYmISOubkvpyIvNiUAlI80quWViRF5ObmEJLmGawiJqEzTuwX6yEipxqQyWEicsyAFI/9aPkSVyF5sV5CUkTC14lJ6ICT+sV6CMhkuhEpIGF9rYSkiNzMHCIyMedgGDE5Qc+6+endTLC5mJP55Xp4jQvIzY0VkOKxf62uSgrJ5YRkO8w72ISYhIY4kV+uh4BMDh+RLmPdjoCcjhZDUkQuJyLbYe7BNsQkjMxJfDURuZxVyM0JyOlpLSRF5HIish3mH+yCmJwol7q2zQl8PT28hl3KuhkBya60FpFJGyEpIjcjImEYMQkH4uS9PhF53hQCMhGR7FZrISkizxOR7TAPYR/E5IRZnRyfE/cwPbxeReRwApJdaykiWwjIRERuSkTCdsTkxAnKw3LC3kwPr1EROZyIZNdE5HktRWQvAZmISNgVMTkDgnK/nKw308trcmoROcVLWQXk9InIq7UUkEk/ETmHgEzMSzgsMTkTgnJ3nKS308vrUEQOJyLZpZYCMhGRy4jItpifMAYxOSOCcjNOzrvRw2tvindmFZH0pqWIFJDn9RKQiYiEQxCTMyMoV3NS3r0eXnNWIocTkexKSwGZiMhleonIuQRkYr5CG8TkDJ0++fQwyd8nJ+L96uH1JSKHc2MddkFAnicgNzeXiDRvoTVicuZOTko9TPq35QR8OD28nqYWkYm7s9I+AXmegNzcXAIyMYehXWKSJNOKSifc8fTw+hGRmxGSbEpAntdaQCb9RKSAhLaISa6y7MTVaiA4ybaj1dfIaSJyM2NEZCIkeycgzxOQ2xGR0CYxyUpjBKYTaR9E5HJCcjtCsj+txWMiIC8iINtk3kOvxCQbcdJDSJ4nIrcnJPvRWkC2EI+JgNzWnAIyMZ+if2ISGEREnneIiEyEJONqLR4TAXmRnuIxEZDQMzEJrK31kHRJ63bGDEna0mI4JuLxMgKybQKSqRKTwEqtR2RiNXJbQnLexONqAnJ7AhKmR0wCFxKRywlJetZqOJ5oJSDF427MLSATEcm8iElgKSG5nJCkJ62HY9JOPCYCchfmGI+JgGS+xCRwFRG53KEiMjlcSLbmpue/0E14tiAchxOPuyEgYb6uHB0dHfQXPvz4E4f9hcDahORyUw3JVlclBeXleojGE+JxPQKyHwKSKbrlhuuvbPqzViaBJELyIlMNyZZZoewrGE8Tj+sRj30RkHAxMQkz10NEJkJybs7G1BTjstdgPNFaOCbicZfEI7AOMQkz1kNIjhGRiZBsTW9x2XsoniUch+kxHhMBCQwnJmGmhOTFDhmSbOayWNt1aE4tDFdpMRwT8bgPc47HREDCLohJmCEhebFDh6RVyd2bW/xtQzhuRjz2S0DCbolJmJEeIjIZLyTn5vkvvLXZO7qyW61G44mW47HXcEzEYyIeYd/EJMyEkFxtjpe3CsrpEY7bEY/9E5BwOGISZkBIrjbHkDwhKPvTejCeEI77JR6PiUcYj5iEiROSrENQtqmXaEzaD8dEPE6FeIR2iEmYsF5CcmxzXpU87SRcROVh9RSMJ4Tj/gnHrxOP0C4xCRPVU0halWzL6bgRltvrMRZP9BCNSf/hmIjHswQk9OHK0dHRQX/hw48/cdhfCDMkJIcZe2Wyp48HEZdX6zkUzxKOhyMczxOPMJ5bbrj+yqY/a2USJqankOTYHbff3E1QXhRPU4vMKUXiWb1EYzKNcEzE4zLiEaZBTAKjaWFVkt1YJ74OHZxTDsJ19BSNiXCcOvEI0yQmYUKsSvbrZOLfywrlJuYed/vQWzCeEI7TJx5hHsQkTERvIWlVcrmeLnnlcETj+ITjxYQjzJeYhAnoLSS53BxWKTmv12BMphWNiXBcRTwCJ8QkMHsv++YbRr+j6zKiclp6jsXThOO8CEfgMj4aBDrX66pki5e5thiUy4jLNk0lFk9MLRoT4bgO8Qjz46NBAGbkdLQIy8OYWiieJhrnSzgC27IyCR3rdVUyaXNlMulndfIyAnN9U47E06YYjCeE43qEI3CRva5MllKem+S9Sf4wyU1JbkzyI7XWr57Z7pNJvnLqW2+stX5x0x0DGEOr758c4rJAmnpoziUOLyIaOSEegUNY5zLXm5N8uNb64SQppXwsyfcm+ZUz232s1vru3e4ewOGdTMh7j8plNo2tfUbo3ANwqCkHYyIaNyEcgbGsjMla6/1J7j/1rWuSLJth3VFKeVuSpyT5XK31A7vZRWCZni9x7cWUo3IowXdYUw/GRDRuSjgCLRl0A55Synck+askH13y8PtqrZ8qpVxJ8qFSytNqre/fxU4C0/PiZ97Y7PsmzxKV7NocYvGEaNyMaAR6sPYNeEop357krUl+qNb6lRXbviXJd9Va33D2MTfggd2YyspkL0F5mqhkHXMKxkQ0bkM4AmPa+0eDlFJekeQ1Sd6c5Eop5Q1Jfi3JN9Vav1BKeVGSl9daf3nxIy9I8sCmOwXQstORICznaW6heJpo3I5wBKZknbu5fluSX0/yB0l+J8fvmbw/yZ8luTfJnUkeS/LqUsptSZ6a47u+vms/uwxMSU+Xuy4jLKdnzqF4QjBuTzQCc+BzJqFTU7nM9UTPQXkZgdkWoXg10bg90Qj0bpvLXMUkdGpqMZlMNyiXEZm7IxAvJxh3RzgCUyQmYYamGJMn5hSVF5ljbIrCzQnG3RKNwJyISZihKcdkIii3te8YFX6HJRb3QzQCiEmYrakHZSIqmQ/BuB+CEeBye/9oEICxnEywRSW9E4v7JRoBDk9MAl04PREXlrRILO6fYARoi8tcoXNzuNT1MsKSQxGLhyMaAQ7HeyZh5uYelCeEJZsSiocnGAHa4D2TADkfBOISkTguwQgwbVYmYSKsTq5HYE6DSGyHYATom8tcgSSCchsic1zisG2CEWC6xCRwFVG5H4JzNVHYL8EIME9iEjhHULaj5QgVf/MhFgFYRkwCFxKVMB+CEYCh3M0VuNCzbn66oISJEIsAtERMwgycTEBFJbRNLALQEzEJMyIqYVxiEYApEZMwQ6cntMISdkMoAjA3YhJmTljCakIRAM4Tk8BfE5bMkVAEgM2ISWApYckUCEUA2B8xCay0bEIuMBmTSASA8YlJYCNnJ/Pikm0JRADoi5gEduKiEBCZ8yYQAWC6xCSwV5fFhNDsizAEAE4Tk8Bo1o0T0bk7ghAA2BUxCTRvmwDqNURFHwDQOjEJTJooAwDYj2vG3gEAAAD6IyYBAAAYTEwCAAAwmJgEAABgMDEJAADAYGISAACAwcQkAAAAg4lJAAAABhOTAAAADCYmAQAAGExMAgAAMJiYBAAAYDAxCQAAwGBiEgAAgMHEJAAAAIOJSQAAAAYTkwAAAAwmJgEAABhMTAIAADCYmAQAAGAwMQkAAMBg163aoJTy3CTvTfKHSW5KcmOSH6m1fvXMdvckuSvJI4s/95211qNd7zAAAADjW2dl8uYkH661/lyt9e1Jnpfke09vUEq5Lcm9SX601vruJC9I8rod7ysAAACNWLkyWWu9P8n9p751TZLHz2z2yiT311qfXHz98SSvTfKRXewkAAAAbRn0nslSynck+askHz3z0K1JHjv19WOL7wEAADBBa8dkKeXbk/yzJG+stX7tzMNfyvF7KU/cuPgeAAAAE7RWTJZSXpHk+5K8OclXSylvKKVcV0p59mKT+5J8aynl2sXXd+X86iUAAAATceXo6PIbrpZSvi3H74H8g8W3rsnxeyh/I8m9tdY7F9vdk+TuJI8utlt6N9eHH3/CHV4BAAAacMsN11/Z9GdXxuSuiUkAAIA2bBOTg27AAwAAAImYBAAAYANiEgAAgMHEJAAAAIOJSQAAAAYTkwAAAAwmJgEAABhMTAIAADCYmAQAAGAwMQkAAMBgYhIAAIDBxCQAAACDiUkAAAAGE5MAAAAMJiYBAAAYTEwCAAAwmJgEAABgMDEJAADAYGISAACAwcQkAAAAg4lJAAAABhOTAAAADCYmAQAAGExMAgAAMJiYBAAAYDAxCQAAwGBiEgAAgMHEJAAAAIOJSQAAAAYTkwAAAAwmJgEAABhMTAIAADCYmAQAAGAwMQkAAMBgYhIAAIDBxCQAAACDiUkAAAAGE5MAAAAMJiYBAAAYTEwCAAAwmJgEAABgMDEJAADAYGISAACAwcQkAAAAg4lJAAAABhOTAAAADCYmAQAAGExMAgAAMNh162xUSnlWkp9O8rJa64sv2OaTSb5y6ltvrLV+cftdBAAAoDVrxWSSu5J8JMnLL9nmY7XWd2+9RwAAADRvrZistf5qKeXuFZvdUUp5W5KnJPlcrfUD2+4cAAAAbVp3ZXId76u1fqqUciXJh0opT6u1vv/sRrfccP2VHf5OAAAARrCzG/DUWj+1+OdRkv+U5B/u6s8GAACgLRvHZCnlulLKsxf//qJSyptPPfyCJA9su3MAAAC0ad27ub4iyRuT/M1Syk8m+cUkfy/JvUnuTPJYkleXUm5L8tQkNyV51172GAAAgNFdOTo62ssfXEq5Nsl9SR6qtb5pyeP35PgusY/kOGrfubhElpGs8Zz5+JeGlFJ+PcnfOPWtH6u1/tGZbRxnDVnzOXOcNaSU8reS/HiS/5Pj/4j6mVrrz53ZxnHWkDWfM8dZI0opz83x26M+v/jWdUm+odb60jPbOc4aMeA5c5w1pJTy+iRvSPJfk7wsyU/WWv/bmW0GH2e7vAHPWe9K8idJnnH2gcUK5r1JSq31yVLKB5O8LscfP8J4LnzOFnz8S1v+6LLnw3HWpEufswXHWVv+Q5K31lo/X0q5JsnfPf2g46xJlz5nC46zdjye5C211vuSpJTyfUm+4fQGjrPmrHzOFhxnbfmlJN9Za/3jUspbk7w7yetPHtz0ONtLTC7K96Ekn01y95JNXpnk/lrrk4uvP57ktXFSGM0az1ni419a8+xSyjuSfC3J/03y72qtXz31uOOsPaues8Rx1oxSyq1J/n6SV5ZSnpHkG5P86zObOc4asuZzljjOmlFr/bMcXxV14h8n+e4zmznOGrLmc5Y4zlrzv5PckuSPF//8gzOPb3Sc7exuridKKSXH1XvuY0FOuTXH77M88djie4xgzecsOf74l59L8rNJvruU8oP73zsu8UtJ/kWt9b1Jnpvkp8887jhrz6rnLHGcteQ5SW5L8mCt9V8leTDJr5zZxnHWlnWes8Rx1qTFPTr+S631K2cecpw16pLnLHGcteYtSX6mlPIvk7wiya+deXyj42znMZnke5L8RSnl7UleneP/KvH2UsrTT23zpSQ3nvr6xsX3GMc6z5mPf2lMrfX3T13H/js5/3w4zhqzxnPmOGvLXy7++YnFP/9zkn9wZhvHWVvWec4cZ+16a5J/u+T7jrN2XfScOc4aUkr520k+nOS1tdafSPK+nF9x3Og423lM1lrfU2v9qcV/ef+tJJ9e/PsTJx8lkuOl8W9d3PAlOX6j50d3vS+sZ53nzMe/tKWUckMp5fQdk1+Q5IHTH9kTx1lT1nnOHGfN+R85vvT/eYuvnxPHWetWPmeOszaVUr4lyZ/WWh9ZfO04a9xlz5njrDk357j7/mLx9Z8meeoujrO93YCnlPKqJK9Jcnsp5Ydy/AK6N8mdtdb/tXjf0C+UUh5N8t9zfqmVA7vsOYuPf2nNE0m+pZTyM0m+muRFSX4ixwe+46xNK5+zOM6aUmv9ainlHyV5RynlM0lekuRNcZw1a53nLI6zVv14kvec+tpx1r4Ln7M4zpqyuOnOLyT5pVLK/8zxe8vfnB0cZ3v7aBAAAACmax/vmQQAAGDixCQAAACDiUkAAAAGE5MAAAAMJiYBAAAYTEwCAAAwmJgEAABgMDEJAADAYP8fdbYL9e8l7ggAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "distribution = chaospy.GaussianMixture(means, covariances)\n", "\n", "xloc, yloc = numpy.mgrid[4:8:100j, 1.5:4.5:100j]\n", "density = distribution.pdf([xloc, yloc])\n", "pyplot.contourf(xloc, yloc, density)\n", "\n", "pyplot.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Like `scikit-learn`, `chaospy` also support higher dimensions, but that would make the visualization harder." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Low discrepancy sequences\n", "\n", "`chaospy` support low-discrepancy sequences through inverse mapping.\n", "This support extends to mixture models, making the following possible:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "execution": { "iopub.execute_input": "2021-05-18T10:57:26.615714Z", "iopub.status.busy": "2021-05-18T10:57:26.615225Z", "iopub.status.idle": "2021-05-18T10:57:26.733413Z", "shell.execute_reply": "2021-05-18T10:57:26.733665Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA4sAAAFvCAYAAADjWV7eAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAuJAAALiQE3ycutAABcLElEQVR4nO3df4yl113n+U/dqlt9u9rd7rbTC3GU0M4KHQQOs9gSm1nUgV0pNiGCiOyA2Emzk10ywk60Y6P5A5RdQmchDIu0wgnZabNohJbYf0AYfmhYYXsEAjuD+RXD7joDDzDYMw4OVtt0u39W1a26d/+ouuXbt54f55znnOc5z3PfLylyuurWfX6d58f3Od/zPSvT6VQAAAAAAMwbtL0CAAAAAID0ECwCAAAAAA4hWAQAAAAAHEKwCAAAAAA4hGARAAAAAHAIwSIAAAAA4JC1phd48eJVr7k6Tp7c0OXLN0KvDtAo2jG6jjaMPqAdo+towwjp9OnjK0W/60zP4krhJgDdQTtG19GG0Qe0Y3QdbRhN6UywCAAAAABoDsEiAAAAAOAQgkUAAAAAwCHWBW6MMauSnpL05SzLPpzz+z+QtDn3o+/Lsuzvaq8hAAAAAKBxLtVQPyHpLyTdVvD7J7MsO197jQAAAAAArbMKFo0xH5T0ZUkvSvq2go+9yxjzw5LWJb2UZdnngqwhAAAAAKBxlcGiMcZI+tYsyx42xny45KM/nWXZHxpjViT9sjFmlGXZzy9+6OTJDa9yv6urA506teH+h0BCaMfoOtow+oB2jK6jDaMpK9PptPQDxpj/WXu9hTcl3SfpnZI+L+nTWZbdLPibhyT9N1mWfc/i7y5evFq+wAKnTm3o0iUmH0W30Y7RdbRh9AHtGF1HG0ZIp08fL+zKq+xZzLLsU7P/v9+z+G1Zlv2UMWbNGPP2LMteNsZ8naRvybLsX+1/9Gsl/WXN9QYAAAAAtMSlGur7JH2npLuNMQ9qLxj8SUnvlnRF0vuNMW+TdETSKe0VxAEAAAAAdJB1sJhl2W9J+q2FH797/3evSPpgwPUC0KLBYEWTiVfGOAAAAHrCZeoMAD23PhpqNBpqvDvRcHWgzc2xtjfHba8WAAAAWkCwCEDSXqD4yrVtPfTYc3rp9Ru6+y3H9C8/dK/uum2dgBEAAGAJDdpeAQBpGI2GeujxL+ql1/eqq7342nV99InnNRoNW14zAAAAtIFgEYAGgxWNdycHgeLMi69d13h3osHAY3JUAAAAdBrBIgBNJlMNVwc6c+etE/ze/ZZjGq4OKHYDAACwhAgWAUiSNjfHunDuvoOAcTZmcZPxio2g9xYAAKSGAjcAJEnbm2Pdddu6nnrkPVRDbRAVaAEAQKoIFgEc2N4PVAaDFd0k9TQ6KtACAICUkYYK4BDGKDaDCrQAACBlBIsA0AIq0AIAgNQRLAJAC6hACwAAUkewCAAtoQItAABIGQVuAKAlVKAFAAApI1gEgBZRgRYAAKSKNFQASABjFAEAQGoIFgEAAAAAhxAsAkgW00cAAAC0hzGLAJKzPhpqNBpS9AUAAKBFBIsAkrI+GuqVa9t66LHn9NLrNw6mk7jrtnUCRgAAgAaRhgogKaPRUA89/kW99PoNSdKLr13XR594XqPRsOU1exPpsQAAYBnQswggGYPBisa7k4NAcebF165rvDvRYLDSatVQ0mMBAMAyIVgEkIzJZKrh6kBn7ty4JWC8+y3HNFwdtDoPIemxAABg2ZCGCiApm5tjXTh3n87cuSFJB0HZZssBWRfSYwEAAEKiZxFAUrY3x7rrtnU99ch7kkn3TD09FgAAIAaCRQDJ2d4PDgeDlVZTT2dSTo8FAACIhTRUAMlKqbcu1fRYAACAWOhZBAALKabHAgAAxESwCACWUkuPBQAAiIk0VABwlFJ6LAAAQCwEiwAAAACAQwgWAQAAAACHECwCAAAAAA4hWAQAAAAAHEKwCAAAAAA4hGARAAAAAHAIwSIAAAAA4BCCRQAAAADAIQSLAAAAAIBDCBaBwAaDlbZXAQAAAKhtre0VAPpifTTUaDTUeHei4epAm5tjbW+OG1+PwWBFk8m08eUCAACgXwgWgQDWR0O9cm1bDz32nF56/Ybufssx/csP3au7bltvLGBMJVgFAABAP5CGCgQwGg310ONf1Euv35AkvfjadX30iec1Gg0bWf4sWH3g0Wd0z/mn9e2fflavXNvWekPLBwAAQP8QLAI1DQYrGu9ODgLFmRdfu67x7qSRMYxtB6sAAADoH4JFoKbJZKrh6kBn7ty45ed3v+WYhquD6OMHUwhW+4D9BAAAcCuCRSCAzc2xLpy77yBgnI1Z3GxgzGDbwaqr1IKy9dFQJ05u6MixIzpxcoPUXQAAgH0UuAEC2N4c667b1vXUI+9ppcDMLFidpaI2GazaSrEATwqFiQAAAFK1Mp022+tw8eJVrwWeOrWhS5duVH8QaFnZ1BUx23GKwdjMQVC2EMy2HZSdOLmhBx595pYU3rvfckxPPnxWVy5373rTxLQpXIvRB7RjdB1tGCGdPn28MO2LNFQgsLbSPrc3x7py+Ya2rm/pyuUbyQSKUvgCPCFSWfs01pNUWgAAEANpqEDPpDhGsSoos13nkL2n82M9F3sWh6sD3UxsPxYhlRYAAMRCzyKAqEIV4Ikxl2SbhYlCYdoUAAAQC8EigOhCBGUxgqL5wkQvnL9fTz58tlM9cn1KpQUAAOkhDRVAdHWrxYZMZc1bt+3NsQaDlc6kns70JZUWAACkyTpYNMasSnpK0pezLPtwzu+/V9J7JL22/70/mmUZTyoAJNULypoIilIb62mrC9OmAACAbnLpWfyEpL+QdNviL4wxb5P0k5JMlmW7xphfkvTdkn41yFoCaF2oaRl8v4OgKF/bc3wCAID+sgoWjTEflPRlSS9K+racjzwg6fksy3b3//2MpA+IYBHovFTmbyQoKtblVFoAAJCuymDRGGMkfWuWZQ8bYz5c8LGvknRl7t9X9n92yMmTG1rxqLmwujrQqVMb1R8EEta1drwzmeql12/owblpGS6cu1dn7tjQWovFU1bX9mpzHTs61LGjVP1sUtfaMJCHdoyuow2jKTY9i/9I0mVjzI9Iuk/SO/f//6ezLLu5/5lXJX3T3N+c2P/ZIZcv38j7caVTpzZ06ZLf3wKp6Fo7PnFyQw8uVCB96PHn9eTDZzu1HQina20YyEM7RtfRhhHS6dPHC39XOXVGlmWfyrLsx7Is+ylJ/7ek/2///4+NMW/f/9hTku7dL4Ij7RW6+Y16qw2gTUzLAAAAsNys51k0xrxP0ndK+kZjzIPaCwg/L0lZlv2tpI9L+lljzE9I+itJvxZ+dQE0Zb4C6bxZBdLQ1UMJPgEAANKyMp02Wwzh4sWrXgukux190LV2vD4a6pVr24cqkIacuD6VAjo4LK8CbtfaMJCHdoyuow0jpNOnjxe+sXeZOgPAkoldgfQgGJ0roBM6GIU7AngAACA5pKECWE7bm2NduXxDW9e3dOXyjaBBw2g0POi1lPbGQ370iec1GlHhtC2zAP6BR5/RPeef1rd/+lm9cm1b6xwTAACWDsEiACsxxihSQCc9BPAAAGCGYBFAK5ouoINqBPAAAGAewSLQEak+qNdZr83NsS6cu+8gYJyNWdxkfFwrCOABAMA8CtwAiUu12EiI9YpdQAfuZgH8YgVcAngAAJYPwSKQsCaqheZNj9Dkem3vB4eDwYpu0nPVOgJ4AAAwQxoqkLCYxUbWR0OdOLmhI8eO6MTJDadql7br5ZKiSopjOmJWwAUAAN1BzyKQKJtiI74BVp2eQZv1WltfSzJ11ked/dx1y7rdAABgDz2LQKJiFhup02NZtV5r62tB5+lrq7BPnZ5XAACAPiBYBBIWo1poiOkRytYrVOpsm8EaE9MDAACQhgokLUaxkfmewfmAcdYzaFNkpmi9drZ3NB6u1k6dbaKwT5nRaHiwbOnNgPfJh892Np0WAADAFT2LQOJiFBsJ0WOZt16hUmdjFvapwsT0AAAAe+hZBDoiZLERmx5L217Axc/UnacvZmGf2feX/X2InlcAAIA+IFgEOs43eCqa33B9NKxVybRu6mysYM1lu5iYHgAAgGAR6Ky6Qd3MZCFQDDFWsCgQtRU6WHPdLiamBwAAkFam02ZTqi5evOq1wFOnNnTp0o3qDwIJC9WOD4KfhWCqbgGYEyc39MCjzxzq0Xvy4bO6cjnc+WfTGxoqGJbqbdcyz7OYh2sx+oB2jK6jDSOk06ePFxZkoMANYCmlwiYxCsA0UdjFZTqMUIV96m4XgSIAAFhWBItAhSbm+3MJxGIFdaEqmRbxnbuw7nJjbxcAAEBfESwCJWJPzu4TiMYMfkJMqVGkzekwYm4XAABAXzFmESgRegzffDuuM+4w1pjF2XcXjRX0Hb83GKzoyLEjuuf804d+98L5+7V1fSt6D1/IMZDLjGsx+oB2jK6jDSOksjGLVEMFCsSe7280Gh5U55x970efeF5PPny2MoiJWa0zr5Jp3UArhbkL61ZoBQAAWDakoQIFYqZ7hhh36FsAxnZM42QuUAyRiptKKihjFAEAAOwQLAIlYgU4IQNR28/6FuoJNdZwvjf0hfP368mHzwZJmwUAAEAcpKECJWKme4aeeL6M66T0M6FTcUkFBQAA6A6CRaBCrACnKhANORm87/jIWGMNi7Yr5DYDAACgHtJQAUsxgpi8cYeh53WsOz6yibGGTcxlCQAAADf0LAIJWCwm45ouWvXdVb2DZT16MVNxpTjbDAAAgProWQRqsq0uaiPWxPVFvYPj8a5Vj55v5VUbsbYZAAAA9dCzCHgKPcm7bbqoTzpsXu/geLyrl9/YPFRgp6xHzzcVt6jnMvZclgAAAPBHzyLgIdTcg/OqptO47cTRWmP6FnsHh8PV6D16VWMRY85lCQCAi5CZQkBfECyit2Je9JtMF73woXt16Ua4wHSyP0axTtEbG7YBdRMFdAAAKEKRNaAYaajondDpoYtipk7mpYtOplO9/zPPOk97UcZ2Sow622I7XUfsAjoAABShyBpQjp5F9EqM9NBFsVMn59NFr125qd3JVH/zWvgewLIevbpvWV17LmMW0AEAoAhF1oByBIvolaYu+lWpk2VBnG2AN5lMgwSmRcub79F74fz9evLhs7rrtnVJqh1w+643YxQBdA3j3LqriSEZQNcRLKI3mrzolwVaRT1y66OhJlpx7q3zHdNn0zuY16MXKuBmLCKAPmOcW/dRZA2oxphF9IbtOLxQtvfH1Q0GK7o5mZaOe5DkPSbCZ0yf6xiMydwYxVDjMVMYi8jUGwBiYJxbf8xebC5OI8WLTWAPPYtoXcgevzZ6s2bBSFmPXN3eOtcxfb7LC/2Wta2xiLzxBxAT49z6oyhTiKAf2EPPIloTo2ppW71ZZT1yO5OJplNV9tbZ9ILZjlG0Scct+q4Yb1nzluXS6+fyWd74A4gpZkVstGMxUwjAmwgW0YqYD/QhLvquN/uyFNi1wV4HflF67O76WtCguSodd3jiaOmyYgfcLi8JfF4o2E7ZAQA+mh7ygOYQ5AOHkYaKVjSRwuNz0a+TvliWAlv0u/F4N8pUH3nLu/Che3Xpht2yYqWPukxt4jMNCpXtADSBAl4AlsXKdNrsW5SLF696LfDUqQ1dunSj+oNI3mCwVxH0nvNPH/rdC+fv19b1rVbe7h30di6kX7r0dpb1hK2PhhodHWq88+bvRqOhHnj0mUNvp598+KyuXL5RK51pcV0m06ne/5lnb5mzcX5ZoZSt84mTG6Xb6/vZEH8HO1yL0Qch2nGMoRSALa7FCOn06eOFb9PpWUTjUi1VHaK3s6xHbntzrFXp4Hc72zulvWB1C7TMr8u1Kze1O5neEijOL2swWKnd61bVK+vS61enh5A3/kBzlrm3vq0CXgDQJMYsohWplaq2LVhg29NX9pnZ70rHOa4O9F2f/YKyV6/VHs85W16dcYxVbMaguozzqTMmKIUpO4C+o1ftTYxzA9Bn9CyiFamVqq7q7VxbX4syFUPR2MLP/8nLyl69JulwD6fvm/y64xjL2PbKuvT61ekh5I0/EI/PeGIAQDcxZhGtS6XMeNGYxbffPtLLb2zWGss4k9eOF9/Q7+xO9K5PPq3FXfLC+fs1XB3U7gEMPY7RdQxq7GqoiItrMfowLph2jK6jDSOksjGLpKGidSkEilJx+uJwuJrbazY/FUOdgHdxqo8TJzf0jjvipKbOL2tL0pFjR0rHMdqm3LqkjLpMbcLcV0BamGMQAJYLaahIXpMFFBbTF6uK0Bw5uh4sPXX2gOWTmuqzrJCFhnxSRl3nsWxCisU6UlwnxNGFY51qgbKu68KxB7Cc6FlEstpMQbQpQjNcHejFSzf14IXfd05PLXv7ntfDubM70Xd85ku3fM7lTX7RZ2wKDdl8f9eLyqSY7priOiGOrh3r1AqUdVnXjj2A5cOYRSQpxJyHsdflzJ0bet+nn3Uat3P89qNaWVmxfjCYBWq+Y4RsHkSKPuP7ENO1NLSU2lrK6zTDtTislI91ma4HOSm0464ee6QhhTaM/igbs0iwiCSlVkBh8cFoa2uswdqqdVGX2Xd85dq2HvR4MPB5qHD9m/kgb5keYlJra6mu0wzX4rBSPtY2uvByKG8dU2jHXT/2OKzJ8yGFNoz+KAsWGbOI5LhOyN7EWI/FsYxbN8fO43ZGo+FBoCjljzss2hafqUZsp7OYmV9n17/tKte2tqzrhDhiH+sm2krKgeL6aBhlyqMQOM/7JeW2BtRFsIjk2BZQaOPiPP9g5FLUperBwKZQjsvcgXUeRGI9xKT48JNisY4U1wlxxDrWPLimPxfkMp7nKd4DQki9rQF1VQaLxpiBMeY3jTGfMMb8uDHmXxtjbs/53B8YY3537n9fHWeVsQyqArGQF2ffG5hLb1/Vg8ErV7est8V1Oou85U0m08LtdnmIsdl3qT+4+lRyXcZ1QhyhjzUPrnu6kB2xLOd56veAurrQ1oA6KscsGmMGkj6eZdlP7P/7lyT9uyzLPrPwufNZlp2vWiBjFmGrrIBCiLEeIQs02IxTWB8N9ZXr23rwc26FcnzHQBSNO3z77SMNh6uVRW/Kxiza7ruujH1MsVhHiuskcS2OIeSxZhzc3vX4yLEjpWPKb7/9aBLtONXzPJSu3AN82bS1WL3EXIsRUtmYxcqpM7Ism0iaBYqrkr5G0s/nfPRdxpgflrQu6aUsyz7nt7rAnqIJ2UNMCn1wA3vsucIbmEuQZvO57c2xztyxccv0EltbY+1OpoXbcuLkhvdDRN50FuPxrl5+Y7Pyxl02FYbNvpsZjYYHn5tt10efeF5PPnw2qQeForbWphTXCXGEOtYhro19UDXlUUrnU9/P867cA3x1qa0BvqznWTTGvFfSxyQ9J+l3cz7y01mW/aExZkXSLxtjRlmWHQoqT57c0IpH1t/q6kCnTm1UfxBLYaKV/Ivz2kCj249a/X3uDeyRszpyZE2DwUDbuxOtrw40mUy0FmisxerqQLu7E62u7WWAb4yGhduytjrQd332C8pevaa733JMF87dqzN3bHity2x5K+truekyTz5yVseO5qfMzP722NGhjh0dlu67xe/Y3C0IhCdTzueO4lqctrrXxr7YmUz12Pffd0smx4Vz92o63bv20I6bsQz3gKq2FgttGE2xDhazLPu3kv6tMebjki5I+qcLv//D/f9OjTG/I+l+5fRAXvZMg6G7HfPWR8P8SaFvVve+zdJGcm9gOxP93ZVNfeT/+pMgPY6L8tpx3rZc+NC9+vyfvKzs1WsH6/bQ43tvY33Pg6rtvmaRLuP6HSdObuQ/uA5WcvdDn9Ox+mLZr8Wpt9M618a+eeuxw9kRV9+4KYl23BSXe0CXlbW1WGjDCOn06eOFv6sMFo0xXy/pnVmW/eb+j/5a0vuMMWuS3ppl2cvGmK+T9C1Zlv2r/c98raS/rLfaQLGyNMkqpWkjawN97Innc1NmJHk/JJYFmHnbsrM70Xd85ku3fM42lazo9yHSZVy/Y1bA4dCDa9EYSYvUViynFFIou9BO61wbY2nr2PU9xbMLbO8BXUdbQ5/ZFLj5zyX9tKQ/2//RN0n6lKTbJf1klmXvNsbcJemz+585IukuSf9TlmXXFr+PAjcIzedBpGjQ/TtOHdU3FAxUr+pxLFrOfIA5nU5L3zbOtsWnSIVNj0dVsQGn8Z6WBQts1ouiHN3R9LU4pZ68rrXTtgPslI7doqp23Pa+65OU20GX8VyMkMoK3FQGi6ERLCIVeTew0WiY/zD4yFl9189+4SAt9ODnFcHbYlB14dy9euux6sDMJyCz/XzedktuvaY+N/+i7W2zmhzcNXktTqmSIu3UTUrHLk9ROyawiYcA3F/evuO5GCGVBYuV8ywCfZU3yX3RvFc7O5NbAkWpepL6vLmXHnp8b+6lqjmnXOZwLFpW0TxPi9styXletrx9V6XoIWEZJ6eGnZTmL6Oduknp2Nlijsq4OEfc9X2OSnQDwSKi8p3wvknzN7CiIG3N8SGxrIT99u5EH/6FP6p8GLENyGzK5Zdtd52HulA3/2WZnBr2fNt1TLRTOykeOxtdDHDRX7y8QCqsq6ECLrqcylM0UN1lkH5ZIRhJ+tOXL0vKn3NqMd2kKiCzLTqTl8aSyrxsKRbliIl0rGopzl+2bO3UV4rHrkoq10Jgpu9zVKI76FnsmRTe2Pq8DWtqvV2WY9PjWHbBzuuFuPChe/WLz/1HbY4nB5+bPYwcObrunW5S1uNRlsaSUmqdT2pr15BS5CbFnrxlaKchpHjsyoS6FqZwD0b3dbV3Hv1EgZueSKknz6VioMt613mzW7Ucl+92+ezicidT6f2feUZ/89rhffPia9f1YI1iEEWFa6qKTKReiKIv+rKfl7kaKtykfOyK5rz1PUdT3lZ0U9WzVKrPxfTCdxPVUHsu1kOozwnvUjHQdr1tb8JF61u2HKm8Cmioi97se47fflQvvn7j0LqcuXND7/v0s0FK8s+vs23gzoNOfF2bdqFIW9fiZX8A6fL226x709sXshpqH14Edbl99VVVu0rtuZjniG4rCxYZs9gDofPa65zwLmNVbNbbZhLsqvUtW86Lr13XgznfPfs734te0bjDtcHKoTFPW1tj7U6mwcbKzD7rMgaHCYXjYjxUfcu6f/rwAFZ27FLbPp9rYZfHlqWy/7kGHtalMdI2z2roLsYsdlzovPYQ1bdsxqrYrndVdbqq9S1dzmSqT/6bL+V+91cq9kHRfrUZk7Y45mnr5rjWWJmyaqeu3xvqZs14ilulNDYU3ZFKNcRY53OI7Yu1bi7DEro6tiyF9sU47nJdGSNNJeF+I1jsuNAPoSFOeJtiMDbrbXMTrlrf8uWsHFQlveW7S4LIshubTeA6b/7Y+BSDsLnJNl1kos0bf8oPZVL3Cn6gfW0/gMU+n+tsXypBRpdfBKXQvtoOVrsi5XbU5RcmsEMaag/MHkJtp3UoEjJVziaVp2q9q1JatySr9S1azmQy1VefGOV8d0EQuTvRxZtjfaQgzaIsFWmW5jPRitZHw0NvB13TTWxTPppMY2krDSWVNKoqXUopQvvaTl2OfT7X2b7UUt5C3YOb1Hb7krqdvos3dXGqHLihZ7EHfKZ1yBPjDWnZ39isd1lvjO36Fi5nayf3u2dB5KHvXBvoY088n/sWtuzGu7070Yd/4Y90z/mn9cCjzxS+OXVJN3F5I1z2vSHf+LXxlrprb6a7klKE9rXdYxX7fK6zfW33iC0KdQ9uUtvti96ofiFzpt8IFnsi1ENo0yd81XpX3YRt1zdvOa5B5M7ORNmr12753tmNTVLhjVfSQU/l4kNN3g3RZoyiz012/ntDp3C1deNP7aHRVsopRcsstQfUth7AmjqffbYv1SCjiy+C2nzAbztYRVhdfGECe6Sh9kzdC2xbqXJVPZBFKa2u67u4nKLvzvvOtYo0i7xUpAsfule/+Nx/1OZ4cvA3s4eaEyc3oleczRMjhauNNJQU0qjQD6mmMrd5PW7ifPbZvrrrFvu60KVrTtup8V1M30Uxqqr3F/MsolDXHrZjrO/8d9rMpbX40DmZTvX+zzyrv3nt1oea33r4rD7w2S8oe/WaV7BWtS5l+yLWfH9tzDXWl7kLu6Rv1+KuzJFXdk7HuPY1vV9ctsFn3RavzdPpVFffuBlyEzqtrft9qi9quqBv12K0q2yeRYJFLJ06N0XbG9tsGXkPNRc+dK/++KW/14/+xpcOPu8T4OSti1Q+P+RgsKIjx47onvNPH/q+F87fr63rW7UeGJq+8bf9oN+1Fyoh9O1a3OUXDrHPN9frXZNctj33OnzuXr31WFovBJbZMl5L6+rbtRjtKgsWSUPF0qh6uLC5WdmmWeQV1xnvTjRcG2hnZ6Lv+MyXbvn8Yuqkz7rYpJfGTi9rOg2lrTQq3ob3Q1upzCG+t4mKoFXnc5vngcu1Jq/q5kOP21XdJIhpBvsYSBcFbtAI26IDVZ/zLV5QVjWzqtiLTxGaefOFD1Ylra0O9I478gf1r62vOReema2LbbGXJooaNHnjb7qwRNcqsKas7YIyTRfZCFlYqsniTnn7IZXzIFZBsFTmcQSAthEsIirbG27V5+reuMserL7iGUS6mj3UFAVr4/Gu98OXywNRX6uWNRWgdrUCa0pCnlt1A86mKkKGDK5SqAjalfPA54VAKoEw3tT2iyVgmZGGimhs06SqPlc33ar0wWoy1Sf/zZdyJwV+8bXrerBGildR+lJR6uRwuJr78GWTKuWaXkrVMj9UYK0vVPpkqBTIplKZQ05A3vYk2D7nQZvnRm6l6nPFLwRCHauuXw9SWH9S/oH20bOIaGzfPFd9ru4b7PI3yysHcyDOlAWRNsu06TVZTJ3c2d6p3VPg00PS9oNA1zA3WH0heqRC9/zETmWO0RPYlTnyUkjnzMumOHPHRmGxnrrHKoVtriOV9aeHF0gDwSKisL3hVn1ubW0Q5CGr6MFqMpnqq0+MbvlsaRBZsUzXm9vsoSpEENLX9NLUtPmQ3nWhgqZYKZBlU1PU/d7QLxnaPt9tzoOUHvYXXwisFRzTuscqpW32kdL6dyXVGeg7gkVEYXvDrfrczs4kyENW4YPV1o5jEFm+zDo3txBBSNPFXpZR2w/pKfANnkIETU2O1wvZw2J7frusf8jz3XW/2ZwHKT7s27SxOtfiFLfZRSrrn8K4XAB7GLOIaPLGieTdcKs+Z/s9VYrG6eWOV9oPIl2WWffmFnLsVJfSIVMYF+NqWcd8ro+Gmmhvrk7f9ln3fG5qvF7oqSmqzu86Y7Pamhu17DxoenxvyO/zvRZ3fUxzSuvf9rhcAG8iWEQ0tjfcqs/VCaLybm6L/3YKIivG/5Xd3IYnjmpzd6oTJzcKv2uZgpA+FC5I+cEvtFDBU4iXIqFeIJUJWZBmpuj8bmLOxDyhlpt3HjQZ1Me4jvhci7se4KS2/k2c5wCqrUynzZ78Fy9e9VrgqVMbunTpRvUHkSTbN5JVn7P9nqoHCJc3pC6fPXj4mq+696F7dXJjqO/7P/+g0QfBlOXtp2XfJ6k7cXJDDzz6zKGHyCcfPqsrl/2uzXV6KmK+bBgM9npP7zn/9KHfvXD+fm1d3wr6oiDGvk1hubHP8zrfH+uZouvXttTWvw8vFWPhuRghnT59vDD9jZ5FNML2warqc7aBYtHbcknONx6Xh8K8XpPJdKr3f+bZoD0UvlJJg4rRa4N4YqWn1WmLMXvhm+xhaSv1r4npL2JPS5LidaSpqVhiSW39lynbBkgVBW56rquDwOusd9kA/a9YVHmru8/mi05cu3JTu5Op/ua1dgfpp1IKXaJwQRc1PWWISxuI9fKjqaq3Lvs25LnR1PQXsYpu+V5HZtsyGxIQ41rY9UJjKa5/Ci85gWVFsNhTKQUH86oeduqud+kDRMXciaH32WQytXogC/EAGHI6j9iYq7CbmgieUrpu2VT7DBW8Ve3bWPulyekvQp/XPteRW7blx56Kfi3s+rWsi+vPy0YgPMYs9lDTYw5sUpNsxh2EWu/CcTiPnNU/+OTT2hxPbvn8C+fv12Rnor+9uuW97LJ9ULRdb799pOFwtVaqj81+bWs8VJnUxsXAzvpoqNHRocY74dPTUrxuFX02xjiqou9sYtxf2ba4XD+aTnN33Td1r4WppPEXSX39YlrGsY08FyOksjGLBIs91FRwYHtxtr2hh1rvouWduXND7/v0s7nfL8lr2S77YDQaajyZajhY0Xi8q5ff2Kz1AGizX5su1OFiGW/ufXDq1IbeeONm8HaT2nWr7O9jBm+LD/xN7ZeiMYo21482z2XbZde5FqZ+rWpy/VIMSJf15SPPxQipLFgkDbVnmhoP5pKaZDPJb8j1Lkwf2587cTHlamtr7D32xXYfzMaAjFZXdOXyDQ2Hq7UnPrbZrymnfKY4LmaZtTlOMMXrVpHYk5YvjlG03S9191HV9Bfz5q8fbae5215HfK+FbW9flabWL6UU8UWxz0lg2REs9kxTwYHtxdn2YSf0euc9QBQFkVs3x17L9r1BhXgwdvmOpgp1+ErtLXVXhAygQj0E+q5TatetIk0XZ7IN1lyOn+s6Vl0/UnlQt2kjPtfCVLavSBPrl3LATME0ID6CxR6KHRy4XJxdHgJjrPfiA0TRW2jXZde5QYV4MHb5DptCHeiOkMFdqIfAEOuU0nWrSOig1maZZfvF5fj5HqOy60dTD+qhvueWbfnkA5XXwtQDkabWL+WAOeXsGaAvmGexh2LPk+Q6B9nsYWdxPMHiQ2Cd9XYdR5EXRLos23YfFK2X7T4p4/IdMeeqSnEMS1+VzSHqc36HmKcu1DrZnIN12lrdc3YmxLnrMsasbL+cOLlhdfzqHqOi60fs+ShjjMWbbYvNeK8m59v00cT6tTUPqIsQ5ySAYhS46blYF3LXAeWuN33b9bb5Xtd94LLson0gKXe95ttxiAehLhSWQDghi534FvxYvBbHqJYZq/qozzlb9D2+61OnGMf8fnE5fjErmsYqLuL7vbbrb/tMkXrxlCbWL8WK2ouW8X7EczFCKitwQ89iz8V64+faE+fas+UUrBW8Lfe9edjus6J9IKlwvRb/vm5vX8wewzKhe7hQLfQb/hC9Erbr5HouLgaKodqazzmbt4w6512d3tyJR69erGM0EyuTxXU/xQoWYmfq1NXE+nWh566teyGwDBizCG8+1Szn34pXqfpM2TgKm7E8IcZz5O0D1/EdIQL6ptOAUh7D0lcxxubUHSfYRLXM0G0txDk74zNG0WaMme21yeb4NXGMQlc2dh2LF7sAS+qVm2OvX5fGvbedEgv0EcEianNNWaoqsmDzmaqHiapAMnQJ8InlW/w+SL3oQ5dV7bvQRWBCPATGrJYZs63ZnrOuAVzVMssCt7X1Nadrk+3xC3mMqop3heD6YqRuG3NZr5TFXL/UA2YA8ZCGimCq0uBs0slsU86qUrDKHv4u3hzrI5HSJ6vWqw9SL/pgK4XCDDO2KXQxUs7qpm+VrVPd1Nkm2lrVMnbX14KmNxal9I3Hu17ptjbHL8QxanpMmG3qo28b25lMdeLkRlKppSldk4qkvn5Aqrpwfhfpx9NrQpaxV8W2p87m7a/LG+Kyt+WFb6XXBvrYE89bzQ/pK/V5DUPo8jamNrm0awpdrDf8RTcxm3OhaJ1CpM6GaGu+PbazAC5kemNRb+BwuGo9d22eqn1Z5xi1Mc+eba+pTxtbHw310us3kpk3MNY1qcnnkWV89gFspPbM4YNqqIEsYyUuyb4Sm03lPknO1RmL9nvRer3j1FF9Q8n3r1n2Itj0olZVQ+26Lrb5FCsbplppMFQbrtrnNm9bfduay9/lfXY0GkY9NrNtt7k+2l6bfFQdI9c2GroCtXXWiuV5ndI5F+Oa1OS1uYv3gb7o0/NEX6X4zFGEaqiR9bEypO3N3rZinW06mWvKWVEKVlHa1VpFuplNmqxtumBRaliXUxHmdbH6XIh5BUNKdQ6zsmuaq7IKpLZpgD5tzfW6vLiMwWBF4+FqkGNT9NnZz2xSYWPeY0KlErsGDraft+k1tU3PTu2cC31NavJ5pI/PPkBIqT1z+CINNYA+VYZ06S53LT5hk07mm3KWd3PPS7sq+/6q4+iTirU4BcBEK51ORcjTlcA3xcI8MSqchhC7Aqkkr7RGl/1Rt8JpiGPjcj2tc20KoW4qsev1MXRqq216dkrnXIxrUpPPI3169gFCS/GZwxfBYk19agyuN2/Xm67NGJQYJbonBT2O89+/s71TeRzr3Bhn+/b+n/m9JMbILKOUHhLnpTb+M2ZF39k+jv2QGarCaZ1j4zMW1ffa5Mq1oqnNfnA9pj5twGZbbcfBPvb97Z9zoa9JTT6P9OnZB4gh1WcOH6Sh1tSXypCSX3e562S9tpX7fNIbbdOHir6/7DhuSbVSl/qSitB1vpNLx0xNS23Sb9uKvr77xDYNsM4+D1XhtM6x8Tnnfa5NLtdI3/FlVfvBNbXT9fOhx8Vtb4515o6NJM65kBPeN/k80qdnHyCWkOd3mwgWA+hDY/Adx+H7MGXzEGj7oOj7ILH4/WXHsc6NMbUxMsvMtb02Vbyh6fGfVW2u7Fw4cmSt1pQDTU1VEWqKCp9j43LO5x0Ll2uTLZfxZXnrVLYfXK+PLp/3HRdX1cbXBiu6dOlG62OuQ78savJ5pA/PPkBMqb0M9kU11ED6UBGsboU4n6CnbqAUutJU2XGss6yUqu9hT+gKi11QtzqoJH3l2rYerLlPivbt228f6eU3NqOezzYVTkO8wKk652MVg/Fdn7rLcD1fbD/veu203YYUnylCvTikGupySLENo1jqHQNl1VAJFgNLvTGUCflw7Du1hCvbB4lQpdzrlPDvW+DRdykF+CGuK75tcH7ZIfdJk1NV2E5RMdmZ6MiRtSAPvmX7W5L1lEOLx933pZzN1Bx1r1GhA2Cb9V5MV7XdhrxnilD371SeA5pcj1S2eZl05bkY3cDUGQ3q8sUyRHe5zcNCqHLbNqlevnOTFR1H3ykxDvbtD71H4x3ewKYuldThkG/tfcfNzqdLhtwnvlNV+MzJZztFxYuXburBC78f5IVO2fX0xMmN0mNRdtx92p1N2qdr+3BNVS3aR2Wfd01v9W3joc6zkOdriGtMk88jXX72AVCOaqgtSrFamG358Ty21f9CVUKsqjQ1e1MeqjT74rJnbMvjb2+OtSp57Vs0K4UqZiGnFghRuTDWPrGdqmJtfa30PLM9D4sqe06m04P0WilMhda862nVsThydD3KdausoqlL+7DZz65twWb8bFXlUt82Huo8C/k9ttOtAEATCBZb0IWbwXxvgi2bIDB0ue225yarO/eirxRfNPRN21NahGy/toGezVQSdaccKFtG0T6fFaUpOs9czsO8KSrednxdu5NpkOtS3mfnz/mqY3HkyFqw4z6/LmXTEsWaSzEE2+mUfF9mhDrPQnxP7P3bxn2DexXQfaShNixUCmZsruk0tilqocttF6V67WzvWKW01dX0lBgUE2hO7CpmZW0wRhpsWeVC23ZVZ8oBm2UU7fPhcDX3QXx2ntmch/P7LC/98cTJYa3rksu5WXQstrbGGqzVv24VrUtZ2qdNZcvY17uyqts26a2u1TlDnWehvifW/m3jvsG9CugPgsWGdWG+PZ+A1iUIDF1uu+7cZLHniwvF9bhQcKC+GFNa2DxExZjDrCgQk+TUrnymHHBpu65jGdfWBpVpnUVFayaBrks227cYrOYdi62b4yBBa9W6+EyDFPN6ZxtYVH2/6wueUOdZiO+JtX/beEHdlZfiAOxUpqEaYwbGmN80xnzCGPPjxph/bYy5Pedz32uM+awx5rwx5ieMMeQeLAidghmLbzqNbdqebVqRq6K5yYrWp246cNPj2myPSxfSnLsm1LF0STNzGatlK28Mne/57rJPqpZRlrpZdZ7t7BeMKvr9K1e3vNNTba9LZdtXdD4WjQ+vm/7scjwX93vZmPVY1zuf1MuyNu867j5Uunnd74m1f5sYjpHCMgHEY9uz+AdZlv2EJBljfknSP5H0mdkvjTFvk/STkkyWZbv7n/luSb8aeH07LUZvQWh13m66vNWN0WPjsj6h3nw2NSmx7XFJ/Y1uG72dKfWwumQWVJ1PddK8XCuc1tmHVcs4cXKjchuqzrOi3xcVrXFJT627fRdvjvURh16+OunPLteJsrZTdKxjXO9czgmXNm/bXkOlm4f4ntD7t42qzqlUkgYQTmWwmGXZRNIsUFyV9DWSfn7hYw9Iej7Lst39fz8j6QMiWDykqeDCV92A1vVhK/ZNo2h9QqUDxx7XNmN7XFJNc2bMjN9DVFH7dXkpUPZwVtWudnOmnnFVtoy11YG+67NfUPbqtcrU1LLzLO/3W1vj0qI1tumpdbZvuDbQx554vvJ8XDxGvi/TbK4TdV4ohb7euZwTMV+EhXp5Wfd7Qu/fNl5Qd+GlOAA31mMWjTHvlfQxSc9J+t2FX3+VpCtz/76y/7NDTp7c0IpHtuXq6kCnTm1Uf7AD7j6ytjff3u5Uw9UVTSdTrQ1WdOxoGikaO5OpHvv++/Tg594MaC+cu1fT6bQ3x2Bzt+AhcuK/jatre1ndx44OC49lnXZsc1xibFddO/sP7A/OPeRdOHevztyxobVIqddtLNPGRCuFQcXo9qNO35P7UuCRswdtb2cy1WAw0PbuROurA00mk9xtL2pXu5OpvnJt+9A+fOfRoXNbKlrG5//kZWWvXivchjxV59ns9xujYen+funSTf3g3JyKddpH0fbt7E4Ptm9m/ny0PUYh1mV2nbBpOzZsrnc2bM+JUOt9sP6JP1OE2r9t3M+X4RkiBam3YfTHynTq9pbHGPNxSXdnWfZP5372P0r69izLvnf/3x+T9M1Zlv2Txb+/ePGq12ulU6c2dOnSjeoPdkjeG/9UekNSWY9Y63Ti5IYeePSZQw8oTz58Vlcux2tnddtx1T5oa7vKtLFOdZYZM03qoHdkIbPApXdkMFjRkWNHdM/5pw/97oXz92vr+tbBHKO2y8lrV6PRMH8fPnJWVzza8OIydnYnetcnn9birn7h/P2a7EwKe/1cl5m3H87cuaH3ffrZ0vbh2g6c9uHDZ7W5Oa7dFlzWZdbjVdV2bNJUQ7I5J2zXO++7i7ajj88URcjs6KdlasOI7/Tp44VvKit7Fo0xXy/pnVmW/eb+j/5a0vuMMWuS3ppl2cuSnpL0cWPM6n4q6nsk/VL9Ve+3vEAxlfFmodJyQj14h943qacDF6k6LqltV5fGzDTxcGOTZla1T2zSvFzTkd0qkNqPYSwbD3ji5IbecUf+Nrx46aYenOv18z3XQ6enVi0r79wsOh9tj5HNvrZNY42dpurD5pzwSW1M5X6awhi9JmoEpLBMAHFUVkOVtCXpfzDG/Kgx5kclfZ+kR7QXEH5ekrIs+1tJH5f0s8aYn5D0V5J+Lcoa91iKFcTqVLgLWZEz9L6JVZG1KTZFhlLYrqYrxvous8nJxouqNbqcM2WVF22rLvtXIK1++C3bltnfFm1DUVEal3N9cTL6+f29dXMcpHpqkcWgLe983NneqTxGNu2h6jN5x6mqamfM+1BRFVObCqau1Ubbvp+mWJW6jaC17UAZQH02BW7+g6T/tuDX75773C9L+uVA67V0ulBBzGVC6JBvdGPtm76++XTdrthtq43eTtdltlEYaH6fu54zVb0xrkVrXCqQTi0CRZtt8en1q2qrZb3DNnMq2lRPtWFTYbXqGFXtQ9/rbFnbiTnXX915FF2Kv7Q9TVUqvZoAEIJ1gRvElXIFMdf0vNAP3i77xudhpu0gPBab3p8mxpQ0VTHWd5lNvKip+g6fc6bspUBRMDQe71o/xBbtw6oCLDbbMtsfedvgOym9ywN6jEB1tg42wapU/kLDZh+6tJmQaaquQgZOti/C2r6fxnj5lMILYwDLySYNFQ0JNTlwSK7pebHe6FbtmxRTflLWZNql5D5RdpPLjJkqa9Mu66SNztZ/UVH643C4WpmaV5bCWXXcqrblyNH13P2R1+vneh20STusk55qEyjanFOzdaiTomrbZmKkqbrySQetuk/YnJNt3U9D3wO5twFoGz2LCWmjB6aK6xvSWG90y/YNKT/u6rz5rvOGO9UxMzFSZW3bpc9chz4FV8qL1oSZd7BqW6oK1wwGK17XQZsg1Wbb6rSDqnOqqNfRNUV19rlYhWpC3odce+1DZju0dT8NeQ9s4t5GjyWAKvQsJiZWD4xPj57vG9JYb3SL9k3bhQy6xve49vkNd4zCQC7tsuicmaWNhii4UtWDWrewS9W2lBWuWWxbkpyug3W3raq3L0SPatlxLEpRLbuGhipUk3e+h7oPufTax8h2aCOjQQp3D4x5b+vz9RxAWPQsJirUm746b2p935DGfqM7v2+6UBgoNV0uQ1+lzvEOWfDItV0WnTNFaaO+Y5/qFHax3bc+4wEv3hzrIzXbls+2SbLu7StTdU5pdaCH9ntUF9chb//aXEPrFqpZs+ixDnHttO2t9cl2sG2TTd8DQtwDY97bunI9B5AGgsUeC3FD8E3LaqrSaNuFDLqqC9VCXYRMXwvxYGnbLsuqZlaljfo8LPoEcospnDsFy6yqAFpYuGZtoI898bx325ovlhMySHXdt0Xn1NbWWIM19/Rfm2uob6Eam2qrodgETm2mq8ZS9x4Y896W+vUcQFpIQ+2xECksddPzmnijm2JhoHmxy7T7cDmubZehr1InfS3mupe1S5t5CEMW3gk57+BLf3/jlvW12Zay/bGzM1H26rVblm3TtvKW67RtBUGqa4pfVQpr3XkdfXvOytqfy73B9RzxSWttO101pjr3wBj3ttSv5wDSszKdNtvzcvHiVa8Fnjq1oUuXblR/EJL2bghHjh3RPeefPvS7F87fr63rW843sZRTOlN805y3TseODpNrxzbH9cTJDT3w6DOH3nA/+fBZXbnc7vb4rFtT7SVvOZL2enUWeqDygvWD7ACLz9bZzqLlnLlzQ+/79LOF+9Zm/ebbV976jEZDr+NXdx++49RRfUON62PZvs3rCfPZv3XlrePO9o7VvcH1HKl7Ttke01DXoq48U8S4VqV8PYe9rrRhdMPp08cL3xSRhtpTMVJYUg0UpebSXm0VpQDffSS9U66taqEh+IzraXK8Tn465oZ1ClidsU9NzDtYls4m2Y8HdG1brnM55s4XWeP6WLVvF9uc7/6t+0LPt9qq6zkS4pyKka4aQ9MvTWPc21K9ngNIE2moPZZ6emYMqQS0RWleKx1N8fFJR24incknVTNmhcGquRB95lT0regYe95BSaXbcvHm2Kr6p2vb8pnLMW8f1rk+uu5bn/3rci1znU8xVCVVl/1hI2S6amhtVw8NuW0xqj8D6K/0ujkQTFvzTC278ofZadLpvGVs33A3nRLs8pY8Vs+E7TbXmVPRdYxirHkHL5zb27el22JRtKaqIE6ROnM52lSgDTFNRtG+tdm/ri/0fHr16lZSrVuRuuo8K/tdG71ifawemlo2DoB0ESz2XKwbQuyAJ8b3NxWklT/Mdv/GXLYP23iocnnoj5Ge7brNRQ+7szkVQ+y7OgGVVJ7COZ1OdfWNm6XbUla0xjZQLRNqCpAY02RU7dt5IV7o2Va2tE1RdT1HXD4f4kVSGy9B+1w9tIsvLgE0iwI3cOJ7s7cN1GL0SrVR/KaoWMPdd24cPGj3UduFE2zaWd2iMYtCFdjxKfTis51lRVVm61FWtGXxWuy6LS++dv0goKsqiFO1ffPLnU1RUVS4ZbIzKe3xi71v54+hT29bHptCZjbzKdpuX+WYxZLPhz7vZtvvG+zYPlOELBbX1awSpInnYoREgRsc4nPT8uk1cgnUYvRKtZU+VFhUo6NjFm3USfEMUcxDsntLbltIw+a7XLa5qTkVm5530KVoTVnPn5RfEKds+w4XDyqYy9Gxx6+Kb8GassDNJ4iINZ+ia++dzedj9M61nylil42QYrVuALBFgZslU2eQvmsRA9f5sGIUHmmjmMmMb2GSrvIpPuHTHn3+pqrYyHzvh8t322xzU3Mqtj3v4GQ/6JXyC2i87fh6ZaBadq0o2o82czkWBam+cyrOttGlYM3afuAWen7AUPMpLnK9fpV93nduv1Tm/KtTDKlr80ICwCKCxSVSd/Jy15u968TPoScKjjX5sE9AsSxcHqp82qPr37hUifQ9P8q22eU7Yz6Q2gRUZeMMq86V2IGqzX4cDFa8glSb64BNwF+2b+sGbmWKKlvubO84Xf+qKvnayvu868uQtiuPLqpTPTTmC0sAaAJpqEukThqQayqOa0pirHkh2y5msmxc0td82qPL37geK9/zo2ybm5pTsc15B13mHPQpiOMzn+OVyzes01Przqk4r2jf7mzvBEszlsLNpzjbviZSJG2rmKZ6jfUphhSr+jJjHwE0iZ7FJRGil82l58MnrS7GvJChv5O3xNVs0td82qPr3zTZs523zU3Nqdj2vIMucw4W9dCslVwrJL/5HG17/KrUna9yFjCHmB8w9HyKTaZI2vbOpX6NdQnSQs8LmVqPK4DlQM/ikgjRy+ba8+E6H1aMkughvzPWW+K+qppHzbU9uvxNWz3bLt8Zak7FNucd9J1z0LYgTt35HKt6U6uuA6Hmq5Tqzw8Yej5FqfkpIap653yvsSlfe0PNC5lqjyuA/qNncYmE6GVz6fnwGecRoyhMqO8M/ZZ42fm0R9u/Sb1nezanYqgenTqFXer0bFbt51eublmNo6y6VviMs7TpTa1SZ/sW1Rn3Jvn3uBVtd6wx3TaKrpVdH9uYp+5xn0m9xxVAfzHP4pJpq4R3ym9+XdSdK6xr7Tj2cfNpj7Z/Y3us5rcxxPmxuM+amFMxbzmx5h3Mm2fRd85Bl/kGXfbjbz18Vh/47BeUvXqtdg9M3TkV8/iMUbSd68/lu2PNi1rnumF73oa8Fjd1f6ozRjHUXI/oj649TyBtzLOIAz6D9EPoy40sRqpsipp6qeDTHm3/pupYFW2j7/lR9n3z3xljTsWi/dLEvIN5+7lqzkGX9M355dikr1740L36/J+8fNDjWFbsJ9b2zb6/aDk+FUar0qR9ztlQKZIzIa4bttfYECm0Tb889b0PxijWtiz68qIaaBPB4pLi4umvrYC7KT5jY+rekKv+Nu/7i/6mKoCS3Cp42rDZZ3lzKsaY5Num+mhReqrLg3ZbgWrVOMud3Ym+4zNfuuVvbILUMq7bVzYe1VdZYOc7ni3ky6+QY+pijW2ctzOZdmoMYOjAvu/ayqIC+og01ETw9stNV/dXzHZsu0+qPueSmuZzQ3Y5di7f7/LZ0Ol3rt9XJ4XO5m/LUmur0lOrUtqO335UKysr3mnAddI3Z8rSV8uOxYuvXT8IkmOlp7799pFefmPTOz2yapl5bTxEe657TY2V0hprebefOqb7f+b3GlvfEAiA7NRNUe6Kvj8Xo1mkoSaMi78b9tdhLmP4qj7n8sbetSfB9di5fL/LZ0NXtXX5vrrVOaXy9LvZ76tSa33nHVwfDfXS6zdKA66ybbRN3yxbflWPasze1HlFx3A4XM0tRBKiwmjesQzVnusEim1Uiq7T0zYYrGi7g5Wt+57VEkrTVX6BviNYbBGlsN0s6/4qe3Cx3Se2n3NJkXS5IfscO5fvd/ls6PE/dcaTuT74VT2UX7w51kcsUmt9H7RHo6EebDFQtWlDMYLUInXGo9YJSHynfYkVBLUxpq7OC5fJZKr1Do8BTDGQTQVTXAHhMXVGiyiF7WbZ9pdNWXjbfeKy72ymkHAtu+967Fy+32cKgNDTZJR9X9XE58Em+S6YczBvmgyfcv42garNNvrue5s2NL9989NFbN0cB5v2Jq895Y1HLVqOy3QPttNXVO3TJqaY8DmudafnKJoWxMZkMgk+VQ7axxRXQHj0LLaEt19ulm1/2fSi2O4T131n88betTfD9di5fL9Pr0aowh42aaUnTm5U9nq6tN+iXsGyOQeLCru49GyW7ueCQDVvG332fVUbKtq+EL2pM7Zp1CEK0bimbJft06YyMlyOa+jhBD7X/rXBSu1rQN/uO31BMSAgLILFlqRUCrsLN7yU9lcTbNIqbfeJbzBVFUjY3pB9j53LDd/2szaVUm3YppX6Bjllih7K10r2cVn1UZdzf3NzrMe+/z49+LmEAlXL6qp1XhC4BFx1XxzUqWyat0+bHL9lc1xTGk7gew1g7HzaQlb5BUA11Fa1XbGraze8tvdXCDbt2GUC5qYmsC7iUlzHZ/lV319W9dNmTkXfbXbZlpjVOXMLDnlUH616YTT/+7xqqKPRMHoF0pDVVV1fkPlW3lwco2hzXoesKpriZO5NV00t4vtM0Yf70DLpwstwX318LkZ7qIaaqDbffqX0dtfWsrwtdOmJs90nsfad7Zt53+UXfb9LwZjQbd21pyZmdc6q9OGqwi5VPZt5+3ltsKJLl24c2s+xK5DWKVyz+MAYq9rtItdCNKHT7VPLyOjDcAIqbXZL6u0J6AKCxZa1VQq7qze8ZSkd7pKC6RKwxdp3NjfkOsufOAR/i+sSsq2HmiYjZHXOqtTasuqjZembRfv57iN7t426gapPYOC6fbvra7V7lG0DLpvtqTqvYwR3KY3fSi14ddWHYBcAXFENNRFN3mB8qkeWfVcb+n5D9qlWabtP2t53dZfvUlk1ZFuX/KtdxqjOWVbh0qb6aFGv32w/Fu3nlYJ9FqsCaVn10bLtG493S6vQuqiqdmtbadTmvLatKmrbdn2uJTGFrkTcJCptAlhG9Cx2QOi3lSHe7nZtvGMX9aEXNXTbdX2z33RPjUuvZ50en7oFV6p6/dbWBiX7+c20zqJKtqG20eYaU5TiPByu5ga7oVJgZ9vgmuJcdV5XpWz7XHtTupZ0fThBSj21ANAECtwkLGZAVmeQPgP8/fWlHVcFgTHbrmuBjBjttWj7bNbNtihP3X2Qd4zmf1b1HYW/f+SsNm+OrdfbZxt9j9l8EBursIvLPqxr8RimdO0N8SLI9jtipHfWuRbzshQp6MvzBNJAgZsOil2Aps7b3a6Od0R9Ng9JsduuzZv9xXF8IedUnH3nYk+N7zQZrj0+dabjcOn1K/r97u7U6fj6bKPvNWa2fTHHxk0sj3eIACfm+FtfIQOlqv2TalCWUk8tAMRGsJioJh4KYjyoMsC/v2yDwNhttyz4K6uS6vtwV/bA6pLyGmquwxBzDlbtx7Lfr60OKtM7q9JTy9heY6quNbHTBV0C0lC9cG1fe5usot2Fit3c6wAsAwrcJCh0UY4qvg+q8xjg76atwkB12BSWidl25/92sZjKfPXOsoImru3T5jvn+RaTcVVnOVX7cd7i73e2d7Rd0atpW+ylSNU1Zm19zWoZIQq7VLXXqmItLsVvqqRw7XUpLtWlZbWhi/cAAMuJnsUEpV5enAH+/nb2J91OLa2qim2vRoy2a9uzF6NH0/U7m5xCIuScilXLn0/vXA/Qq1ml6Bozq3AaMwVWql9cp2zqkTo9Y21ee5vs2UyhFzWWVFNrAaAIwWKiUg7Iul7Nri3ro6Feev3GQc9PimlVRVyCwJBt1/aBO8bDpW86pOtcgC5BTKw5FV1MJpPc41vUqxmy+qhvhVOXY+8a5BUFpDFeXthee2MEU02+xEz9hamvLqTWAsAi0lATldrcWIuq0tdw2Gg0DJqKaCNkqpPt/Ggh265tKppLip7tPqmbDmkzF6BtAB17TkUXa4OVQ8f3bcfXK3tPJbf2mJcC20R6vm/642Ibi7WuZdfekGmveZqcI7HL8zEW6XtqLYB+omcxYV2ouNbVVKCmNZ1WFSPVyaVHOUTbdd1nVT2aPvskZDqkb2987DkVfdNgXXo1d9fXvNuja4XTOudSqPO0iZ6xxfVooteqyaySvmWw9Dm1FkC/0bPYAdxAuq/J4hSuRVlcuPYo19ku131W1qPpu0+KvrMoHXKxh8ClmEwR28JCRcvZujkO0vbyesNsejVngbXNvq9TUCZEr1rI87TpnrGmeq2azCrpUwZLCgWKAMDHynTa7AXq4sWrXgtk8lF03fpoqK9c39aDn4s7ZjH2ROFFYrwZt5mEvGryecl9n5R9p82E72s1etIWl1m2rMnOpLBozbw6k7kv9shOp1NdfeOm1Wc3N8cajYaV+96l1zfvs5KCTVYfcuL7poqZ2LRJgpFbtfFMEbJtATwXI6TTp48Xvq0lDRUIqCxg2t4c68wdG8HSqooCmqZTnWI+EPvMqSi5jx+bfd7mO6tSDHfX14KlAzY1p2KRvNTGC+fu1VuP2RV7GQxWNB6ulu77Ncf9lZ8CuxGsmIzNvrI9j5oaSpBqQRhSK2/Vt9RaAMuBNFQgANsUuLXBSu20qqpiJ02mOsVMeZ3xnVNxxnafuHxnWYph6HTAUEVrfFL68rblocfti73Y7Pu6BWViFJMp2le+qa5NBEwpFYSJXWiny/qUWgtgOdCzCNTkU1jC9+HRZlmhp10p6x2IMT1A0TLrzKlos09cvrOoh2Bne6eyJ62JORWLpvVY3I9lQvVSl+37EMuI2as2v+zUpz1Ipdcq9f2UCnpcAXQFwSJQU6yAyXdZoR4aq9JLY6S82qS0+iy3ap/4zKlYlGIYI3BpsvroTKjqo1X7PsT+amJe2ibPc18pVNDuwn4CANgjDRWoIeZ8anWWVTfVySYlM3TKq20aqO9yy/ZJnTkVF5cXKh2wieqjVfKWceGce/XRsn1vs7+qzqPY89I2eZ6H0FavVdf2EwCgGj2LQA1NFpbwWZbvQ6Nt70DIHh2XHgnb5bqkYaYwp6JkXzCoaDlF03qEKvYynU61tbXjlWqYt+99ixjlfU+sXrU2C8h0qUhMqoV2AAD+mDoDqMmlHHrddtxE6XXXMvwhqqH6lP4vW67vOvlO/VC2XS4P+r7H12VajzoT1k8mU506taHdqaJMz7JYmTalaQaaXp+YVYZjBqCpHbciPFOg62jDCImpM7B0mnwb32RhiSaW5do7EKJHx6dHomi5dQps+Ez9UNbOXNug73gv22k9QhV7iTU9S50iRrE1eZ7HKhLTxLyPqRTaAQCEURksGmPOSPopSX8q6ZSkE5L+WZZlOwuf+wNJm3M/+r4sy/4u3KoC1ZqaBHtRk4UlmliWT3ppUYAQompmmcXvDhFk+ARfdV9Q+BTZydNEsZfYqYZtzBdq851NnecxAuUmq5SmUGgHABCGTc/iHZJ+JcuyX5EkY8yTkv6xpF9c+NyTWZadD7t6gL0USrY3ObYo5rJC9A64Bu4hJkOPEWRUBV+hXlBUBaa2FU6b6NmxCUjr9jA2NfbN5/jFPPdiBcpt9NR2ZawlAKBYZbCYZdnzkp6f+9FA0tWcj77LGPPDktYlvZRl2efCrCJgJ+bDUJeKTIRSp3fAN3AvSy21eaCPEWRUFWAJ+YIiZJGduj07ZW0+VFGaMk30kIY+fiGuEzHacBs9tQCAfnAqcGOM+YeSfljSB7Msmyz87r/MsuwPjTErkn5Z0tNZlv384neMx7vTFY/q2aurA+3uTqo/iKW1uTvVPT/21KGfv/DJBzRa9SvZvjOZajAYaHt3ovXVgSaTidZqlH9flnY80Yru/5nfO1wA5ZGzWnX8rp39yecfnAsaLpy7V2fu2Mg9FjuTqV76+xt68HN2n68j5HbO7EymWhmsaLw71XB1RdP9Nhh6OWXLL2vzZW3Y9VjZrMvivgh5DEMdv9DXiRhtOEZb7bJluRajv2jDCGltrfhB2brAjTHmmyX9oPbGIh5qnVmW/eH+f6fGmN+RdL+kQ8HiZc9KeVR9elNb4/JSd+LkRv7b+MGKV9uJUdWvD+24qv3NKnLm9mLsTHRtoSJnVa/GiZMbB8HH7Hseenyvx7hoX7712OFer6tv3Kyz2Ye4bqfP98+K7MRczjybNl/Whn2OlY0YY99C7ddY1T9Dt+H10TC/p/bmct4/+nAtxnKjDSOk06ePF/5uYPMFxpj3Sjon6Qck7RhjvscYs2aMefv+77/OGPMDc3/ytZL+0n+VUcR24vJlFGoy9JnRaJg7X91oife1TfurmuB+fmqEqkndfSf5LpsE3lbVBOK22+k7EXlekZ2y5YRQp83HnJA9RopkqP0a6zoRog0vft8sdfiF8/fryYfPJjedBQAgPZXBojHmPkm/LukbJf22pN+R9C2S3iPp8/sfuyLp/caYTxhjPqW9qqn/IsYKLzsCmGIhH4ZiPvjGFHu9bNtfVeBu+9Kj7gO9T5BhE8TabKfL91SxeRFS99jXbfNNBrXz61xH3RdMTVwnQu630AFoH6V6bQeAttgUuPmipGMFv373/mdekfTBgOuFHBQpqBaqZHuT1RhDaCI12aX9VVXkdClGFLvy5jzXgidF2ykpaOGUJgrKhGjzTRSlkcK197qVY7t2nZhZ9vtEHoZ3AEA+6zGLaF9XH0zaEOJhqKkH37qaqujo2v6KAnfXlx5NBEozPhV187bzxMmN4JV585bTVDVW2zbfxLQdobe57gumrlwnUCyFaZcAIFVO1VBDuHjxqtcCGci7J1YxBeQLHYzEaMcnTm7ogUefOVzl8OGzuuJQUMpmW0O1P991ng8kQ58Ls4In95x/+tDvXjh/v7YsC56E+h4boY79vKp2YNuGY2U6xNjmuvrcK9XXjJX5dpximwKq8FyMkE6fPl6/GirS0MSbe7wpVFprLKFSk23frIdqf769MfPbEnpezVA9901lAMRKSw+Zyh1aqqn4qV8nfPQ5AJ6XapsCgFRYVUNFWihS0LxUHxbaqOgYov3VLUYUq7BIqEIyoSvz5oldUGbx72cFezZ3p7UL9tRZp6aL6Lhoe/mhLFPV7dTbFAC0jTRUoEEx2nHddMwm0yaLlp/3/TbzL8ZIHSvqUXHtaWmiZ6aptPSU0t9TWpe+Woa0zPlrMW0KXcRzMUIqS0MlWAQaMAt8fNtxVeBUNzAJ9XAYImXLdltiP+CFGh8ZO42tiaA0teBhWVIk29D2y6OmLF6LaVPoGp6LERJjFoGWLD6A7Dg+ZNk+wLRd0THUg5ZLVcLY43dDjY+M/WAde7xcimO6+jhGMBXLWnWbNgUA+QgWgUjyAp8L5+7VW4/Z9Xz5lHP3fWivE3iFLDvvGpQ18YDXZrDk8t1Fn6u7fikHD33o4UrRMk8HQpsCgFtR4AaIJK9ozEOP5xeNsf37oqIzIfgWrgm1nnWK1sR8wGujAMasmMyRY0e8i8mE+I6ZJgr2IB11C1ABAPqDnkUkratly+v2RrXZm+XyvSHXM+UerCZ7WkL01MaYuP6g53ky1XCwwpiuniMtEwAg0bOIRIXsFWlD3d6orpRzD72eqfZgNdnTEqKnNkav9KznebS6slRT9vhOw9KEJtYtlWsNAKAdBItITl/m+MoLfC6csw98Ug2cFoVczzbS32wfuJuY3zTE/JGx5qBsSirrl/ILq5TXzVUqxxsAkI80VCSnTuXJlOQVjZlOp7r6xk3vv/dJ/QuRslr2HaGrkjaV/uZbwbUqfbjtYjIpp/OWSWnqgtBpvMuybi5SOt4AgGL0LCIpXe8VWbTYG7XmuP51erOaLJISo9ct9tyEIXuvUysm05Ve6ZnUsgmaLi7lIuV1s5Xa8e66rt0XAXQLPYtISld7RarUDXxc/76tIilNjG8K0VMasvc6ajEZz14X1+9ou5BUStkEKc4r2YV1c5HS8e4yemcBNIGeRdwihTeUXesVSVGqRVLqCNV7F7r3OmYxmTo9tTbfkcLYt9SyCVIuLpXyutlK7Xh3Fb2zAJpCsAhJaTw0zjDHVz19LJIS8sEo5AN37P0U4uG/6DtSedhMMQBK+YVVyutmI8Xj3UWpvcwD0F8Ei0jmoXFeE5Un+yrEw1hqD3ShH4xCPXCntp9cpPSwmVoAlPILq5TXzVZqx7trUnuZB6DfGLOIpMePpPywLbU/1qtIiEnkm5yIvkyMcVohK7imsp9cpDb2LXRF3VDrlOqk9Cmvm40Uj3eX9HVsP4A0ESwuudQeGrsi9cICbRRJiSXWg1GoB+5U9pOL0Ps0xHUiZAAU8rqV8vUv5XWr0vWAt21dfEkFoJsIFpccbyjddWWesxAPY6k80MV8MArxwJ3KfnIRYp/GeGlS53ik/hIHh3U54G1TF19SAeimlem02Qv1xYtXvRZ46tSGLl26Uf1BODsIfhYeGlMLfkKp2+tw4uSGHnj0mUPB9ZMPn9WVy+VtNLV23KWeYwKBN4U6bj77dNaGU7tupLY+SFtq1+I6unQdRzh9asNo3+nTxwsHO9OziKV5Qxki2OhL2m4XA68u9t6FFvq41dmnqY11Tm19gKZ04Z4DoLsIFiGp/w/ioVJH+5C225U02iLL+mAU87i57tPUXpqktj4AAPQFU2fgFn19oAo5TUBKZd99SqSnNGXCTKql3lNar5SOW2pThqS2PgAA9AU9i+i90L0OKaTt+qYjptYDk2o6bGrrldpxk9Krxpja+gAA0AcEi+i9GKmjbabt1klHTCmNNtV02BTXK6XjNpPCS5OU1wcAgD4gDRVLIVbqaBvpbXXTEVNJo00prbIL65XKcZu3vTnWlcs3tHV9S1cu32g9MEttfQAA6Dp6FrEU+tLrECIdMeS+8E1/TDGtMuX1ktJuw6mNCUxtfQAA6CqCRSyNPlR8DZWOWHdf1B3Tl2JaZcrrNdOHNgwAALqDNFQsna73OoRMR/TZF7MxfQ88+ozuOf+0vv3Tz+qVa9tad0zTjJFWGaJ6aYrpnou63oYBAEA3rEynzT50XLx41WuBp05t6NKlG9UfBBIWqh23Wa3zxMkNPfDoM4d63p58+KyuXHbbtlDbEXp/pFYNNSVci9EHtGN0HW0YIZ0+fbzwbTtpqEAHtZWOGGMakrrbEaN6KemeAAAApKECndaXyc/rbEfM6qWkewIAgGVGsAjASUpj+mx6OgEAAOCHNFQATlKawiH16qUAAABdRs8iAGcpTX6eUk8nAABAn9CzCMBbCmP6UurpBAAA6BOCRQCdR/VSAACA8EhDBdAbKfR0AgAA9AXBIgAAAADgEIJFAAAAAMAhBIvAAubmAwAAAChwAxxYHw01Gg2pqAkAAACIYBGQtBcovnJtWw899pxeev3GwVx9d922TsAIAACApUQaKiBpNBrqoce/qJdevyFJevG16/roE89rNBq2vGYAAABAOwgWsfQGgxWNdycHgeLMi69d13h3whhGAAAALCWCRSy9yWSq4epAZ+7cuOXnd7/lmIarA+buAwAAwFIiWAQkbW6OdeHcfQcB42zM4ibjFQEAALCkKHADSNreHOuu29b11CPvoRoqAAAAIIJF4MD2fnA4GKzoJqmnAAAAWHKkoQILGKMIAAAAECwCSaICKwAAANpWmYZqjDkj6ack/amkU5JOSPpnWZbtLHzueyW9R9Jr+9/7o1mW0UUDOFgfDTUaDRk3CQAAgNbZjFm8Q9KvZFn2K5JkjHlS0j+W9IuzDxhj3ibpJyWZLMt2jTG/JOm7Jf1q+FUG+ml9NNQr17b10GPP6aXXbxxUZL3rtnUCRgAAADSuMg01y7LnZ4Hi3N9cXfjYA5Kez7Jsd//fz0j6QJhVBJbDaDTUQ49/US+9fkOS9OJr1/XRJ57XaDRsec0AAACwjJyqoRpj/qGkG5J+Y+FXXyXpyty/r+z/7JCTJze04jEca3V1oFOnNqo/CCSsrB1v7k4PAsWZF1+7rvFkSttHMrgWow9ox+g62jCaYh0sGmO+WdIPSvq+LMsmC79+VdI3zf37xP7PDrl8+UbejyudOrWhS5f8/hZIRVk7PnFyQ2fu3LglYLz7Lcc0HKzQ9pEMrsXoA9oxuo42jJBOnz5e+DuraqjGmPdKOifpByTtGGO+xxizZox5+/5HnpJ0rzFmdf/f79Hh3kcguD5VDd3cHOvCuft05s69N4WzMYubjFcEAABAC2yqod4n6dcl/bGk39ZegPm8pNe1V9Tm3VmW/a0x5uOSftYY8/eS/krSr8VaaaCPVUO3N8e667Z1PfXIe3q1XQAAAOimlem02dktLl686rVAutsxc1A1dL8YTJeqhtq248FgRZMJM88gPVyL0Qe0Y3QdbRghnT59vDBVzyoNFUjJMlQNJVAEAABA2wgW0SmDwYrGu5P8qqG7k16NYQQAAADaRLCITplMphquDg6KwMzc/ZZjGq4O6JEDAAAAAiFYnEOvVDdQNRQAAACIz3qexT7rY2XNPqNqKHxROAgAAMDe0geLB5U1H3uuc5U1l9n2fnA4GKzoJg//qMALIQAAAHdLn4a6DJU1+4xeIlSZvRB64NFndM/5p/Xtn35Wr1zb1jrnOAAAQKmlDhaprAn0Hy+EAAAA/Cx1sEhlTaDfeCEEAADgb6mDRYnKmkCf8UIIAADA39IXuKGyJtBvsxdCs1RUXggBAADYWfpgUaKyJtBnvBACAADwQ7A4h5Q0oJ94IQQAAOBu6ccsAlgevBACAACwR7AIAAAAADiEYBEAAAAAcAjBIgAAAADgEIJFAAAAAMAhBIsAAAAAgEMIFgEAAAAAhxAsAgAAAAAOIVgEAAAAABxCsAgAAAAAOIRgEQAAAABwyMp0Om17HQAAAAAAiaFnEQAAAABwCMEiAAAAAOAQgkUAAAAAwCFrba+ADWPMqqSnJH05y7IPt7w6gBdjzK9LOjn3o0eyLPuzVlYG8GCMOS3phyRdlPQPJP15lmX/W7trBdgxxpyR9DuS/tP+j9YkHcuy7JtaWynAgzHmg5K+R9L/I+m/kvS/ZFn2/7a7VuirTgSLkj4h6S8k3db2igA1/FmWZefbXgmghl+Q9NEsy/6TMWYg6RvaXiHAwVVJD2VZ9pQkGWPOSTrW7ioBXn5O0rdmWfbvjTEflXRe0gfbXSX0VfLB4v7bky9LelHSt7W7NkAtbzfGfFzSRNI1SY9lWbbT8joBVowxXyXpv5D0gDHmNkm3S3q0zXUCXGRZ9rr2spRm/ntJ39XS6gB1fEXSfybp3+//94/bXR30WdJjFo0xRntvTn6+7XUBAvg5Sf8iy7KfknRG0ifbXR3AyddIepukv86y7Gck/bWkX2x3lQA/xpj3Svp3WZZttr0ugIeHJP24MeZ/l/ReSb/W8vqgx5IOFiX9I0mXjTE/Iun9kt5ljPkRY8zRltcLcJZl2R9lWTab2PS3Jd3f5voAjt7Y/+/v7//3C5L+65bWBajro5L+j7ZXAnBljPlqSb8i6QNZlv1zST8t6VfbXSv0WdJpqFmWfWr2/40xH5b0bfu9MkCnGGOOS/qhLMv+1/0ffa2kv2xxlQBX/0F7wwHeKenPtdfTSBtG5xhjvlHSK1mWvdb2ugAe7tBeZ8/l/X+/IulIa2uD3ks6WJwxxrxP0ndKutsY82CWZY+1vU6Ao7GkbzTG/LikHUlfJ+mft7tKgL0sy3aMMf+dpI8bY/5c0tdL+nC7awV4+SFJn6r8FJCg/aI2Pyvp54wxf6O9seQ/0O5aoc9WptNp9acAAAAAAEsl9TGLAAAAAIAWECwCAAAAAA4hWAQAAAAAHEKwCAAAAAA4hGARAAAAAHAIwSIAAAAA4BCCRQAAAADAIQSLAAAAAIBD/n+Fm2RCV+be1QAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "pseudo_samples = distribution.sample(500, rule=\"additive_recursion\")\n", "\n", "pyplot.scatter(*pseudo_samples)\n", "pyplot.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Chaos expansion\n", "\n", "To be able to do point collocation method it requires the user to have access to sampler from the input distribution and orthogonal polynomials with respect to the input distribution.\n", "The former is available above, while the latter is available as follows:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "execution": { "iopub.execute_input": "2021-05-18T10:57:26.735997Z", "iopub.status.busy": "2021-05-18T10:57:26.735682Z", "iopub.status.idle": "2021-05-18T10:57:26.757446Z", "shell.execute_reply": "2021-05-18T10:57:26.757112Z" } }, "outputs": [ { "data": { "text/plain": [ "polynomial([1.0, q1-3.0518, 0.2121*q1+q0-6.4705])" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "expansion = chaospy.generate_expansion(1, distribution, rule=\"cholesky\")\n", "\n", "expansion.round(4)" ] } ], "metadata": { "jupytext": { "formats": "ipynb,py:percent" }, "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.9.5" } }, "nbformat": 4, "nbformat_minor": 4 }