hw2_phi_X_to_Y.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "a0129eae",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAEKCAYAAAASByJ7AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAh+0lEQVR4nO3deXzV9Z3v8dfnJCFhSUCWsCUh7FsIi8EFlM0NlUFG0erUjk7rcG+vTpd7u9h2ep2xnT5sO63tnU7HoR0f6HTUcUMZa0VEUBBBQBbZdySsIWwhZD/f+8c5cQBZEjgn39855/18PHiEk/M753xOzvL+/b7bz5xziIiIhHwXICIiwaBAEBERQIEgIiJRCgQREQEUCCIiEqVAEBERIACBYGZpZrbKzN7wXYuISCrzHgjA14GNvosQEUl1XgPBzPKA24Hf+6xDREQg3fPj/wr4DpB9vg3MbAYwA6Bt27ZXDho0qGUqExFJEitXrjzsnOtyse28BYKZTQEOOedWmtmE823nnJsJzAQoKSlxK1asaJkCRUSShJntbsp2PpuMxgJTzWwX8AIwycz+4LEeEZGU5i0QnHPfc87lOecKgXuBd51z9/uqR0Qk1QVhlJGIiASA705lAJxzC4GFnssQkQRVV1dHaWkp1dXVvkvxKisri7y8PDIyMi7p9oEIBBGRy1FaWkp2djaFhYWYme9yvHDOUV5eTmlpKb17976k+1CTkYgkvOrqajp16pSyYQBgZnTq1OmyjpIUCCKSFFI5DBpd7t9AgSAiIoACQUQkriZMmECiTKhVIIiICKBAEBGJiV27djFo0CAeeOABiouLmT59OqdOnTpjm3bt2n32/5dffpkHH3wQgJdeeomioiKGDx/OuHHjWrLsMygQRCQlhcNw8CA4F7v73Lx5MzNmzGDt2rXk5OTw29/+tkm3e/zxx5k7dy5r1qxhzpw5sSuomRQIIpJywmGYOBHy8mDChMjlWMjPz2fs2LEA3H///SxevLhJtxs7diwPPvggv/vd72hoaIhNMZdAgSAiKaesDJYsgfr6yM+ystjc79nDPi90+fT5Ak899RQ//vGP2bNnDyNGjKC8vDw2BTWTAkFEUk5uLowZA+npkZ+5ubG5308//ZQPP/wQgOeff57rrrvujOu7du3Kxo0bCYfDzJ49+7Pfb9++nauvvprHH3+czp07s2fPntgU1EwKBBFJOWawYAGUlsLChZHLsTB48GCeeeYZiouLOXLkCF/96lfPuP6JJ55gypQpTJo0ie7du3/2+29/+9sMGzaMoqIixo0bx/Dhw2NTUDNpLSMRSUmhEHTtGuv7DPHUU0+d8buFCxd+9v/p06czffr0z93u1VdfjW0hl0hHCCIiAigQRERiorCwkHXr1vku47IoEEREBFAgiIhIlAJBREQABYKIiEQpEEREWsicOXN44oknfJdxXpqHICLSQqZOncrUqVN9l3FeOkIQEYmBxuWvH3roIYqKivjiF7/IO++8w9ixY+nfvz8fffQRs2bN4pFHHgHgwQcf5Gtf+xpjxoyhT58+vPzyy0BkItuUKVM+u99HHnmEWbNmAfDoo48yZMgQiouL+da3vhXz56AjBBFJTeFwZFW73NyYrV2xbds2XnrpJWbOnMno0aN57rnnWLx4MXPmzOEnP/kJ06ZNO2P7/fv3s3jxYjZt2sTUqVPPOYu50ZEjR5g9ezabNm3CzDh27FhMaj6djhBEJPXEaf3r3r17M2zYMEKhEEOHDuWGG27AzBg2bBi7du363PbTpk0jFAoxZMgQDh48eMH7zsnJISsri4ceeohXX32VNm3axKTm0ykQRCT1xGn968zMzM/+HwqFPrscCoWor6+/4PYueqae9PR0wqcFVOMy2enp6Xz00UfcddddvPbaa0yePDkmNZ9OgSAiqSde61/HQK9evdiwYQM1NTUcP36c+fPnA3Dy5EmOHz/Obbfdxq9+9StWr14d88dWH4KIpJ7G9a9j3IcQC/n5+dxzzz0UFxfTv39/Ro4cCUBFRQV33HEH1dXVOOd48sknY/7Y5mJ5QtHmPLBZFvA+kEkkmF52zj12oduUlJS4FStWtER5IpJANm7cyODBg32XEQjn+luY2UrnXMnFbuvzCKEGmOScO2lmGcBiM/uTc26px5pERFKWt0BwkUOTk9GLGdF/fg5XRETEb6eymaWZ2WrgEDDPObfMZz0ikrh8NX8HyeX+DbwGgnOuwTk3AsgDrjKzorO3MbMZZrbCzFaUxWhomIgkl6ysLMrLy1M6FJxzlJeXk5WVdcn3EYhRRs65Y2a2EJgMrDvrupnATIh0Krd8dSISdHl5eZSWlpLqO41ZWVnk5eVd8u29BYKZdQHqomHQGrgR+KmvekQkcWVkZNC7d2/fZSQ8n0cI3YFnzCyNSNPVi865NzzWIyKS0nyOMloLjPT1+CIiciYtXSEiIoACQQIkHIaDByGFB4qIeKVAkECI02rEItIMCgQJhDitRiwizaBAkEAI8GrEIikjEBPTRAK8GrFIylAgSGCEQtC1q+8qRFKXmoxERARQIIiISJQCQcQTzbuQoFEgiHigeRcSRAoEEQ8070KCSIEg4oHmXUgQadipiAeadyFBpEAQ8UTzLiRo1GQkIiKAAkFERKIUCCIiAigQREQkSoEgIiKAAkFERKIUCCIiAigQREQkSoEgIiKAZiqLnJNzjuM11RypquJETQ0na2upqa+nNtyAcw4w0kNGVnoGrTPSyW6VyRVZrenYujUZaWm+yxe5JAoESVnOOUpPnGBzeRlbj5Sz89hR9hw/TumJExyqPEndJa5J3TGrNT2ys8nLaU9hhyvo27EjAzp1ZkDHTmSm6yMnwaV3p6SMqro6Vu7fx/J9pazav5+1hw5woqbms+u7tGlLQfv2lPToSbd27ejSpi0dW7cmJzOLdq1akZmeTqtQiLRQiLBz1IfD1DTUc6qujoqaGo5WV3P4VCWHKivZV3GCzeWHmb9z+2fBkmbGwE6dGd6tO6N79OTqnvl0z8729ecQ+RwFgiS1HUePMH/ndhbu2sXKfXupDTcQin4x39ZvAEW5XRncuQv9OnYiOzMz5o9fHw6z+9hRNpcfZkNZGWsPHuC/tmzi+XVrAehzxRWM69WbSYV9uLpnnpqbPAiHtepsI3Oezt9nZvnAs0A3IAzMdM79+kK3KSkpcStWrGiJ8iSBfXr8GHM2b+KNrZvZUn4YgAGdOjOuoBdj8ntxZfcecfnyb6qGcJhNh8v4sHQPH+zZzdLSUmoa6mmfmcVNffsydcBgrs3LJy2kMR/x1njmuiVLIuelWLAgsgptsjGzlc65kotu5zEQugPdnXMfm1k2sBKY5pzbcL7bKBBS18X24mrq63lr+1b+c90nLN27B4CSHj25td8Abu7Tj545OS1ccdNV1dWx+NPdvLV9K/N2bONkbS3d27Vj+pAivjB0GD2yg1t7ojt4MHIa0/r6yMmKSkuTc0nywAfC2czsdeA3zrl559tGgZCaLrQXd6TqFP++djV/WLuG8qpT9GrfgelDhjJt0BB6JuAXaXV9HfN37ODljet4f/cuzIxb+vbnr0eVMKJbd9/lxV1LN984FzmndeN7a+HC5Gw2SqhAMLNC4H2gyDl34qzrZgAzAAoKCq7cvXt3yxcoXp1rLy7UrpJ/Xbmc59atobq+nomFfXhw+EjGFvQilCSf6NITx/nD2tW8sP4TTtTUcG1ePl+76lquzsv3XVpc+Gq+SYU+hIQJBDNrB7wH/INz7tULbasjhNR0+l7cNeNqmPqj5cxavZLahgbuGDiY/1lyFf06dvJdZtxU1tby/Lq1/O7jFZSdquT6gl58d+w4hnRJrhMxp0rzjQ8JEQhmlgG8Acx1zv3yYtsrEFJXfYPjmY/W89S6RZRXnWLKgIF885qx9O5whe/SWkx1fR1/WLuG365YxvHqau4ZOoxvj7mOjq3b+C4tJlKl+caHwAeCmRnwDHDEOfeNptxGgZCath8p53vvzmPFvr2M6tadH46fxPCu3XyX5c2Jmmr+6aOlPLNmFe1ateL7143nrsFDsST49kyF5hsfEiEQrgMWAZ8QGXYK8H3n3Jvnu40CIbWEnePpVSv5xw8X0zo9g+9dN47pQ4qSpo/gcm0pP8wP3p3Hyv37GN+rNz+98WZy27bzXZYEUOAD4VIoEFJHWWUl//vtN/lgz6fc1KcvP550E13atPVdVuCEnePZNav42ZJFtEnP4Oc3T2ZiYR/fZUnANDUQknAKhiS6Ffv2MuX5f2fl/n38ZNJNPHX7HQqD8wiZ8eCIUcy59366tmvHV+bM5pcffkA4gXb0JDgUCBIoL29YxxdffZE2GRm8es9fcG9RcVK0jcdbv46deOWe+7h7SBG/Wb6Ur/7xdU7V1fkuSxKMAkECwTnHr5ct4TvvzOWqnnm89oUvMqhzF99lJZSs9AyeuOFmHhs/kfk7d3DfK/9J+alTvsuSBKJAEO+cc/xo0UJ+vexD7ho8lKen3kn7rCzfZSUkM+OB4aP419vvYOuRcr7wygscOFnhuyxJEAoE8co5x4/eX8Cs1R/z4IhR/OzGW7TiZ6NwODJb6xL6A27o05dnpt3FoZOV3PfKixw8eTIOBUqyUSCIV08uXcKsNav48ogr+eH1E9Rf0KhxHYe8vMhsrUs4Wc/oHnnMmnYXh09V8sBrL3Osuir2dUpSUSCIN8+vW8tvli/lC0OH8YPrxysMTldWFpmyW18f+VlWdkl3M6p7D2ZOmcauY8f4H2+8Tk19fYwLlWSiQBAvlpXu4bGF8xnfqzc/mnhjbMLgMppYAic3N7J+Q3p65Gfupa9bdG1+AT+76RaW79vLYwvnk0hzjy5VMr0VWpICQVrcocqTPPKnNyho357/N/l20mOxpGUMmlgCxSyy3GdpaUwW9Zk6cDAPj76aFzes48UN62JTY0Al21uhJSkQpEWFneNbb79FZV0t/3LbHbE7c1mMmljiqrm7raFQZLnPGDWlfePqMYzNL+Dv33uX7UfKY3KfQZQIb4WgUiBIi/qPT9aweM9u/vb6CfTvFMMlq2PYxBIXAdhtTQuF+MXNt9I6PZ1vzXuLhiTddQ76WyHIFAjSYg6crODnHyziuvxe3FdUHNs7j3ETS8wFZLc1t207Hhs/iTUHD/Ds2tVeaoi3oL8VgkyBIC3mpx8sojbcELtO5LPFuIklpgK02/pnAwYxrqCQJ5d+wOEknckc5LdCkCkQpEV8cuggr2/eyEMjS+jVoYPvclpegHZbzYwfjptAVV0dv/noQ291SPAoEKRF/HrpEjpkZTHjytG+S/EnQLutfTt24u4hRbyw7hP2VwRjaQsNFfVPgSBxt+lwGe/u2sGXR1xJTqxGFclle3j0NTS4MP+2aqXvUoLQ5y4oEKQFPLtmFVnp6dxfPNx3KXKanjk53Np/AC9u+MT7UtkB6XNPeQoEiavK2lrmbNnElAED6ZDV2nc5cpYvFY/gZG0tf9y62WsdAepzT2kKBImrd3Zu51RdHdMHF/kuRc6hpHtPerXvwOubN3qtI0B97ilNgSBxNXfbVnLbtqWkR0/fpcg5mBm39x/IstI93ldDDVCfe8pSIEjc1IfDLN6zm0mFfQjpUx5YN/TuQ4NzLP50t+9SxDMFgsTN+rJDnKytZWx+L9+lyAUUd+1Gu4xWLN1b6rsU8UyBIHGz+sA+ILImf2Bp8DtpoRAjunVn9f59vktpEr1k8aNAkLjZWFZGx6zWdGvXzncp56bB758ZkpvL1iPl1Af8b6CXLL4UCBI3O48dpU/HjsE9E5oGv3+m7xUdqQuHKT1x3HcpF6SXLL4UCBI3+yoqyMvO8V3G+cVy8HuCt2M0vk77ArKMxflovkJ8pfsuQJLXkapTdGrTxncZ59c4+L2sLPLNcqlHMo3tGEuWRL6lFiyIjKFMII2v05GqYK9+GquXTM7N67vWzJ42s0Nmltzn9EtBDeEwVfX1tGvVyncpFxaLwe9J0I7R+DpV1tZ6ruTiNF8hfnzvxswCJnuuQeKgIdp0kh5K81xJC0iCdozG81rXJ2iTl8TGRQPBzL5pZnnxeHDn3PvAkXjct/iVFt19qw83eK6kBSTBuguNo4vSz1F7gnePSDM05QghB5hrZovM7GEz6xrvok5nZjPMbIWZrShLwEPxVJUWCpGVnu59Fc0Wk+DtGI2vU5uzmvg0zDO1XDQQnHN/75wbCjwM9ADeM7N34l7Zfz/+TOdciXOupEuXLi31sBIDV2S1pjxJT9HoXYx32xtfpyvOWpE2CbpHpBma04dwCDgAlAOJ10gqLa5HdjZ7K074LiP5xGG3vfF16pGdfcbvk6B7RJqhKX0IXzWzhcB8oDPw18654ngXJomvsMMV7Dh61HcZyScOu+07jh4lzYy8nPZn/D4JukekGZpyhNAL+IZzbqhz7jHn3IZYPbiZPQ98CAw0s1Iz+0qs7lv8G9S5C2WnKimrrPRdSnKJw277+rJD9O3YiVZpnx8VluDdI9IMF52Y5px7NF4P7py7L173Lf6N6NYNgJX79zG5X3/P1SSRGM/OCjvH6gP7uaVvvxgVKInK9zwESWLDcrvROj2dJXu0zn5MhcMxnaq7vuwQx2uquSYvPwbFSSJTIEjctEpLY0x+AQt27cRpEHtsNKVDuZkjkN7duR0DrisojGWlkoAUCBJXN/ftz96KE6w5eMB3KcnhYh3KzRyB5Jzjza1buLJHTzoHed0paREKBImrW/r2JzMtnVc2rvddSnK4WIdyM0cgfXLoIFuPlHPHwMFxK1kznROHAiGFtcQHNSczk9v6D+C1TRuoqKlpscdNWhcbB9rMEUh/+GQ1rdPT+bMBg+JSrmY6JxYFQopqyQ/qA8NHUllXxwvr1+oLIhYuNA60GRMHDlWeZM6mTdw5eCg5mZlxKVUznROLAiFFteQHtbhrN67NK+B3H69gz/661PmC8HUo1MSJA0+tXE6DC/PQyJK4lRK0mc46Or0wBUKKaukP6jeuuZbDp07xxt6PA/UFETcBPxQqPXGc59au4c7BQ+nVoUPcHidIM50D/pIEgs6YlqJa+sxTo3vkcXOffvzLymW8/V9DSKvKTu4zXp3rEKxriy4UfEE/WfQeoZDxjavHxP2xGg9YfAv4SxIIOkJIYS29JMH3rx9PQ9jx+KJ3k38phKC1lZxm/s7tvLV9Kw+PvobuZy1ml8wC/JIEhgJBWkxB+w58/epreXv7NuZs3ui7nPgKUlvJaY5WVfGDd+cxoFNn/npU/PoOgiigL0mgKBCkRT00qoRR3brztwve4dPjx3yXE18BWxXOOcf35r/N0aoqfnHT5HMuZJfsAvaSBI4CQVpUeijEk7fcjmE88uZ/UV2fImdUC4Dfr1rB2zu28e0x1zM0V43nGnH0eQoEaXH57dvzy1tuZX3ZIb4zby5hfSLj7t2dO/jpB4u4td8AvjLySt/leKcRR+emQBAvbujdl++MvZ43tm7mZx+836zbas+uedYcPMDX3nqDIZ278PObJmNqL9GEufNQIIg3M0aN5kvFI5j58Qr+efnSJt1Ge3bNs+lwGX/1+it0at2G30/9c9pkZPguKRA04ujcNA9BvDEzHhs/iZO1tfziww9wDh4effUF92A1lrzpNpQd4i9nv0xmWjrPTptObtt2vksKjJaeh9NSLvdUGTpCEK9CZvzsxlu4c9AQfrn0A/5h0XsX7FPQnl3TfLS3lPteeZHM9HSeu+ueuM5GTlTJNuIoFkfPOkIQ79JCIX5202SyMzN5evVK9p+s4B9vmkzrczRvJOueXSy9vnkj3503l/z27Zk17S56Zuf4LklaQCyOnnWEIIEQMuP/jpvI968bz1vbtnD3S8+z5/jxc2+bZHt2sVIfDvPE4vf45tw3Gdm9Oy/dfa/CIIXE4ujZEunUhiUlJW7FihW+y5A4W7BrB9+c+yYAT9xwC5P79fdcUfDtqzjBN+e+yfJ9e/nisOH8cNzElJx4lurO14dgZiudcxedmq4jBAmciYV9mHPv/fRq34H/9eYcvvvOXE5ET64jZ3LO8frmjdz23LNsKDvEL2++lR9NvFFhkKIu9+hZfQgSSAXtO/DS3ffx62VL+NeVy1n86S7+bvwN3NS3n+/SAmNvxQn+buF85u/cwchu3fnFzbdS2OEK32VJAlOTkQTe6gP7eXT+22wpP8zEwj787bgJ9E7hL77q+jr+bdXH/DY6d+Ob14zlr0aMIi0Uuuxhh5KcmtpkpECQhFDX0MCsNR/zT8uWUt1Qz31FxTwy+hq6tG3ru7QW0xAO89rmjTy59AP2VVRwc99+/PD6ifTMiXQcNw47XLIk0qm4YEGkCUFEgSBJqayykl8tW8KL6z8hIy2NvygazkOjrqRbu+Rd17+uoYE5Wzbx2+XL2HnsKMNyu/Lo2HFcm19wxnYHD0bGoNfXR0aalJZq0p5EKBAkqe06dpTffLSU1zdvJGTGnw0YxAMjRjEsiVbxPFpVxYsbPuHZNavYf/Ikgzt34W+uupZb+vY752xu5yITkhqPELTmvzRKiEAws8nAr4E04PfOuScutL0CQc625/hxnl69kpc2rONUXR3DcrvyhaHDuL3/QNpnZfkur9kawmGW7t3DyxvW86dtW6htaODavHweGlXChF69L7ownfoQ5FwCHwhmlgZsAW4CSoHlwH3OuQ3nu40CQc7nRE0Nszet5/l1n7Cl/DCtQmmMLyzk1n4DmNS7DzmZwQ2HhnCYjw/s461tW/nT1i0cqDxJdqtMpg0azH1FxQzq3MV3iZLgmhoIPoedXgVsc87tADCzF4A7gPMGgsj55GRm8sDwUfxl8UjWlR3itU0beHPrFubt2E56KMSV3XtwfUEhY/MLGJrblXTPva37Kk6wtHQPiz7dzaLduzhSXUWrUBrjehXy/YHjubFPX7LStTKptCyfRwjTgcnOuYeil78EXO2ce+Ss7WYAMwAKCgqu3L17d4vXKokp7ByrD+znnR3beW/3TjYejix63zYjg+HdujO8azeKcrsyqHMXCnLakxankDhSdYrNhw+zvuwQaw8dYNX+/eytOAFAp9atub6gkEm9+zC+V2+yMzPjUoOktkRoMrobuOWsQLjKOfc357uNmozkcpSdqmRZ6R6W79vLqv372FR+mProkpCt0tIobN+BgvYdyMvJoVu7bLq0acsVrVvTPjOTtq1a0To9g4y0ECEznIN6F6amvp7Kujoqamo4Vl3F4VOnOFh5kn0VFew5foxdx45xpLrqsxq6t8tmRLfulPToyTU98xjYuQshNfZLnCVCk1EpkH/a5Txgn6daJAV0adOWKQMGMWXAIABq6uvZXH6YTYfL2H70CDuOHuHT48f4sPRTKusu/VzP6aEQXdu2Iz+nPTf37UefKzoysFNnBnfJpXObNrF6OiIx5zMQlgP9zaw3sBe4F/gLj/VIislMT6e4azeKu3b73HUVNTUcrjrF0aoqTtTUUFlbS1V9HXXhMGHnMCLLdmempdM2I4PszEw6ZGXRqU0bOrVuo73+S6RRUn55CwTnXL2ZPQLMJTLs9Gnn3Hpf9YicLjszk+zMzJReIqOlaaa1f14Xt3POvQm86bMGEQkGnR7VP+WviASCTo/qn5a/FpFA0OlR/dMRgog0WTgcWUQvXqPVdXpUvxQIItIkjZ2+eXmRRfSiUzgkiSgQRKRJztXpK8lFgSAiTaJO3+SnTmURaRJ1+iY/HSGINFO8O1aDTJ2+yU2BINIM6liVZKZAEGkGdaxKMlMgiDSDOlYlmalTWaQZ1LEqyUyBINJMjR2rIslGTUYiIgIoEEREJEqBICIigAJBRESiFAgiIgIoEEREJEqBICIigAJBRESiFAgiIgIoEEREJEqBICIigAJBRESiFAgiIgIoEEREJMpLIJjZ3Wa23szCZlbiowYRETmTryOEdcCdwPueHl9ERM7iJRCccxudc5t9PHaiC4fh4EFwznclIpJs1IeQQMJhmDgR8vJgwoTIZRGRWInbKTTN7B2g2zmu+oFz7vVm3M8MYAZAQUFBjKpLTGVlsGQJ1NdHfpaV6VSOIhI7cQsE59yNMbqfmcBMgJKSkpRuKMnNhTFjImEwZkzksohIrMQtECT2zGDBgsiRQW5u5LKISKz4Gnb652ZWClwL/NHM5vqoIxGFQpFmIoWBiMSalyME59xsYLaPxxYRkXPTKCMREQEUCCIiEqVAEBERQIEgIiJRCgQREQEUCCIiEqVAEBERQIEgIiJRCgQREQEUCCIiEqVAEBERQIEgIiJRCgQREQEUCCIiEqVAEBERQIEgIiJRCgQREQEUCCIiEqVAEBERQIEgIiJRCgQREQEUCCIiEqVAEBERQIEgIiJRCgQREQEUCCIiEqVAEBERQIEgIiJRXgLBzH5uZpvMbK2ZzTazDj7qEBGR/+brCGEeUOScKwa2AN/zVIeIiER5CQTn3NvOufroxaVAno86RETkv6X7LgD4MvCf57vSzGYAM6IXa8xsXYtU5Udn4LDvIuIomZ9fMj830PNLdAObspE55+Ly6Gb2DtDtHFf9wDn3enSbHwAlwJ2uCYWY2QrnXElsKw0OPb/ElczPDfT8El1Tn1/cjhCcczde6HozewCYAtzQlDAQEZH48tJkZGaTge8C451zp3zUICIiZ/I1yug3QDYwz8xWm9lTTbzdzDjWFAR6fokrmZ8b6PkluiY9v7j1IYiISGLRTGUREQEUCCIiEpVwgZDMy16Y2d1mtt7MwmaWNEPgzGyymW02s21m9qjvemLJzJ42s0PJOj/GzPLNbIGZbYy+N7/uu6ZYMbMsM/vIzNZEn9vf+64pHswszcxWmdkbF9s24QKB5F72Yh1wJ/C+70JixczSgH8GbgWGAPeZ2RC/VcXULGCy7yLiqB74P865wcA1wMNJ9PrVAJOcc8OBEcBkM7vGb0lx8XVgY1M2TLhASOZlL5xzG51zm33XEWNXAducczucc7XAC8AdnmuKGefc+8AR33XEi3Nuv3Pu4+j/K4h8sfT0W1VsuIiT0YsZ0X9JNcrGzPKA24HfN2X7hAuEs3wZ+JPvIuSCegJ7TrtcSpJ8oaQaMysERgLLPJcSM9HmlNXAIWCecy5pnlvUr4DvAOGmbByEtYw+pxnLXtQD/9GStV2upjy3JGPn+F1S7YWlAjNrB7wCfMM5d8J3PbHinGsARkT7ImebWZFzLin6g8xsCnDIObfSzCY05TaBDIRkXvbiYs8tCZUC+addzgP2eapFLoGZZRAJg/9wzr3qu554cM4dM7OFRPqDkiIQgLHAVDO7DcgCcszsD865+893g4RrMjpt2YupWvYiISwH+ptZbzNrBdwLzPFckzSRmRnwb8BG59wvfdcTS2bWpXGUopm1Bm4ENnktKoacc99zzuU55wqJfO7evVAYQAIGApe+7EXgmdmfm1kpcC3wRzOb67umyxUdAPAIMJdIh+SLzrn1fquKHTN7HvgQGGhmpWb2Fd81xdhY4EvApOjnbXV0jzMZdAcWmNlaIjsu85xzFx2amcy0dIWIiACJeYQgIiJxoEAQERFAgSAiIlEKBBERARQIIiISpUAQERFAgSAiIlEKBJHLYGaFp58Lwcy+ZWZ/57EkkUumQBAREUCBICIiUQoEkctTz5mfoyxfhYhcLgWCyOU5COSaWSczyySyLLtIQgrk+RBEEoVzrs7MHidyFrGdJNHyyZJ6tNqpiIgAajISEZEoBYKIiAAKBBERiVIgiIgIoEAQEZEoBYKIiAAKBBERifr/wGnen2iQmuMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "np.seterr(invalid='ignore', over='ignore')  \n",
    "# suppress warning caused by division by inf\n",
    "\n",
    "N =30\n",
    "np.random.seed(0)\n",
    "X = np.random.randn(2,N)\n",
    "y = np.sign( X[0,:]**2+ X[1,:]**2 -0.7)\n",
    "theta = 0.5\n",
    "c, s = np.cos(theta), np.sin( theta )\n",
    "X = np.array ([[c,-s],[s,c]])@X\n",
    "X = X + np.array([[1] ,[1]])\n",
    "\n",
    "def phi(X):\n",
    "    temp = np.ones( (5, X.shape[1]))\n",
    "    for idx in range(X.shape[1]):\n",
    "        temp[1, idx] = X[0, idx]\n",
    "        temp[2, idx] = X[0, idx]**2\n",
    "        temp[3, idx] = X[1, idx]\n",
    "        temp[4, idx] = X[1, idx]**2\n",
    "    return temp.T\n",
    "\n",
    "class LR():\n",
    "    def __init__(self, X, Y):\n",
    "        self.X = X\n",
    "        self.Y = Y\n",
    "        self.N = X.shape[0]\n",
    "        self.p = X.shape[1]\n",
    "        self.w = np.zeros(self.p)\n",
    "        self.learning_rate = 0.01\n",
    "\n",
    "    def loss(self):\n",
    "        temp = np.mean(np.log(1 + np.exp(-self.Y*(self.X@self.w))))\n",
    "        return temp\n",
    "\n",
    "    def train(self):\n",
    "        idx = np.random.randint(self.N)\n",
    "        self.w -= self.learning_rate * \n",
    "        ( (np.exp(-self.Y[idx]*(self.X[idx,:]@self.w)))\n",
    "         /(1 + np.exp(-self.Y[idx]*(self.X[idx,:]@self.w))) ) \n",
    "        * (-self.Y[idx]*self.X[idx,:])\n",
    "\n",
    "\n",
    "epochs = 10000\n",
    "\n",
    "model = LR(phi(X), y)\n",
    "loss_list = []\n",
    "\n",
    "for epoch in range(epochs):\n",
    "    model.train()\n",
    "    loss_list.append(model.loss())\n",
    "\n",
    "xx = np.linspace(-4, 4, 1024)\n",
    "yy = np.linspace(-4, 4, 1024)\n",
    "xx, yy = np.meshgrid(xx, yy)\n",
    "Z = model.w[0] + ( model.w[1]*xx + model.w[2]* (xx**2)) \n",
    "+ ( model.w[3]*yy + model.w[4]*(yy**2))\n",
    "\n",
    "plus_mask, minus_mask = [], []\n",
    "for idx in range(len(y)) :\n",
    "    if y[idx] > 0 :\n",
    "        plus_mask.append(idx)\n",
    "    else :\n",
    "        minus_mask.append(idx)\n",
    "\n",
    "plt.scatter(X[0,plus_mask], X[1,plus_mask], color = 'b', \n",
    "            s=5, label = 'plus')\n",
    "plt.scatter(X[0,minus_mask], X[1,minus_mask], color = 'r', \n",
    "            s=5, label = 'minus')\n",
    "\n",
    "plt.contour(xx, yy, Z, 0)\n",
    "plt.xlim((-2,4))\n",
    "plt.ylim((-2,4))\n",
    "plt.xlabel(\"u\")\n",
    "plt.ylabel(\"v\")\n",
    "plt.legend()\n",
    "plt.show()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e5a4ce95",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}