491 lines
192 KiB
Plaintext
491 lines
192 KiB
Plaintext
|
{
|
||
|
|
"cells": [
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"id": "1e9456",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"# Generowanie zmiennych losowych o dowolnym rozkładzie w Pythonie. Dyskretne zmienne losowe wielowymiarowe w Pythonie."
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"id": "382d67",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"## Metoda odwracania dystrybuanty\n",
|
||
|
|
"\n",
|
||
|
|
"Istnieje wiele różnych metod generowania liczb losowych o zadanym rozkładzie prawdopodobieństwa. Dzisiaj omówimy jedną z nich - będzie to metoda odwracania dystrybuanty."
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"id": "dcfe80",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"### Dystrybuanta i dystrybuanta empiryczna\n",
|
||
|
|
"\n",
|
||
|
|
"Przypomnijmy najpierw, że dystrybuanta jest jedną z funkcji, za pomocą której można zdefiniować rozkład prawdopodobieństwa zmiennej losowej. Zgodnie z definicją dystrybuantę zmiennej losowej $X$ liczymy ze wzoru:\n",
|
||
|
|
"\n",
|
||
|
|
"$$ F_{X}(x)=\\mathbb{P}(X\\le x) = \\mathbb{P}(X^{-1}\\left((-\\infty,x]\\right)),$$\n",
|
||
|
|
"\n",
|
||
|
|
"dla dowolnego $ x\\in\\mathbb{R}$. Ponadto, jeśli zmienna losowa $X$ jest dyskretna, a jej zbiór atomów to $A=\\{a_1,a_2,\\ldots\\}$, to \n",
|
||
|
|
"\n",
|
||
|
|
"$$F_{X}(x)=\\sum_{a\\in A, a\\le x}\\mathbb{P}(X=a). $$\n",
|
||
|
|
"\n",
|
||
|
|
"Chcąc zbadać, jak wygląda rozkład pewnej zmiennej losowej na podstawie próbki, możemy spróbować oszacować właśnie dystrybuantę. W tym celu definiujemy tzw. **dystrbuantę empiryczną**.\n",
|
||
|
|
"\n",
|
||
|
|
"**Definicja (dystrybuanta empiryczna)**\n",
|
||
|
|
"\n",
|
||
|
|
"Niech $x_1,x_2,\\ldots, x_n$ będzie próbką wylosowaną z rozkładu o dystrybuancie $F$. Funkcję $\\hat{F}:\\mathbb{R}\\rightarrow [0,1]$ zdefiniowaną wzorem\n",
|
||
|
|
"\n",
|
||
|
|
"$$\\hat{F}(x)=\\frac{1}{n}\\sum_{i=1}^n\\pmb{1}_{(-\\infty,x]}(x_i)=\\frac{\\#\\{i:x_i\\leq x\\}}{n}$$\n",
|
||
|
|
"\n",
|
||
|
|
"nazywamy **dystrybuantą empiryczną**.\n",
|
||
|
|
"\n",
|
||
|
|
"Innymi słowy, chcąc wyznaczyć wartość dystrybuanty empirycznej w punkcie $x$ sprawdzamy, ile wyników z naszej próbki jest mniejszych lub równych od $x$ i porównujemy tę liczbę z całkowitą liczbą otrzymanych wyników.\n",
|
||
|
|
"\n",
|
||
|
|
"Chcąc wyznaczyć dystrybuantę empiryczną w Pythonie, możemy posłużyć się funkcją `ecdf(x)` z pakietu `stats` z biblioteki `scipy`. Argumentem tej funkcji jest wektor próbek. Załóżmy, że `f` jest dystrybuantą empiryczną wygenerowaną za pomocą polecenia `ecdf(x)`. Wówczas możemy badać jej właściwości za pomocą następujących poleceń:\n",
|
||
|
|
"\n",
|
||
|
|
"* `f.cdf.quantiles` - zwraca wektor unikalnych wartości w próbce,\n",
|
||
|
|
"* `f.cdf.probabilities` - zwraca wektor skumulowanych prawdopodobieństw związanych z wartościami z poprzedniego podpunktu,\n",
|
||
|
|
"* `f.cdf.evaluate(t)` - zwraca wartość dystrybuanty empirycznej w punkcie $t$,\n",
|
||
|
|
"* `f.cdf.plot(ax)` - rysuje wykres dystrybuanty empirycznej na zadanych osiach.\n",
|
||
|
|
"\n",
|
||
|
|
"**Przykład 1**\n",
|
||
|
|
"\n",
|
||
|
|
"Załóżmy, że w $10$ rzutach czworościenną kostką do gry otrzymaliśmy następujące wyniki: $1, 1, 4, 1, 3, 2, 2, 1, 4, 4$. Na podstawie otrzymanych wyników wyznacz dystrybuantę empiryczną zmiennej losowej $X$ oznaczającej liczbę wyrzuconych oczek w rzucie czworościenną kostką. Następnie zbadaj własności tej dystrybuanty i zwizualizuj ją na wykresie."
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 3,
|
||
|
|
"id": "9e518c",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"name": "stdout",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"Wartości otzymane w eksperymencie: [1. 2. 3. 4.]\n",
|
||
|
|
"Skumulowane prawdopodobieństwa dla otrzymanych wartości: [0.4 0.6 0.7 1. ]\n",
|
||
|
|
"Wartości dystrybuany empirycznej w punktach 0, 2, 3.5: 0.0 0.7 1.0\n"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAABaAAAANwCAYAAADZTs3dAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuNSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/xnp5ZAAAACXBIWXMAABYlAAAWJQFJUiTwAABAiUlEQVR4nO3deZxlZ13n8e/PhEAIEAgQGEWJICEZQSFBCIuQsIRNTWQRRQFRnEFwIgiOOA4CjoxRGAygjixCBBTBhcVBMYIgS0ARghgNJCBRkD2QhSQkkPzmj3MuqRRV3dXd9eRWdb/fr1e9Ttc555771K3b/Uo+9dRzqrsDAAAAAACb7ZuWPQAAAAAAAPZOAjQAAAAAAEMI0AAAAAAADCFAAwAAAAAwhAANAAAAAMAQAjQAAAAAAEMI0AAAAAAADCFAAwAAAAAwhAANAAAAAMAQAjQAAAAAAEMI0AAAAAAADCFAAwAAAAAwhAANAAAAAMAQAjQAAFtWVfX8cdiyx8L2VFXPnN9Dpy57LAAA+6L9lz0AAIC93Ry+HrNq99eSXJjkS0nOTPL3SV7d3R+/Zke3tqo6NsmxST7Y3a9f5lj2JVV1YpI7JHl7d799qYMBAIBNYAY0AMA156tJPjt/fDHJdZPcOskJSZ6d5GNV9dqqusnyhvh1xyZ5RpITlzuMfc6JmV73Y5c7jL3KF5J8JMmnlz0QAIB9kQANAHDNOb27bz5/3Ky7D0xyoyQPTPKaJJ3k4Uk+WFW3WOZAYW/R3b/V3Ud09y8ueywAAPsiARoAYIm6+/zufnN3/3CSByf5SpJvSfInyx0ZAADAnhOgAQC2iO5+c5Knzp/epaq+f3Gsql4230hth2G6qp41n3f6qv2HVtVzqurMqrq4qr5SVZ+oqtOr6leq6pbzeYdVVWdaBiJJHrPiRoBXuyHg4tz5/FTVMVX1J1X16aq6oqpOqapfns/5h52M+7HzeZ+oqjX/G7WqbldVf1RVn5nH/+GqenpVXXud88+dr3nsDp53zZscVtW1q+rhVfWKqvrHqvrC/Jz/VlV/UFVH7+CaX3/eqjqkqp5XVR+vqsuq6j+q6iVV9Z9WPebY+XVcrBX+jNWv+6rzj6qqk6vqXVX17/O1z6uqt1fV46pqv/XGtxHza/2yedxfqarzq+rdVfX4qrrWGuevfi/cuareUFWfr6qL5vfZg1acf0BV/cL8frykqj5bVS+qqkM28Jp+W1W9dH6vfGUe43Or6uB1HrvuTQhXfv+r6siq+v35ul+tqtdX1aPn45+pqnXvn1NVx83nXbLWOOZr/25VnT2fc35V/VNVvWDle2nFWHf28fYdvPbr/T05YJ2x36KqnlpVb66qc+bxXVhVZ9T078kN1/u6AQA2wk0IAQC2lpck+eUkhyZ5ZJI/n/e/NMljk3x/Vd24u89b/cCawu0iYL5sxf5bJnlPkkX0vCLTDRC/Jcktktw1yaeS/O587LNJrpfkoEwzsi9Y9VRXrPHcj0jyqkz/fXnBinNelilmH11Vt+/uf1rn6/6Jefv73X3lGsfvluTF85guTFJJbpvkV5I8qKru191fXufau+N+SV47/7mTnD9vvy3T9+WHquonuvuVO7jGLZKcmuSWSS6ZH//NSR6X5L5VdVR3f2k+9/JMr/vBSa6T5OIkO/p6Tkty4/nPl8wfhyS51/zxg1V1Qnd/bYNf79dV1c8keX6umqzy5Uzvh7vNH4+oqgd39yXrPP6EJH+c6b1w4fzYuyb586r64Uzv6b/MtM71VzK9Locm+S9Jvqeqjunuy9cZ3ndk+r7cdB5XJzksyVOSnFBV9+zu3Vnr+Xszvf+vm+SiTDcJzfx1vCDJzTItlfPnaz76qvfvn3X31f6+VNV/S/KbSRY/FLh4Hvft5o/vylVrfn850/tgPTdKsmZInp/r+CSvT3Jgpr+H18pVf0+Oztprup+S5KHzny+fx3DDTDfDvEOSH62qY7v7kzsYFwDAusyABgDYQubw9tb50+9dsf/0JP+SKT796DoPv3em2HlxpjWlF56RKT5/NMk9kxzQ3YdkilS3T/KrST4zP88nuvvmSZ47P/Y1K9atXnx8Yo3nfmmSNyT59u6+YaaQd8ocrf5qPuexaw26qm6T5B6ZotzL1/nafmf++r+ruw9Ocv35epcmOSbJ89Z53O76cqbweM8k1+vuQ+Y1u2+ZKdjtn+TFVfVtO7jGC5N8KcnduvugTCH2hEwx+7AkX1+TuLtPn1/3xfftuatf91XXPi3JjyT5T919UHffaL7+ozJ9Lx+U5Mm7+kVX1YnzuC9O8t+T3LS7r5/p+/mAJOdkiqW/uYPL/H6SV8xju2GmuPyGTP/v8ZuZ3ltHJPm+eczXz/S6XJTkjpkC/Xqemymsfu88roMyRdUvZIrTv79LX/BVfifJ+5LcvrtvkOnrfUp3X5rkD+dz1nv/3iDJQ+ZPX7bq2MMzvY/2y7Sszn/u7uvN368bJ/mxJO9fnN/d3/B9X/H9f2Cu+sHOX67zdbwmUyRf/D28Qab3WWcK9A9a4zFnJTkpyeFJDuzuG2f6Icix82ty6yQvWuf5AAB2SoAGANh6FrOEv2XVcgcvnbdrhrBcNQvzT7r7ohX7j5m3/7O737mYYdzdl3X3md399O5+/R6O+R+T/FB3nztf+2uLP2ea1Z0kP7bW8g256ut5R3d/bJ3rX5bkAYsZ1N19eXefmuQJ8/Gf3EkM3iXd/fbu/tn59bpkxf5/7+4nZwqN18n634vFmO/b3e+ZH/u17n5jpuCfJA/bg/E9srv/qLs/s2Lfxd39qiQ/NO96wtqPXtu8bMcp86cP7+7ndPcX5mtf3t1/lSmCXpLkJ2rVMiIrfKC7H9fdn50f+/lMPzRZzLp/YpJHdPebuvuK+eONSZ4zP35Hr8u1kzywu981X/vK7n7Diq/5flV1j135umefm6975nzdXvFeXLx/v6+qbrrGY384U7D+eJK3LXbO7/VFqH91dz+8u89aHO/uL3b3H3T3U3Y2uKo6NFfNbH5td//6Oqe+L8kPr/h7eHF3n5zkTfPxb3ht57//L+zuc1b82/DV7v7bTD90+HySB9aqZWoAADZKgAYA2Hq+tOLPK9fEfUWmX5G/Q1XdceUD5nVaf3D+9GqzMDOFv+SqJThG+D/rLJ2RTDMyP5tp2YTvW3lgXjbk0fOnq8e90u929xfX2P+KJJ/M9N+1D1nj+CiLpRjuvoNzXrzWUimZQmKSfHtVHbSpo0rS3e/MPMu6qr55Fx56bKYZ3mfOsXmta38syXszzQA/dp3rnLzG4y6eH5ckp3f3O9Z43GLm/+12MMbXdvdH17j+25Is1j3fnbD/W/Ns52/Q3Wck+UCm5Sx+bI1TFj+EOLW7V67VfZ9Mwf2KJD+/G2NK8vWQ/aeZln85Izv+ocfJq8aw8Pp5u6PX9hvMf+dOz7Tkzd125bEAAAsCNADANjHHzNfPn66OUD+SaUbuOWvEvb+Yt79eVb893zDtwE0e3nvWO9DTOsSLpRFWj/v+mSLdhZmWKFjP29e59pVJ3jl/etRGBrpRNd1A8OnzDfTOq6qvrbjZ2+vm03YUeN+3zv7/WPHnG+7B+B4+3yjv36vq0rr6DQsX192VAL0IjLeZb2C35seK8751neust8735+btmescX6x9fKMdjPHtOzj2t/N2d94H675/Z2v+9kFVHZnpNwyuzLTe90qL3zz4x+7+j+y+38q0RM3nkpzY66y9PdvZe27N17amm0a+bL5h4ZdXvZdOmE/blfcSAMDXCdAAAFvPyki0etbvIoQ9sqpW3oxssfzGWmso/3qSN2ZaP/oJSf4myYVzWP35efb0nvr8To4vxv3AqrrZiv2Lcf/RTsLajgLe4thayyPslqr6z5nWnP6VTDfQOyTT0hOfyxRKF7PUdzSD+aK1dnb3V1Z8utaSJDsb2/5V9WeZbsZ3QqYQXJnWQf7s/LGYjb4rM6wXM+Svnemme+t9XGc+77prXaTXvwngYv3inR3f0Y3SR70Pdvb+/cNM3//bV9XRK/Yv3r9v6e5/X/WYxft89f4Nq6onZLo54+VJHrrGc1z
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 864x504 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"execution_count": 3,
|
||
|
|
"metadata": {
|
||
|
|
"image/png": {
|
||
|
|
"height": 440,
|
||
|
|
"width": 720
|
||
|
|
},
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "execute_result"
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"from scipy import stats\n",
|
||
|
|
"import matplotlib.pyplot as plt\n",
|
||
|
|
"\n",
|
||
|
|
"# Wprowadzamy dane z zadania i generujemy dystrybuantę empiryczną\n",
|
||
|
|
"x = [1, 1, 4, 1, 2, 2, 3, 1, 4, 4]\n",
|
||
|
|
"f = stats.ecdf(x)\n",
|
||
|
|
"\n",
|
||
|
|
"# Badamy własności otrzymanej funkcji\n",
|
||
|
|
"x = f.cdf.quantiles\n",
|
||
|
|
"print('Wartości otzymane w eksperymencie:', x)\n",
|
||
|
|
"p = f.cdf.probabilities\n",
|
||
|
|
"print('Skumulowane prawdopodobieństwa dla otrzymanych wartości:', p)\n",
|
||
|
|
"print('Wartości dystrybuany empirycznej w punktach 0, 2, 3.5:', f.cdf.evaluate(0),f.cdf.evaluate(3),f.cdf.evaluate(5))\n",
|
||
|
|
"\n",
|
||
|
|
"# Rysujemy wykres dystrybuanty empirycznej\n",
|
||
|
|
"os = plt.subplot()\n",
|
||
|
|
"f.cdf.plot(os)\n",
|
||
|
|
"plt.title('Dystrybuanta empiryczna')\n",
|
||
|
|
"plt.xlabel('x')\n",
|
||
|
|
"plt.ylabel('F(x)')\n",
|
||
|
|
"plt.xticks(x)\n",
|
||
|
|
"plt.yticks(p)\n",
|
||
|
|
"plt.show()"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"id": "39c337",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"### Dystrybuanta pseudoodwrotna\n",
|
||
|
|
"\n",
|
||
|
|
"Metoda odwracania dystrybuanty opiera się na następującym twierdzeniu:\n",
|
||
|
|
"\n",
|
||
|
|
"**Twierdzenie (metoda odwracania dystrybuanty)**\n",
|
||
|
|
"\n",
|
||
|
|
"Niech $F$ będzie dystrybuantą pewnej zmiennej losowej, a $U$ będzie zmienną losową o rozkładzie jednostajnym na odcinku $(0,1)$. Jeśli $F^{-1}$ jest funkcją odwrotną do $F$, to funkcja $X=F^{-1}(U)$ jest zmienną losową o dystrybuancie $F$.\n",
|
||
|
|
"\n",
|
||
|
|
"Niestety, nie wszystkie zmienne losowe mają dystrybuanty, które są funkcjami odwracalnymi. Dlatego w praktyce będziemy korzystać z tzw. dystrybuanty pseudoodwrotnej, którą możemy zdefiniować dla dowolnej zmiennej losowej.\n",
|
||
|
|
"\n",
|
||
|
|
"**Definicja (dystrybuanta pseudoodwrotna)**\n",
|
||
|
|
"\n",
|
||
|
|
"Niech $F:\\mathbb{R}\\rightarrow [0,1]$ będzie dystrybuantą zmiennej losowej $X$. Wtedy funkcję $F^-:(0,1)\\rightarrow \\mathbb{R}$ definiowaną wzorem\n",
|
||
|
|
"\n",
|
||
|
|
"$$F^-(u)=\\inf\\{x:F(x)\\geq u\\}$$\n",
|
||
|
|
"\n",
|
||
|
|
"nazywamy **dystrybuantą pseudoodwrotną**.\n",
|
||
|
|
"\n",
|
||
|
|
"Jeśli $F$ jest dystrybuantą, która posiada funkcję odwrotną $F^{-1}$, to $F^-=F^{-1}$.\n",
|
||
|
|
"\n",
|
||
|
|
"**Twierdzenie (ogólna metoda odwracania dystrybuanty)**\n",
|
||
|
|
"\n",
|
||
|
|
"Niech $F$ będzie dystrybuantą pewnej zmiennej losowej, a $U$ będzie zmienną losową o rozkładzie jednostajnym na odcinku $(0,1)$. Jeśli $F^-$ jest dystrybuantą pseudoodwrotną do $F$, to funkcja $X=F^-(U)$ jest zmienną losową o dystrybuancie $F$.\n",
|
||
|
|
"\n",
|
||
|
|
"Powyższe twierdzenia stanowią podstawę następującego algorytmu, który pozwala generować zmienne losowe o dowolnym rozkładzie.\n",
|
||
|
|
"\n",
|
||
|
|
"**Algorytm (Generowanie zmiennych losowych przez odwracanie dystrybuanty)**\n",
|
||
|
|
"\n",
|
||
|
|
"1. Generujemy $n$ liczb pseudolosowych $u_1,u_2,\\ldots,u_n$ z rozkładu jednostajnego na odcinku $(0,1)$.\n",
|
||
|
|
"2. Obliczamy $x_i=F^-(u_i)$ dla $i=1,2,\\ldots,n$.\n",
|
||
|
|
"\n",
|
||
|
|
"Przypomnijmy, że do generowania liczb pseudolosowych o rozkładzie jednostajnym $U(0,1)$ możemy użyć polecenia `uniform.rvs(size=n)` z pakietu `scipy.stats`, gdzie $n$ jest rozmiarem próbki, którą chcemy otrzymać.\n",
|
||
|
|
"\n",
|
||
|
|
"**Przykład 2**\n",
|
||
|
|
"\n",
|
||
|
|
" Niech $X$ będzie zmienną losową o dystrybuancie $F:\\mathbb{R}\\rightarrow [0,1]$ danej wzorem:\n",
|
||
|
|
" \n",
|
||
|
|
" $$F(x)=\\begin{cases} 0, &\\text{ dla }x<0,\\\\\n",
|
||
|
|
"x^2, &\\text{ dla } x\\in[0,1],\\\\\n",
|
||
|
|
"1, &\\text{ dla } x>1.\\end{cases}$$\n",
|
||
|
|
"\n",
|
||
|
|
"Wygeneruj $100$ liczb zgodnie z rozkładem zmiennej losowej $X$. Następnie wyznacz dystrybuantę empiryczną na podstawie otrzymanej próbki i porównaj ją z rzeczywistą dystrybuantą zmiennej losowej $X$.\n",
|
||
|
|
"\n",
|
||
|
|
"Wyznaczmy najpierw dystrybuantę pseudoodwrotną do dystrybuanty $F$. Niech $u\\in (0,1)$. Wówczas:\n",
|
||
|
|
"\n",
|
||
|
|
"$$F^-(u)=\\inf\\{x:F(x)\\geq u\\}=\\inf\\{x:x>0\\wedge x^2\\geq u\\}=\\sqrt{u}.$$\n",
|
||
|
|
"\n",
|
||
|
|
"Zobaczmy, jak możemy wykorzystać dystrybuantę pseudoodwrotną do generowania liczb o rozkładzie zgodnym z rozkładem zmiennej losowej $X$."
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 4,
|
||
|
|
"id": "ad3784",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAABYYAAANyCAYAAADWxR9BAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuNSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/xnp5ZAAAACXBIWXMAABYlAAAWJQFJUiTwAAD9E0lEQVR4nOzdd3gUVdvH8e8JKXQIUkWaiIKoIKCggjRBmoAC6gMiKIqKDeyor5RHBQsKyCOiiKjYQUGa9GpFsCF2AiodQiDUkOS8f8xsdjfZTU825fe5rr2yc2bmzL1nZybJvWfPMdZaRERERERERERERKT4CAt1ACIiIiIiIiIiIiKSv5QYFhERERERERERESlmlBgWERERERERERERKWaUGBYREREREREREREpZpQYFhERERERERERESlmlBgWERERERERERERKWaUGBYREREREREREREpZpQYFhERERERERERESlmlBgWERERERERERERKWaUGBYREREREREREREpZpQYFhERERERERERESlmlBgWERERERERERERKWaUGBYREREREREREREpZsJDHYCIiEhBY4y5COgOHAQmW2ttiEMSERERERERyVVKDIuIiPgwxlQCPgaqA12UFBYREREREZGiSENJiIiI+HsdOAMYZq1dEepgCipjzExjjDXGjA51LL6MMdvcuNqFOpZQMsbUddtBH2zkA2PMYLe9VxemugubgnrfkcJD51DeMMa0c9t1W6hjERGRrFFiWERE8pzPP2KpH4eNMd8bY54zxpxRAOIcBvQGJlhrXwtxOFJIGWN6G2NGF/fkdH5zE6ijjTFNQx1LUWOMqei27ehQx1KcuB/wjDbGDA91LCIiIlI0aSgJERHJT6eAWPe5AaoATdzHLcaYq6y160MRmDHmPGAC8CnwUChiKGR2Ab8B+0MdSAHUGxjkPl8dujCKncFAW2Ab8H0+H/sQzvXwdz4fN79UBEa5z0eHLoxipy5Ou28HJoY0kqJDv7vyxjGcdt0R6kBERCRrlBgWEZH89IW1tp1nwRhTGugDTMZJPHxkjDnTWns8vwOz1m4GSuX3cQsra+1IYGSo4xApCKy1nwCfhDoOEUmffnflDWvtN0DDUMchIiJZp6EkREQkZKy1x6y1bwP3uEXVcXpbioiIiIiIiEgeUmJYREQKgg+BZPd5c98Vxpjy7hiLPxhjjriPH40xY4wxFQJV5m5v3bGNw4wxdxljvjHGxLnlTX22jTLG3GeM+doYc8gYc9wY85sx5gVjTPUAdQ9069gQYF1lY0yyu35xgPXnuOtOGGNK+pSvdssHG2NKufH/5say1xjzvjGmQZDXGmWM6WeMectto/1u/duNMe8YY5oH2s/dN2WiNmNMJfc1xxhjThpjdhhjXjPG1Aiyb4YT+BhjrjLGzDPG7DbGJLivZb4x5spg+2SGMWaAMeYr91yINcasNMZ0T2f7lW6sz2dQ75vudu+mKq9njJlqjPndfU+Oue272hgz0hhT2d2unXEme/MMIzHKpBpX26dOv4l6jDFdjTGL3TZKNsYMN8bMcLeZnUHcY9ztvkhnm8uMMQuMMfvc+L93r4uAfwv6xFw3yPqgk9sZY8q55/KHxpjN7nV33BjzpzHm1WDncurjGmNqu+fgv+45GWOMed4YUz7VPoPdONq6RW+kavdt7na51p4B9snRBHHGmNPdttnhXr9b3euxYpDtL3ePd9IYc1o69Z5pvPekc3zKw9yYVxljDhhjTrnnxs9uO3Xx2XY1EOOznHqs+NG+2xrvvayiMeYZY8yv7jkX576vnnjOSyfussa5vq0xpnOQbUq679WvxnuvfM8Yc3aQ7VN+L6Rz3KD3NWNMM2PMeGPMemPM327bH3Bf8y3GmBKZOa4xZpBxft/EG2ec/VXGmE4B9tsGrHIX6wRo98G52Z7pMca0Ns7voX99XvdyY8x/jDEmwPap729XutvHuufBMmPMJT7bVzDGPGW899h/3HMn4Ld4jP994jw3tt3utfOrMeb/jDFRQfYN+B6bVPc0Y0wrY8xsY8wuY0ySMWaiMeYJd5tvM2ivm9zt/jEB7rFu3W8Z53fwCeP83t5kjBln/K/TYPMzpH7MTKftPff+/W7b/mCce3+a9y1VjFn6/Z36uCIiUohYa/XQQw899NAjTx/ATMACq9PZZo+7zas+ZWfhjBdq3cdR9+FZ3g40CFDXaHf9m8Bc93kicNB93tTdrgqwyae+E8Bhn+VYoFWqumv71Fcu1bprfPY9BJRItX6ou25NqvLVbvk9PvGcwBmzz1PfAaB+gNfaw2ebZDfm4z5lp4CBQdrc07Y3+Dw/6h7bs38MEJ3Oezo6wLoIYJZPHZ728F1+Jpvn0hSfOpLc9zTZp/08r6Odzz793bLdQHiQesv5nFtX+JQ3S3VOJPicR55HF3fbS91jeNr/iLuc8vCpt527zTbgfp/376B7bg1367PASeC0IHGH+bzmW3zK6/rE18c9D6xb/ymfdZ8EahOf9XWDHDel/gDr7vLZPxHn3D3pU3bEt42DHLeXu59129835g1AhM8+17ntm+Bzrvm2+waf9ydb7ZmJ83IwGdzj0tm3EbA3Vft4rv0/gPsC1Y0znqcF7k6n7v+626xPVf6Oz/EsEJfqPfrKZ9uPgX0+63anejwQ4F72IPAX/vfVOHebpW75hHTiHoL3Hh8W4L4zDvjS5/30vb8cBS4PUOdod/3MdI7rqT/QfW1/qmOkvg8sJPC1lHJcYDre68I35iSgT6r9NuDczz3rU7f7dTltz0yen8+kep2H8N5zLfBe6jrxv78Nc7dPSvWajwOX4fwe/gnvue97Hi7I4D7R393HE5fvvl8CZTP7HuN/z7wO7z0nDufeMhE4w30dFjg/nTZb527zZKpyE6Q9fX/HzPTZflKA9933kRhgH9+2H+xuk+y+Dt/jTgwSe7Z+f/seN6v3QD300EMPPUL7CHkAeuihhx56FP0HGSSGccb29fyj+axbFgn84Jb9DXRy/6kyQEecf3AtsBmISlXfaHddPE5S4g6gtLuuKlDefb7Y3S4W6IebyAVaAD/iTYJUTlX/NnwSgj7lk/AmsizQItV6TzJmbKry1XiTdjHAlUAJnCRVG+Afd/2HAdqunXvcNp7X6JbXBl7E+w947QD7bvM57nfAJW55ONATb+Lj2XTe09EB1nmO+4fbrmXc8nLue+Fpn/9k8Twa4POP6XNARbe8Gs6HAAl4k7vtfPYriTfB0jNI3bfg/Wfa+JSvdMu/Ai70KS/tnicvetotM22T6n3zvDeJwP+Aaj7xnuE+/9nd7p4g9VyBN6FSzqe8rk9bxeGc6/XcdWVwEneeBMejAer17Fs3yHFT6g+w7nrgSeAiINItMzjjT3oSDns950WQ4x4EVgDnueVRwM14P7QYFmDf1e66wem0e7baMxPn5mCykRjGScJ4YvoLN6GJc+1f5bZTXKC6cSbJtMCmIHWH4dw7LXCzT/nleJOTwz2v032PauD0eH8+s+93kPcg3j12F9ykIXCW+/Nad5s9BP+gZj2B75Uzfc7po8BA3A8JgKbARrz37ehU+44mZ4nhd3HO7eo+ZWVwPljb5e73YID9PMc9iHO9347391E9YI27fmfq9iATybactGcmzs97fdrzVqCCW14KJ3nqed0jg8R9FCdZ+xTe+3Vd4At3/TfAHOBXoLV7DkbiJLI9idlu6dwn4tw6znfLI3GuRc8HK68G2Dfge4z/PTMemI17/8P5neh5vsjd5oUgbdbAXZ9Mqg9zce67nmP8D6jjs64GcBvwWCbfm2547+HXp9P2L+H93VIRZ04HT3yNA9T7Itn4/Y0Sw3rooYcehfYR8gD00EMPPfQo+g8yTgz79jDs45YNdJcTcJNDqfZpjLeH4M2p1o32qW9okGO28dnmygDrq+FNJqZOTrzplo9LVf69W/60+/P+VOv/dcs7pipf7ZYfw02epFrfx11/AjfRloW2f93dd1SAddvw/tOfpgcl3p6sW9N5T0enKm/g/sO5F6gVJKbr3X03Z+F1GOBPgiR23PXLfN7TdqnWe5L2nwS
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 864x504 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"execution_count": 4,
|
||
|
|
"metadata": {
|
||
|
|
"image/png": {
|
||
|
|
"height": 441,
|
||
|
|
"width": 707
|
||
|
|
},
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "execute_result"
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"import numpy as np\n",
|
||
|
|
"\n",
|
||
|
|
"# Definiujemy dystrybuantę pseudoodwrotną\n",
|
||
|
|
"def F_pseudoinverse(u):\n",
|
||
|
|
" x = [i**(1/2) for i in u]\n",
|
||
|
|
" return x\n",
|
||
|
|
"\n",
|
||
|
|
"# Generujemy 100 liczb o rozkładzie zgodnym z rozkładem zmiennej losowej X, używając dystrybuanty pseudoodwrotnej\n",
|
||
|
|
"u = stats.uniform.rvs(size=100)\n",
|
||
|
|
"probki = F_pseudoinverse(u)\n",
|
||
|
|
"\n",
|
||
|
|
"# Wyznaczamy dystrybuantę empiryczną na podstawie wygenerowanej próbki\n",
|
||
|
|
"F_empirical = stats.ecdf(probki)\n",
|
||
|
|
"\n",
|
||
|
|
"# Narysujemy wykresy dystrybuanty zmiennej losowej i dystrybuanty empirycznej\n",
|
||
|
|
"x = np.linspace(0, 1, 200)\n",
|
||
|
|
"y1 = [i**2 for i in x]\n",
|
||
|
|
"y2 = F_empirical.cdf.evaluate(x)\n",
|
||
|
|
"\n",
|
||
|
|
"plt.plot(x, y1, label='dystrybuanta zmiennej losowej X', color='red')\n",
|
||
|
|
"plt.plot(x, y2, label='dystrybuanta empiryczna wygenerowanej próbki', color='blue')\n",
|
||
|
|
"plt.title('Porównanie dystrybuanty i dystrybuanty empirycznej')\n",
|
||
|
|
"plt.xlabel('x')\n",
|
||
|
|
"plt.legend()\n",
|
||
|
|
"plt.grid()\n",
|
||
|
|
"plt.show()"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"id": "ba7094",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"**Przykład 3**\n",
|
||
|
|
"\n",
|
||
|
|
"Załóżmy, że $X$ jest zmienną losową o następującym rozkładzie:\n",
|
||
|
|
"\n",
|
||
|
|
"\n",
|
||
|
|
"| $k$ | $0$ | $1$ |\n",
|
||
|
|
"| --- | --- | --- |\n",
|
||
|
|
"| $\\mathbb{P}(X=k)$ | $0{,}8$ | $0{,}2$ | \n",
|
||
|
|
"\n",
|
||
|
|
"Znajdź dystrybuantę i dystrybuantę pseudoodwrotną dla tej zmiennej losowej. Następnie wygeneruj próbkę $100$ liczb zgodnie z rozkładem zmiennej losowej $X$, używając metody odwracania dystrybuanty i porównaj rozkład otrzymanych liczb z rozkładem tej zmiennej losowej.\n",
|
||
|
|
"\n",
|
||
|
|
"Zacznijmy od wyznaczenia wzoru dystrybuanty:\n",
|
||
|
|
"\n",
|
||
|
|
"$$F(x)=\\begin{cases}\n",
|
||
|
|
"0, &\\text{ dla } x<0,\\\\\n",
|
||
|
|
"0{,}8, &\\text{ dla } 0\\leq x<1,\\\\\n",
|
||
|
|
"1, &\\text{ dla } 1\\leq x.\n",
|
||
|
|
"\\end{cases}$$\n",
|
||
|
|
"\n",
|
||
|
|
"Przejdziemy teraz do dystrybuanty pseudoodwrotnej. Spróbujemy najpierw wyznaczyć jej wartość w przykładowym punkcie $u=0{,}1$. Zgodnie z definicją dystrybuanty pseudoodwrotnej musimy znaleźć możliwie najmniejszą liczbę rzeczywistą $x$, dla której zachodzi nierówność $F(x)\\geq u$. Patrząc na wzór dystrybuanty, możemy łatwo zauważyć, że taką liczbą będzie $x=0$. Możemy też łatwo stwierdzić, że dystrybuanta pseudoodwrotna będzie przyjmować wartość $0$ dla wszystkich argumentów $0< u\\leq 0{,}8$. W podobny sposób możemy wyznaczyć pozostałą część wzoru na dystrybuantę pseudoodwrotną:\n",
|
||
|
|
"\n",
|
||
|
|
"$$F^-(u)=\\begin{cases}\n",
|
||
|
|
"0, &\\text{ dla } 0<u\\leq 0{,}8,\\\\\n",
|
||
|
|
"1, &\\text{ dla } 0{,}8<u\\leq 1.\n",
|
||
|
|
"\\end{cases}$$"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 5,
|
||
|
|
"id": "63f4ab",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"data": {
|
||
|
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAABaIAAANyCAYAAACQc7zrAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuNSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/xnp5ZAAAACXBIWXMAABYlAAAWJQFJUiTwAADG+klEQVR4nOzdd5gUVdbH8d8hpyGJCAICihkVQVByNKKImRUVDOual33VXeMCuoi6ZlExLZjXsCrKqqDIoIABwbAKKCpBFEElg6SZ8/5RVTM9Pd0Tu5kZ+H6ep5+eubfq1q2qW9VVp2/fMncXAAAAAAAAAADpUqmsKwAAAAAAAAAA2LERiAYAAAAAAAAApBWBaAAAAAAAAABAWhGIBgAAAAAAAACkFYFoAAAAAAAAAEBaEYgGAAAAAAAAAKQVgWgAAAAAAAAAQFoRiAYAAAAAAAAApBWBaAAAAAAAAABAWhGIBgAAAAAAAACkFYFoAAAAAAAAAEBaEYgGAAAAAAAAAKQVgWgAAAAAAAAAQFpVKesKAEB5ZGYdJfWXtErSfe7uZVwlAAAAIK3M7I+Smkma6e6Ty7o+2DFxrwXsvOgRDQBxzKyhpJclXS/pSy6MAAA7KzOrbWbzzMzN7OKyrg+A9DGz0yU9IukMSR+XcXWwg+JeC9i5EYgGgPwel9Rc0iXuPqWsK1Nemdn4MDAxYjsvt5qZ3RgGRjaFdSjTC1gza7U96xEty8xaxaUPDdMzt0c9sPMysxFhWxtf1nXZHtKxvqU9XpOdB9LgAUn7SbrL3R9K87KQBOf38qegfWJmmWHe0O1fs5Ixs5YKgtC/SOrv7qtTVG6JzlXJtuH2vuZKxMwWhXXoVVZ1qOAKvdfinAfsuAhEA9guYoKW8a+1ZvaZmf3TzJqXg3peImmgpDvd/dEyrg4Se0DSTQoCI9skLQ9fAIAUMrOzJQ2R9Iqkq9NQ/ojwVT/VZW9PZlbLzL4Pr2ueKmTa/uF02WbWY3vVcUdiZvWjtpPGZQwLl9EqXcsoT8ysiqRnJVWXdKK7f1/GVcIOinstAIwRDWB72yppZfi3SdpV0iHh6wIzO8Hdp5dFxcysraQ7Jb0m6a9lUYcKZpmkryX9ur0WaGb1JA0N/z3F3V/eXssuZ74O37fGpa8J85Zs3+oAKAPJzgMpYWb7SnpQ0ixJg909Ow2LGR6+j5e0Og3lbxfuvjEcV/cdSWeZ2bPu/mb8dGZWV9LD4b8Puft7xVgM5/dc9ZXbdkakaRnDJLWUlClpUQnmX6Jgf61JWY3Sa7ikzpIGufsHKS471eeqrTFlogLhXguARCAawPY30917Rf+YWS1Jp0i6T8GNxYtmtqe7/769K+buX0qqub2XW1G5+7WSrt3Oi91XwWfXbztxEFruvl+S9FcU9FwEsINLdh5IoYMk3SFpbFl8Jlc07j7FzB6XdL6ksWbW1t3XxU12p4KHwC2RdE0xy+f8XoG4+zllXYficPcbJd2YprJTeq5y9x8V/CoOFQz3WgAkhuYAUMbcfaO7PyXpijCpiYKfawGJRBev68u0FgCwg3P3l9x9pLsz9FHRXaXg10J7SLo1NsPM+kq6IPz3ogRBagAAgB0egWgA5cULkqKf/XaIzTCzuuE4fZ+b2frw9YWZjQyHasgn9sFSZlbJzC4zs4/NbHWY3i5m2upm9n9m9pGZrTGz383sazO7y8yaJCj77LCMWQnyGoXjPrqZJfpZ7r5h3iYzqxGTnvNAFjOrGdb/67AuK8zs32a2d5J1rW5mp5nZk+E2+jUsf7GZPWNmHRLNF86b87AVM2sYrvNCM9tsZj+a2aNm1jTJvIU+rNDMTjCzCWb2s5ltCdfldTM7Otk8ScoZasFDaTLDpJZxY40PDacr9IFiyeodbgM3s0Xh/13NbGK4PX8Pt+1lZmbFqXtY1iFmtjws/+lwLMbY/KpmdqGZTTGzX8Ltv9jMJofpteOmT+nDChOs+9Fm9o6ZrQyPmbfNrHPM9PXMbJSZfRNumx/M7DYzS9jLxcyam9lVZvaWmS0ws40WjA//qQXHcf0C6tbazB6KWdbGcNtkmtm1ZtYonK52WKab2fEFlGdhG3czu7AY2yj2WNnDzB4L13tTWN4dlvx8lNPmLDher7fgHLYuTK8fTpcR7sMXzOzLcNv/bmbfmtkjluAcYMHDOzeG5RyYIP/1mPayW4L8D2KPobi8fc3sOQuO29/NbL6ZDTez6kXYXsU6r4bz5Gm/ZjbEzD4M9+ua8Pg4JtXLLav1TVBGsdbXkpwHYvLrmNl1ZjYrLG+TBcfffWbWIsk8Jf4sCuc/0czesOB8t9WCc8jX4XY9I2a68Zb3QWMLLe85fXw43ePh//9MsKxTYqb/W4L8P1mC86GV4nxUmPDhbpeE/15sZt3DZdaWFI2F+mSiYTsKE398lGB+M7MzzOy/FnwmR5/z75nZX8xsl5hpo/NdYa+hCZYT3+42mtlXZnazBUOTJKpbpXD9pprZb2Hb+SWc71+xx0G4/gtj/o+v04gE5fc2s5ct91rkZzN7xcz6JJh2RNg2W4ZJU+PKzyzi9s73oD0z6x6mrUiyDaJr1HkJ8uuE2yXhMW9mh1pwffFDuG9/NbNJZnZKAXWsZmZ/NrOZ4bK3hsfu52b2gMV87sfNl5JrlpKyIjys0Mx2CY/n2eG6bbTgOuLfZjYwZrpeVrS2XtCyinVNUEA5afs8L2CZsef8BmZ2twXj3W8ys6UWXHskuw9I271WkuUV+5ognO8YM3spXJ/NFhz/H5rZDZb8s7CtBeeeheG2WG1mM8zsIjOrWpT6AiiAu/PixYtX2l8Kxn90SZkFTLM8nOaRmLQ2Csbm8/C1IXxF/y+WtHeCskaE+U9IejX8e5ukVeHf7cLpdpU0J6a8TZLWxvy/UtIRcWXvEVNeRlzeyTHzrpFUOS7/wjBvWlx6Zph+RUx9NknaGFPeb5L2SrCux8dMkx3W+feYtK2Szk6yzaNte1bM3xvCZUfzL5TUoIB9OiJBXlVJT8eUEW2P2P9vK0b7OUPSz+G6uaSs8P/odUbcfh9fhLY4Ii69V5i+SME41NvC7bk6rt73JCizVZSfIK9LTLt7UJLF5TeT9GlM+Vnhvt4ck9Yrbp4ovVVc+lAVcpwl2Sax635JuN5Zcfvsd0ldFRwz/wvT1sfVc2KS8l+KmWZzuH5ZMWnfSmqeYL72yns8bonZltHrmJjpHwnT/lPAuvZVbjuvW4xttCic7wJJK8K/1ynvsbZAUtMC2tytkj6KWZeobdUPp7sspqxtCdrBekn9EpT/bph/cVx6pbjtdVpcfm0F5weX1Dour4fynmvXxNRlpqRblORYUwnOq/HtV9Ldyj0eVilok9H8VyXZRyVabgVe34TngTBvf+X97Nwatp/YOnVNMF+mSv5ZNCpmGg+3Qezx8XPMtPcqOHdHeb8o7zn93nC6c8L8jxIs776Y+f+bIP/ZMG9kKs5HxTynPh+W9bWkGuH6erhuDUtYZk57KcG89SS9HbOOia4VhsZMPytuf8S+ViSaJ0m726i8n6HfStojQf2eiWs7q5X33PdhzLQvh+0lp13Fva6KK/sfcesdf3yNjpv+qrCcqE2sjCv/5aLsE+UeS7HbtXrMNt8/bvr2cdugcVz+UWH64gTLujCuDa9S8BkS/f+U8l+PVompY+y2iZ3v3wmWlbJrliK023zbMExvFZWZZL7uCp5fEn+cb4ufT8E1WrK2/rOCz/l8y1IprgkKWeeUf54XcRtfqeD4jI7b2M+LFfHtNZx3RJif8nut+ONLJfuMrKag7cefW2LXbUSC+S5T3uNpnfIeF1Ml1SrOfuXFi1feV5lXgBcvXjvHS4UEohUMuRBdVNweplWT9HmYtkTSkQoecGgKgkmLw7wvJVWPKy+6OFoXXvBcHF00SGqsMAAl6c2Yi6DTFF6oSzpM0hfKvclpFFf+ojDvmLj06GYzusA6LC4/utm6KS49U7k3DwslHS2psoILz+6SfgjzX0iw7XqFy+0ee2GkIGAeXbj9rsQ3f4tilvuppM5hehVJA5R7MXl7Aft0RIK8aLkLwu1aO0zPCPdFtH3
|
||
|
|
"text/plain": [
|
||
|
|
"<Figure size 864x504 with 1 Axes>"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
"execution_count": 5,
|
||
|
|
"metadata": {
|
||
|
|
"image/png": {
|
||
|
|
"height": 441,
|
||
|
|
"width": 721
|
||
|
|
},
|
||
|
|
"needs_background": "light"
|
||
|
|
},
|
||
|
|
"output_type": "execute_result"
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"# Definiujemy dystrybuantę zmiennej losowej X\n",
|
||
|
|
"def F(x):\n",
|
||
|
|
" n = len(x)\n",
|
||
|
|
" y = [0]*n\n",
|
||
|
|
" for i in range(n):\n",
|
||
|
|
" if x[i] < 0:\n",
|
||
|
|
" y[i] = 0\n",
|
||
|
|
" elif x[i] < 1:\n",
|
||
|
|
" y[i] = 0.8\n",
|
||
|
|
" else:\n",
|
||
|
|
" y[i] = 1\n",
|
||
|
|
" return y\n",
|
||
|
|
"\n",
|
||
|
|
"# Definiujemy dystrybuantę pseudoodwrotną zmiennej losowej X\n",
|
||
|
|
"def F_inv(u):\n",
|
||
|
|
" n = len(u)\n",
|
||
|
|
" y = [0]*n\n",
|
||
|
|
" for i in range(n):\n",
|
||
|
|
" if u[i] < 0.8:\n",
|
||
|
|
" y[i] = 0\n",
|
||
|
|
" else:\n",
|
||
|
|
" y[i] = 1\n",
|
||
|
|
" return y\n",
|
||
|
|
"\n",
|
||
|
|
"# Generujemy próbkę 100 liczb używając metody odwracania dystrybuanty\n",
|
||
|
|
"u = stats.uniform.rvs(size=100)\n",
|
||
|
|
"probki = F_inv(u)\n",
|
||
|
|
"\n",
|
||
|
|
"# Aby porównać rozkład zmiennej losowej X z rozkładem wylosowanych liczb możemy porównać funkcję masy prawdopodobieństwa zmiennej losowej X z częstotliwością wylosowanych liczb\n",
|
||
|
|
"atomy = [0, 1]\n",
|
||
|
|
"pr = [0.8, 0.2]\n",
|
||
|
|
"\n",
|
||
|
|
"jedynki = np.count_nonzero(probki)\n",
|
||
|
|
"zera = 100 - jedynki\n",
|
||
|
|
"czestotliwosc = [zera/100, jedynki/100]\n",
|
||
|
|
"\n",
|
||
|
|
"plt.bar(atomy, pr,label='Funkcja masy prawdopodobieństwa zmiennej losowej X', width=0.3, color='blue', align='center')\n",
|
||
|
|
"plt.bar(atomy, czestotliwosc, label='Częstotliwość liczb w próbce', width=0.3, color='red', align='edge')\n",
|
||
|
|
"plt.legend()\n",
|
||
|
|
"plt.xlabel('Atomy')\n",
|
||
|
|
"plt.ylabel('Prawdopodobieństwa / częstotliwości')\n",
|
||
|
|
"plt.title('Porównanie funkcji masy prawdopodobieństwa X i częstotliwości liczb w próbce')\n",
|
||
|
|
"plt.grid()\n",
|
||
|
|
"plt.show()"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"id": "90ade4",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"Ogólnie, w przypadku zmiennej losowej dyskretnej o rozkładzie $\\mathbb{P}(X=x_i)=p_i$ dla $i=1,2,\\ldots,n$ dystrybuantę pseudoodwrotną możemy przedstawić graficznie w następujący sposób:\n",
|
||
|
|
"\n",
|
||
|
|
"\n",
|
||
|
|
"\n",
|
||
|
|
"Zatem jeśli chcemy znaleźć wartość $F^-(u)$, możemy podzielić przedział $(0,1)$ na mniejsze przedziały, których długości są równe kolejnym wartościom $p_i$ i które odpowiadają kolejnym atomom $x_i$. Następnie wystarczy sprawdzić, do którego z przedziałów należy $u$."
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "markdown",
|
||
|
|
"id": "3f5655",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"source": [
|
||
|
|
"## Dyskretne zmienne losowe wielowymiarowe\n",
|
||
|
|
"\n",
|
||
|
|
"W celu badania dyskretnych zmiennych losowych wielowymiarowych w Pythonie, będziemy posługiwać się macierzami, czyli obiektami klasy `ndarray`. Sposób ich generowania został już omówiony na początku poprzedniego przewodnika *03LRAP*. Warto jednak znać klika innych przydatnych funkcji z biblioteki NumPy. \n",
|
||
|
|
"\n",
|
||
|
|
"* `A.sum(axis=k)` - jeśli *A* jest obiektem klasy `ndarray`, to ta funkcja zwróci sumy elementów tej macierzy. Opcjonalny argument `axis` pozwala sumować elementy z poszczególnych kolumn (jeśli $k=0$) lub wierszy (jeśli $k=1$);\n",
|
||
|
|
"* `numpy(A, weights=None)` - oblicza średnią wartość elementów z macierzy *A*. W opcjonalnym argumencie `weights` możemy podać wagi, aby wyznaczyć średnią ważoną. Wagi powinny być podane jako macierz liczb o tych samych wymiarach, co macierz *A*.\n",
|
||
|
|
"\n",
|
||
|
|
"**Przykład 4**\n",
|
||
|
|
"\n",
|
||
|
|
"Dany jest dyskretny wektor losowy $(X,Y)$ o rozkładzie podanym w poniższej tabeli:\n",
|
||
|
|
"\n",
|
||
|
|
"| $X$ \\ $Y$ | $1$ | $2$ | $3$ | $4$ |\n",
|
||
|
|
"| --- | --- | --- | --- | --- |\n",
|
||
|
|
"| $-2$ | $0{,}1$ | $0{,}2$ | $0$ | $0{,}1$ |\n",
|
||
|
|
"| $2$ | $0{,}3$ | $0$ | $0{,}2$ | $0{,}1$ |\n",
|
||
|
|
"\n",
|
||
|
|
"* Znajdź rozkład brzegowy zmiennej losowej $Y$.\n",
|
||
|
|
"* Oblicz $\\mathbb{E}Y$ oraz $Var(Y)$.\n",
|
||
|
|
"* Oblicz $\\mathbb{E}(\\sqrt{|XY|})$.\n"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 10,
|
||
|
|
"id": "4f9d77",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
{
|
||
|
|
"name": "stdout",
|
||
|
|
"output_type": "stream",
|
||
|
|
"text": [
|
||
|
|
"Rozkład brzegowy zmiennej losowej Y:\n",
|
||
|
|
"P(Y=1)=0.4\n",
|
||
|
|
"P(Y=2)=0.2\n",
|
||
|
|
"P(Y=3)=0.2\n",
|
||
|
|
"P(Y=4)=0.2\n",
|
||
|
|
"EY=2.2\n",
|
||
|
|
"Var(Y)=1.3599999999999994\n",
|
||
|
|
"E(sqrt(|XY|))=2.021268798455112\n"
|
||
|
|
]
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
"# Definiujemy wektory atomów zmiennych losowych X i Y oraz tabelkę z prawdopodobieństwem łącznym\n",
|
||
|
|
"atomy_X = np.array([[-2], [2]])\n",
|
||
|
|
"atomy_Y = np.array([[1, 2, 3, 4]])\n",
|
||
|
|
"pr_laczne = np.array([[0.1, 0.2, 0, 0.1], [0.3, 0, 0.2, 0.1]]) # Zakładamy, że wiersze tej tabeli odpowiadają atomom zmiennej losowej X, a kolumny atomom zmiennej losowej Y\n",
|
||
|
|
"\n",
|
||
|
|
"# Szukamy rozkładu brzegowego zmiennej losowej Y\n",
|
||
|
|
"pr_Y = np.array([pr_laczne.sum(axis=0)])\n",
|
||
|
|
"print('Rozkład brzegowy zmiennej losowej Y:')\n",
|
||
|
|
"for k in range(4):\n",
|
||
|
|
" print('P(Y=', atomy_Y[0, k], ')=', pr_Y[0, k], sep=\"\")\n",
|
||
|
|
" \n",
|
||
|
|
"# Szukamy EY oraz Var(Y)\n",
|
||
|
|
"EY = np.average(atomy_Y, weights=pr_Y)\n",
|
||
|
|
"print('EY=', EY, sep=\"\")\n",
|
||
|
|
"EY2 = np.average(atomy_Y**2, weights=pr_Y)\n",
|
||
|
|
"VarY = EY2 - EY**2\n",
|
||
|
|
"print('Var(Y)=', VarY, sep=\"\")\n",
|
||
|
|
"\n",
|
||
|
|
"# Szukamy E(sqrt(|XY|)). Skorzystamy z ,,prawa leniwego statystyka''\n",
|
||
|
|
"XY = np.dot(abs(atomy_X)**(1/2), atomy_Y**(1/2)) # pomocnicza tabelka, która zawiera wszystkie iloczyny postaci (pierwiastek z |x_i|)*(pierwiastek z y_j), gdzie x_i, y_j to odpowiednio atomy zmiennych losowych X i Y\n",
|
||
|
|
"E = np.average(XY, weights=pr_laczne)\n",
|
||
|
|
"print('E(sqrt(|XY|))=', E, sep=\"\")"
|
||
|
|
]
|
||
|
|
},
|
||
|
|
{
|
||
|
|
"cell_type": "code",
|
||
|
|
"execution_count": 0,
|
||
|
|
"id": "a1029e",
|
||
|
|
"metadata": {
|
||
|
|
"collapsed": false
|
||
|
|
},
|
||
|
|
"outputs": [
|
||
|
|
],
|
||
|
|
"source": [
|
||
|
|
]
|
||
|
|
}
|
||
|
|
],
|
||
|
|
"metadata": {
|
||
|
|
"kernelspec": {
|
||
|
|
"argv": [
|
||
|
|
"/usr/bin/python3",
|
||
|
|
"-m",
|
||
|
|
"ipykernel",
|
||
|
|
"--HistoryManager.enabled=False",
|
||
|
|
"--matplotlib=inline",
|
||
|
|
"-c",
|
||
|
|
"%config InlineBackend.figure_formats = set(['retina'])\nimport matplotlib; matplotlib.rcParams['figure.figsize'] = (12, 7)",
|
||
|
|
"-f",
|
||
|
|
"{connection_file}"
|
||
|
|
],
|
||
|
|
"display_name": "Python 3 (system-wide)",
|
||
|
|
"env": {
|
||
|
|
},
|
||
|
|
"language": "python",
|
||
|
|
"metadata": {
|
||
|
|
"cocalc": {
|
||
|
|
"description": "Python 3 programming language",
|
||
|
|
"priority": 100,
|
||
|
|
"url": "https://www.python.org/"
|
||
|
|
}
|
||
|
|
},
|
||
|
|
"name": "python3",
|
||
|
|
"resource_dir": "/ext/jupyter/kernels/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.10.12"
|
||
|
|
}
|
||
|
|
},
|
||
|
|
"nbformat": 4,
|
||
|
|
"nbformat_minor": 4
|
||
|
|
}
|