{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "$\\newcommand{\\si}{\\sigma}\n", "\\newcommand{\\al}{\\alpha}\n", "\\newcommand{\\tta}{\\theta}\n", "\\newcommand{\\Tta}{\\Theta}\n", "\\newcommand{\\Si}{\\Sigma}\n", "\\newcommand{\\ld}{\\ldots}\n", "\\newcommand{\\cd}{\\cdots}\n", "\\newcommand{\\Ga}{\\Gamma} \n", "\\newcommand{\\bet}{\\beta}\n", "\\newcommand{\\cU}{\\mathcal{U}}\n", "\\newcommand{\\cN}{\\mathcal{N}}\n", "\\newcommand{\\R}{\\mathbb{R}}\n", "\\newcommand{\\p}{\\mathbb{P}}\n", "\\newcommand{\\f}{\\frac}\n", "\\newcommand{\\ff}{\\frac{1}}\n", "\\newcommand{\\ds}{\\displaystyle}\n", "\\newcommand{\\bE}{\\mathbf{E}}\n", "\\newcommand{\\E}{\\mathbb{E}}\n", "\\newcommand{\\bF}{\\mathbf{F}}\n", "\\newcommand{\\ii}{\\mathrm{i}}\n", "\\newcommand{\\me}{\\mathrm{e}}\n", "\\newcommand{\\hsi}{\\hat{\\sigma}}\n", "\\newcommand{\\hmu}{\\hat{\\mu}}\n", "\\newcommand{\\ste}{\\, ;\\, }\n", "\\newcommand{\\op}{\\operatorname} \n", "\\newcommand{\\argmax}{\\op{argmax}}\n", "\\newcommand{\\lfl}{\\lfloor}\n", "\\newcommand{\\ri}{\\right}\n", "\\newcommand{\\supp}{\\operatorname{supp}}$\n", "\n", "\"Open\n", "\n", "# Echantillonnage préférentiel (*importance sampling*)\n", "\n", "\n", "## Exercice 1\n", "\n", "On souhaite calculer numériquement l'intégrale \n", "$$\\int_0^{10} e^{-2|x-5|}dx.$$\n", "Une première manière de faire est d'écrire cette intégrale comme \n", "$$10\\E[g(X)]$$\n", "avec $g(x) = e^{-2|x-5|}$ et avec $X$ de loi uniforme sur $[0,10]$. \n", "\n", "1. Utilisez la méthode de Monte-Carlo pour approcher l'intégrale précédente : tirez $N$ réalisations indépendantes $x_1\\dots x_N$ d'une loi uniforme sur $[0,10]$ et l'approximation est donnée par \n", "$$10 \\frac 1 N \\sum_{i=1}^N e^{-2|x_i-5|}.$$" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3X10XPV95/H3d2b0LEuyLMmWbMs2tgx+AAIISEIgEEgKCYE0JSkkbWlOG9I2NFmSbpfubtIesg9tsts0PWWTcGieaAohpLtxgZY2hEBCArEcwGCDbdlYlixZGj1Yz9JoZn77x8zIsixbY3mkOzP38zpHZ2bu3Ln365H1md/87u93rznnEBGR/BLwugAREck8hbuISB5SuIuI5CGFu4hIHlK4i4jkIYW7iEgeUriLiOQhhbuISB5SuIuI5KGQVzuuqalx69ev92r3IiI5adeuXb3Oudr51vMs3NevX09LS4tXuxcRyUlm1pbOeuqWERHJQwp3EZE8pHAXEclDCncRkTykcBcRyUPzhruZfcPMeszstdM8b2b2t2bWama7zezSzJcpIiJnI52W+7eAG8/w/E1AU/LnLuCr516WiIici3nD3Tn3HNB/hlVuBb7jEl4AqsysPlMFiqTLOcejO9uZmIp5XYqI5zLR574aaJ/xuCO57BRmdpeZtZhZSzgczsCuRU549eggf/qD3fzo9W6vSxHxXCbC3eZYNudVt51zDzjnmp1zzbW1886eFTkr3UOTAPQkb0X8LBPh3gGsnfF4DdCZge2KnJXwcCLUwyMKd5FMhPsO4HeSo2beCgw657oysF2RszId7sMKd5F5TxxmZg8D1wI1ZtYB/DlQAOCc+xrwJPBeoBUYAz62WMWKnEl4ZCJxq3AXmT/cnXN3zPO8Az6ZsYpEFkgtd5ETNENV8ob63EVOULhL3kiFet/IJLH4nAO2RHxD4S55wTlHeHiSkoIgcQf9oxGvSxLxlMJd8sLIZJSJqTgX1C8D1O8uonCXvJAK8631FYnH6ncXn1O4S16YDveGipMei/iVwl3yQqqlvqVe4S4CCnfJE6kwX7+ijLLCoMJdfE/hLnkhPDxJKGBUlRRQu6xIfe7iewp3yQvh4UlqyosIBCwR7sMTXpck4imFu+SF8MgktcuKAJLhrpa7+JvCXfJCeHhGuJcr3EUU7pIXekcmqS0/0XIfmogyGdXl9sS/FO6S8+JxR+9I5KRuGYDeEZ2CQPxL4S45b2AsQizuTgl3dc2InyncJeelhj2e6HMvTixXuIuPKdwl56VCXC13kRMU7pLzpsM9eUB1RXnhSctF/EjhLjlvdsu9IBiguqxw+pqqIn6kcJecl7pIR1nRiUsC15QXquUuvqZwl5w3c3Zqimapit8p3CXnzZydmlJbrpOHib8p3CXnhYdPzE5NSbXcndOFssWfFO6S807XLTMxFWdkMupRVSLeUrhLTpuMxjg+NjVnuIOGQ4p/Kdwlp/Ulzx9zap+7ZqmKvyncJafNnsCUMt1y10FV8SmFu+S07qHERKW6ipPDvS4Z7scGNZFJ/EnhLjmtKxne9ZUlJy2vKi2guCCgcBffUrhLTuscHKcwGGBFWeFJy82MhsqS6fAX8Zu0wt3MbjSzfWbWamb3zvF8o5k9Y2YvmdluM3tv5ksVOVXX8QlWViYujD3bqspiOgfHPahKxHvzhruZBYH7gZuArcAdZrZ11mr/FXjUOXcJcDvwfzJdqMhcjg1OnNIlk1JfWaJuGfGtdFruVwCtzrlDzrkI8Ahw66x1HFCRvF8JdGauRJHT6xwcp6GyeM7nGqqK6R6aIBqLL3FVIt5LJ9xXA+0zHnckl830F8BvmVkH8CTwxxmpTuQM4nFH99AE9VVzt9xXVRYTd9Cjse7iQ+mE+6mdmYmW+kx3AN9yzq0B3gs8ZGanbNvM7jKzFjNrCYfDZ1+tyAy9o5NMxRz1p2u5J7trdFBV/CidcO8A1s54vIZTu11+D3gUwDn3C6AYqJm9IefcA865Zudcc21t7cIqFknqOj73MMiU+qpE6HfpoKr4UDrhvhNoMrMNZlZI4oDpjlnrHAGuBzCzLSTCXU1zWVSp0D5dyz0V+qkPARE/mTfcnXNR4G7gKeB1EqNi9pjZfWZ2S3K1zwIfN7NXgIeB33U616osss5kaDecps+9ojhEaWFQwyHFl0LzrwLOuSdJHCiduezzM+7vBa7KbGkiZ3ZsaIKiUIDlpQVzPm9m1FcWazik+JJmqErO6jw+Tn1lMWZzHfNPaKgqoVPhLj6kcJec1XWGCUwpqyqK6TqubhnxH4W75Kyu4+PTI2JOp76qhPDIJJGoJjKJvyjcJSfF4o7u4cnTjpRJaagsxjnoGVbXjPiLwl1yUnh4kljczd8tU5ka665wF39RuEtOSg1vbJinWyY1TLJT/e7iMwp3yUnzzU5NqVfLXXxK4S45ab7ZqSnLigtYVhTSWHfxHYW75KSuwQlKCoJUlsw9gWmmVZXF6pYR31G4S07qGkwMgzzTBKaU+ipdbk/8R+EuOanz+MT0KX3n01BZrDNDiu8o3CUnHRucmB7mOJ9VlcX0jkSYjMYWuSqR7KFwl5wTjcXpGZ447eX1Zku18LsHdUUm8Q+Fu+Sc7uFJ4o7TXl5vttQpCnTqX/EThbvknGPJkE63WyY1Fl7DIcVPFO6Sc46mLtKR5gHV1Fj4oxoOKT6icJec094/BsDa6vTCvawoxIqywunXifiBwl1yzpG+MWqXFVFamNaFxABoXFHKEYW7+IjCXXJOW/8ojdWlZ/WaxupS2voU7uIfCnfJOe3942cd7uuqS+kaHNdFO8Q3FO6SUyajMToHzz7c11aXEnc6qCr+oXCXnHJ0YBznOPuW+4oyAPW7i28o3CWntCXDed2Ks+9zBzjSN5rxmkSykcJdckpqOOPZttzrlhVRGAqo5S6+oXCXnHKkb4ziggC1y4rO6nWBgNFYreGQ4h8Kd8kpbf1jNFaXpnUe99k0HFL8ROEuOaU9Ge4L0VhdSnv/GM65DFclkn0U7pIznHMc6R+jsbpsQa9vrC5lNBKjbzSS4cpEso/CXXJG70iEsUiMxjTPKTNbaoSN+t3FDxTukjOO9CeGMabGrJ+tE8MhFe6S/xTukjOOTJ8NcmF97qnXqeUufpBWuJvZjWa2z8xazeze06zzYTPba2Z7zOwfM1umCLT1jWEGa5YvrFumuCDIyooijZgRX5j3nKlmFgTuB94NdAA7zWyHc27vjHWagD8DrnLODZhZ3WIVLP51pH+MVRXFFBcEF7yNddVlOq+7+EI6LfcrgFbn3CHnXAR4BLh11jofB+53zg0AOOd6MlumSGIY5EK7ZFLWaiKT+EQ64b4aaJ/xuCO5bKbNwGYze97MXjCzG+fakJndZWYtZtYSDocXVrH4VlvfGOvOMdzXrSjl2NAEE1OxDFUlkp3SCfe5pgLOngUSApqAa4E7gAfNrOqUFzn3gHOu2TnXXFtbe7a1io+NR2L0DE8ueAJTSur1HQNqvUt+SyfcO4C1Mx6vATrnWOeHzrkp59ybwD4SYS+SEe3JMG48y7NBzpZ6vQ6qSr5LJ9x3Ak1mtsHMCoHbgR2z1vl/wHUAZlZDopvmUCYLFX9LjU3PVMtd/e6S7+YNd+dcFLgbeAp4HXjUObfHzO4zs1uSqz0F9JnZXuAZ4D865/oWq2jxn8N95zaBKWVFWSHlRSEO9+q87pLf0rp8vHPuSeDJWcs+P+O+Az6T/BHJuNaeEarLCqkuKzyn7ZgZG2vLaA2PZKgykeykGaqSE1p7RthUW56RbW2sK6e1R+Eu+U3hLlnPOceBnhE2rcxMuDfVLaN7aJKhiamMbE8kGyncJev1jkQYHJ/KWMt9U11iO2q9Sz5TuEvWS4VwKpTPlcJd/EDhLlkvdfCzKUPdMmuXl1AYCnBQ4S55TOEuWa+1e5jyohCrKoozsr1QMMB5NWUcULhLHlO4S9ZrDY+wsbZsQRfFPh2NmJF8p3CXrNfaM8KmumUZ3WZTXTntA2M6gZjkLYW7ZLWhiSm6hyYzdjA1ZVNdOc7BQU1mkjylcJeslumRMikaMSP5TuEuWS0Vvk0ZDvcNNWUEDI2YkbylcJes1tozQmEocM5XYJqtKBRk3QqNmJH8pXCXrNbaM8J5NWUEA5kbKZOysVYjZiR/KdwlqyVGymS2SyalaWU5h/tGicbii7J9ES8p3CVrTUzFaB8YW7Rw31RbzlTM0aYLd0geUrhL1joYHsG5zI+USUlt90C3umYk/yjcJWudGCmT2QlMKRuT4a6x7pKPFO6StfZ3DxMMGOtrMjtSJqW8KMTqqhL2HRtelO2LeEnhLllrT+cQTXXlFIWCi7aPLfUV7OkcXLTti3hF4S5Za0/nENsaKhd1H9tXV3Cod5SxSHRR9yOy1BTukpV6hiYID0+yfXXFou5nW0MlzsHrXUOLuh+RpaZwl6y0pzMRtkvRcp+5P5F8oXCXrPTa0UQ/+NaGxW25r6ooprqscHp/IvlC4S5Z6bXOQTbUlFFeFFrU/ZgZ2xoqeO2oWu6SXxTukpUSB1MXt9Wesn11JQd6hpmM6sIdkj8U7pJ1jo9F6BgYZ/vqxe1vT9nWUMFUzGmmquQVhbtknb3TB1OXqOWePGir8e6STxTuknVeS4bsYo+USWmsLqW8KKR+d8krCnfJOns6h2ioTIxiWQqBgLG1QTNVJb8o3CXrvHZ0kG1L1N+esq2hgr1dQ8Tibkn3K7JY0gp3M7vRzPaZWauZ3XuG9W4zM2dmzZkrUfxkdDLKod7R6X7wpbK9oZKJqTiHdIZIyRPzhruZBYH7gZuArcAdZrZ1jvWWAZ8CXsx0keIfbxwbwrmlO5iask0zVSXPpNNyvwJodc4dcs5FgEeAW+dY7wvAF4GJDNYnPpM6qLlUwyBTNtWWUxQKaKaq5I10wn010D7jcUdy2TQzuwRY65x7PIO1iQ+90nGcmvIiVlYULel+Q8EAW+or2N2hcJf8kE64z3XZ+emjTmYWAL4MfHbeDZndZWYtZtYSDofTr1J8Y1fbAJetq8Jsrv92i+uydct5peO4ZqpKXkgn3DuAtTMerwE6ZzxeBmwHfmJmh4G3AjvmOqjqnHvAOdfsnGuura1deNWSl3qGJ2jrG+Py9dWe7P/y9cuZjMY13l3yQjrhvhNoMrMNZlYI3A7sSD3pnBt0ztU459Y759YDLwC3OOdaFqViyVu7Dg8AiRa0Fy5bl/hQ2dXW78n+RTJp3nB3zkWBu4GngNeBR51ze8zsPjO7ZbELFP/YeXiAolBgyWamzla7rIj1K0rZmfyQEcllaZ1P1Tn3JPDkrGWfP8261557WeJHu9r6ecvaKgpD3s2tu2xdNc/s68E550m/v0imaIaqZIWxSJTXOodoXu9Nl0zK5euX0z8a4VDvqKd1iJwrhbtkhZfbjxOLO5o9Opiakvpw2aWuGclxCnfJCrsOD2AGlzZ623LfWFvO8tICWnRQVXKcwl2yws62ATbXLaOypMDTOsyMy9Ytp0Utd8lxCnfxXCzueKltwPP+9pTm9dUc6h2lb2TS61JEFkzhLp7bd2yY4clo9oR7cpx9S5ta75K7FO7iudSkoeZ13h5MTdm+upLCYIBdCnfJYQp38dyLb/azsqKINctLvC4FgOKCIBetqeTFQ31elyKyYAp38VQs7vhZay/v2FSbVZOG3tFUw+6jgwyMRrwuRWRBFO7iqVePDnJ8bIprNtd4XcpJrtlci3Pws9Zer0sRWRCFu3jq2X1hzODqpuw6S+jFa6qoLCng2f06NbXkJoW7eOq5A2EuWl1JdVmh16WcJBgw3tFUw08PhHFOF82W3KNwF88Mjk3x0pEBrtmcXa32lHc21dI9NMm+7mGvSxE5awp38czzB3uJO7I23K9OHgd4dp+6ZiT3KNzFM8/tD7OsOMQla6u8LmVO9ZUlnL9yGc8dULhL7lG4iyecczy3P8xVG2sIBbP3v+E1m2vY+eYAY5Go16WInJXs/auSvNbaM0Ln4ETWdsmkXLO5lkgszgua0CQ5RuEunkgNMcy28e2zXb6+muKCAM/t13h3yS0Kd/HET/aF2VhbxprlpV6XckbFBUHedt4KfvxGj4ZESk5RuMuS6x+N8ItDfdy4fZXXpaTl17at4kj/GHs6h7wuRSRtCndZck/tOUYs7njvhfVel5KW92xbRTBgPPlql9eliKRN4S5L7ondXWyoKWNrfYXXpaSluqyQt29cwROvdqlrRnKGwl2WVN/IJL841Md7L1yVVWeBnM/7LqynrU9dM5I7FO6ypJ7a000s7njfhQ1el3JWfi3ZNfOEumYkRyjcZUk98WonG2rK2FK/zOtSzsryVNfMbnXNSG5QuMuS6RuZ5BcH+3jfhfU51SWTcvNF9Ro1IzlD4S5L5l/3HCPuyJlRMrO9Z2uia+bx3eqakeyncJcl8/grXTnZJZMy3TXzaifxuLpmJLsp3GVJHO4d5ReH+vjgJatzsksm5TcuXUN7/zi/0LlmJMsp3GVJPLzzCMGA8eHL13pdyjm5cfsqqkoL+McXj3hdisgZKdxl0UWicR5r6eD6C+pYWVHsdTnnpLggyG9cuoan9hwjPDzpdTkip5VWuJvZjWa2z8xazezeOZ7/jJntNbPdZva0ma3LfKmSq/5t7zH6RiN85MpGr0vJiDuuaCQadzy2q8PrUkROa95wN7MgcD9wE7AVuMPMts5a7SWg2Tl3EfAY8MVMFyq56+FfHmF1VQnXNGX3udvTtamunCs2VPPIziM6sCpZK52W+xVAq3PukHMuAjwC3DpzBefcM865seTDF4A1mS1TctXh3lGeb+3jjivWEgjk7oHU2T5yRSNtfWP8/KAOrEp2SifcVwPtMx53JJedzu8B/zLXE2Z2l5m1mFlLOKzrUvrB9IHU5tw+kDpb6sDqw7/UgVXJTumE+1zNrTm/i5rZbwHNwJfmet4594Bzrtk511xbmx9f0eX0RiejPLqznRu21FGX4wdSZysuCHJb8sBqx8DY/C8QWWLphHsHMLPZtQbonL2Smd0A/BfgFuechhEI//jiEQbGpvjEOzd6Xcqi+Ng7NgDwwHOHPK5E5FTphPtOoMnMNphZIXA7sGPmCmZ2CfB1EsHek/kyJddMTMV44KeHePvGFVzauNzrchbF6qoSPnjpah7Z2U7P8ITX5YicZN5wd85FgbuBp4DXgUedc3vM7D4zuyW52peAcuD7Zvayme04zebEJ76/q4Pw8CR3X7fJ61IW1R9eu4loLM7f//RNr0sROUkonZWcc08CT85a9vkZ92/IcF2Sw6Zicb7+7EEuaazibRtXeF3OotpQU8b7LmrgH15o4w+v3UhVaaHXJYkAmqEqi2DHy510DIxz93Wbcvo8Mun65HUbGY3E+NbPD3tdisg0hbtk1FQszv0/aWVLfQXvuqDO63KWxAWrKrhhy0q++fxhBsenvC5HBFC4S4Z994U2DoVH+cy7N/ui1Z5yz7ubGJqY4m+fPuB1KSKAwl0yaGA0wpd/dICrm2q4YYs/Wu0p2xoquf3ytXz754dp7RnxuhwRhbtkzpd/tJ+RySifu3mrr1rtKZ99z/mUFAT570/s9boUEYW7ZMYbx4b4hxfa+OiVjWxemZtXWjpXNeVFfOr6Jp7ZF+aZNzTdQ7ylcJdz5pzjC4/vZVlxAffcsNnrcjx159vXs6GmjC88sZdINO51OeJjCnc5Z4/sbOf51j7+5D2bWV7m73HehaEAn3//Vg6FR/mbH+33uhzxMYW7nJO2vlG+8Phertq0go9eqWu0AFx3fh0fbl7D1549yK62fq/LEZ9SuMuCxeKOzz76CsGA8aXbLs6r87Wfq8/dvJWGqhLu+d4rjE5GvS5HfEjhLgv29ecO0tI2wH23bqOhqsTrcrLKsuIC/veHLqZ9YIz/9sTrXpcjPqRwlwVpOdzPl/99PzdtX8UH3nKma7f415XnreDjV5/Hw788wj+/cspZskUWlcJdzlrHwBifeGgXq6tK+J8fvNCXY9rT9dn3bKZ53XL+5PuvsLvjuNfliI8o3OWsjE5G+f1vtxCJxXnwzst1FsR5FIWCfO23L6OmvIiPf6eF7iGd912WhsJd0haLO+753svs7x7m/o9cyqa6cq9Lygk15UU8eGczwxNR7vpOC+ORmNcliQ8o3CUtsbjjTx/bzb/t7eZzN2/lms26Bu7Z2FJfwd/85lvYfXSQj3+nhYkpBbwsLoW7zCsWd/ynH+zmB7/q4DPv3szHrtrgdUk56T3bVvGl2y7m+YO9CnhZdAp3OaNY3HHvD3bz2K4O7rlhM5+6vsnrknLabZet4Yu/cRE/a+3lrod2KeBl0Sjc5bRGJqN84qEWvr+rg09d38Snb1CwZ8KHmtfyVx+8iJ8eCHP7Ay/o4tqyKBTuMqeOgTFu++rPeWZfmPtu3cZn3u3vE4Jl2ocvX8tXP3op+44N84G/e569nUNelyR5RuEup3i+tZcP3P88R4+P883fvZzfedt6r0vKSzdur+f7f/A24g5u+9rP+eHLR70uSfKIwl2mTUzFuO+f9/LRB1+ksqSA//tHb9eomEW2fXUlO+6+ii31FXz6kZf59CMvMTim67DKuQt5XYBkh18dGeDeH+xmf/cId75tHffetIWSwqDXZflCXUUx37vrrXz1Jwf5ytMH+OWb/fyPX7+Q63xygXFZHOac82THzc3NrqWlxZN9ywndQxP81b+8wT+9dJSVFUV88baLeada657Z3XGce773MgfDo1x3fi2fu3kr59VqspicYGa7nHPN866ncPengdEI33z+Tf7+Z28yFXP8/tUb+OR1mygr0pc5r0Wicb7988N85ekDTEZjfOSKRj7xzo0686YACnc5je6hCb7xszd56IU2xiIxbtq+intvuoB1K8q8Lk1mCQ9P8tf/vo/vt3Rglhgjf9c1G9lQo9+VnyncZVo87vhZay/ffbGNH73eg3OO91/cwCev2+Tbi1nnko6BMb7+7CG+19JOJBrnHZtq+OiVjdywdSUFQY2J8BuFu88553ip/TiPv9LFk692cWxoguqyQj7UvIaPXNGolnoO6hme4NGd7Tz8y3aOHh9neWkBN25fxc0XNXDlhmpCCnpfULj7UP9ohOdbe3l2f5jn9ofpGZ6kMBjgnefX8v6LG/i1bSspCmkETK6LxR3P7u/hhy938qO93YxGYlSWFHB1Uw3v3FzL1U21rKos9rpMWSTphruOnuWoqVicA90jvNY5yEtHBvjlm/0cDI8CTP+hv+uCOm7YupKK4gKPq5VMCgaMd12wknddsJKJqRjPvNHD02/08Oz+MI/v7gJgbXUJl6+r5tJ1y7lwdSXnr1pGcYE+2P0krZa7md0IfAUIAg865/5y1vNFwHeAy4A+4Dedc4fPtE213NMzFolypH+Mtr4xDoZHaO0eYX/PMPu7R4hE4wBUFIdoXl9N8/rlXLlhBW9ZW0VQF6v2Hecce7uGeOFQPzvf7KelrZ/ekQiQ+EDYVFtO08pymuqWsamunHUrSmlcUaoP/xyTsW4ZMwsC+4F3Ax3ATuAO59zeGev8EXCRc+4PzOx24Nedc795pu36Odydc4xFYvSPRqZ/wiOThIcTP12D43Qen6BrcHz6jzOlvrKYTXXlbKmvYFtDBdtXV7JhRRkBhbnM4pyjY2CcPZ2DvHZ0iL1dQ7T2jNA+MMbMP/uq0gIaKktoqCqmvrKEumVF1CZ/qssKqS4rZHlZIcuKQrqkYhbIZLfMFUCrc+5QcsOPALcCe2escyvwF8n7jwF/Z2bmvOrQn4dzjrhL9F3GnSMad4n78RP3p2JxYnFHNB5nKpZ4PBWLE4m65G2cSCzOZDTG5FSciakY49O3MUYno4xHYoxMRhmNRBmZiDI8EWVwfIrB8Smi8bnfmmVFIVZVFlNfVcK2hgrWVpfSWF3KuhWlbKgpY5laWZImM2NtdSlrq0u5cXv99PLxSIw3e0c50j9KW98YR/rH6BqcoGNgnJ2HBxgcn/v0B8GAUVEcoqKkgIriAsqLQpQXhygrDFJalLgtKQhSnLotCFIUClAUClIYClAYClAQNIpCAQqCAUKBAIUhIxgIEAoYoaARDBihQIBgIHE/aEYgAEFLPNaHS/rSCffVQPuMxx3AladbxzkXNbNBYAXQm4kiZ3rohTb+7scHcA4cJG6dwwFx53AucUvy+bhLBHg8nrgfS66zmIoLApQVhigtClJWGKK8KMTyskLWVpdSWVJAZUkBFSUFiVZRaaJVVLesiJryIk35l0VXUhhka0MFWxsq5nx+MhqjdyRCeHiSgdEIfaMRBkYj0w2TwfEpRiYTDZb2/rFkYybGWCTK+FRsUf++zBJBHzDDDAJmBJK3JG9Ty43EB1zyqeTtzMd20nZPuk2ul7g/c/8zXnNSYXPePe2H0aevb+L9Fzec3T/+LKUT7nNVN/vXl846mNldwF0AjY2Naez6VGuXl3Dt5rrEL8gSuw7M+qUm9nXiF2/J/wzBQOKXFki2CIIBCAYCBAOJdUMBIxgMELREKyLRmghQEDAKggEKZrQ8CoOJ1khxQaJlUhQKUFKYuFXrQnJZUSjI6qoSVi9gRqxzjslo4hvsZDSe+FYbjRGJxpmMnvjGG019E445YvE40VjiW3M0+Q16KhZPNMZSjbLUt+tUY80l5m84ErczG20zG3mJL8gu2QgEl7o/XW9iWeLBiZtUp8PMEJv5oXXycjfn8lMT8ITKksX/Bp5OuHcAa2c8XgN0nmadDjMLAZVA/+wNOeceAB6ARJ/7Qgq+9vw6rj1fJ1QSyUZmRnGyS0a8lc6sh51Ak5ltMLNC4HZgx6x1dgB3Ju/fBvw4W/vbRUT8YN6We7IP/W7gKRJDIb/hnNtjZvcBLc65HcDfAw+ZWSuJFvvti1m0iIicWVqTmJxzTwJPzlr2+Rn3J4APZbY0ERFZKJ2MQkQkDyncRUTykMJdRCQPKdxFRPKQwl1EJA95dj53MwsDbZ7sPHNqWIRTLOQwvR8n6L04md6Pk53L+7HOOTfvVew9C/d8YGYt6ZydzS/0fpyg9+Jkej9OthTvh7plRETykMJdRCQPKdzPzQNeF5Cy1FCTAAACrUlEQVRl9H6coPfiZHo/Trbo74f63EVE8pBa7iIieUjhvgBmttbMnjGz181sj5l92uuavGZmQTN7ycwe97oWr5lZlZk9ZmZvJP+PvM3rmrxkZvck/05eM7OHzazY65qWipl9w8x6zOy1GcuqzezfzexA8nb5Yuxb4b4wUeCzzrktwFuBT5rZVo9r8tqngde9LiJLfAX4V+fcBcDF+Ph9MbPVwKeAZufcdhKnDffTKcG/Bdw4a9m9wNPOuSbg6eTjjFO4L4Bzrss596vk/WESf7yrva3KO2a2Bngf8KDXtXjNzCqAa0hc4wDnXMQ5d9zbqjwXAkqSV2kr5dQrueUt59xznHpVuluBbyfvfxv4wGLsW+F+jsxsPXAJ8KK3lXjqb4A/BeJeF5IFzgPCwDeT3VQPmlmZ10V5xTl3FPhfwBGgCxh0zv2bt1V5bqVzrgsSDUVgUa4bqnA/B2ZWDvwA+A/OuSGv6/GCmd0M9DjndnldS5YIAZcCX3XOXQKMskhfu3NBsj/5VmAD0ACUmdlveVuVPyjcF8jMCkgE+3edc//kdT0eugq4xcwOA48A7zKzf/C2JE91AB3OudQ3ucdIhL1f3QC86ZwLO+emgH8C3u5xTV7rNrN6gORtz2LsROG+AGZmJPpUX3fO/bXX9XjJOfdnzrk1zrn1JA6U/dg559uWmXPuGNBuZucnF10P7PWwJK8dAd5qZqXJv5vr8fEB5qQdwJ3J+3cCP1yMnaR1DVU5xVXAbwOvmtnLyWX/OXmtWZE/Br5rZoXAIeBjHtfjGefci2b2GPArEqPMXsJHs1XN7GHgWqDGzDqAPwf+EnjUzH6PxIffolx/WjNURUTykLplRETykMJdRCQPKdxFRPKQwl1EJA8p3EVE8pDCXUQkDyncRUTykMJdRCQP/X9ViTYCUsUdgQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "x = np.linspace(1,10,100)\n", "def g(x):\n", " return np.exp(-2*np.abs(x-5))\n", "plt.plot(x,g(x))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Cependant, la fonction $g$ atteint son maximum en $x=5$ et décroît\n", "rapidement après, il est\n", "donc sans doute plus malin d'utiliser une fonction d'échantillonnage\n", "préférentiel gaussienne $f_Y$ centrée en 5 et de variance 1 (par exemple). On réécrit donc \n", "$$\\int_0^{10} e^{-2|x-5|}dx = \\int \\mathbf{1}_{[0,10]}(x)\n", "\\frac{e^{-2|x-5|}}{\\frac{1}{\\sqrt{2\\pi}}e^{-(x-5)^2/2}}\\frac{1}{\\sqrt{2\\pi}}e^{-(x-5)^2/2}dx\n", "= \\E[\\frac{g(Y)}{f_Y(Y)} \\mathbf{1}_{[0,10]}(Y)],$$\n", "avec $f_Y(x) = \\frac 1 {\\sqrt{2\\pi} }e^{-(x-5)^2/2} $ et $Y$ de loi gaussienne $\\mathcal{N}(5,1)$.\n", "\n", "Pour calculer l'intégrale précédente, on tire donc $N$ réalisations indépendantes $x_1\\dots x_N$ d'une loi $\\mathcal{N}(5,1)$ et l'approximation est donnée par \n", "$$ \\frac 1 N \\sum_{i=1}^N \\sqrt{2\\pi}e^{+(y_i-5)^2/2}e^{-2|y_i-5|} \\times \\mathbf{1}_{[0,10]}(y_i).$$\n", "\n", "2. Utilisez cette nouvelle approche pour calculer l'intégrale et comparer les variances et vitesses de convergence des deux approches. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercice 2\n", "\n", "On souhaite calculer \n", "$$\\int_{-\\infty}^{+\\infty} x^2 \\frac 1 2 e^{-|x|}dx,$$\n", "mais on ne sait pas échantillonner suivant la densité $p(x) = \\frac 1 2 e^{-|x|}dx$.\n", "\n", "Réécrire cette intégrale à l'aide d'une fonction d'échantillonnage préférentiel gaussienne de $\\cN(0,4)$ et estimer sa valeur par une méthode de MOnte-Carlo. " ] } ], "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.5" } }, "nbformat": 4, "nbformat_minor": 2 }