diff --git a/NM4_ex_solution.ipynb b/NM4_ex_solution.ipynb new file mode 100644 index 0000000000000000000000000000000000000000..ceb9f72e4d5882b5a0caa07fef5bf2db951c49df --- /dev/null +++ b/NM4_ex_solution.ipynb @@ -0,0 +1,140 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "broken-executive", + "metadata": {}, + "source": [ + "# Exercice : méthode de la transformée inverse\n", + "\n", + "But : générer des points selon une distribution de forme $f(x) = \\cos x$ dans l'intervalle $\\left[0.1, 1.5\\right]$\n", + "\n", + "Premièrtement, on ajoute une **constante de normalisation** $\\alpha$ de façon à ce que l'intégrale entre les bornes de notre distribution soit bien égale à $1$ (la probabilité de tirer un nombre entre les deux bornes. On a $f(x) = \\alpha \\cos x$\n", + "\n", + "On calcule $F(x) = \\int_a^x f(x') dx'$\n", + "\n", + "$$\n", + "\\int_a^x f(x') dx' = \\big[ \\alpha \\sin x \\big]_a^x = \\alpha ( \\sin x - \\sin a )\n", + "$$\n", + "\n", + "Notez que je n'utilise pas la constante $\\beta$ comme aux le cours, qui constitue une étape inutile. On fixe la valeur de $\\alpha$ pour avoir $F(b) = 1$:\n", + "\n", + "$$\n", + "F(b) = \\alpha ( \\sin b - \\sin a ) = 1 \\Rightarrow \\alpha = \\frac{1}{\\sin b - \\sin a}\n", + "$$\n", + "\n", + "La fonction $F(x)$ est donc après normalisation:\n", + "\n", + "$$\n", + "F(x) = \\frac{\\sin x - \\sin a}{\\sin b - \\sin a}\n", + "$$\n", + "\n", + "Et la fonction inverse:\n", + "\n", + "$$\n", + "T(u) = \\arcsin((\\sin b - \\sin a) u + \\sin(a))\n", + "$$" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "id": "french-candle", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "inf: 0.1 sup: 1.5000000000000018\n" + ] + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAD4CAYAAADy46FuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAUCUlEQVR4nO3df4zcdX7f8eerOEe4u8Dx41ghm9SkuEkA9dTD5WjSRts6Cr5LVFMJJKfN4Z5cWaXkeq0q9SB/lKoREkhVSVAKkXVQfjQ6cMkpuG0gRabbaxUwgeRyHFCKe6Tg4kIIhGAqCCbv/jEfW+O93c+ud3Z2dvHzIY125j3fz2feY+/otZ/vd+Y7qSokSZrPn5t0A5Kk1c2gkCR1GRSSpC6DQpLUZVBIkrrWTbqB5XbOOefUxo0bJ93GMe+++y6f+MQnJt3GotnveNnveNnv0j399NNvVNWn57rvIxcUGzdu5Kmnnpp0G8fMzMwwPT096TYWzX7Hy37Hy36XLsn/nu8+dz1JkroMCklSl0EhSeoyKCRJXQaFJKnLoJAkdRkUkqQug0KS1GVQSJK6PnKfzB7Vxuv/08hz/MHNP70MnUjS6uCKQpLUZVBIkroMCklSl0EhSeryYPYK8SC5pLXKFYUkqcugkCR1GRSSpC6DQpLU5cHsNc6D5JLGzRWFJKnLoJAkdRkUkqQug0KS1GVQSJK6FgyKJHcleT3Jd4ZqZyV5NMmL7eeZQ/fdkORAkheSXDFUvzTJM+2+25Kk1U9N8kCr70+ycWjMjvYYLybZsVxPWpK0eItZUdwNbJ1Vux7YV1WbgH3tNkkuArYDF7cxtyc5pY25A9gFbGqXo3PuBN6qqguBW4Fb2lxnATcCnwMuA24cDiRJ0spYMCiq6pvAm7PK24B72vV7gCuH6vdX1ftV9RJwALgsyXnA6VX1eFUVcO+sMUfnehDY0lYbVwCPVtWbVfUW8CjfG1iSpDFb6gfupqrqEEBVHUpybquvB54Y2u5gq33Qrs+uHx3zSpvrSJK3gbOH63OMOU6SXQxWK0xNTTEzM7PEp7U8hh//8OHDy9bPuJ7XXP3+vUfeHXneu7d+YuQ5FrKc/74rwX7Hy37HY7k/mZ05atWpL3XM8cWq3cBugM2bN9f09PSCjc7rkdE/6Tz8+DMzM4PbyzzvMWut32V2rN81wn7Hy37HY6nvenqt7U6i/Xy91Q8C5w9ttwF4tdU3zFE/bkySdcAZDHZ1zTeXJGkFLTUo9gJH34W0A3hoqL69vZPpAgYHrZ9su6neSXJ5O/5wzawxR+e6CnisHcf4LeCnkpzZDmL/VKtJklbQgrueknwdmAbOSXKQwTuRbgb2JNkJvAxcDVBVzybZAzwHHAGuq6oP21TXMngH1WnAw+0CcCdwX5IDDFYS29tcbyb5ReB32nb/sqpmH1SXJI3ZgkFRVT87z11b5tn+JuCmOepPAZfMUX+PFjRz3HcXcNdCPUqSxsdPZkuSugwKSVKXQSFJ6jIoJEldfhWqVtRyfHUr+PWt0kpyRSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXQaFJKnLcz3pI8FzSEnj44pCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklS10hBkeSfJHk2yXeSfD3J9yc5K8mjSV5sP88c2v6GJAeSvJDkiqH6pUmeaffdliStfmqSB1p9f5KNo/QrSTpxSw6KJOuBfwRsrqpLgFOA7cD1wL6q2gTsa7dJclG7/2JgK3B7klPadHcAu4BN7bK11XcCb1XVhcCtwC1L7VeStDSjnj12HXBakg+AjwOvAjcA0+3+e4AZ4KvANuD+qnofeCnJAeCyJH8AnF5VjwMkuRe4Eni4jfkXba4HgV9JkqqqEfuWFsWz0kojrCiq6v8A/wp4GTgEvF1V/xmYqqpDbZtDwLltyHrglaEpDrba+nZ9dv24MVV1BHgbOHupPUuSTtySVxTt2MM24ALgj4F/n+TnekPmqFWn3hszu5ddDHZdMTU1xczMTKeN8Rt+/MOHDy9bP+N6Xmut3+G5l7Pf4XmX27j6HTf7Ha+10u8ou55+Enipqv4QIMk3gB8DXktyXlUdSnIe8Hrb/iBw/tD4DQx2VR1s12fXh8ccTLIOOAN4c3YjVbUb2A2wefPmmp6eXvqzemT0XQ3Djz8zMzO4vczzHnMS9js893L2OzzvMePqd42w3/FaK/2O8q6nl4HLk3y8vUtpC/A8sBfY0bbZATzUru8Ftrd3Ml3A4KD1k2331DtJLm/zXDNrzNG5rgIe8/iEJK2sJa8oqmp/kgeB3wWOAL/H4K/6TwJ7kuxkECZXt+2fTbIHeK5tf11Vfdimuxa4GziNwUHsh1v9TuC+duD7TQbvmpIkraCR3vVUVTcCN84qv89gdTHX9jcBN81Rfwq4ZI76e7SgkSRNhp/MliR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXaOePVbSEnhWWq0lrigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl99HIX2E+D0XGgdXFJKkLoNCktRlUEiSugwKSVKXQSFJ6hopKJJ8KsmDSf5HkueT/NUkZyV5NMmL7eeZQ9vfkORAkheSXDFUvzTJM+2+25Kk1U9N8kCr70+ycZR+JUknbtQVxS8Dj1TVjwCfAZ4Hrgf2VdUmYF+7TZKLgO3AxcBW4PYkp7R57gB2AZvaZWur7wTeqqoLgVuBW0bsV5J0gpYcFElOB34CuBOgqv60qv4Y2Abc0za7B7iyXd8G3F9V71fVS8AB4LIk5wGnV9XjVVXAvbPGHJ3rQWDL0dWGJGlljPKBux8C/hD4t0k+AzwNfAWYqqpDAFV1KMm5bfv1wBND4w+22gft+uz60TGvtLmOJHkbOBt4Y7iRJLsYrEiYmppiZmZmhKc1uuHHP3z48LL1M67ntdb6HZ57Ofsdnne52e/KsN/xGCUo1gGfBb5cVfuT/DJtN9M85loJVKfeG3N8oWo3sBtg8+bNNT093WljAY+M/snW4cefmZkZ3F7meY85Cfsdnns5+x2e9xj7nfv/cZWy3/EY5RjFQeBgVe1vtx9kEByvtd1JtJ+vD21//tD4DcCrrb5hjvpxY5KsA84A3hyhZ0nSCVpyUFTV/wVeSfLDrbQFeA7YC+xotR3AQ+36XmB7eyfTBQwOWj/ZdlO9k+Tydvzhmlljjs51FfBYO44hSVoho54U8MvAryX5GPBd4EsMwmdPkp3Ay8DVAFX1bJI9DMLkCHBdVX3Y5rkWuBs4DXi4XWBwoPy+JAcYrCS2j9ivJOkEjRQUVfUtYPMcd22ZZ/ubgJvmqD8FXDJH/T1a0EiSJsNPZkuSugwKSVKXQSFJ6vIb7iQtyG/OO7m5opAkdRkUkqQug0KS1GVQSJK6DApJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1+Q13kibGb85bG1xRSJK6DApJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkrpGDIskpSX4vyX9st89K8miSF9vPM4e2vSHJgSQvJLliqH5pkmfafbclSaufmuSBVt+fZOOo/UqSTsxyrCi+Ajw/dPt6YF9VbQL2tdskuQjYDlwMbAVuT3JKG3MHsAvY1C5bW30n8FZVXQjcCtyyDP1Kkk7ASEGRZAPw08DXhsrbgHva9XuAK4fq91fV+1X1EnAAuCzJecDpVfV4VRVw76wxR+d6ENhydLUhSVoZo57C45eAfwb8wFBtqqoOAVTVoSTntvp64Imh7Q622gft+uz60TGvtLmOJHkbOBt4Y7iJJLsYrEiYmppiZmZmxKc1muHHP3z48LL1M67ntdb6HZ57Ofsdnne52e/x8y63cfU7bmul3yUHRZKfAV6vqqeTTC9myBy16tR7Y44vVO0GdgNs3ry5pqcX0848Hhn93DPDjz8zMzO4vczzHnMS9js893L2OzzvMfa7NvtdI9ZKv6OsKH4c+FtJvgB8P3B6kn8HvJbkvLaaOA94vW1/EDh/aPwG4NVW3zBHfXjMwSTrgDOAN0foWZJ0gpZ8jKKqbqiqDVW1kcFB6seq6ueAvcCOttkO4KF2fS+wvb2T6QIGB62fbLup3klyeTv+cM2sMUfnuqo9xvesKCRJ4zOO04zfDOxJshN4GbgaoKqeTbIHeA44AlxXVR+2MdcCdwOnAQ+3C8CdwH1JDjBYSWwfQ7+SpI5lCYqqmgFm2vU/ArbMs91NwE1z1J8CLpmj/h4taCRJk+EnsyVJXQaFJKnLoJAkdfmd2ZI+cvwu7uXlikKS1GVQSJK6DApJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkLoNCktRlUEiSugwKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl9+ZLUmLdLJ+F7crCklSl0EhSeoyKCRJXQaFJKlryUGR5Pwk/yXJ80meTfKVVj8ryaNJXmw/zxwac0OSA0leSHLFUP3SJM+0+25LklY/NckDrb4/ycalP1VJ0lKMsqI4AvzTqvpR4HLguiQXAdcD+6pqE7Cv3abdtx24GNgK3J7klDbXHcAuYFO7bG31ncBbVXUhcCtwywj9SpKWYMlBUVWHqup32/V3gOeB9cA24J622T3Ale36NuD+qnq/ql4CDgCXJTkPOL2qHq+qAu6dNeboXA8CW46uNiRJK2NZPkfRdgn9ZWA/MFVVh2AQJknObZutB54YGnaw1T5o12fXj455pc11JMnbwNnAG7MefxeDFQlTU1PMzMwsx9NasuHHP3z48LL1M67ntdb6HZ57Ofsdnne52e/x8y43+x2vkYMiySeBXwf+cVX9SecP/rnuqE69N+b4QtVuYDfA5s2ba3p6eoGuOx4Z/QM1w48/MzMzuL3M8x5zEvY7PPdy9js87zH2a7/DxtXvKjfSu56SfB+DkPi1qvpGK7/WdifRfr7e6geB84eGbwBebfUNc9SPG5NkHXAG8OYoPUuSTswo73oKcCfwfFX966G79gI72vUdwEND9e3tnUwXMDho/WTbTfVOksvbnNfMGnN0rquAx9pxDEnSChll19OPA18EnknyrVb7BeBmYE+SncDLwNUAVfVskj3AcwzeMXVdVX3Yxl0L3A2cBjzcLjAIovuSHGCwktg+Qr+SpCVYclBU1X9n7mMIAFvmGXMTcNMc9aeAS+aov0cLGknSZPjJbElSl0EhSeoyKCRJXQaFJKnLoJAkdRkUkqQuvzNbkiZstX8XtysKSVKXQSFJ6jIoJEldBoUkqcugkCR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXQaFJKnLoJAkdRkUkqQug0KS1GVQSJK6DApJUpdBIUnqMigkSV0GhSSpy6CQJHUZFJKkrjURFEm2JnkhyYEk10+6H0k6maz6oEhyCvBvgM8DFwE/m+SiyXYlSSePVR8UwGXAgar6blX9KXA/sG3CPUnSSSNVNekeupJcBWytqr/fbn8R+FxV/fzQNruAXe3mDwMvrHij8zsHeGPSTZwA+x0v+x0v+126P19Vn57rjnUr3ckSZI7acelWVbuB3SvTzolJ8lRVbZ50H4tlv+Nlv+Nlv+OxFnY9HQTOH7q9AXh1Qr1I0klnLQTF7wCbklyQ5GPAdmDvhHuSpJPGqt/1VFVHkvw88FvAKcBdVfXshNs6Eatyl1iH/Y6X/Y6X/Y7Bqj+YLUmarLWw60mSNEEGhSSpy6BYBgudYiTJ303y7Xb57SSfmUSfQ/0s6pQoSf5Kkg/bZ1kmZjH9JplO8q0kzyb5ryvd46xeFvp9OCPJf0jy+63fL02iz6F+7kryepLvzHN/ktzWns+3k3x2pXuc1c9C/a6211u336HtVsXrbU5V5WWEC4MD7P8L+CHgY8DvAxfN2ubHgDPb9c8D+1dzv0PbPQb8JnDVau4X+BTwHPCD7fa5q7zfXwBuadc/DbwJfGyCPf8E8FngO/Pc/wXgYQafabp8kr+/i+x31bzeFtPv0O/NxF9v811cUYxuwVOMVNVvV9Vb7eYTDD4LMimLPSXKl4FfB15fyebmsJh+/w7wjap6GaCqJtnzYvot4AeSBPgkg6A4srJtDjVT9c3Ww3y2AffWwBPAp5KctzLdfa+F+l1lr7fF/PvC6nm9zcmgGN164JWh2wdbbT47Gfx1NikL9ptkPfC3gV9dwb7ms5h/378InJlkJsnTSa5Zse6+12L6/RXgRxl8cPQZ4CtV9Wcr096SnOjv+Goy6dfbglbZ621Oq/5zFGvAgqcYObZh8jcY/OL+tbF21LeYfn8J+GpVfTj4o3eiFtPvOuBSYAtwGvB4kieq6n+Ou7k5LKbfK4BvAX8T+AvAo0n+W1X9ybibW6JF/46vJqvk9bYYq+n1NieDYnSLOsVIkr8EfA34fFX90Qr1NpfF9LsZuL/90p4DfCHJkar6jZVp8TiL6fcg8EZVvQu8m+SbwGeASQTFYvr9EnBzDXZOH0jyEvAjwJMr0+IJW3On0VlFr7fFWE2vtzm562l0C55iJMkPAt8Avjihv3KHLdhvVV1QVRuraiPwIPAPJ/hLu5hTuDwE/PUk65J8HPgc8PwK93nUYvp9mcHqhyRTDM54/N0V7fLE7AWuae9+uhx4u6oOTbqp+ayy19uCVtnrbU6uKEZU85xiJMk/aPf/KvDPgbOB29tfDUdqQmeMXGS/q8Zi+q2q55M8Anwb+DPga1XVfSviJPsFfhG4O8kzDHbrfLWqJnaq6SRfB6aBc5IcBG4Evg+O9fubDN75dAD4fwxWRBOziH5XzesNFtXvqucpPCRJXe56kiR1GRSSpC6DQpLUZVBIkroMCklSl0EhSeoyKCRJXf8fYcOeX7mQY9cAAAAASUVORK5CYII=\n", + "text/plain": [ + "<Figure size 432x288 with 1 Axes>" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "import math\n", + "import random\n", + "import matplotlib.pyplot as plt\n", + "def fm1(y, inf, sup):\n", + " return math.asin((math.sin(sup) - math.sin(inf)) * y + math.sin(inf))\n", + "\n", + "inf, sup = 0.1, 1.5\n", + "print(\"inf:\", fm1(0, inf, sup), \"sup:\", fm1(1, inf, sup)) # should be inf, sup\n", + "uni = [random.random() for i in range(1000000)]\n", + "new = [fm1(x, inf, sup) for x in uni]\n", + "plt.hist(new, 14, rwidth=0.8, zorder=2)\n", + "plt.grid()\n", + "plt.show()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "suspected-namibia", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.10" + }, + "varInspector": { + "cols": { + "lenName": 16, + "lenType": 16, + "lenVar": 40 + }, + "kernels_config": { + "python": { + "delete_cmd_postfix": "", + "delete_cmd_prefix": "del ", + "library": "var_list.py", + "varRefreshCmd": "print(var_dic_list())" + }, + "r": { + "delete_cmd_postfix": ") ", + "delete_cmd_prefix": "rm(", + "library": "var_list.r", + "varRefreshCmd": "cat(var_dic_list()) " + } + }, + "types_to_exclude": [ + "module", + "function", + "builtin_function_or_method", + "instance", + "_Feature" + ], + "window_display": false + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}