Skip to content

Instantly share code, notes, and snippets.

@oxinabox
Created March 11, 2020 11:52
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save oxinabox/1674e5b8e0e86cfa80535ba615ac545f to your computer and use it in GitHub Desktop.
Save oxinabox/1674e5b8e0e86cfa80535ba615ac545f to your computer and use it in GitHub Desktop.
NonNegative Matrix Factorization
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[32m\u001b[1mActivating\u001b[22m\u001b[39m environment at `~/Documents/misccode/NNMF/Convex/Project.toml`\n"
]
}
],
"source": [
"using Pkg: @pkg_str\n",
"pkg\"activate .\"\n",
"# use EricPHanson's branch as it has updated packages and Compat setup right\n",
"#pkg\"add Plots Distributions Convex ECOS SCS https://github.com/ericphanson/SequentialConvexRelaxation.jl.git\""
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"┌ Info: Precompiling Convex [f65535da-76fb-5f13-bab9-19810c17039a]\n",
"└ @ Base loading.jl:1273\n",
"┌ Info: Precompiling SequentialConvexRelaxation [154aca32-cc25-11e9-1929-3b69ed34445f]\n",
"└ @ Base loading.jl:1273\n",
"┌ Info: Precompiling ECOS [e2685f51-7e38-5353-a97d-a921fd2c8199]\n",
"└ @ Base loading.jl:1273\n",
"┌ Info: Precompiling SCS [c946c3f1-0d1f-5ce8-9dea-7daa1f7e2d13]\n",
"└ @ Base loading.jl:1273\n"
]
},
{
"data": {
"text/plain": [
"Plots.GRBackend()"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using Convex\n",
"using SequentialConvexRelaxation\n",
"using ECOS\n",
"using SCS\n",
"using Distributions\n",
"using LinearAlgebra\n",
"using Plots\n",
"gr(fmt=:png)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"WARNING: redefining constant true_cluster_assignment\n",
"WARNING: redefining constant cluster_distributions\n",
"WARNING: redefining constant data\n"
]
},
{
"data": {
"text/plain": [
"5×10 Array{Float64,2}:\n",
" 101.024 101.004 101.03 101.023 101.024 … 98.9384 99.9862 99.9652\n",
" 101.058 100.993 101.013 100.972 101.038 98.9286 99.9543 100.056 \n",
" 100.949 101.003 100.992 101.022 101.065 98.9416 99.9211 99.9866\n",
" 101.005 100.968 101.0 101.014 100.985 99.0079 100.023 100.029 \n",
" 101.027 100.97 100.978 100.954 100.993 99.0447 100.02 99.9276"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"const true_cluster_assignment = [fill(1, 5); fill(2, 3); fill(3, 2)]\n",
"const ndims = 5\n",
"const cluster_distributions = [ # Use IsoNormal because only allow circular distributions\n",
" MultivariateNormal(fill(101, ndims), 0.001I),\n",
" MultivariateNormal(fill(99, ndims), 0.002I),\n",
" MultivariateNormal(fill(100, ndims), 0.003I), \n",
"]\n",
"\n",
"# features × observations\n",
"const data = hcat((rand(cluster_distributions[cl]) for cl in true_cluster_assignment)...)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3df3hU9Z3o8e8kmckPQgSDGQihChjEAg1GCgEWBWqhStckvd2UEMCio1KhEJaYFUN6qZRSSLhK3Qa8YCksWIZaulK7NuYp9nkSpVVWHbqYln2kIINhMRMw45DA5Jzv/ePU3DRMkjkzJ5zMzPv1nD9mvufHfDw99JPvz2ORUgoAAGJVnNkBAABgJhIhACCmkQgBADGNRAgAiGkkQgBATCMRAgBiWoLZAQAA0Ju33377hRde8Hg86enp3/nOd6ZMmaKVe73eLVu2nDx5cuLEieXl5YMHD76+JJjrUyMEAAxciqL88Ic/XLFixeHDhx966KHq6urOXU6n0263O53OjIyMQ4cOBSwJBokQADBwKYqybt26u+66q7293Wq1Dho0qHNXQ0NDfn6+zWbLz8+vr68PWBIMmkYBAAOXzWabNm1aW1tbfn6+EOK5557r3OXxeOx2uxDCbre3tLQELAnGwE2EH3/8cVpaWmpqqlEX7OjoSEgYuP+9A5+qqkKIuDhaEULHQxgmRVHi4uIsFovZgUQqKaWqqvHx8f35I61SKLpOOHv23NYtO7sVPv3001lZWZ1fk5OTjxw58qtf/aqmpuYnP/mJViil1B4G7b8rYEkwBu6/yWefffaBBx6YM2eOURe8fPnysGHDjLpaDGpvbxdCpKSkmB1IBOMhDJPP50tKSrLZbGYHEqlUVW1tbR06dGg//oRskOKSrlMuXPjw4sWL5eXlXQtvvvlm7UNTU9Orr7766KOPJicn33///T//+c87j0lPT7948WJWVlZzc7P2L+v6kmDw1z0AwDCqVFTZoXNTMjIypv69zr+509PTf/Ob3/zpT3+SUv7+97+//fbbhRAul0sIkZeXV1tbK6Wsra2dMWNGwJJgkAgBAIaRokNKnVuvTak2m+373/9+TU3NN77xjTfeeGPNmjVCiLKyMiHEkiVLTp8+XVxcfObMmZKSkoAlwRi4TaMAAAghcnJyduzY0bWkrq5OCJGamrpp06au5deXBINECAAwjFQVKTv0ntJPwQSJRAgAMIyUuhOhkCRCAMANpyhKQ0PD+fPnc3JyJkyYYHY4ZiIRAkDMuXDhwj333OvxeD77zHfTTYNnzJj5q18dNmaCpuzQ3TRqdo2QUaMAEHMWLSr58MPTLS2Xr13zf/JJy9Gjb3QbjRIyrWlU56hRnU2pRiMRAkBskVI2Njaqquws8Xo/+8UvfmHQxRWpdujc6CMEANxoUoiuK5BZLJbYrReRCAEgtlgslvj4BCEsQlg+L5HjxmUbcvFQR42auX4siRAAYo6i/F3ukdLy17/+1ZAryxAGywjF3GQUu3VhAIhZUspuJZ988okpkQwEJEIAiDnDh9u7frVYxKeffvrZZ58ZcGmpCLVD58b0CQDAjbVv376hQ28SQmqblOq5c+cfeughAy6txmQiVBRl2bJl2mev17t+/frCwsLKykqv19ttb0/HAABupC996Uvz58/7PBEKIeL8/o633nor+JfZ9kgqQnbo2yJ9HuHhw4dXr17tdru1r06n0263O53OjIyMQ4cOddsb8JgwAwAABOPKlSvf+c53vvjFLw4fPmL69Onnz58XwiJEXOeoGSnllStXzA3SFOEO1BkzZkxmZmZlZaX2taGh4ZlnnrHZbPn5+d/73vdKS0u77g14zCOPPBLwylLKTz755Ny5c50lqampaWlpIYeqKIqimFwBj2ja3eMehoOHMEzK58wOJCItXLiwtrbu2jW/EOJ//ufiTTcNTkpKbG+/1nlAcnJKcnJyL7c3Li6u72XYZIdQddbwzG4aDTcRTp48uetXj8djt9uFEHa7vaWlpdvegMf0dGWfz1ddXZ2ent5Z8rWvfW3p0qUhh+r1em02W8ino62tTQjR0WFyI0ZE4yEMk9frvXbtGvcwBF6v99ixP2hZUAghhMXna7vllvRPP229cqXNYrHcfPOQ7dufa21t7eUiaWlp8fHxffyS1keoS6Qnwm6klNrfC1LKntqagzlGCJGamrply5Y5c+YYFZuiKEOHDjXqajEoMTFRCJGSkmJ2IBGMhzBM8fHxSUlJNywRLlq06NNPP70xv9Xf2traEhNtdvst2tdPP21tb29PS0tLT0+32WyzZ89euXLl6NGjzQ3SLAYnwvT09IsXL2ZlZTU3Nw8bNizkYwDAdH/+859fe+01s6Mw3nvvvbdkydKrV9v+8pf/FkKkpKSoqtyyZYsxV1c7hOrv+7Cu9K5EYzSDE2FeXl5tbe3DDz9cW1s7Y8aMbntdLldOTk7vxwDAABEXF6f140SZm2++WQgp5d8GS1650nbq1Kn9+/d/+9vfNuDq2qhRvaeYyuB5hEuWLDl9+nRxcfGZM2dKSkq67S0rK+vzGABAf+vW1XflStsbb7xhzKWl3kmEHaYnQmNqhHV1ddqH1NTUTZs29bRX+xDwGADADWO1Wrt+tVgsMdtBKFh0GwBi0JAhQy5cuNjR8beq2LBhNy9ZssSYSzNqFAAw8OXl5Y0de/uxY29ZLJaUlJQXXnhh7NixhlzZIlWLzqZOvccbjkQIADHHarX++7//qqOj48qVK+EsVBIdSIQAEKMSEhKMz4JS0d3UKcNe4DQ8JEIAgGEsqmrRmwjDX+k7PCRCAOh3b7755pEjR9xu9xe+8IWCgoJp06aZHVG/UfXXCBksAwBR7OOPP/7617/+0UcfeTyXtJJdu3aPHTv2178+kpGRYW5s0PBiXgDoL21tbffcc8/775/weC4LYdE2j+fS8eP/ee+99/r9OpciC+T6d76aTRVS0bcJk5tGSYQA0F+effbZ8+ebpOxerqryo4/cO3bsCPP6Ad/5ai6tj1DnRiIEgCi1f//+9varAXddudL2s5/9TNfVnn322V/+8pedn51O55gxYxYvXhxmkCARAkB/6f31fs3NzbquNmvWrIaGBiGE3+9vaGiYM2fO5MmT8/LywgrRcNpgGb2bqUiEANBfen+fe1ycvv8Hnjx58kcffXTp0qV33333tttuG5hjbSxSf9Mo8wgBIFplZ2e73R8LETAdypycHF1XS0hImD59+ltvvXXixIn77rvPkAiNF4HTJ6gRAkB/Wb9+/dChQwLuuvnmoevWrdN7wVmzZh09evSdd96ZNWtW2NHhb0iEANBf5s6dW1T0TzfdNLhb+ZAhg5ct+3YI3Xt33XXXhx9+mJOTk5qaalCMRlOlRVX1bdcPq72xaBoFgH60c+fOadOmff/7329vb1dV1WKJGzQo5Qc/+MGiRYtCuJrNZrv11lu7tYt2vvN1QAhhrdG+pk/U19fv3bu3ubl59OjRa9euzcrK0sq/+tWvdjvyzjvvbGxs1D4vWLCgtLQ0mN8nEQJA/1q2bNmyZcuampqamppGjhxpt9tDu05HR8eZM2c+/vjjL3/5y8ZGOJA1NTVVVVVt3bp19OjRr7zySlVV1fbt27VdR44c6TzM6XT6/f7XXnvN6XQmJycLIeLj44P8CZpGAeBGGDFiRG5ubshZUAhx7Nixp5566rvf/a7NZjMwMIOpqu65E72+j7CpqWnu3Lnjx49PTEycN2/euXPnOnclf+7ChQsffPBBfn6+oigVFRVFRUWbN2/2+XxBhkyNEAAiw6xZswb+GBmL0U2jubm5ubm5QghFUfbu3Tt79uxuB/j9/m3btj355JOXL1/Ozs5evnx5RkbGjh07ampqKioqgvl9EiEAwDiq1P1aJVWtr69fsWJF17Knn3565MiRnV+PHz++e/fuKVOmXL+w6ssvvzx+/Phbb71VCFFdXa0VOhwOh8MR5O+TCAEAJrvjjjuWLl3atWTo0KHaBynlrl27Ghsb169f3zlMppOiKK+++urWrVuFEKdOnfL7/RMmTBBCWK1Wq9Ua5K+TCAEAhrHon1BvkWpGRkZP72g8ceLEsWPHnn/++fj4+La2NiFEcnKyy+XSliN4//33b7nlFq3u2N7evnHjxurq6szMzP3798+cOTPIAEiEAADjSNXYlWVcLpfb7S4sLOwsqaurKysr0yaNvP7663fffbdWPmnSpJKSksrKSp/PN3Xq1JUrVwb5+yRCAMDAtXTp0m6tpqLL1Mmuq/NYLJaCgoKCggK9P0EiBAAYJ4S1Rll0GwAGpqtXr/7oRz8yOwqDSSnPnz/ff9e36B81alFZYg0ABqSDBw92rtcVHVRVbWtry87OHjt2rNmxDCAkQgAIbMKECdpY/KihKEpra2vnzIR+IVXd8whpGgUARA81hFGjJEIAQNQIoUZodiJk0W0AQEyjRggAMI5KHyEAIIZF4vQJmkYBADGNGiEAwDghvI+w1xfz3gAkQgCAcUJ5HyErywAAokYEDpahjxAAENOoEQIAjBPKhHqaRgEAUUOVuhOb2YmQplEAQEyjRggAMI6U+ge/sLIMACBq0EcIAIhp9BECABBZqBECAIwTwsoyZk+oJxECAIwTgX2ENI0CAGIaNUIAgHFCGCwjGTUKAIgakrdPAABimap/frzJeZA+QgBAbKNGCAAwTgg1QpNnT5AIAQDGCWGpUbOnEdI0CgCIbdQIAQDGicDBMiRCAIBx6CMEAMQ0GXk1QvoIAQAxjRohAMA4Un8Nz+waIYkQAGAYqVqkatF5Sh+ZsL6+fu/evc3NzaNHj167dm1WVlbnrlWrVjU2NmqfFyxYUFpa6vV6t2zZcvLkyYkTJ5aXlw8ePLjPAGgaBQAMXE1NTVVVVWVlZU6nc/r06VVVVZ27pJRut9vpdB45cuTIkSNPPPGEEMLpdNrtdqfTmZGRcejQoWB+gkQIADCO/HzgaPBbrxXCpqamuXPnjh8/PjExcd68eefOnevc5fF4FEWpqKgoKiravHmzz+cTQjQ0NOTn59tstvz8/Pr6+mBCJhECAIyjCqFadG89y83NLS0tFUIoirJ3797Zs2d37mppacnOzl67du2BAwcGDRpUU1MjhPB4PHa7XQhht9tbWlqCCZk+QgCAYaS0SKmzj1BaGhoavvvd73YtfOqpp0aOHNn59fjx47t3754yZcqyZcs6C8eNG1ddXa19djgcDodDCCGltFgs2gc1uBdCkQgBACbLzs5etGhR15KhQ4dqH6SUu3btamxsXL9+fddhMkKIU6dO+f3+CRMmCCGsVqvVahVCpKenX7x4MSsrq7m5ediwYcH8OokQAGAc2UdTZ8BT7Hb79OnTA+48ceLEsWPHnn/++fj4+La2NiFEcnKyy+XKyclpb2/fuHFjdXV1Zmbm/v37Z86cKYTIy8urra19+OGHa2trZ8yYEczvkwgBAMbpq88v0Cm97XS5XG63u7CwsLOkrq6urKysrq5u0qRJJSUllZWVPp9v6tSpK1euFEIsWbJk8+bNxcXF2dnZTz31VDC/TyIEAAxcS5cuXbp0abfCuro6IYTFYikoKCgoKOi6KzU1ddOmTbp+gkQIADBMCINlhN7jjRZ6IlQUxeFw7NmzRwhx/Uz+Xub2X78QQJj/DQCAgUKbPqGL2YkwxHmEhw8fXr16tdvt1r5eP5O/p7n9ARcCAADALCEmwjFjxixevLjz6/Uz+Xua2x9wIQAAQJSQcULVv5kqxKbRyZMnd/16/Uz+nub2awsBLF++PCMjY8eOHTU1NRUVFQF/oq2t7ac//enRo0c7S+6+++777rsvtICFED6fLykpKeTTceXKFSFEkBNUERAPYZh8Pl9HR4fNZjM7kEilKIrP59Pm24UgJSUlLq6vpBXKhPrQwjGMMYNlrp/J39Pc/oALAQSOLCFh7Nix2kxJTVZWVkJC6AEnJCSEczq0u8c9DAcPYZji4+O5h+GwWCz9fgNDmD5hdh+hMbfj+pn8Pc3tD7gQQEBWq/Xee++dM2eOIREKIRITE/ljPBzaHzTcw3DwEIbp2rVrSUlJ1AhDpiiKdg/NDmRgMaZlVpvJL6XsnMl/fYnL5RJCtLe3b9iw4ezZs36/v3MhAABAdJDyb68k1LWZG7MxiXDJkiWnT58uLi4+c+ZMSUlJwJKysjIhROdCAAsXLvR6vb00jQIAIo+M072Z/R6ksJpGtbn9ItBM/utLelkIAAAAs9DnDAAwTAhNndLsoegkQgCAYaTU3+dndh8hiRAAYBzVonuCvDS5j9DknwcAwFzUCAEAhgnh7RNRsrIMAABChPSGerP7CGkaBQDENGqEAADjqHFS52AZafZgGRIhAMAwUgr9fYRMnwAARI1Q+gj7J5Kg0UcIAIhp1AgBAIaR9BECAGJZCPMITX8xL02jAICYRo0QAGCYUN4+YXaNkEQIADCQRW+fH4kQABBFVJZYAwAgolAjBAAYhj5CAEBMkyKE1zDRNAoAgHmoEQIADCNVi96VZUwfLEMiBAAYJqQ31JMIAQDRIhIHy9BHCACIadQIAQAGCqFptJ8iCRaJEABgmFAGy/S1JFt9ff3evXubm5tHjx69du3arKys3netWrWqsbFRO2DBggWlpaW9X59ECAAYuJqamqqqqrZu3Tp69OhXXnmlqqpq+/btveySUrrdbqfTmZycLISIj4/v8yfoIwQAGEYbNap36+WCTU1Nc+fOHT9+fGJi4rx5886dO9f7Lo/HoyhKRUVFUVHR5s2bfT5fnzFTIwQAGEZK/aNGez0+Nzc3NzdXCKEoyt69e2fPnt37rpaWluzs7OXLl2dkZOzYsaOmpqaioqL3AEiEAADDhDSPUDQ0NKxevbprYXl5+ciRIzu/Hj9+fPfu3VOmTFm2bFm307vtGjduXHV1tbbL4XA4HI4+AyARAgBMdvvttxcVFXUtGTJkiPZBSrlr167Gxsb169d3HSbT065Tp075/f4JEyYIIaxWq9Vq7fPXSYQAAOPoHzUqZdzw4cNnzpwZcO+JEyeOHTv2/PPPx8fHt7W1CSGSk5NdLldOTk7AXe3t7Rs3bqyurs7MzNy/f39Pl+2KRAgAMEwITaOi1+NdLpfb7S4sLOwsqaurKysrq6urC7hr0qRJJSUllZWVPp9v6tSpK1eu7PP3SYQAgIFr6dKlS5cu7VZYV1fX0y6LxVJQUFBQUBD8T5AIAQCGicT3EZIIAQCGCW3UqLmYUA8AiGnUCAEAhonE1zCRCAEAhpHSova1iPZ1p5jcNkkiBAAYJhLfUE8fIQAgplEjBAAYJpRFt80eNUoiBAAYhqZRAAAiDDVCAIBhIrFGSCIEABhH6k9sJEIAQNQIZa3RfgolaPQRAgBiGjVCAIBhVGlRddYI9R5vOBIhAMAwhr+Y9wagaRQAENOoEQIADMP0CQBAbNM/fYJECACIHqHUCPsplKDRRwgAiGnUCAEAhgllQj1NowCAqBHCPELTEyFNowCAmEaNEABgHKZPAABiGfMIAQAxTZUWVerrdJOCRAhEEVVVDxz4t9dff9VqtT744D8VFBSaHRGAPpAIASM99FDxiEz/uooZfr+y/bkdf/zjm5s3V5sdFHDjMH0CiGnvvPNOh3LxR1se1b7u2v3QV+ZUnz9/fuTIkeYGBtwwMgJfw8T0CcAw77//3uzZt3d+tVgss+4Zd+LECRNDAtAnEiFgmMzMkW53a9cSt/tyZmamWfEAN542alTnZnLMJELAMPfee+9vX/vA5fpI+/r7Nz44/eGnEydONDcq4EZS5d8Wlwl+Y9QoED1SU1P37395zZonPv30E0WRo0aNfumlX8bHx5sdF3DjMI8QiHV33HHHf/zH765duxYXF5eQwL8vIALwDxUwns1mMzsEwCQRWCMMq49QUZRly5Zpn71e7/r16wsLCysrK71eb8CSno4EAEQHVejrIFSDSJz19fUOh6OgoGDNmjVut7uzPPi807vQE+Hhw4dXr17dGZPT6bTb7U6nMyMj49ChQwFLejoSAICAmpqaqqqqysrKnE7n9OnTq6qqOncFn3d6F3oiHDNmzOLFizu/NjQ05Ofn22y2/Pz8+vr6gCU9HQkAiA5SWPRvvWlqapo7d+748eMTExPnzZt37ty5zl3B553ehd5HOHny5K5fPR6P3W4XQtjt9paWloAlPR0ZUGtr66OPPpqamtpZUlhYuGLFipAD7uW3EIy2tjYhxJUrV8wOJILxEIbJ6/UmJibSBRsyRVG8Xq+iKKGdPmTIkD6HgBm+skxubm5ubq4QQlGUvXv3zp49u3NX8Hmnd4YNlpFSWiwW7YOqqgFLejoyoLS0tF27ds2ZM8eoCIUQw4YNM/BqsUZLgSkpKWYHEtl4CMNhs9mSkpJIhCFTFMVmsw0dOrT/fkLKEAa/WN588801a9Z0LXryySe7LkZx/Pjx3bt3T5kypXNgitCTd3pnWCJMT0+/ePFiVlZWc3Oz9k/9+pKejgQAxLIxY8Z84xvf6Fpy0003aR+klLt27WpsbFy/fn1WVlbXY4LPO70zbGWZvLy82tpaKWVtbe2MGTMClrhcroDlAIDooMpQthEjRsz6e4MGDdIueOLEiWPHjm3cuDE9Pb2trU3ro+kpm4SWXwxLhEuWLDl9+nRxcfGZM2dKSkoClpSVlQUsB6JVa2vrv/7rj9esWbFzZw3dq4gFIa012ltTqsvlcrvdhYWFD35O9JxNQssvFmn6cqc9ePLJJx944AED+whpiQ0TfYR6nTt3Lj9/vuPRGRMnjXz3P8/u/7f/PHDg5TvuuMPsuCJYa2srfYThUBSltbW1X/sI31xU7v3vs7pO+Yuv5S+zs2tqavoppD6xsgzQX556as3254tmzhwnhPiHf7jj9tszqqo27d69z+y4gH4khe5Ro5G9sgyAXpw69WctC2q+Om/Sn/7kMjEe4AZQhUX/ZjISIdBfEhOTrl71d371etu6zosFMECQCIH+UlDwzY3P/FrrhldV+b+/98o//uM3+jwLiGiGD5a5AUiEQH/5538uv3Z1ZN7UTUsW/3Tal38wLD1n8eKHzA4K6F+SF/MC6BQXF1ddvf3KlSsfffTRbbfdlpSU1NzcbHZQQP/ixbwAuktJSRk/frzZUQDoEYkQAGAYVQi9o0BNHzVKIgQAGCaUwS9mN40yWAYAENOoEQIADBPCyjKq2Qt9kggBAIbRXjqv95R+CiZIJEIAgGFU2ccb569n+vQJ+ggBADGNGiEAwDBS6u7zY/oEACB6RGIfIU2jAICYRo0QAGCYkNYa7adYgkUiBAAYJoQl1szOgyRCAICBIvDtE/QRAgBiGjVCAIBhaBoFAMS0SHwxL02jAICYRo0QAGAYNYS3T5g9oZ5ECAAwkt4+P/oIAQDRg7dPAAAQYagRAgAMI2kaBQDEMin1D5YxOxPSNAoAiGnUCAEAhgmhadR0JEIAgGHUUJpGmUcIAIgWqtTd52f6+wjpIwQAxDRqhAAAA1mkziXT9B5vOBIhAMAwqv7pEH2+tklRFIfDsWfPns6Sr371q92OufPOOxsbG7XPCxYsKC0tDT4AEiEAYOA6fPjw0aNH3W5318IjR450fnY6nX6//7XXXnM6ncnJyUKI+Ph4XT9BHyEAwDCqsOjdem8aHTNmzOLFi7sVJn/uwoULH3zwQX5+vqIoFRUVRUVFmzdv9vl8umImEQIAjCOF1L/1YvLkyXl5eQF3+f3+bdu2rVix4vLly9nZ2WvXrj1w4MCgQYNqamp0hUzTKADAMCG8fUKVlrfeemvt2rVdC9euXZuZmdn7iS+//PL48eNvvfVWIUR1dbVW6HA4HA6HrgBIhAAAk912220PPvhg15Kbbrqp91MURXn11Ve3bt0qhDh16pTf758wYYIQwmq1Wq1WXb9OIgQAGCa0t09kZmbee++9wRzscrlycnKEEO+///4tt9wycuRIIUR7e/vGjRurq6szMzP3798/c+ZMXQHQRwgAMIy2soyuTdfKMmVlZdqH119//e6779Y+T5o0qaSkpLKycuHChV6vl6ZRAEC0qaur6/Zh3bp1nXstFktBQUFBQUFoFycRAgAMI4VF1blSjGp24ySJEABgmD6nQwQ8xVz0EQIAYho1QgCAYdQg1g7txuwKIYkQAGAcGcKLeXn7BAAgasgQ+vzMrhLSRwgAiGnUCAEAhgmhj1Dv8YYjEQIAjMP0CQAAIgs1QgCAYVT9K8v0/mLeG4BECAAwjLaOtt5TzEUiBAAYJrTXMJmLPkIAQEyjRggAMIzU3zRq+qhREiEAwDCh9BH2TyTBo2kUABDTqBECAAwjhUXvdAimTwAAoofUPx2CPkIAQPRg+gQAABGGGiEAwDCsLAMAiGmROI+QplEAQEyjRggAMEwkDpYhEQIADBOJK8sYkAgvXbq0c+fO9957Ly4ubtq0aY8//nhKSkrAws5TVq1a1djYqH1esGBBaWlp+GEAAEwn9Se2aKgRbtu2bezYsQcOHJBS7t27d9++fcuXLw9YqB0vpXS73U6nMzk5WQgRHx8ffgwAAITGgETocrnWrVtntVqFEAsXLly+fPny5csDFmrHezweRVEqKircbvddd921Zs0am812/WU7Ojrefvvtq1evdpZkZWVlZ2eHHOe1a9e6Xg16aXePP1zCwUMYpqtXr1osFmn6KMOIpSjK1atXQ34IbTabxdLHcmhS6h4Favr/ngYkwnHjxh08ePBb3/qW3+9/6aWXWlpaeirUtLS0ZGdnL1++PCMjY8eOHTU1NRUVFddf9tq1a++8886FCxc6S/Ly8m699daQ47x27dq1a9dCPh1+v18IwT0MBw9hmPx+f1wcY91DpyiK3+8P+SFMSEjo80/hGJ1HWF5e/uMf/7i4uHjIkCGFhYWDBw/uqVAzbty46upq7bPD4XA4HAEvm5KSsmLFijlz5oQfoebq1atdw4Be2j+Arn290IuHMExSyqSkpIBtSAiGoihSSh7CbgxIhElJSRs2bNBaQV0u16hRo3oq1Jw6dcrv90+YMEEIYbVatWMAAFEghOkTpjOgkWHXrl3PPfecz+fzeDwvvvhiQUFBT4Uul0sI0d7evmHDhrNnz/r9/uZvZgwAAA+/SURBVP3798+cOTP8GAAAA4H8vHVU12YuAxLhY4895vV6Fy1aVF5ePn/+/FmzZvVUWFZWJoSYNGlSSUlJZWXlwoULvV5vT02jAICII/82XEbfZm7MBjSNpqWlPfPMM8EU1tXVCSEsFktBQYFWRwQAwFysLAMAMEwITZ3RMH0CAABNf6w1qiiKw+HYs2dP18LrVyjzer1btmw5efLkxIkTy8vLgx8cSyIEAAxchw8fPnr0qNvt7loYcIUyp9Npt9u/973vvfDCC4cOHXrkkUeC/AmmpgIADCP1DxntfW3SMWPGLF68uFth5wplRUVFmzdv9vl8QoiGhob8/HybzZafn19fXx98zCRCAIBhpP6td5MnT87Ly+tWqK1Qtnbt2gMHDgwaNKimpkYI4fF47Ha7EMJut3ddzqxPNI0CAAwT2hJrx44de/LJJ7sWrlmzJjMzs6dTAq5QJqXUlkKVUqqqjndgkAgBACa79dZbH3jgga4laWlpvRwfcIWy9PT0ixcvZmVlNTc3Dxs2LPhfJxECAIwT0tsnMjMzg1xZ2uVy5eTktLe3b9y4sbq6OjMzs3OFsry8vNra2ocffri2tnbGjBnBB0AfIQDAMKqQejddK8v0skLZkiVLTp8+XVxcfObMmZKSkuCvSY0QADDQaQuTiV5XKEtNTd20aVMIFycRAgAMI/WvFMPKMgCA6KG1duo9pZ+CCRJ9hACAmEaNEABgGCn6WCkm4CnmIhECAAwjpZQ6O/3oIwQARA9Vf41Q7/GGo48QABDTqBECAAwjQ3gxb/9EEjwSIQDAMFIIXSvFCP3HG46mUQBATKNGCAAwTGivYTIXiRAAYBipf2UZ05tGSYQAAMNInW+TEAMgEdJHCACIadQIAQCGicRFt0mEAADDhNRHaDKaRgEAMY0aIQDAMKoUqs5VtPUebzgSIQDAMCGMGjW9cZRECAAwTCQOlqGPEAAQ06gRAgAME8KoUdNrhCRCAIBhQmgaZfoEAABmokYIADCMFFIKVecp+o43HIkQAGAY+ggBADFNtaiqRWeN0ML0CQAAzEONEABgGCmkqrPPT+/xhiMRAgAMowpVb2LjxbwAAJiJGiEAwDBSqEyfAADELtUiVJ2jQPUebzgSIQDAMJI+QgAAIgs1QgCAYUIYNcr0CQBA9OiPtUYVRXE4HHv27OlaWF9fv3fv3ubm5tGjR69duzYrK0sIsWrVqsbGRu2ABQsWlJaWBhMAiRAAMHAdPnz46NGjbre7a2FTU1NVVdXWrVtHjx79yiuvVFVVbd++XUrpdrudTmdycrIQIj4+PsifoI8QAGAYKRRV59Z7jXDMmDGLFy/uVtjU1DR37tzx48cnJibOmzfv3LlzQgiPx6MoSkVFRVFR0ebNm30+X5AxUyMEABhGDWWJtd5GjU6ePPn6wtzc3NzcXCGEoih79+6dPXu2EKKlpSU7O3v58uUZGRk7duyoqampqKgIJgASIQDAMFIoUih6T/nDH/7wL//yL10LS0tLR4wY0fuJx48f371795QpU5YtWyaEGDduXHV1tbbL4XA4HI4gAyARAgBMNmrUqHnz5nUtGTx4cC/HSyl37drV2Ni4fv16bZiMEOLUqVN+v3/ChAlCCKvVarVag/x1EiEAwDChTagfOXLkV77ylWAOdrlcOTk5J06cOHbs2PPPPx8fH9/W1iaESE5Obm9v37hxY3V1dWZm5v79+2fOnBlkACRCAIBhQuoj1HF8WVlZXV2dy+Vyu92FhYWd5XV1dZMmTSopKamsrPT5fFOnTl25cmWQ1yQRAgAGurq6uq4fli5dunTp0m7HWCyWgoKCgoICvRcnEQIADBPSYBlWlgEARIsQ3lBveiJkQj0AIKZRIwQAGIamUQBATOuPRbf7G4kQAGAYVaqq1JkIJS/mBQDAPNQIAQCGkUKlaRQAEMOklLqbRkmEN1xra2ttbW2LxzN9xowvfelLZocDADBTzCXC4++88+1vfmuGbehgGf9zdXvOffdu/787zQ4KAKJESE2jJg+WiblE+Njih55OvzPDliKEuF+IZ39X/3pt7bz5882OCwCigRSq3nmEQvfxBoutUaMXLlxI7ZBaFtTcaxtW9+vfmBgSAEQTVajaDAo9G9MnbqDExMSOv++V9UslKSXZrHgAAKaLrUQ4dOhQS1rqX698qn2VQr7e/sn9BfnmRgUAUUSV+jdzIzamj/DSpUs7d+5877334uLipk2b9vjjj6ekpAQs1I73er1btmw5efLkxIkTy8vLBw8ebEgYwfjZLw4Wff3BMVcvpon499svlzzmmDFjxg37dQCIbjKE6RNmD5Yxpka4bdu24cOHHzhwYN++fampqfv27eupUON0Ou12u9PpzMjIOHTokCExBCk7O/udk38q3/dC0XM/+O3bb5VVrLuRvw4AGGiMSYQul6uoqMhqtdpstoULF9bX1/dUqGloaMjPz7fZbPn5+V3Lb4yEhIRp06bdf//9drv9Bv80AEQ37e0TujbTR40a0zQ6bty4gwcPfutb3/L7/S+99FJLS0tPhRqPx6MlIbvd3rW8K5/Pt2HDhp/85CedJV/5ylcWLlwYcpCXL1+Oj48P+XS0tbUJIa5evWp2IBGMhzBMXq83MTHRZrOZHUikUhTF6/WGfHpaWlqfD3AoTaNmjxo1JhGWl5f/+Mc/Li4uHjJkSGFhodbnF7BQI6W0WCzaB1UNfMuSk5MffvjhqVOndpYMHTo0LS0t5CCvXbsWzulISEgQQnR29CIEPIThS0pKIhGGTFEUIUTID2FcXN+NiLG71mhSUtKGDRusVqsQwuVyjRo1qqdCTXp6+sWLF7Oyspqbm4cNGxbwmnFxcV/4whfuvPNOQyIUQsTHx/PHeDi0u8c9DAcPYZjiP2d2IBGMG3g9Y/oId+3a9dxzz/l8Po/H8+KLLxYUFPRU6HK5hBB5eXm1tbVSytraWgZtAkD0kKrUvUXFqNHHHnvM6/UuWrSovLx8/vz5s2bN6qmwrKxMCLFkyZLTp08XFxefOXOmpKTEkBgAAKbTlliLxcEyaWlpzzzzTDCFdXV1QojU1NRNmzYZ8tMAAIQj5hbdBgD0H621U+cpUTFqFAAAIYQUMkZHjZroww8//Pjjj7/4xS+mp6ebHQsAxDopFSn19fnxhvrQXbp0qejrDyoXPPaEpMa2S/mLir//ox/quoLP59tete3thjczhtuX/3Npbm5uP4UKABiwIvjtE098+5F7Lomnb5nwyNCxVZl3v3PolX8/fDj4030+3z1TprYc/I9veuInvOd+LP+buk4HAFxPaxqNrLdPRGoilFKefN81PW249tUiLAWpIw8fOBj8FXY+/6+zlEFfGzJqeOKgCYPTK4d/6Xtl/9I/wQJA7JBC6tyi4+0TN56iKN1CT4qL9332WfBXOP7WsUnJQzu/Doq3JktLa2urQQECACJDpCbChISEofaMM23/P2/VX/lkzv3zg7/CqNtu++RaW9eSK4r/Rr4ZEQCikap/Y/pEqHb+28/+1/wHZvmG3hJnPan6rmQN27FyRfCnL33MsfRXD96eMmSINVEK8cvm07Pnz9OWAgcAhEj/2yeE2X2EEZwI77jjjmN/ev/w4cNNH517bOqX583Tl8YmTpy4eVdNRela0X71itIx/x8X/OjZ/9N/0QJAbFD1JzZqhGEYPHjwQw89FPLp87/2tfl//lpbW1tycrKBUQEAIkhkJ0JDkAUBwCghTYdg+sSN8otf/MLsECLbyZMn/+u//svsKCIbD2GY3nzzTbfbbXYUEay1tfW1117r39+QUkhV59ZH06iiKMuWLeta4vV6169fX1hYWFlZ6fV6A5YEL4YS4QsvvGB2CJHt7bff/uMf/2h2FJGNhzBMv/3tbz/88EOzo4hgly5dOnhQx3zrgeDw4cOrV6/u9geQ0+m02+1OpzMjI+PQoUMBS4IXQ4kQANDfQlhWpvem1DFjxixevLhbYUNDQ35+vs1my8/Pr6+vD1gSPPoIAQAGMnjU6OTJk68v9Hg8drtdCGG321taWgKWBG/gJsLk5OQ9e/bU1tYadUGv1/vUU08ZdbUYdPLkSSnl+fPnzQ4kgvEQhundd9/1eDy/+93vzA4kUnm93rNnz4b8EK5evXrEiBG9H/P+++/qvezbb7+9fPnyblH1/ltSSm2+nJRSVdWAJcEbuIlw3bp1b775poEXvO+++wy8WgziBoaPexgmbmD4CgsLQz63n9beGj9+/NatW3X9Vnp6+sWLF7Oyspqbm4cNGxawJHgDNxEmJyfz0ANA1EtLSwv+/+1dLldOTk5eXl5tbe3DDz9cW1s7Y8YMIcT1JcFjsAwAIGKUlZUJIZYsWXL69Oni4uIzZ86UlJQELAmeRfY1gQMAgChGjRAAENMGbh+hXpcuXdq5c+d7770XFxc3bdq0xx9/PCUlJWChdrzX692yZcvJkycnTpxYXl7OC5j03kAhxKpVqxobG7XPCxYsKC0tNSn2AUFRFIfDsWfPHhHo6erpeeM51IR29wQP4XV6v5Nd92p4AkU01Qi3bds2fPjwAwcO7Nu3LzU1dd++fT0VasJZhiAq6b2BUkq32+10Oo8cOXLkyJEnnnjCvNjN123xi+CXveA5FGHcPR7Cbnq/k0Eu0RKDoicRulyuoqIiq9Vqs9kWLlyorSwQsFATzjIEUUnvDfR4PIqiVFRUFBUVbd682efzmRe7+botfhH8shc8hyKMu8dD2E3vdzLIJVpiUPQkwnHjxh08ePCzzz67dOnSvn37tJUFAhZqwlmGICrpvYEtLS3Z2dlr1649cODAoEGDampqzIvdfJMnT87Ly+v8GvyyFzyHIoy7x0PYTe93stvegMfcyGgHjuhJhOXl5drY2VWrVo0YMUJr6Q5YqAlnGYKopPcGjhs3rrq6+vbbb09LS3M4HMePHzcv9gEn+GUveA6vF/zd4yHsXTBPF0+giKbBMklJSRs2bLBarUIIl8s1atSongo14SxDEJX03sBTp075/f4JEyYIIaxWq3YMNMEve8FzeL3g7x4PYe+Cebp4AkU01Qh37dr13HPP+Xw+j8fz4osvFhQU9FTocrnE58sQSClDWIYgKum9ge3t7Rs2bDh79qzf79+/f//MmTNN/g8YSK5/uq4v4TnsSfB3j4ewd70/XTyBnaInET722GNer3fRokXl5eXz58+fNWtWT4U9LUwQ4/TewEmTJpWUlFRWVi5cuNDr9TocDpP/AwaSYJa94DnsSfB3j4ewd70/XTyBnVhZBgAQ06KnRggAQAhIhACAmEYiBADENBIhACCmkQgBADGNRAgAiGkkQgBATCMRAgBiGokQABDTSIQAgJhGIgQAxDQSIQAgppEIAQAxjUQIAIhp/w9IIKsz50FcYgAAAABJRU5ErkJggg=="
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"scatter(data[1, :], data[2, :]; zcolor = true_cluster_assignment, legend=:bottomleft)"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"nnmf"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"\"\"\"\n",
" nnmf(V, latent_size)\n",
"\n",
"### Args\n",
" - `V` the data, a matrix of features × observations\n",
" - `latent_size` the number of latent variables. For clustering this is the number of clusters.\n",
"\n",
"Returns `W,H`, which are nonnegative matrixes such that ``W*H ≈ V``\n",
"\"\"\"\n",
"function nnmf(V, latent_size)\n",
" all(V.>=0) || throw(ArgumentError(\"V contains negative entries\"))\n",
" \n",
" W = Variable(size(V, 1), latent_size)\n",
" H = Variable(latent_size, size(V, 2))\n",
" \n",
" # The main goal `opnorm(W*H - V))` will be handled by the BilinearConstraint\n",
" # thus we set the problem to be a dummy value\n",
" problem = minimize(\n",
" 1,\n",
" [\n",
" W>=0,\n",
" H>=0,\n",
" ]\n",
" )\n",
" \n",
" # Constructing the bilinear equality constraint\n",
" # In general this is of the form A*P*B=C.\n",
" # Here A=W, P=I, B=H, C=V .\n",
" eye = Matrix(I, latent_size, latent_size)\n",
" reconstruct_constraint= BilinearConstraint(W, eye, H, V, λ=0.1)\n",
" eye_orth_inner = Matrix(I, size(H,2), size(H,2))\n",
" orthogonality_constraint= BilinearConstraint(H, eye_orth_inner, H', eye, λ=0.1)\n",
"\n",
" # The bilinear problem is the convex problem with the bilinear constraint\n",
" bilinear_problem = BilinearProblem(\n",
" problem,\n",
" [\n",
" reconstruct_constraint,\n",
" # orthogonality_constraint\n",
" ],\n",
" )\n",
" res = solve!(bilinear_problem, () -> SCS.Optimizer(verbose=0), iterations=10)\n",
" (W=W, H=H), bilinear_problem\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"scrolled": true
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"┌ Warning: Bilinear constraint 1 not satisfied (gap = 0.11234117682821705)\n",
"└ @ SequentialConvexRelaxation /Users/oxinabox/.julia/packages/SequentialConvexRelaxation/9jc2O/src/SequentialConvexRelaxation.jl:408\n"
]
},
{
"data": {
"text/plain": [
"((W = Variable\n",
"size: (5, 3)\n",
"sign: real\n",
"vexity: affine\n",
"id: 140…734\n",
"value: [0.1570462845229664 11.417119780254275 1.9790771743623983; 0.021141235108979976 11.299729854032721 2.297251144320985; 0.023640789088924436 11.484763948519115 1.844308740524739; 0.2620201397676304 11.303159303534848 2.2257761979127215; 0.47901684406437545 11.362428001936705 2.0262235624224973], H = Variable\n",
"size: (3, 10)\n",
"sign: real\n",
"vexity: affine\n",
"id: 114…020\n",
"value: [0.9090761727022365 0.8018998057589781 0.8156148196829522 0.793127657006549 0.7303882764054443 0.9560068400380448 0.8835317223972944 1.108916235520905 1.0571211581975082 0.689856014912743; 8.22145778191618 8.255554482564627 8.245501346939792 8.260128763462827 8.265130415875387 8.132590595705121 8.0724488049853 8.075654836102757 8.147332052116685 8.13525181295928; 3.5356276815757015 3.3466106276014056 3.4074959247868417 3.3267055113350907 3.316543654963961 3.048534318704198 3.386482639601667 3.333416861023218 3.4313193606284322 3.536654868595232]), BilinearProblem(minimize\n",
"└─ 1\n",
"subject to\n",
"├─ >= constraint (affine)\n",
"│ ├─ 5×3 real variable (id: 140…734)\n",
"│ └─ 0\n",
"└─ >= constraint (affine)\n",
" ├─ 3×10 real variable (id: 114…020)\n",
" └─ 0\n",
"\n",
"termination status: OPTIMAL\n",
"primal status: FEASIBLE_POINT\n",
"dual status: FEASIBLE_POINT, BilinearConstraint[BilinearConstraint(Variable\n",
"size: (5, 3)\n",
"sign: real\n",
"vexity: affine\n",
"id: 140…734\n",
"value: [0.1570462845229664 11.417119780254275 1.9790771743623983; 0.021141235108979976 11.299729854032721 2.297251144320985; 0.023640789088924436 11.484763948519115 1.844308740524739; 0.2620201397676304 11.303159303534848 2.2257761979127215; 0.47901684406437545 11.362428001936705 2.0262235624224973], Bool[1 0 0; 0 1 0; 0 0 1], Variable\n",
"size: (3, 10)\n",
"sign: real\n",
"vexity: affine\n",
"id: 114…020\n",
"value: [0.9090761727022365 0.8018998057589781 0.8156148196829522 0.793127657006549 0.7303882764054443 0.9560068400380448 0.8835317223972944 1.108916235520905 1.0571211581975082 0.689856014912743; 8.22145778191618 8.255554482564627 8.245501346939792 8.260128763462827 8.265130415875387 8.132590595705121 8.0724488049853 8.075654836102757 8.147332052116685 8.13525181295928; 3.5356276815757015 3.3466106276014056 3.4074959247868417 3.3267055113350907 3.316543654963961 3.048534318704198 3.386482639601667 3.333416861023218 3.4313193606284322 3.536654868595232], [101.02364745133809 101.00418121769547 … 99.98623365604394 99.96517226835647; 101.05839042646836 100.99314613005879 … 99.95431490867371 100.05551695615708; … ; 101.00457311317994 100.96845799149115 … 100.02296898927675 100.0289186610699; 101.02662908998043 100.9704462880352 … 100.01970963605815 99.92763650541038], Variable\n",
"size: (5, 3)\n",
"sign: real\n",
"vexity: constant\n",
"id: 988…862\n",
"value: [-0.1570462845229664 -11.417119780254275 -1.9790771743623983; -0.021141235108979976 -11.299729854032721 -2.297251144320985; -0.023640789088924436 -11.484763948519115 -1.844308740524739; -0.2620201397676304 -11.303159303534848 -2.2257761979127215; -0.47901684406437545 -11.362428001936705 -2.0262235624224973], Variable\n",
"size: (3, 10)\n",
"sign: real\n",
"vexity: constant\n",
"id: 114…974\n",
"value: [-0.9090761727022365 -0.8018998057589781 -0.8156148196829522 -0.793127657006549 -0.7303882764054443 -0.9560068400380448 -0.8835317223972944 -1.108916235520905 -1.0571211581975082 -0.689856014912743; -8.22145778191618 -8.255554482564627 -8.245501346939792 -8.260128763462827 -8.265130415875387 -8.132590595705121 -8.0724488049853 -8.075654836102757 -8.147332052116685 -8.13525181295928; -3.5356276815757015 -3.3466106276014056 -3.4074959247868417 -3.3267055113350907 -3.316543654963961 -3.048534318704198 -3.386482639601667 -3.333416861023218 -3.4313193606284322 -3.536654868595232], Variable\n",
"size: (5, 5)\n",
"sign: real\n",
"vexity: constant\n",
"id: 169…881\n",
"value: [1.0 0.0 0.0 0.0 0.0; 0.0 1.0 0.0 0.0 0.0; 0.0 0.0 1.0 0.0 0.0; 0.0 0.0 0.0 1.0 0.0; 0.0 0.0 0.0 0.0 1.0], Variable\n",
"size: (10, 10)\n",
"sign: real\n",
"vexity: constant\n",
"id: 624…272\n",
"value: [1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0; 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0; 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0], 0.1)], SequentialConvexRelaxation.Result(10, false, AbstractFloat[1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0], Array{AbstractFloat,1}[[0.47002954935838154], [0.13459561941770679], [0.12111592265281247], [0.11550643127095846], [0.11349153619789144], [0.11283465517897813], [0.11252425804889135], [0.11240977173745709], [0.11236063100153294], [0.11234117682821705]], Any[])))"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res, bilinear_problem = nnmf(data, 3)"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.11234117682821705"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"norm(res.W.value * res.H.value .- data)"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"5×10 Array{Float64,2}:\n",
" -0.0182321 -0.000390748 -0.0182582 … -0.010305 0.0236244 \n",
" -0.0166948 -0.00265244 0.00409316 0.0132875 0.00979968\n",
" 0.0147802 0.000936427 0.00941258 0.0024681 -0.016155 \n",
" 0.0315861 0.00430964 -0.00136609 -0.0180409 -0.022313 \n",
" -0.0114721 -0.00219804 0.00604102 0.0127632 0.00508257"
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"(res.W.value * res.H.value .- data)"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"5×3 Array{Float64,2}:\n",
" 0.157046 11.4171 1.97908\n",
" 0.0211412 11.2997 2.29725\n",
" 0.0236408 11.4848 1.84431\n",
" 0.26202 11.3032 2.22578\n",
" 0.479017 11.3624 2.02622"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.W.value"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"3×10 Array{Float64,2}:\n",
" 0.909076 0.8019 0.815615 0.793128 … 1.10892 1.05712 0.689856\n",
" 8.22146 8.25555 8.2455 8.26013 8.07565 8.14733 8.13525 \n",
" 3.53563 3.34661 3.4075 3.32671 3.33342 3.43132 3.53665 "
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.H.value"
]
},
{
"cell_type": "code",
"execution_count": 39,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.3"
]
},
"execution_count": 39,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"predicted_cluster_assignments = vec(first.(Tuple.(findmax(res.H.value, dims=1)[end])))\n",
"mean(predicted_cluster_assignments .== true_cluster_assignment)"
]
},
{
"cell_type": "code",
"execution_count": 40,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3df3hU9Z3o8e8kmckPQgSDGQihChjEAg1GCgEWBWqhStckvd2UEMCio1KhEJaYFUN6qZRSSLhK3Qa8YCksWIZaulK7NuYp9nkSpVVWHbqYln2kIINhMRMw45DA5Jzv/ePU3DRMkjkzJ5zMzPv1nD9mvufHfDw99JPvz2ORUgoAAGJVnNkBAABgJhIhACCmkQgBADGNRAgAiGkkQgBATCMRAgBiWoLZAQAA0Ju33377hRde8Hg86enp3/nOd6ZMmaKVe73eLVu2nDx5cuLEieXl5YMHD76+JJjrUyMEAAxciqL88Ic/XLFixeHDhx966KHq6urOXU6n0263O53OjIyMQ4cOBSwJBokQADBwKYqybt26u+66q7293Wq1Dho0qHNXQ0NDfn6+zWbLz8+vr68PWBIMmkYBAAOXzWabNm1aW1tbfn6+EOK5557r3OXxeOx2uxDCbre3tLQELAnGwE2EH3/8cVpaWmpqqlEX7OjoSEgYuP+9A5+qqkKIuDhaEULHQxgmRVHi4uIsFovZgUQqKaWqqvHx8f35I61SKLpOOHv23NYtO7sVPv3001lZWZ1fk5OTjxw58qtf/aqmpuYnP/mJViil1B4G7b8rYEkwBu6/yWefffaBBx6YM2eOURe8fPnysGHDjLpaDGpvbxdCpKSkmB1IBOMhDJPP50tKSrLZbGYHEqlUVW1tbR06dGg//oRskOKSrlMuXPjw4sWL5eXlXQtvvvlm7UNTU9Orr7766KOPJicn33///T//+c87j0lPT7948WJWVlZzc7P2L+v6kmDw1z0AwDCqVFTZoXNTMjIypv69zr+509PTf/Ob3/zpT3+SUv7+97+//fbbhRAul0sIkZeXV1tbK6Wsra2dMWNGwJJgkAgBAIaRokNKnVuvTak2m+373/9+TU3NN77xjTfeeGPNmjVCiLKyMiHEkiVLTp8+XVxcfObMmZKSkoAlwRi4TaMAAAghcnJyduzY0bWkrq5OCJGamrpp06au5deXBINECAAwjFQVKTv0ntJPwQSJRAgAMIyUuhOhkCRCAMANpyhKQ0PD+fPnc3JyJkyYYHY4ZiIRAkDMuXDhwj333OvxeD77zHfTTYNnzJj5q18dNmaCpuzQ3TRqdo2QUaMAEHMWLSr58MPTLS2Xr13zf/JJy9Gjb3QbjRIyrWlU56hRnU2pRiMRAkBskVI2Njaqquws8Xo/+8UvfmHQxRWpdujc6CMEANxoUoiuK5BZLJbYrReRCAEgtlgslvj4BCEsQlg+L5HjxmUbcvFQR42auX4siRAAYo6i/F3ukdLy17/+1ZAryxAGywjF3GQUu3VhAIhZUspuJZ988okpkQwEJEIAiDnDh9u7frVYxKeffvrZZ58ZcGmpCLVD58b0CQDAjbVv376hQ28SQmqblOq5c+cfeughAy6txmQiVBRl2bJl2mev17t+/frCwsLKykqv19ttb0/HAABupC996Uvz58/7PBEKIeL8/o633nor+JfZ9kgqQnbo2yJ9HuHhw4dXr17tdru1r06n0263O53OjIyMQ4cOddsb8JgwAwAABOPKlSvf+c53vvjFLw4fPmL69Onnz58XwiJEXOeoGSnllStXzA3SFOEO1BkzZkxmZmZlZaX2taGh4ZlnnrHZbPn5+d/73vdKS0u77g14zCOPPBLwylLKTz755Ny5c50lqampaWlpIYeqKIqimFwBj2ja3eMehoOHMEzK58wOJCItXLiwtrbu2jW/EOJ//ufiTTcNTkpKbG+/1nlAcnJKcnJyL7c3Li6u72XYZIdQddbwzG4aDTcRTp48uetXj8djt9uFEHa7vaWlpdvegMf0dGWfz1ddXZ2ent5Z8rWvfW3p0qUhh+r1em02W8ino62tTQjR0WFyI0ZE4yEMk9frvXbtGvcwBF6v99ixP2hZUAghhMXna7vllvRPP229cqXNYrHcfPOQ7dufa21t7eUiaWlp8fHxffyS1keoS6Qnwm6klNrfC1LKntqagzlGCJGamrply5Y5c+YYFZuiKEOHDjXqajEoMTFRCJGSkmJ2IBGMhzBM8fHxSUlJNywRLlq06NNPP70xv9Xf2traEhNtdvst2tdPP21tb29PS0tLT0+32WyzZ89euXLl6NGjzQ3SLAYnwvT09IsXL2ZlZTU3Nw8bNizkYwDAdH/+859fe+01s6Mw3nvvvbdkydKrV9v+8pf/FkKkpKSoqtyyZYsxV1c7hOrv+7Cu9K5EYzSDE2FeXl5tbe3DDz9cW1s7Y8aMbntdLldOTk7vxwDAABEXF6f140SZm2++WQgp5d8GS1650nbq1Kn9+/d/+9vfNuDq2qhRvaeYyuB5hEuWLDl9+nRxcfGZM2dKSkq67S0rK+vzGABAf+vW1XflStsbb7xhzKWl3kmEHaYnQmNqhHV1ddqH1NTUTZs29bRX+xDwGADADWO1Wrt+tVgsMdtBKFh0GwBi0JAhQy5cuNjR8beq2LBhNy9ZssSYSzNqFAAw8OXl5Y0de/uxY29ZLJaUlJQXXnhh7NixhlzZIlWLzqZOvccbjkQIADHHarX++7//qqOj48qVK+EsVBIdSIQAEKMSEhKMz4JS0d3UKcNe4DQ8JEIAgGEsqmrRmwjDX+k7PCRCAOh3b7755pEjR9xu9xe+8IWCgoJp06aZHVG/UfXXCBksAwBR7OOPP/7617/+0UcfeTyXtJJdu3aPHTv2178+kpGRYW5s0PBiXgDoL21tbffcc8/775/weC4LYdE2j+fS8eP/ee+99/r9OpciC+T6d76aTRVS0bcJk5tGSYQA0F+effbZ8+ebpOxerqryo4/cO3bsCPP6Ad/5ai6tj1DnRiIEgCi1f//+9varAXddudL2s5/9TNfVnn322V/+8pedn51O55gxYxYvXhxmkCARAkB/6f31fs3NzbquNmvWrIaGBiGE3+9vaGiYM2fO5MmT8/LywgrRcNpgGb2bqUiEANBfen+fe1ycvv8Hnjx58kcffXTp0qV33333tttuG5hjbSxSf9Mo8wgBIFplZ2e73R8LETAdypycHF1XS0hImD59+ltvvXXixIn77rvPkAiNF4HTJ6gRAkB/Wb9+/dChQwLuuvnmoevWrdN7wVmzZh09evSdd96ZNWtW2NHhb0iEANBf5s6dW1T0TzfdNLhb+ZAhg5ct+3YI3Xt33XXXhx9+mJOTk5qaalCMRlOlRVX1bdcPq72xaBoFgH60c+fOadOmff/7329vb1dV1WKJGzQo5Qc/+MGiRYtCuJrNZrv11lu7tYt2vvN1QAhhrdG+pk/U19fv3bu3ubl59OjRa9euzcrK0sq/+tWvdjvyzjvvbGxs1D4vWLCgtLQ0mN8nEQJA/1q2bNmyZcuampqamppGjhxpt9tDu05HR8eZM2c+/vjjL3/5y8ZGOJA1NTVVVVVt3bp19OjRr7zySlVV1fbt27VdR44c6TzM6XT6/f7XXnvN6XQmJycLIeLj44P8CZpGAeBGGDFiRG5ubshZUAhx7Nixp5566rvf/a7NZjMwMIOpqu65E72+j7CpqWnu3Lnjx49PTEycN2/euXPnOnclf+7ChQsffPBBfn6+oigVFRVFRUWbN2/2+XxBhkyNEAAiw6xZswb+GBmL0U2jubm5ubm5QghFUfbu3Tt79uxuB/j9/m3btj355JOXL1/Ozs5evnx5RkbGjh07ampqKioqgvl9EiEAwDiq1P1aJVWtr69fsWJF17Knn3565MiRnV+PHz++e/fuKVOmXL+w6ssvvzx+/Phbb71VCFFdXa0VOhwOh8MR5O+TCAEAJrvjjjuWLl3atWTo0KHaBynlrl27Ghsb169f3zlMppOiKK+++urWrVuFEKdOnfL7/RMmTBBCWK1Wq9Ua5K+TCAEAhrHon1BvkWpGRkZP72g8ceLEsWPHnn/++fj4+La2NiFEcnKyy+XSliN4//33b7nlFq3u2N7evnHjxurq6szMzP3798+cOTPIAEiEAADjSNXYlWVcLpfb7S4sLOwsqaurKysr0yaNvP7663fffbdWPmnSpJKSksrKSp/PN3Xq1JUrVwb5+yRCAMDAtXTp0m6tpqLL1Mmuq/NYLJaCgoKCggK9P0EiBAAYJ4S1Rll0GwAGpqtXr/7oRz8yOwqDSSnPnz/ff9e36B81alFZYg0ABqSDBw92rtcVHVRVbWtry87OHjt2rNmxDCAkQgAIbMKECdpY/KihKEpra2vnzIR+IVXd8whpGgUARA81hFGjJEIAQNQIoUZodiJk0W0AQEyjRggAMI5KHyEAIIZF4vQJmkYBADGNGiEAwDghvI+w1xfz3gAkQgCAcUJ5HyErywAAokYEDpahjxAAENOoEQIAjBPKhHqaRgEAUUOVuhOb2YmQplEAQEyjRggAMI6U+ge/sLIMACBq0EcIAIhp9BECABBZqBECAIwTwsoyZk+oJxECAIwTgX2ENI0CAGIaNUIAgHFCGCwjGTUKAIgakrdPAABimap/frzJeZA+QgBAbKNGCAAwTgg1QpNnT5AIAQDGCWGpUbOnEdI0CgCIbdQIAQDGicDBMiRCAIBx6CMEAMQ0GXk1QvoIAQAxjRohAMA4Un8Nz+waIYkQAGAYqVqkatF5Sh+ZsL6+fu/evc3NzaNHj167dm1WVlbnrlWrVjU2NmqfFyxYUFpa6vV6t2zZcvLkyYkTJ5aXlw8ePLjPAGgaBQAMXE1NTVVVVWVlZU6nc/r06VVVVZ27pJRut9vpdB45cuTIkSNPPPGEEMLpdNrtdqfTmZGRcejQoWB+gkQIADCO/HzgaPBbrxXCpqamuXPnjh8/PjExcd68eefOnevc5fF4FEWpqKgoKiravHmzz+cTQjQ0NOTn59tstvz8/Pr6+mBCJhECAIyjCqFadG89y83NLS0tFUIoirJ3797Zs2d37mppacnOzl67du2BAwcGDRpUU1MjhPB4PHa7XQhht9tbWlqCCZk+QgCAYaS0SKmzj1BaGhoavvvd73YtfOqpp0aOHNn59fjx47t3754yZcqyZcs6C8eNG1ddXa19djgcDodDCCGltFgs2gc1uBdCkQgBACbLzs5etGhR15KhQ4dqH6SUu3btamxsXL9+fddhMkKIU6dO+f3+CRMmCCGsVqvVahVCpKenX7x4MSsrq7m5ediwYcH8OokQAGAc2UdTZ8BT7Hb79OnTA+48ceLEsWPHnn/++fj4+La2NiFEcnKyy+XKyclpb2/fuHFjdXV1Zmbm/v37Z86cKYTIy8urra19+OGHa2trZ8yYEczvkwgBAMbpq88v0Cm97XS5XG63u7CwsLOkrq6urKysrq5u0qRJJSUllZWVPp9v6tSpK1euFEIsWbJk8+bNxcXF2dnZTz31VDC/TyIEAAxcS5cuXbp0abfCuro6IYTFYikoKCgoKOi6KzU1ddOmTbp+gkQIADBMCINlhN7jjRZ6IlQUxeFw7NmzRwhx/Uz+Xub2X78QQJj/DQCAgUKbPqGL2YkwxHmEhw8fXr16tdvt1r5eP5O/p7n9ARcCAADALCEmwjFjxixevLjz6/Uz+Xua2x9wIQAAQJSQcULVv5kqxKbRyZMnd/16/Uz+nub2awsBLF++PCMjY8eOHTU1NRUVFQF/oq2t7ac//enRo0c7S+6+++777rsvtICFED6fLykpKeTTceXKFSFEkBNUERAPYZh8Pl9HR4fNZjM7kEilKIrP59Pm24UgJSUlLq6vpBXKhPrQwjGMMYNlrp/J39Pc/oALAQSOLCFh7Nix2kxJTVZWVkJC6AEnJCSEczq0u8c9DAcPYZji4+O5h+GwWCz9fgNDmD5hdh+hMbfj+pn8Pc3tD7gQQEBWq/Xee++dM2eOIREKIRITE/ljPBzaHzTcw3DwEIbp2rVrSUlJ1AhDpiiKdg/NDmRgMaZlVpvJL6XsnMl/fYnL5RJCtLe3b9iw4ezZs36/v3MhAABAdJDyb68k1LWZG7MxiXDJkiWnT58uLi4+c+ZMSUlJwJKysjIhROdCAAsXLvR6vb00jQIAIo+M072Z/R6ksJpGtbn9ItBM/utLelkIAAAAs9DnDAAwTAhNndLsoegkQgCAYaTU3+dndh8hiRAAYBzVonuCvDS5j9DknwcAwFzUCAEAhgnh7RNRsrIMAABChPSGerP7CGkaBQDENGqEAADjqHFS52AZafZgGRIhAMAwUgr9fYRMnwAARI1Q+gj7J5Kg0UcIAIhp1AgBAIaR9BECAGJZCPMITX8xL02jAICYRo0QAGCYUN4+YXaNkEQIADCQRW+fH4kQABBFVJZYAwAgolAjBAAYhj5CAEBMkyKE1zDRNAoAgHmoEQIADCNVi96VZUwfLEMiBAAYJqQ31JMIAQDRIhIHy9BHCACIadQIAQAGCqFptJ8iCRaJEABgmFAGy/S1JFt9ff3evXubm5tHjx69du3arKys3netWrWqsbFRO2DBggWlpaW9X59ECAAYuJqamqqqqrZu3Tp69OhXXnmlqqpq+/btveySUrrdbqfTmZycLISIj4/v8yfoIwQAGEYbNap36+WCTU1Nc+fOHT9+fGJi4rx5886dO9f7Lo/HoyhKRUVFUVHR5s2bfT5fnzFTIwQAGEZK/aNGez0+Nzc3NzdXCKEoyt69e2fPnt37rpaWluzs7OXLl2dkZOzYsaOmpqaioqL3AEiEAADDhDSPUDQ0NKxevbprYXl5+ciRIzu/Hj9+fPfu3VOmTFm2bFm307vtGjduXHV1tbbL4XA4HI4+AyARAgBMdvvttxcVFXUtGTJkiPZBSrlr167Gxsb169d3HSbT065Tp075/f4JEyYIIaxWq9Vq7fPXSYQAAOPoHzUqZdzw4cNnzpwZcO+JEyeOHTv2/PPPx8fHt7W1CSGSk5NdLldOTk7AXe3t7Rs3bqyurs7MzNy/f39Pl+2KRAgAMEwITaOi1+NdLpfb7S4sLOwsqaurKysrq6urC7hr0qRJJSUllZWVPp9v6tSpK1eu7PP3SYQAgIFr6dKlS5cu7VZYV1fX0y6LxVJQUFBQUBD8T5AIAQCGicT3EZIIAQCGCW3UqLmYUA8AiGnUCAEAhonE1zCRCAEAhpHSova1iPZ1p5jcNkkiBAAYJhLfUE8fIQAgplEjBAAYJpRFt80eNUoiBAAYhqZRAAAiDDVCAIBhIrFGSCIEABhH6k9sJEIAQNQIZa3RfgolaPQRAgBiGjVCAIBhVGlRddYI9R5vOBIhAMAwhr+Y9wagaRQAENOoEQIADMP0CQBAbNM/fYJECACIHqHUCPsplKDRRwgAiGnUCAEAhgllQj1NowCAqBHCPELTEyFNowCAmEaNEABgHKZPAABiGfMIAQAxTZUWVerrdJOCRAhEEVVVDxz4t9dff9VqtT744D8VFBSaHRGAPpAIASM99FDxiEz/uooZfr+y/bkdf/zjm5s3V5sdFHDjMH0CiGnvvPNOh3LxR1se1b7u2v3QV+ZUnz9/fuTIkeYGBtwwMgJfw8T0CcAw77//3uzZt3d+tVgss+4Zd+LECRNDAtAnEiFgmMzMkW53a9cSt/tyZmamWfEAN542alTnZnLMJELAMPfee+9vX/vA5fpI+/r7Nz44/eGnEydONDcq4EZS5d8Wlwl+Y9QoED1SU1P37395zZonPv30E0WRo0aNfumlX8bHx5sdF3DjMI8QiHV33HHHf/zH765duxYXF5eQwL8vIALwDxUwns1mMzsEwCQRWCMMq49QUZRly5Zpn71e7/r16wsLCysrK71eb8CSno4EAEQHVejrIFSDSJz19fUOh6OgoGDNmjVut7uzPPi807vQE+Hhw4dXr17dGZPT6bTb7U6nMyMj49ChQwFLejoSAICAmpqaqqqqysrKnE7n9OnTq6qqOncFn3d6F3oiHDNmzOLFizu/NjQ05Ofn22y2/Pz8+vr6gCU9HQkAiA5SWPRvvWlqapo7d+748eMTExPnzZt37ty5zl3B553ehd5HOHny5K5fPR6P3W4XQtjt9paWloAlPR0ZUGtr66OPPpqamtpZUlhYuGLFipAD7uW3EIy2tjYhxJUrV8wOJILxEIbJ6/UmJibSBRsyRVG8Xq+iKKGdPmTIkD6HgBm+skxubm5ubq4QQlGUvXv3zp49u3NX8Hmnd4YNlpFSWiwW7YOqqgFLejoyoLS0tF27ds2ZM8eoCIUQw4YNM/BqsUZLgSkpKWYHEtl4CMNhs9mSkpJIhCFTFMVmsw0dOrT/fkLKEAa/WN588801a9Z0LXryySe7LkZx/Pjx3bt3T5kypXNgitCTd3pnWCJMT0+/ePFiVlZWc3Oz9k/9+pKejgQAxLIxY8Z84xvf6Fpy0003aR+klLt27WpsbFy/fn1WVlbXY4LPO70zbGWZvLy82tpaKWVtbe2MGTMClrhcroDlAIDooMpQthEjRsz6e4MGDdIueOLEiWPHjm3cuDE9Pb2trU3ro+kpm4SWXwxLhEuWLDl9+nRxcfGZM2dKSkoClpSVlQUsB6JVa2vrv/7rj9esWbFzZw3dq4gFIa012ltTqsvlcrvdhYWFD35O9JxNQssvFmn6cqc9ePLJJx944AED+whpiQ0TfYR6nTt3Lj9/vuPRGRMnjXz3P8/u/7f/PHDg5TvuuMPsuCJYa2srfYThUBSltbW1X/sI31xU7v3vs7pO+Yuv5S+zs2tqavoppD6xsgzQX556as3254tmzhwnhPiHf7jj9tszqqo27d69z+y4gH4khe5Ro5G9sgyAXpw69WctC2q+Om/Sn/7kMjEe4AZQhUX/ZjISIdBfEhOTrl71d371etu6zosFMECQCIH+UlDwzY3P/FrrhldV+b+/98o//uM3+jwLiGiGD5a5AUiEQH/5538uv3Z1ZN7UTUsW/3Tal38wLD1n8eKHzA4K6F+SF/MC6BQXF1ddvf3KlSsfffTRbbfdlpSU1NzcbHZQQP/ixbwAuktJSRk/frzZUQDoEYkQAGAYVQi9o0BNHzVKIgQAGCaUwS9mN40yWAYAENOoEQIADBPCyjKq2Qt9kggBAIbRXjqv95R+CiZIJEIAgGFU2ccb569n+vQJ+ggBADGNGiEAwDBS6u7zY/oEACB6RGIfIU2jAICYRo0QAGCYkNYa7adYgkUiBAAYJoQl1szOgyRCAICBIvDtE/QRAgBiGjVCAIBhaBoFAMS0SHwxL02jAICYRo0QAGAYNYS3T5g9oZ5ECAAwkt4+P/oIAQDRg7dPAAAQYagRAgAMI2kaBQDEMin1D5YxOxPSNAoAiGnUCAEAhgmhadR0JEIAgGHUUJpGmUcIAIgWqtTd52f6+wjpIwQAxDRqhAAAA1mkziXT9B5vOBIhAMAwqv7pEH2+tklRFIfDsWfPns6Sr371q92OufPOOxsbG7XPCxYsKC0tDT4AEiEAYOA6fPjw0aNH3W5318IjR450fnY6nX6//7XXXnM6ncnJyUKI+Ph4XT9BHyEAwDCqsOjdem8aHTNmzOLFi7sVJn/uwoULH3zwQX5+vqIoFRUVRUVFmzdv9vl8umImEQIAjCOF1L/1YvLkyXl5eQF3+f3+bdu2rVix4vLly9nZ2WvXrj1w4MCgQYNqamp0hUzTKADAMCG8fUKVlrfeemvt2rVdC9euXZuZmdn7iS+//PL48eNvvfVWIUR1dbVW6HA4HA6HrgBIhAAAk912220PPvhg15Kbbrqp91MURXn11Ve3bt0qhDh16pTf758wYYIQwmq1Wq1WXb9OIgQAGCa0t09kZmbee++9wRzscrlycnKEEO+///4tt9wycuRIIUR7e/vGjRurq6szMzP3798/c+ZMXQHQRwgAMIy2soyuTdfKMmVlZdqH119//e6779Y+T5o0qaSkpLKycuHChV6vl6ZRAEC0qaur6/Zh3bp1nXstFktBQUFBQUFoFycRAgAMI4VF1blSjGp24ySJEABgmD6nQwQ8xVz0EQIAYho1QgCAYdQg1g7txuwKIYkQAGAcGcKLeXn7BAAgasgQ+vzMrhLSRwgAiGnUCAEAhgmhj1Dv8YYjEQIAjMP0CQAAIgs1QgCAYVT9K8v0/mLeG4BECAAwjLaOtt5TzEUiBAAYJrTXMJmLPkIAQEyjRggAMIzU3zRq+qhREiEAwDCh9BH2TyTBo2kUABDTqBECAAwjhUXvdAimTwAAoofUPx2CPkIAQPRg+gQAABGGGiEAwDCsLAMAiGmROI+QplEAQEyjRggAMEwkDpYhEQIADBOJK8sYkAgvXbq0c+fO9957Ly4ubtq0aY8//nhKSkrAws5TVq1a1djYqH1esGBBaWlp+GEAAEwn9Se2aKgRbtu2bezYsQcOHJBS7t27d9++fcuXLw9YqB0vpXS73U6nMzk5WQgRHx8ffgwAAITGgETocrnWrVtntVqFEAsXLly+fPny5csDFmrHezweRVEqKircbvddd921Zs0am812/WU7Ojrefvvtq1evdpZkZWVlZ2eHHOe1a9e6Xg16aXePP1zCwUMYpqtXr1osFmn6KMOIpSjK1atXQ34IbTabxdLHcmhS6h4Favr/ngYkwnHjxh08ePBb3/qW3+9/6aWXWlpaeirUtLS0ZGdnL1++PCMjY8eOHTU1NRUVFddf9tq1a++8886FCxc6S/Ly8m699daQ47x27dq1a9dCPh1+v18IwT0MBw9hmPx+f1wcY91DpyiK3+8P+SFMSEjo80/hGJ1HWF5e/uMf/7i4uHjIkCGFhYWDBw/uqVAzbty46upq7bPD4XA4HAEvm5KSsmLFijlz5oQfoebq1atdw4Be2j+Arn290IuHMExSyqSkpIBtSAiGoihSSh7CbgxIhElJSRs2bNBaQV0u16hRo3oq1Jw6dcrv90+YMEEIYbVatWMAAFEghOkTpjOgkWHXrl3PPfecz+fzeDwvvvhiQUFBT4Uul0sI0d7evmHDhrNnz/r9/uZvZgwAAA+/SURBVP3798+cOTP8GAAAA4H8vHVU12YuAxLhY4895vV6Fy1aVF5ePn/+/FmzZvVUWFZWJoSYNGlSSUlJZWXlwoULvV5vT02jAICII/82XEbfZm7MBjSNpqWlPfPMM8EU1tXVCSEsFktBQYFWRwQAwFysLAMAMEwITZ3RMH0CAABNf6w1qiiKw+HYs2dP18LrVyjzer1btmw5efLkxIkTy8vLgx8cSyIEAAxchw8fPnr0qNvt7loYcIUyp9Npt9u/973vvfDCC4cOHXrkkUeC/AmmpgIADCP1DxntfW3SMWPGLF68uFth5wplRUVFmzdv9vl8QoiGhob8/HybzZafn19fXx98zCRCAIBhpP6td5MnT87Ly+tWqK1Qtnbt2gMHDgwaNKimpkYI4fF47Ha7EMJut3ddzqxPNI0CAAwT2hJrx44de/LJJ7sWrlmzJjMzs6dTAq5QJqXUlkKVUqqqjndgkAgBACa79dZbH3jgga4laWlpvRwfcIWy9PT0ixcvZmVlNTc3Dxs2LPhfJxECAIwT0tsnMjMzg1xZ2uVy5eTktLe3b9y4sbq6OjMzs3OFsry8vNra2ocffri2tnbGjBnBB0AfIQDAMKqQejddK8v0skLZkiVLTp8+XVxcfObMmZKSkuCvSY0QADDQaQuTiV5XKEtNTd20aVMIFycRAgAMI/WvFMPKMgCA6KG1duo9pZ+CCRJ9hACAmEaNEABgGCn6WCkm4CnmIhECAAwjpZQ6O/3oIwQARA9Vf41Q7/GGo48QABDTqBECAAwjQ3gxb/9EEjwSIQDAMFIIXSvFCP3HG46mUQBATKNGCAAwTGivYTIXiRAAYBipf2UZ05tGSYQAAMNInW+TEAMgEdJHCACIadQIAQCGicRFt0mEAADDhNRHaDKaRgEAMY0aIQDAMKoUqs5VtPUebzgSIQDAMCGMGjW9cZRECAAwTCQOlqGPEAAQ06gRAgAME8KoUdNrhCRCAIBhQmgaZfoEAABmokYIADCMFFIKVecp+o43HIkQAGAY+ggBADFNtaiqRWeN0ML0CQAAzEONEABgGCmkqrPPT+/xhiMRAgAMowpVb2LjxbwAAJiJGiEAwDBSqEyfAADELtUiVJ2jQPUebzgSIQDAMJI+QgAAIgs1QgCAYUIYNcr0CQBA9OiPtUYVRXE4HHv27OlaWF9fv3fv3ubm5tGjR69duzYrK0sIsWrVqsbGRu2ABQsWlJaWBhMAiRAAMHAdPnz46NGjbre7a2FTU1NVVdXWrVtHjx79yiuvVFVVbd++XUrpdrudTmdycrIQIj4+PsifoI8QAGAYKRRV59Z7jXDMmDGLFy/uVtjU1DR37tzx48cnJibOmzfv3LlzQgiPx6MoSkVFRVFR0ebNm30+X5AxUyMEABhGDWWJtd5GjU6ePPn6wtzc3NzcXCGEoih79+6dPXu2EKKlpSU7O3v58uUZGRk7duyoqampqKgIJgASIQDAMFIoUih6T/nDH/7wL//yL10LS0tLR4wY0fuJx48f371795QpU5YtWyaEGDduXHV1tbbL4XA4HI4gAyARAgBMNmrUqHnz5nUtGTx4cC/HSyl37drV2Ni4fv16bZiMEOLUqVN+v3/ChAlCCKvVarVag/x1EiEAwDChTagfOXLkV77ylWAOdrlcOTk5J06cOHbs2PPPPx8fH9/W1iaESE5Obm9v37hxY3V1dWZm5v79+2fOnBlkACRCAIBhQuoj1HF8WVlZXV2dy+Vyu92FhYWd5XV1dZMmTSopKamsrPT5fFOnTl25cmWQ1yQRAgAGurq6uq4fli5dunTp0m7HWCyWgoKCgoICvRcnEQIADBPSYBlWlgEARIsQ3lBveiJkQj0AIKZRIwQAGIamUQBATOuPRbf7G4kQAGAYVaqq1JkIJS/mBQDAPNQIAQCGkUKlaRQAEMOklLqbRkmEN1xra2ttbW2LxzN9xowvfelLZocDADBTzCXC4++88+1vfmuGbehgGf9zdXvOffdu/787zQ4KAKJESE2jJg+WiblE+Njih55OvzPDliKEuF+IZ39X/3pt7bz5882OCwCigRSq3nmEQvfxBoutUaMXLlxI7ZBaFtTcaxtW9+vfmBgSAEQTVajaDAo9G9MnbqDExMSOv++V9UslKSXZrHgAAKaLrUQ4dOhQS1rqX698qn2VQr7e/sn9BfnmRgUAUUSV+jdzIzamj/DSpUs7d+5877334uLipk2b9vjjj6ekpAQs1I73er1btmw5efLkxIkTy8vLBw8ebEgYwfjZLw4Wff3BMVcvpon499svlzzmmDFjxg37dQCIbjKE6RNmD5Yxpka4bdu24cOHHzhwYN++fampqfv27eupUON0Ou12u9PpzMjIOHTokCExBCk7O/udk38q3/dC0XM/+O3bb5VVrLuRvw4AGGiMSYQul6uoqMhqtdpstoULF9bX1/dUqGloaMjPz7fZbPn5+V3Lb4yEhIRp06bdf//9drv9Bv80AEQ37e0TujbTR40a0zQ6bty4gwcPfutb3/L7/S+99FJLS0tPhRqPx6MlIbvd3rW8K5/Pt2HDhp/85CedJV/5ylcWLlwYcpCXL1+Oj48P+XS0tbUJIa5evWp2IBGMhzBMXq83MTHRZrOZHUikUhTF6/WGfHpaWlqfD3AoTaNmjxo1JhGWl5f/+Mc/Li4uHjJkSGFhodbnF7BQI6W0WCzaB1UNfMuSk5MffvjhqVOndpYMHTo0LS0t5CCvXbsWzulISEgQQnR29CIEPIThS0pKIhGGTFEUIUTID2FcXN+NiLG71mhSUtKGDRusVqsQwuVyjRo1qqdCTXp6+sWLF7Oyspqbm4cNGxbwmnFxcV/4whfuvPNOQyIUQsTHx/PHeDi0u8c9DAcPYZjiP2d2IBGMG3g9Y/oId+3a9dxzz/l8Po/H8+KLLxYUFPRU6HK5hBB5eXm1tbVSytraWgZtAkD0kKrUvUXFqNHHHnvM6/UuWrSovLx8/vz5s2bN6qmwrKxMCLFkyZLTp08XFxefOXOmpKTEkBgAAKbTlliLxcEyaWlpzzzzTDCFdXV1QojU1NRNmzYZ8tMAAIQj5hbdBgD0H621U+cpUTFqFAAAIYQUMkZHjZroww8//Pjjj7/4xS+mp6ebHQsAxDopFSn19fnxhvrQXbp0qejrDyoXPPaEpMa2S/mLir//ox/quoLP59tete3thjczhtuX/3Npbm5uP4UKABiwIvjtE098+5F7Lomnb5nwyNCxVZl3v3PolX8/fDj4030+3z1TprYc/I9veuInvOd+LP+buk4HAFxPaxqNrLdPRGoilFKefN81PW249tUiLAWpIw8fOBj8FXY+/6+zlEFfGzJqeOKgCYPTK4d/6Xtl/9I/wQJA7JBC6tyi4+0TN56iKN1CT4qL9332WfBXOP7WsUnJQzu/Doq3JktLa2urQQECACJDpCbChISEofaMM23/P2/VX/lkzv3zg7/CqNtu++RaW9eSK4r/Rr4ZEQCikap/Y/pEqHb+28/+1/wHZvmG3hJnPan6rmQN27FyRfCnL33MsfRXD96eMmSINVEK8cvm07Pnz9OWAgcAhEj/2yeE2X2EEZwI77jjjmN/ev/w4cNNH517bOqX583Tl8YmTpy4eVdNRela0X71itIx/x8X/OjZ/9N/0QJAbFD1JzZqhGEYPHjwQw89FPLp87/2tfl//lpbW1tycrKBUQEAIkhkJ0JDkAUBwCghTYdg+sSN8otf/MLsECLbyZMn/+u//svsKCIbD2GY3nzzTbfbbXYUEay1tfW1117r39+QUkhV59ZH06iiKMuWLeta4vV6169fX1hYWFlZ6fV6A5YEL4YS4QsvvGB2CJHt7bff/uMf/2h2FJGNhzBMv/3tbz/88EOzo4hgly5dOnhQx3zrgeDw4cOrV6/u9geQ0+m02+1OpzMjI+PQoUMBS4IXQ4kQANDfQlhWpvem1DFjxixevLhbYUNDQ35+vs1my8/Pr6+vD1gSPPoIAQAGMnjU6OTJk68v9Hg8drtdCGG321taWgKWBG/gJsLk5OQ9e/bU1tYadUGv1/vUU08ZdbUYdPLkSSnl+fPnzQ4kgvEQhundd9/1eDy/+93vzA4kUnm93rNnz4b8EK5evXrEiBG9H/P+++/qvezbb7+9fPnyblH1/ltSSm2+nJRSVdWAJcEbuIlw3bp1b775poEXvO+++wy8WgziBoaPexgmbmD4CgsLQz63n9beGj9+/NatW3X9Vnp6+sWLF7Oyspqbm4cNGxawJHgDNxEmJyfz0ANA1EtLSwv+/+1dLldOTk5eXl5tbe3DDz9cW1s7Y8YMIcT1JcFjsAwAIGKUlZUJIZYsWXL69Oni4uIzZ86UlJQELAmeRfY1gQMAgChGjRAAENMGbh+hXpcuXdq5c+d7770XFxc3bdq0xx9/PCUlJWChdrzX692yZcvJkycnTpxYXl7OC5j03kAhxKpVqxobG7XPCxYsKC0tNSn2AUFRFIfDsWfPHhHo6erpeeM51IR29wQP4XV6v5Nd92p4AkU01Qi3bds2fPjwAwcO7Nu3LzU1dd++fT0VasJZhiAq6b2BUkq32+10Oo8cOXLkyJEnnnjCvNjN123xi+CXveA5FGHcPR7Cbnq/k0Eu0RKDoicRulyuoqIiq9Vqs9kWLlyorSwQsFATzjIEUUnvDfR4PIqiVFRUFBUVbd682efzmRe7+botfhH8shc8hyKMu8dD2E3vdzLIJVpiUPQkwnHjxh08ePCzzz67dOnSvn37tJUFAhZqwlmGICrpvYEtLS3Z2dlr1649cODAoEGDampqzIvdfJMnT87Ly+v8GvyyFzyHIoy7x0PYTe93stvegMfcyGgHjuhJhOXl5drY2VWrVo0YMUJr6Q5YqAlnGYKopPcGjhs3rrq6+vbbb09LS3M4HMePHzcv9gEn+GUveA6vF/zd4yHsXTBPF0+giKbBMklJSRs2bLBarUIIl8s1atSongo14SxDEJX03sBTp075/f4JEyYIIaxWq3YMNMEve8FzeL3g7x4PYe+Cebp4AkU01Qh37dr13HPP+Xw+j8fz4osvFhQU9FTocrnE58sQSClDWIYgKum9ge3t7Rs2bDh79qzf79+/f//MmTNN/g8YSK5/uq4v4TnsSfB3j4ewd70/XTyBnaInET722GNer3fRokXl5eXz58+fNWtWT4U9LUwQ4/TewEmTJpWUlFRWVi5cuNDr9TocDpP/AwaSYJa94DnsSfB3j4ewd70/XTyBnVhZBgAQ06KnRggAQAhIhACAmEYiBADENBIhACCmkQgBADGNRAgAiGkkQgBATCMRAgBiGokQABDTSIQAgJhGIgQAxDQSIQAgppEIAQAxjUQIAIhp/w9IIKsz50FcYgAAAABJRU5ErkJggg=="
},
"execution_count": 40,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"scatter(data[1, :], data[2, :]; zcolor = true_cluster_assignment)"
]
},
{
"cell_type": "code",
"execution_count": 41,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3dfVxUZf438GtmmBkHRuRBGUV8gA2XVl2ITFHXVq3sl9QCtZGI0qKTWbaKQWQirWnmKrhWbqg3toY/LceMzKwWKWtvSHaLHkYjytJ8mBVFBtBhZGDmnOv+42zcLAzjnMOBw8x83q957WvmOg/z9ezk1+tZRiklAAAAvkoudQAAAABSQiIEAACfhkQIAAA+DYkQAAB8GhIhAAD4NCRCAADwaX5SBwAAAODKZ599tnPnTrPZHBoa+thjj02aNIkrt1gsmzZtqqmpmTBhQm5u7uDBg7uXuHN/1AgBAGDgYhjmhRdeWLZsWWlp6cMPP1xYWNhxyGAw6HQ6g8EQFhZ24MABpyXuQCIEAICBi2GYZ5555pZbbrHZbEqlMiAgoONQZWVlUlKSSqVKSkqqqKhwWuIONI0CAMDApVKppkyZ0trampSURAh58cUXOw6ZzWadTkcI0el0jY2NTkvcMXAT4cWLFwMDA7VarVg3dDgcfn4D9887wFFKWZZVKBRSB+KpWJYlhMjlaIMRiGEYuVwuk8mkDsRT9eNfgNcoYXhdcO7chc2bdnQpXL16dURERMdHjUZz+PDht99+u6io6JVXXuEKKaXcT4L7C8ppiTsGbmLYunXr3LlzZ82aJdYNm5ubhw4dKtbdfA3LsteuXQsODpY6EE/V3t7ucDhE/Iedr7l+/bpKpVKr1VIH4qmuXr0aEhLSD/+SYGklJU28Lrl06XR9fX1ubm7nwpCQEO5NXV3dkSNHHnnkEY1Gc88997zxxhsd54SGhtbX10dERDQ0NHB/vXcvcQf+fQoAAKJhKcNSB88XExYWNvm/+fv7czcMDQ197733Tp48SSn95JNPbrrpJkKI0WgkhCQkJJSVlVFKy8rKpk2b5rTEHUiEAAAgGkoclPJ8uWxKValUzz33XFFR0f333//xxx+vXLmSEJKTk0MIWbhw4ZkzZ9LS0s6ePZuenu60xB0Dt2kUAACAEBIbG7t9+/bOJeXl5YQQrVa7YcOGzuXdS9yBRAgAAKKhLEOpg+8lfRSMm5AIAQBANJTyToSEIhECAEA/+v7777/88suwsLDbb79dqVRKHY70kAgBAHzChQsXRo8e/fDDD3/wwQfXrlk0Gk1QUNCxYx9FRkaK+TXUwbtpFDVCAADoB1sKCkaNGfPOO+9evXqNENLWZm9uvvrggw9WV1eL+C0CmkYp4dmUKjYkQgAA78eyrGHf68OGDrt69SohHdPqZSaT6fLly9yyZKKglKGshw2WwTxCAADvV1FRMSl4RMOlS4RQQtiOF6V07Z/+dO3aNakDlBISIQCAd6KU5ufny+Xy0UPDXsh66n80ujtCR4epNWq5XCEjhMgIkbW1te3fu6/04EHxvpThO6Fe8lGjSIQAAN5JJpOtX7/+8DvvaP39M4aMjdGGxGlCRqgD1HI/hsq5REgojQ0c9r//51WxvpT3sjI3WlmmHyARAgB4s3vvu+/Pf335w+v1LY72jxrOnb7efM3hIIQNUMiDlUrS2v5Q6C+spku/jrzppvCI26f/pqamRuqQ+xsGywAAeDOLxVL27nvXLZasls9ar1379eCw4011hMgpkYWp/Nf/clqwctDzAfEfXav7yP/q30pe41a1Fo4yhOdgGcIy0iYjJEIAAG/27bff/u213Wq5In7k2F8Hjb7Scu1K+/UWh/1yu/WH6837L36/dMyvG9pb32oxnTj13aBBg3r7faywRCglEZpGGYbJzMzk3lssljVr1qSkpOTn51ssli5HezoHAAD6SGRk5MtFr0y67bagdjo3aFRtS6NO7W9jHZTKAhTKQD/VvK/ef7f5gsLBMIwYCYkyhDr4vaSeR9jbRFhaWrpixQqTycR9NBgMOp3OYDCEhYUdOHCgy1Gn5/QyAAAAcCEsLMzf3/+bH099dPmn+vbrNtZxqqWpyd5Oiayddbx16YdrrION/+VPDZfLjx6VOlhp9LZpNCoqKjw8PD8/n/tYWVm5bt06lUqVlJT07LPPZmVldT7q9JzFixc7vTPLshcvXjx9+nRHSWBgYMeexQIwjEj/3vFJzM+kDsRT4QH2Eh6gYO+///6SJUvkbQ61QvH4Nx9lRPzqR2tzVVMdS2kb67A4HEOHhuT96dnVz+YfO3bM9ROWy+U33uOeOjyuabS3iTAuLq7zR7PZzK1QoNPpGhsbuxx1ek5Pd75+/XpRUZHBYOgomTNnjvsbLXZnsVhUKpXgy30cwzAWi0WhUEgdiKdqa2tzOBwsy0odiKdqaWlRKpVqtVrqQDxJfX39o48sMX5ePYihrFymVSjX/2rmCHUAIeR7a9PaU8dbWSYoKDAxMTEiIoIQ8vDDD7ueWR8YGHjjvwQ8sI9Q5MEylFLu3wuU0p7+m3fnHEKIVqt9/vnnZ82aJVZsDMMEBweLdTdfwzCMQqHAAxTMZrM5HA6tVit1IJ7Kz89PpVL1cyKcP3/+1atX+/MbRTfIXxMQHKSRKwIUSrlM9qf6EyMHDSaEfHPZFDxy+J1Tpuj1+jvuuEPqMCUmciIMDQ2tr6+PiIhoaGgYOnSo4HMAACT33XffffDBB1JHIb6vvvrqiczFZ86ea2lpyc3NFfnurIOwdn6X8N2/UGwiJ8KEhISysrJFixaVlZVNmzaty1Gj0RgbG+v6HACAAUIul4u4GvXAERISYmHsarniypXGRx55RNzdJ/4zapTvJZISeWWZhQsXnjlzJi0t7ezZs93783Jycm54DgAA9DWWEhvjIIRcvFgn8q25wTK8XlInQnFqhOXl5dwbrVa7YcOGno5yb5yeAwAA/UtOCNFoej2D3vNhZRkAAB81eHDAY489JvJNMWoUAAAGPrVaPXz4sMceeyw7O1vcO8soK+PZ1Mn3fNEhEQIA+JzExLk7d+6UOoqBAokQAMDn9OHiGJTh3dRJJV5oAokQAABEI2NZGd9EKPWKS0iEAAB97tNPPz18+LDJZBo9enRycvKUKVOkjqjPsPxrhBgsAwDgxS5evHjvvfeeP3/ebG7iSoqLd/3iF794993DYWFh0sYGHJEn1AMAQIfW1tbbb7/9669PmM3NhMi4l9ncVF39xW9/+1u7nedSZM503/NVaiyhDL8XkbhpFIkQAKCvbN269d//rqO0aznL0vPnTdu3b+/l/Z3u+Sotro+Q5wuJEADAS+3du9dma3N66Pr11tdee43X3bZu3frWW291vDcYDFFRUQsWLOhlkIBECADQV1xv79fQ0MDrbjNmzKisrCSE2O32ysrKWbNmxcXFJSQk9CpE0XGDZfi+JIVECADQV1zv5y6X8/sbOC4u7vz5801NTV9++eXYsWMH5lgbGeXfNIp5hAAA3io6OtpkukiI03RIY2Njed3Nz89v6tSpx48fP3HixJ133ilKhOLzwOkTqBECAPSVNWvWBAcHOT0UEhL8zDPP8L3hjBkzjh079vnnn8+YMaPX0cF/IBECAPSV2bNnp6Y+OGTI4C7lQUGDMzP/IKB775Zbbjl9+nRsbKxWqxUpRrGxVMay/F7dh9X2LzSNAgD0oR07dkyZMuW5556z2Wwsy8pk8oAA/+eff37+/PkC7qZSqcaMGdOlXbRjz9cBQcBaozeaPlFRUVFSUtLQ0BAZGZmdnR0REcGV33XXXV3OvPnmm2tra7n3iYmJWVlZ7nw/EiEAQN/KzMzMzMysq6urq6sbOXKkTqcTdh+Hw3H27NmLFy/edttt4kY4kNXV1RUUFGzevDkyMvKdd94pKCh46aWXuEOHDx/uOM1gMNjt9g8++MBgMGg0GsJnYXE0jQIA9IcRI0bEx8cLzoKEkKqqqlWrVv3xj39UqVQiBiYyluU9d8LlfoR1dXWzZ8+OiYlRq9Vz5sy5cOFCxyHNzy5duvTtt98mJSUxDJOXl5eamrpx40ar1epmyKgRAgD0k7a2NrVaLfjyGTNmDPwxMjKxm0bj4+Pj4+MJIQzDlJSUzJw5s8sJdrt9y5YtTz31VHNzc3R09NKlS8PCwrZv315UVJSXl+fO9yMRAgD0k9WrV2/ZskXqKPoYS3lvq8SyFRUVy5Yt61y2evXqkSNHdnysrq7etWvXpEmTui+sevDgwZiYmDFjxhBCCgsLuUK9Xq/X6938fiRCAID+8MMPP7y4devKlSs7xnpAh1/+8pcZGRmdS4KDg7k3lNLi4uLa2to1a9Z0f3QMwxw5cmTz5s2EkFOnTtnt9vHjxxNClEqlUql089uRCAEA+sO+10p+PSTsjT17n1q9SupY+pCM/4R6GWXDwsJ62qPxxIkTVVVV27ZtUygUra2thBCNRmM0GrnlCL7++uthw4ZxdUebzbZ+/frCwsLw8PC9e/dOnz7dzQCQCAEA+sPbhjefjrz1L/8rciLsaWqBZCgr7soyRqPRZDKlpKR0lJSXl+fk5HCTRo4ePXrrrbdy5RMnTkxPT8/Pz7darZMnT37iiSfc/H4kQgAA8V29evXo0aPc+6ampsbGxjCqCFP5+zc7/vKXv2i12o6mv3vuuUfw7HgXUwu8RkZGRpdWU9Jp6mTn1XlkMllycnJycjLfr0AiBAAQ35AhQ6yWlscfWzpz6JhfqgcTQh4aPIoQskA7uvav+wghNbbm400Xi1991f0suHXr1tGjRz/wwAPc+/Dw8OjoaG5qASFkzpw5+/fv76s/j/sErDUq9aLbmEcIANAn/rAo80ujsU6rUCiViWFRYzWBhJAo/yGJYVFqlfqyv9+/qqvT+Kwv030bpvj4eG7xlJ6mFvQ/mYAl1liJl1hDIgQA6CsxMTEVX3x+Oiq0tPFsR+Eb5jNnbx5edeIrbnyj+3rahqm6unrZsmUBAQFdZiCAm9A0CgDQh7Ra7eSpCc2n3yOEtLOMSq4IJIqJ06f7+/vzvVX3bZhcTy2QBmV5zyNE0ygAgHc7fPCtSYOHvdn006IfPjnUdG7y4GGHDAeE3arLNkzc1IL169eHhoa2trZyswskJmCJNb6JU2xIhAAAfai+vv7fJtPm+pqw5DvP/NsUMPc3LzeeOn369NWrVwXcrcs2TB1TC373M7HD54+rEfJ9SQpNowAAfei9I0dUAf6vHHxj8pQphJBNW/9SeX/KgofmlZWVpaam8r1bl22YnE4tAL6QCAEA+lCYTnfi++869wj+ZsaMk99/x43/5MUztmESUMOTuo8QiRAAoA8lJiZ2Lxw8ePA999zD91ZVVVUvvfTS8uXLB/I2TDL+i25LPn0CiRAAwDN4xDZMngiJEAAAxCNgP0KXG/P2AyRCAAAQj5D9CNE0CgAAXsMDB8tgHiEAAPg01AgBAEA8ApZYQ9MoAAB4D5byTmxIhAAAA9P169ffeOMNqaMQ3+nTp6UOYWBBIgQAcK6oqOiLL76QOgpxXL9+XaPRyGQyQohGo+F29+0TlPIf/IKVZQAABqSZM2cOhK1uRWE2m0NCQrhE2LfQRwgAAD7NA/sIMX0CAAB8GmqEAAAgHgEry0g9oR6JEAAAxOOBfYRoGgUAAJ+GGiEAAIhHwGAZilGjAADgNSh2nwAAAF/G8p8fL3EeRB8hAAD4NtQIAQBAPAJqhBLPnkAiBAAA8QhYalTqaYRoGgUAAN+GGiEAAIjHAwfLIBECAIB40EcIAAA+jXpejRB9hAAA4NNQIwQAAPFQ/jU8qWuESIQAACAaysooK+N5yQ0yYUVFRUlJSUNDQ2RkZHZ2dkRERMeh5cuX19bWcu8TExOzsrIsFsumTZtqamomTJiQm5s7ePDgGwaAplEAABi46urqCgoKcnJyDAbD1KlTCwoKOg5RSk0mk8FgOHz48OHDhx9//HFCiMFg0Ol0BoMhLCzswIED7nwFEiEAAIiH/jxw1P2XywphXV3d7NmzY2Ji1Gr1nDlzLly40HHIbDYzDJOXl5eamrpx40ar1UoIqaysTEpKUqlUSUlJFRUV7oSMRAgAAOJhCWFlvF89i4+Pz8rKIoQwDFNSUjJz5syOQ42NjdHR0dnZ2fv27QsICCgqKiKEmM1mnU5HCNHpdI2Nje6EjD5CAAAQDaUySnn2EVJZZWXlH//4x86Fq1atGjlyZMfH6urqXbt2TZo0KTMzs6Nw3LhxhYWF3Hu9Xq/X6wkhlFKZTMa9Yd3bEAqJEAAAJBYdHT1//vzOJcHBwdwbSmlxcXFtbe2aNWs6D5MhhJw6dcput48fP54QolQqlUolISQ0NLS+vj4iIqKhoWHo0KHufDsSIQAAiIfeoKnT6SU6nW7q1KlOD544caKqqmrbtm0KhaK1tZUQotFojEZjbGyszWZbv359YWFheHj43r17p0+fTghJSEgoKytbtGhRWVnZtGnT3Pl+JEIAABDPjfr8nF3i6qDRaDSZTCkpKR0l5eXlOTk55eXlEydOTE9Pz8/Pt1qtkydPfuKJJwghCxcu3LhxY1paWnR09KpVq9z5fiRCAAAYuDIyMjIyMroUlpeXE0JkMllycnJycnLnQ1qtdsOGDby+AokQAABEI2CwDOF7vtiEJ0KGYfR6/e7duwkh3Wfyu5jb330hgF7+GQAAYKDgpk/wInUiFDiPsLS0dMWKFSaTifvYfSZ/T3P7nS4EAAAAIBWBiTAqKmrBggUdH7vP5O9pbr/ThQAAAMBLUDlh+b8kJbBpNC4urvPH7jP5e5rbzy0EsHTp0rCwsO3btxcVFeXl5Tn9CovFkp2d3TGVhBAyd+7czlMp+WpqapLLsZKOQAzDWCwWSqVeJd5jtbW1ORyO9vZ2qQPxVC0tLUqlUq1WSx2Ip2pqaiKEcDPNBRsyZIhCobjBSUIm1AsPSRTiDJbpPpO/p7n9ThcCcEqr1ebn53PzQjgajSYgIEBwkA6HIygoSPDlPo5hGLlcjgcomM1mczgcWq1W6kA8lUKhUKlUSISCMQwTFBTUy0ToVl1CwPQJqfsIxUmE3Wfy9zS33+lCAE7JZLKgoKCwsDBRIiSEyOVy1AgFo5TiAfaG/GdSB+Kp8AB7iXt6vUyE3kqcXxU3k59S2jGTv3uJ0WgkhNhstrVr1547d85ut3csBAAAAN6B0v9sScjrJW3M4iTChQsXnjlzJi0t7ezZs+np6U5LcnJyCCEdCwHMmzfPYrG4aBoFAADPQ+W8X1Lvg9SrplFubj9xNpO/e4mLhQAAAACkgpVlAABANAKaOqlbeyX1ISRCAAAQDaX8+/yk7iNEIgQAAPGwMt4T5KnEfYQYiwwAAD4NNUIAABCNgN0nvGRlGQAAAEIE7VAvdR8hmkYBAMCnoUYIAADiYeWU52AZKvVgGSRCAAAQDaWEfx8hpk8AAIDXENJH2DeRuA19hAAA4NNQIwQAANFQ9BECAIAvEzCPUPKNedE0CgAAPg01QgAAEI2Q3SekrhEiEQIAgIhkfPv8kAgBAMCLsFhiDQAAwKOgRggAAKJBHyEAAPg0SgRsw4SmUQAAAOmgRggAAKKhrIzvyjKSD5ZBIgQAANEI2qEeiRAAALyFJw6WQR8hAAD4NNQIAQBARAKaRvsoEnchEQIAgGiEDJa50ZJsFRUVJSUlDQ0NkZGR2dnZERERrg8tX768traWOyExMTErK8v1/ZEIAQBg4KqrqysoKNi8eXNkZOQ777xTUFDw0ksvuThEKTWZTAaDQaPREEIUCsUNvwJ9hAAAIBpu1Cjfl4sb1tXVzZ49OyYmRq1Wz5kz58KFC64Pmc1mhmHy8vJSU1M3btxotVpvGDNqhAAAIBpK+Y8adXl+fHx8fHw8IYRhmJKSkpkzZ7o+1NjYGB0dvXTp0rCwsO3btxcVFeXl5bkOAIkQAABEI2geIamsrFyxYkXnwtzc3JEjR3Z8rK6u3rVr16RJkzIzM7tc3uXQuHHjCgsLuUN6vV6v198wACRCAACQ2E033ZSamtq5JCgoiHtDKS0uLq6trV2zZk3nYTI9HTp16pTdbh8/fjwhRKlUKpXKG347EiEAAIiH/6hRSuXDhw+fPn2606MnTpyoqqratm2bQqFobW0lhGg0GqPRGBsb6/SQzWZbv359YWFheHj43r17e7ptZ0iEAAAgGgFNo8Tl+Uaj0WQypaSkdJSUl5fn5OSUl5c7PTRx4sT09PT8/Hyr1Tp58uQnnnjiht+PRAgAAANXRkZGRkZGl8Ly8vKeDslksuTk5OTkZPe/AokQAABE44n7ESIRAgCAaISNGpUWJtQDAIBPQ40QAABE44nbMCERAgCAaCiVsTdaRLvbJRK3TSIRAgCAaDxxh3r0EQIAgE9DjRAAAEQjZNFtqUeNIhECAIBo0DQKAADgYVAjBAAA0XhijRCJEAAAxEP5JzYkQgAA8BpC1hrto1Dchj5CAADwaagRAgCAaFgqY3nWCPmeLzokQgAAEI3oG/P2AzSNAgCAT0ONEAAARIPpEwAA4Nv4T59AIgQAAO8hpEbYR6G4DX2EAADg01AjBAAA0QiZUI+mUQAA8BoC5hFKngjRNArQty5fvix1CADgChIhQN/KfvLJ5uZmqaMA6C9URvm/pA0ZiRCgD1ksltKDB996802pAwHoJwKyIBIhgLe5cuUKwzDc+7dLS+/QRe4t/pu0IQH0G5bKWCrn9aIEiRDAW+zdW3LLLTcvWJCoCwsNUA+KGTl68zP5D4ZEqq5cnTTu5qiRoyaOH3/y5EmpwwSA/4JRowDi+PjjY6+/UfR/K3MDAtQsy6Ykbj31mWnV6FuHqwOy1b96w3zmdPSwA0cODx8+XOpIAfqQJ06fQI0QQBz79u1+fkNyQICaECKXyw+9/6QtQH7s2sU2lrnSfr1a1nLk2IfIguD1KJWx/F/SxoxECCCOy5cvDR8+pOOjTCYj19uvt7cvPVN54NqF61ctra2tEoYHAD1BIgQQx623TvmwvKbj408/Xbl89ZoqYWLNmR8fWPPklZZrx44dkzA8gP4haNSoxDGjjxBAHE8+mXvnnTMslrYZt4/76Ux93uq3n3v++aeffpoQkpGZmTB9+scffyx1jAB9jqW8d5yXfNQoEiGAOAIDAz/++Hhx8Y6tW6pHjhx95MhHY8eO7Tg6evTo9PR06aID6CfYjxDApwUEBGRlZUsdBQDwg0QIAADi8cAaYa8GyzAMk5mZyb23WCxr1qxJSUnJz8+3WCxOS3o6EwAAvANLeM+duGEirKio0Ov1ycnJK1euNJlMHeXu5x3XhCfC0tLSFStWdMRkMBh0Op3BYAgLCztw4IDTkp7OBAAAcKqurq6goCAnJ8dgMEydOrWgoKDjkPt5xzXhTaNRUVHh4eH5+fncx8rKynXr1qlUqqSkpGeffXbx4sXdS3o60+n97Xb7J5980tDQ0FESGRk5YcIEwQG3tbXZbDbBl/s4hmHwAHvDZrM5HA4/P3RGCGSz2ViWpZIPtPdY3H+/MlmvGiFVKpVcfoPqEyUyvqNAXf+fWldXN3v27JiYGELInDlz9u/f33HI/bzjmvD/LOPi4jp/NJvNOp2OEKLT6RobG52W9HSmUw6H46effur8u1epVNyzEMbhcDgcDsGX+ziGYfAAeamqqho5cuTo0aO5jw6Hg3uG0kblubinp1AopA7EU3H//fY+Ed7wHMp/pRjX58fHx8fHxxNCGIYpKSmZOXNmxyH3845rov37lFLKPWJKKcuyTkt6OtMpjUbzwAMPzJo1S6wIbTabVqsV626+hmEYlmXxAN331r43oqKjc/Oe4T76+fk5HA48QMEopSqVSq1WSx2Ip2pra9Nqtb1MhO6gVMDgF9mnn366cuXKzkVPPfVUeHh4x8fq6updu3ZNmjSpY2AK4ZN3XBMtEYaGhtbX10dERDQ0NAwdOtRpSU9nAngZh8Px0dHyz//1r45ECAAuREVF3X///Z1Lhgz5z4KFlNLi4uLa2to1a9ZERER0Psf9vOOaaIkwISGhrKxs0aJFZWVl06ZNc1piNBpjY2O7lwN4h0OHDtXX17e1tn5f+11cQOhlm+3Z1XnDdDq1ZlBERMTtt98udYAAfY6lhOXZk8tSMmLEiBkzZjg9euLEiaqqqm3btikUCm7BXo1G01M2EZZfRFtrdOHChWfOnElLSzt79iy3gkb3kpycHKflAN7htttu27195yt/ej74o6/vHzwyM3Cs/cCHG55e/UJ+fkBAgNTRAfQH0XeoNxqNJpMpJSXldz8jPWcTYflFNmBHYT311FNz584VsY8QLbG9wTDMtWvXgoODpQ5koHM4HLkrsmrfP5atG08IefrC8bTlM//nvtiMBa8dOvT3jrEzwJfFYkEfYW+YzeaQkJB+6CP8dH6u5YdzvC753tr4/czooqKiPgrphrD7BICY/Pz8UuY9pFIoPjJf+MpSr1b5PTg/4Ve/ish4ePLf//53qaMD6HO0DybU9zUkQgCRlb6x//NL598jFz8a1niuoelw6ReEkNGjQ/797/NShwbQ51gi4/+SGBIhgMg++/LLwh1FqiEBhz/MzXvhgY+P1RBCKipOT5wYd8NrAaD/YZ0LADG1t7eXvnNIp9OdN/30h4f/9szqe2Injd1S+Peab64+t3aO1NEB9DlswwTg61QqFbewxapVa957L65g0x6r1Tpz5l0HD/71hmtTAXgBio15AaBDYuK9iYn3cu+5tUaljQegH3hijRD/RAUAAJ+GGiEAAIiGJYTvKFDJR40iEQIAgGgENI0SNI0CAABICDVCAAAQDbeyDK9L+C7SLTokQgAAEI2gHeoxfQIAALwFK2AeIfoIAQAAJIQaIQAAiIYK2Ji3byJxHxIhAACIxhP7CNE0CgAAPg01QgAAEI2gtUb7KBZ3IRECAIBoBCyxJnUeRCIEAAARYfcJAAAAz4IaIQAAiAZNowAA4NOwMS8AAICHQY0QAABEwwrYfULqCfVIhAAAICa+fX7oIwQAAO+B3ScAAAA8DGqEAAAgGoqmUQAA8GWU8h8sI3UmRNMoAAD4NNQIAQBANAKaRiWHRAgAAKJhhTSNYh4hAAB4C5by7scsb8cAABaKSURBVPOTfD9C9BECAIBPQ40QAABEJKM8l0zje77okAgBAEA0LP/pEDfctolhGL1ev3v37o6Su+66q8s5N998c21tLfc+MTExKyvL/QCQCAEAYOAqLS09duyYyWTqXHj48OGO9waDwW63f/DBBwaDQaPREEIUCgWvr0AfIQAAiIYlMr4v102jUVFRCxYs6FKo+dmlS5e+/fbbpKQkhmHy8vJSU1M3btxotVp5xYxECAAA4qGE8n+5EBcXl5CQ4PSQ3W7fsmXLsmXLmpubo6Ojs7Oz9+3bFxAQUFRUxCtkNI0CAIBoBOw+wVLZ8ePHs7OzOxdmZ2eHh4e7vvDgwYMxMTFjxowhhBQWFnKFer1er9fzCgCJEAAAJDZ27Njf/e53nUuGDBni+hKGYY4cObJ582ZCyKlTp+x2+/jx4wkhSqVSqVTy+nYkQgAAEI2w3SfCw8N/+9vfunOy0WiMjY0lhHz99dfDhg0bOXIkIcRms61fv76wsDA8PHzv3r3Tp0/nFQD6CAEAQDTcyjK8XrxWlsnJyeHeHD169NZbb+XeT5w4MT09PT8/f968eRaLBU2jAADgbcrLy7u8eeaZZzqOymSy5OTk5ORkYTdHIgQAANFQImN5rhTDSt04iUQIAACiueF0CKeXSAt9hAAA4NNQIwQAANGwbqwd2oXUFUIkQgAAEA8VsDEvdp8AAACvQQX0+UldJUQfIQAA+DTUCAEAQDQC+gj5ni86JEIAABAPpk8AAAB4FtQIAQBANCz/lWVcb8zbD5AIAQBANNw62nwvkRYSIQAAiEbYNkzSQh8hAAD4NNQIAQBANJR/06jko0aRCAEAQDRC+gj7JhL3oWkUAAB8GmqEAAAgGkpkfKdDYPoEAAB4D8p/OgT6CAEAwHtg+gQAAICHQY0QAABEg5VlAADAp3niPEI0jQIAgE9DjRAAAETjiYNlREiETU1NO3bs+Oqrr+Ry+ZQpUx599FF/f3+nhR2XLF++vLa2lnufmJiYlZXV+zAAAEBynriyjAiJcMuWLb/4xS/27dtHKS0pKdmzZ8/SpUudFnLnU0pNJpPBYNBoNIQQhULR+xgAAGAgoPwTm+Q1QhH6CI1GY2pqqlKpVKlU8+bNq6io6KmQYzabGYbJy8tLTU3duHGj1WrtfQwAAADCiFAjHDdu3P79+x966CG73f766683Njb2VMhpbGyMjo5eunRpWFjY9u3bi4qK8vLyut/22rVrjzzyiFar7ShJSUlZtmyZ4Dg7xwB8MQxjsVgYhpE6EE/V1tbmcDhsNpvUgXiqlpYWlUqlUqmkDsRTNTY2siwrk/VqMbOgoCA/vxtkDUp5jwKVfNSoCIkwNzf35ZdfTktLCwoKSklJGTx4cE+FnHHjxhUWFnLv9Xq9Xq93etvAwMDi4uJZs2b1PsIOQ4cOFfFuPoVhGJVKFRwcLHUgnspmszkcjs7/sANe1Gq1SqVSq9VSB+KpZDJZSEhILxOhO3x0HuGgQYPWrl2rVCoJIUajcdSoUT0Vck6dOmW328ePH08IUSqV3DkAAACSEKGPsLi4+MUXX7RarWaz+dVXX01OTu6p0Gg0EkJsNtvatWvPnTtnt9v37t07ffr03scAAAADAeX/kpwIiXDJkiUWi2X+/Pm5ubl33333jBkzeirMyckhhEycODE9PT0/P3/evHkWi6WnplEAAPA49OfWUV4vaYnQNBoYGLhu3Tp3CsvLywkhMpksOTmZqyMCAIA3oYRSntU8vueLDkusAQCAT8MSawAAIBoBTZ3eMH0CAACA0xdrjTIMo9frd+/e3bmw+1KdFotl06ZNNTU1EyZMyM3N7TxtzzUkQgAAGLhKS0uPHTtmMpk6FzpdqtNgMOh0umeffXbnzp0HDhxYvHixm1+BPkIAABAN5T9k1PXapFFRUQsWLOhS6HSpzsrKyqSkJJVKlZSU1HldzxtCIgQAANGIPo8wLi4uISGhSyG3VGd2dva+ffsCAgKKiooIIWazWafTEUJ0Oh2vNTXRNAoAAKIRtsRaVVXVU0891blw5cqV4eHhPV3idKlOSim3hhyllGV57IGBRAgAABIbM2bM3LlzO5cEBga6ON/pUp2hoaH19fURERENDQ28VpZGIgQAAPEI2n0iPDzczS0WjEZjbGyszWZbv359YWFheHh4x1KdCQkJZWVlixYtKisrmzZtmvsBoI8QAABEwxLK98VrZRkXS3UuXLjwzJkzaWlpZ8+eTU9Pd/+eqBECAMBAx63QSVwu1anVajds2CDg5kiEAAAgGsp/pRisLAMAAN6Da+3ke0kfBeMm9BECAIBPQ40QAABEQ8kNVopxeom0kAgBAEA0lFLKs9MPfYQAAOA9WP41Qr7niw59hAAA4NNQIwQAANFQARvz9k0k7kMiBAAA0VBCeK0UQ/ifLzo0jQIAgE9DjRAAAEQjbBsmaSERAgCAaCj/lWUkbxpFIgQAANFQnrtJkAGQCNFHCAAAPg01QgAAEI0nLrqNRAgAAKIR1EcoMTSNAgCAT0ONEAAARMNSwvJcRZvv+aJDIgQAANEIGDUqeeMoEiEAAIjGEwfLoI8QAAB8GmqEAAAgGgGjRiWvESIRAgCAaAQ0jWL6BAAAgJRQIwQAANFQQilheV7C73zRIRECAIBo0EcIAAA+jZWxrIxnjVCG6RMAAADSQY0QAABEQwllefb58T1fdEiEAAAgGpawfBMbNuYFAACQEmqEAAAgGkpYTJ8AAADfxcoIy3MUKN/zRYdECAAAoqHoIwQAAPAsqBECAIBoBIwaxfQJAADwHn2x1ijDMHq9fvfu3Z0LKyoqSkpKGhoaIiMjs7OzIyIiCCHLly+vra3lTkhMTMzKynInACRCAAAYuEpLS48dO2YymToX1tXVFRQUbN68OTIy8p133ikoKHjppZcopSaTyWAwaDQaQohCoXDzK9BHCAAAoqGEYXm+XNcIo6KiFixY0KWwrq5u9uzZMTExarV6zpw5Fy5cIISYzWaGYfLy8lJTUzdu3Gi1Wt2MGTVCAAAQDStkiTVXo0bj4uK6F8bHx8fHxxNCGIYpKSmZOXMmIaSxsTE6Onrp0qVhYWHbt28vKirKy8tzJwAkQgAAEA0lDCUM30v++c9/Pv30050Ls7KyRowY4frC6urqXbt2TZo0KTMzkxAybty4wsJC7pBer9fr9W4GgEQIAAASGzVq1Jw5czqXDB482MX5lNLi4uLa2to1a9Zww2QIIadOnbLb7ePHjyeEKJVKpVLp5rcjEQIAgGiETagfOXLkHXfc4c7JRqMxNjb2xIkTVVVV27ZtUygUra2thBCNRmOz2davX19YWBgeHr53797p06e7GQASIQAAiEZQHyGP83NycsrLy41Go8lkSklJ6SgvLy+fOHFienp6fn6+1WqdPHnyE0884eY9kQgBAGCgKy8v7/wmIyMjIyOjyzkymSw5OTk5OZnvzZEIAQBANIIGy2BlGQAA8BYCdqiXPBFiQj0AAPg01AgBAEA0aBoFAACf1heLbvc1JEIAABANS1mW8kyEFBvzAgAASAc1QgAAEA0lLJpGPYDD4Th58mRbW9vEiRMDAgKkDgcAwItQSnk3jSIR9q8TRmPG7x8aIx+kIrLa1uY/bdr4YNo8qYMCAADJ+FYirKmpyfj9Q08GRIYP0hJC2oYwa55ePSlhSmRkpNShAQB4A0FNoxgs04/mz0sLJyouCxJC1HLFLHXoh0fLpY0KAMBrUMJyUwndfxGe8w5F50M1wpMnT9Z8W+M/fAwJ/f+FKqKwtrRIFxQAgFdhCe/pE6zU0ye8PBG2t7dbrVbuf/cUv5oefvPRpgtnQq8O8VOp5QqVXFHlaNo6a6bUYQIAgGTESYRNTU07duz46quv5HL5lClTHn30UX9/f6eF3PkWi2XTpk01NTUTJkzIzc11vRNxbzgcjiezVpb8757pY6LVLMkeMcEscxSdP1Fz7coY/8BhoaF3PHh/fHx8H307AIDv8bzpE+L0EW7ZsmX48OH79u3bs2ePVqvds2dPT4Ucg8Gg0+kMBkNYWNiBAwdEicEpf3//3SWvFe/Y0exoWxAcqZTLM4ZFK/wH3XfP3GVr1+w4/NbzhZv77tsBAHwNpZRSlt/LOwbLGI3G1NRUpVKpUqnmzZtXUVHRUyGnsrIyKSlJpVIlJSV1Lu8ji5csuSvpvp9ar11otZjbW0dNjCl9792VOdmxsbF9/dUAADDAidM0Om7cuP379z/00EN2u/31119vbGzsqZBjNpt1Oh0hRKfTdS7vrLW1ddeuXR27EhNCJk2aNGfOHGERlr33wUhGftBxZZhD9sPFxubmZj8/L+8fFRfDMFarValUSh2Ip7LZbAwj8dA4j8Z19tvtdqkD8VRWq1WlUslkst7cxN/fXy6/QfVJwO4TXjJqNDc39+WXX05LSwsKCkpJSeH6/JwWciil3P8flFKWdd467OfnFxMTM2HChI6SsWPHCvuL+Pz58xfqL+nXrn18+R8/eP/9JYv11dXVM2bMEHArnyWXy/38/JAIBWMYRiaT4QEKxv388AAF455eLxOhO5dTISvLeMWo0UGDBq1du5b7jRqNxlGjRvVUyAkNDa2vr4+IiGhoaBg6dKjTeyqVyt/85jezZs3qfXhNTU0nar7RarVqtTo5JeU3M2acOHFCrVb3/s6+g2EYtVqNhyYYpVShUOABCtbe3q5SqfAABeOeXi8ToTs8ca1RcfoIi4uLX3zxRavVajabX3311eTk5J4KjUYjISQhIaGsrIxSWlZWNm3aNFFicGHy5Mljx47t+Dh06NDZs2f39ZcCAIBHECcRLlmyxGKxzJ8/Pzc39+677+ZaHZ0W5uTkEEIWLlx45syZtLS0s2fPpqenixIDAABIj++QUcp6SdNoYGDgunXr3CnkBr9otdoNGzaI8tUAADBwcEus8bxI4sEyXrvWaGVlpdQhAACAB/DORNje3p6xYEFbW5vUgQAA+Bb+7aLSN416fCJsbW09f/58lzkYR8vKLl+se//99294+cWLF/ssNAAAn0MJ5QaO8npJG7MHTypvbW19PHPxF59WDRvkX2ezPv5k1uTp065du0YI+WvhXx6NmPjKlq3c5MXAwMDo6Ojud7DZbH9YuHBoUMikKZP/oF8cEhLS338GAADvQilDKb8+P+xQL1z2408EfvnjX0ZOIoS0sczGF4v+9fnnbxj2/8/oX45QqOeG3XT18umXHlnxwfnvcp7Mfn7jC83NzZ0vv3r16qyEaWdOn54aPKLu5PkZrxQd+vCo03wJAABezIObRis//uS+kDHce7VckRY4Wk1l77777mm7dar/MELIVP9hp+0thw4d+nPB5u4Lqj2/5tlfWWVThgw/1dI0J2T0iiHRyxc/0t9/BgAA74Km0f5jt9sV/929GuSnulJff8/cuW8kTLlee6Wq6WKQclBsfPy9993X+TSWZdetW/fnjRsHy/2CFaonxt5S2fjvZd98VNdmZWTkn//8Z0JCQr/+SQAAvAolvAe/eMU8wv6nVCrVgdor7deHqf6zx2H1dXNCyr0sy3748bFvNIMjY8f/dKLmYsU1hmEUCkXHhXK5fO3atXG//vWyjMzlY+NjtCE3a0MPXfrhH40mNngwsiAAgK/x4KbRl4p3brhS848mU21L45uNPx1X2bKeyvnss8+0gYE73z7w1ntHig8dHBISfPz48e7XJt9//zz9oncun7Yy9naWKb3041TdmJTf/77//xQAAN6F5f9CjVCoKQkJf6+q/NvO//PN2XOTpj+4bdEilUql0Wi+rvnG39+fEHLrrbd+/c3J77//3unl1G5XyRWLTx4dOWiwWuHXMmFs0Z9f6N8/AQCA1+G/+wRBH2FvREREPLv+v1Zx67LXrkajiYuLc3rtkfffH39LzMndlUfLjq7OyZk0LQE7vAAA9BrLP7FhQr0ULl26lL/uubf//v6IESMe/sPDNT/+EBgUJHVQAAAgAc+uEQo2fPjwhQsXdnwMCgp65BHMnQAA6C1B0yG8Yj9Cj/Dmm29KHYIHa2lpcWfJOujJjz/+WF1dLXUUHuz48eMXLlyQOgoPdujQofb29v74JkoJZXm+btA0yjBMZmZm5xKLxbJmzZqUlJT8/HyLxeK0xH0+lAh37twpdQgerLm5+fXXX5c6Cg9WU1PzySefSB2FBysrK/vhhx+kjsKDvfbaazabTeoohCgtLV2xYoXJZOpcaDAYdDqdwWAICws7cOCA0xL3+VAiBACAviZgWRnXTalRUVELFizoUlhZWZmUlKRSqZKSkioqKpyWuM9H+wgBAKBviDxq1OnIf7PZrNPpCCE6na6xsdFpifsGbiLUaDS7d+8uKysT64YWi2XVqlVi3c3XtLS0XLhwAQ9QsNOnTzc1NeEBCvbFF19cvnz5H//4h9SBeKqGhobnnntOpVL15iYrVqwYMWKE63O+/vpLvrf97LPPli5d2uW/DtffRSmVyWTcG24bvu4l7hu4ifCZZ5759NNPRbzhnXfeKeLdfFBycrLUIXgw/Px6CQ+wl0R5gNzGdqKLiYnZvHkzr+8KDQ2tr6+PiIhoaGgYOnSo0xL3DdxEqNFo8NMHAPB6gYGB7v9tbzQaY2NjExISysrKFi1aVFZWNm3aNEJI9xL3YbAMAAB4jJycHELIwoULz5w5k5aWdvbs2fT0dKcl7pNR3vtlAAAAeA/UCAEAwKcN3D5Cvpqamnbs2PHVV1/J5fIpU6Y8+uij/v7+Tgu58y0Wy6ZNm2pqaiZMmJCbm9tHncCegu/TI4QsX768traWe5+YmJiVlSVR7AMCwzB6vX737t3E2U+rpx8bfoRE6KMj+AX+N9ePsfNRDn57nXlPjXDLli3Dhw/ft2/fnj17tFrtnj17eirk9GYZAu/D9+lRSk0mk8FgOHz48OHDhx9//HHpYpdel5Uv3F/zAj9CwY8Ov8DOXD9GN1dm8WXekwiNRmNqaqpSqVSpVPPmzeNWFnBayOnNMgTeh+/TM5vNDMPk5eWlpqZu3LjRarVKF7v0uqx84f6aF/gRCn50+AV25voxurkyiy/znkQ4bty4/fv3t7S0NDU17dmzh1tZwGkhpzfLEHgfvk+vsbExOjo6Ozt73759AQEBRUVF0sUuvbi4uISEhI6P7q95gR+h4EeHX2Bnrh9jl6NOz+nPaAcg70mEubm53NjZ5cuXjxgxgmvydlrI6c0yBN6H79MbN25cYWHhTTfdFBgYqNfrsa9CZ+6veYEfYRfuPzr8Al1w53eF315n3jNYZtCgQWvXruV2mTcajaNGjeqpkNObZQi8D9+nd+rUKbvdPn78eEKIUqnkzgGO+2te4EfYhfuPDr9AF9z5XeG315n31AiLi4tffPFFq9VqNptfffVVbj0wp4VGo5H8vAwBpVTAMgTeh+/Ts9lsa9euPXfunN1u37t37/Tp0yX+Awwk3X9a3UvwI3TK/UeHX6ALrn9X+O115z2JcMmSJRaLZf78+bm5uXffffeMGTN6KuxpYQJfxvfpTZw4MT09PT8/f968eRaLRa/XS/wHGEjcWfMCP0Kn3H90+AW64Pp3hd9ed1hZBgAAfJr31AgBAAAEQCIEAACfhkQIAAA+DYkQAAB8GhIhAAD4NCRCAADwaUiEAADg05AIAQDApyERAgCAT0MiBAAAn4ZECAAAPg2JEAAAfBoSIQAA+DQkQgAA8Gn/D5RUvLd5KXJrAAAAAElFTkSuQmCC"
},
"execution_count": 41,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"recon = res.W.value * res.H.value \n",
"scatter!(recon[1, : ], recon[2, :]; zcolor=predicted_cluster_assignments, marker=:star)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"@webio": {
"lastCommId": null,
"lastKernelId": null
},
"kernelspec": {
"display_name": "Julia 1.3.2-pre",
"language": "julia",
"name": "julia-1.3"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.3.2"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment