{ "cells": [ { "cell_type": "markdown", "id": "ff30f263", "metadata": {}, "source": [ "# Математическая статистика (ФБМФ, ФМХФ)\n", "## Домашнее задание 1 — часть A\n", "\n", "**Правила, прочитайте внимательно:**\n", "\n", "* Выполненную работу нужно отправить телеграм-боту `@thetahat_st_bot`. Для начала работы с ботом каждый раз отправляйте `/start`. Дождитесь подтверждения от бота, что он принял файл. Если подтверждения нет, то что-то не так. **Работы, присланные иным способом, не принимаются.**\n", "* Дедлайн см. в боте. После дедлайна работы не принимаются **вообще никак**, кроме случаев наличия уважительной причины.\n", "* До дедлайна можно поменять решение любое количество раз. Начинайте точно так же сдавать задание, бот подскажет.\n", "* Любую уважительную причину нужно подтвердить документально, отправив скан или фото боту. При этом работу можно сдать позже на столько дней, на сколько время ее действия пересекается с временем выполнения задания.\n", "* Прислать нужно **ноутбук в формате ipynb**. Другие форматы не принимаются.\n", "* Выполнять задание необходимо **полностью самостоятельно**. При обнаружении списывания **все участники списывания будут сдавать устный зачет.**\n", "* Решения, размещенные на каких-либо интернет-ресурсах не принимаются. Кроме того, публикация решения в открытом доступе может быть приравнена к предоставлении возможности списать.\n", "* Простой или основной уровень вы выбираете самостоятельно, выполняя или не выполняя задания типа B. При выборе **простого уровня** достаточно выполнить задания *типа A*. При выборе **основного уровня** нужно выполнять *как задания типа A, так и задания типа B*.\n", "* Для выполнения задания используйте этот ноутбук в качествие основы, ничего не удаляя из него. Можно добавлять необходимое количество ячеек. Ячейки с assert'ами **удалять и изменять нельзя**, в противном случае соответствующее задание не будет оценено.\n", "* Комментарии к решению пишите в markdown-ячейках.\n", "* **Если код студента недописан и т.д., то он не оценивается.**\n", "* Каждая задача стоит **5 баллов**.\n", "\n", "_Замечание: перед выполнением задания можно ознакомиться с ноутбуком с семинара._\n" ] }, { "cell_type": "markdown", "id": "43d95d55", "metadata": {}, "source": [ "Данная часть задания проверяется автоматически. Для выполнения задания используйте этот ноутбук в качестве основы, ничего не удаляя из него. Можно добавлять необходимое количество ячеек. Ячейки с assert'ами удалять и изменять нельзя, в противном случае соответствующее задание *не будет оценено*.\n", "\n", "\n", "> *Примечание.* Рекомендуется работать с данным ноутбуком **локально** в Jupyter Notebook (например, используя Anaconda или альтернативные варианты). Если вы используете софт по типу Google Colaboratory, то перед отправкой боту данного ноутбука необходимо проверить, что **в ячейках с assert'ами и \"# Ваше решение тут\"** в метаданных присутствуют поля *metadata* с `nbgrader`. Можно открыть ноутбук с решением в текстовом редакторе (MS Word, Блокнот) и выполнить поиск по документу слова `nbgrader`. Если поиск показал ровно **25 совпадений** — можете отправлять файл боту. Если совпадений меньше, решение может быть не оценено. В таком случае попробуйте скачать файл в форме `ipynb` еще раз и перенесите решения в новый файл. **Внимание! Бот не проверяет решение и не проверяет наличие метаданных.** " ] }, { "cell_type": "code", "execution_count": null, "id": "2e5011a1", "metadata": {}, "outputs": [], "source": [ "# Bot check\n", "\n", "# HW_ID: st_1a\n", "# Бот проверит этот ID и предупредит, если случайно сдать что-то не то.\n", "\n", "# Status: not final\n", "# Перед отправкой в финальном решении удали \"not\" в строчке выше.\n", "# Так бот проверит, что ты отправляешь финальную версию, а не промежуточную.\n", "# Никакие значения в этой ячейке не влияют на факт сдачи работы." ] }, { "cell_type": "code", "execution_count": null, "id": "690a6e0d", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import scipy.stats as sps\n", "import matplotlib.pyplot as plt\n", "from scipy.special import factorial, gammaln" ] }, { "cell_type": "markdown", "id": "0a3e104d", "metadata": {}, "source": [ "### Задача 1.1\n", "\n", "Пусть $X_1, ..., X_n$ — выборка пуассоновского распределения $Pois(\\theta)$, то есть $\\mathsf{P}(X_i = k) = \\frac{\\theta^k}{k!} e^{-\\theta}$ при $k \\in \\{0, 1, 2, ...\\}$.\n", "Реализуйте функции, вычисляющие:\n", "* логарифм правдоподобия;\n", "* градиент логарифма правдоподобия;\n", "* оценку $\\theta$ по методу максимального правдоподобия.\n", "\n", "_Замечание: функция вычисления логарифма факториала уже реализована ниже_\n", "\n", "*Примечание:* вам могут пригодиться методы библиотеки NumPy для работы с массивами и numpy.random или scipy.stats — со случайными величинами. Задача с семинара о распределениях. Задача с семинара о работе с массивами." ] }, { "cell_type": "code", "execution_count": null, "id": "8b76a230", "metadata": {}, "outputs": [], "source": [ "def logfactorial(x):\n", " return gammaln(x + 1)" ] }, { "cell_type": "code", "execution_count": null, "id": "1d9e7c25", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "eabbae1ff0eabca8e9fd95e4d213874c", "grade": false, "grade_id": "cell-a6eee546e2ec3a6c", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def poiss_loglikelihood(x, theta):\n", " ... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "1fdcdfe7", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "16dc0ac9eb870d1795c5c5625ef49889", "grade": true, "grade_id": "cell-64887c82a45dea49", "locked": true, "points": 0.5, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "x = np.array([1, 2, 3])\n", "assert round(poiss_loglikelihood(x, 1.9), 2) == -4.33\n", "assert poiss_loglikelihood(x, 1.9) > poiss_loglikelihood(x, 4.2)" ] }, { "cell_type": "code", "execution_count": null, "id": "ecb052ef", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "8f6b72348ed34ce277e81d0234a251fd", "grade": false, "grade_id": "cell-db5e85f6a0200a00", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def poiss_loglikelihood_grad(x, theta):\n", " ... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "0e97e5fb", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "780a9e644fb78ae1281f0e697a99e90e", "grade": true, "grade_id": "cell-09d1a11f5f5dd198", "locked": true, "points": 1, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "x = np.array([1, 2, 3])\n", "assert round(poiss_loglikelihood_grad(x, 1.9), 2) == 0.16" ] }, { "cell_type": "code", "execution_count": null, "id": "93922628", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "13320cc0f10f8381af4f3fc04e81310b", "grade": false, "grade_id": "cell-c3a633633724fa95", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def poiss_maxlikelihood_estimator(x):\n", " ... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "61f71f24", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "3dc7a9590fbaf6e2a3824a0e9c23a441", "grade": true, "grade_id": "cell-631e62c697b22613", "locked": true, "points": 1, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "x = np.array([1, 2, 3])\n", "theta = poiss_maxlikelihood_estimator(x)\n", "assert np.allclose(poiss_loglikelihood_grad(x, theta), 0)\n", "assert np.allclose(poiss_maxlikelihood_estimator(x), 2)" ] }, { "cell_type": "markdown", "id": "47f63aab", "metadata": {}, "source": [ "### Задача 1.2\n", "\n", "Дана выборка $X_1, ..., X_n$ из нормального распределения $\\mathcal{N}(a, \\sigma^2)$. Реализуйте функции, вычисляющие:\n", "* логарифм правдоподобия для параметра $\\theta = (a, \\sigma)$;\n", "* градиент логарифма правдоподобия $\\theta$;\n", "* оценку $\\theta$ по методу максимального правдоподобия.\n", "\n", "*Примечание:* вам могут пригодиться методы библиотеки NumPy для работы с массивами и numpy.random или scipy.stats — со случайными величинами. Задача с семинара о распределениях. Задача с семинара о работе с массивами." ] }, { "cell_type": "code", "execution_count": null, "id": "5081c5f7", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "b684507ea0d85b59ce8d3d1ca91f8afe", "grade": false, "grade_id": "cell-1aea7dcb70b4f5cb", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def norm_loglikelihood(x, a, sigma):\n", " ... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "8c09822f", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "4df4de5d8c21b521ed8d1ddff546980e", "grade": true, "grade_id": "cell-cc15bd20d5e0e6d4", "locked": true, "points": 0.5, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "x = np.array([1, 2, 3])\n", "assert round(norm_loglikelihood(x, 0, 1), 2) == -9.76\n", "assert norm_loglikelihood(x, 0, 1) < norm_loglikelihood(x, 1, 1)\n", "assert norm_loglikelihood(x, 0, 1) > norm_loglikelihood(x, 0, 100)" ] }, { "cell_type": "code", "execution_count": null, "id": "a465b735", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "ac9bb5c1a4faf0de8e3068a067679a38", "grade": false, "grade_id": "cell-459261bfe20abe10", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def norm_loglikelihood_grad(x, a, sigma):\n", " ... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "2c7937b2", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "14c2988ea6b57b39a8454b1c9ed30132", "grade": true, "grade_id": "cell-eb271e14e658f13d", "locked": true, "points": 1, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "x = np.array([1, 2, 3])\n", "assert np.allclose(norm_loglikelihood_grad(x, 1.9, 3), [0.03, -0.92], atol=0.01)" ] }, { "cell_type": "code", "execution_count": null, "id": "1be9bb0b", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "89b521cf7a90ea62fe5da8634aedb4ba", "grade": false, "grade_id": "cell-78bc50d2c7882218", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def norm_maxlikelihood_estimator(x):\n", " ... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "798deb87", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "875df9b54c4e36fa10951a70965c22e6", "grade": true, "grade_id": "cell-c076cb0e133a2eb6", "locked": true, "points": 1, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "x = np.array([1, 2, 3])\n", "a, sigma = norm_maxlikelihood_estimator(x)\n", "assert np.allclose(norm_loglikelihood_grad(x, a, sigma), 0)\n", "assert np.allclose([a, sigma], [2.00, 0.82], atol=0.01)" ] }, { "cell_type": "markdown", "id": "383aef70", "metadata": {}, "source": [ "---" ] }, { "cell_type": "markdown", "id": "bbfc715e", "metadata": {}, "source": [ "__________________\n", "### Задача 2.1\n", "Даны 6 параметрических моделей — для каждой нужно реализовать подсчёт функции правдоподобия по небольшой выборке, чтобы затем построить график функции правдоподобия.\n", "\n", "*a).* Параметрическая модель $\\mathcal{N}(\\theta, 1)$.\n", "\n", "*b).* Параметрическая модель $Exp(\\theta)$.\n", "\n", "*c).* Параметрическая модель $U[0, \\theta]$.\n", "\n", "*d).* Параметрическая модель $Bin(5, \\theta)$.\n", "\n", "*e).* Параметрическая модель $Pois(\\theta)$.\n", "\n", "*f).* Параметрическая модель $Сauchy(\\theta)$, где $\\theta$ — параметр сдвига.\n", "\n", "*Примечание:* вам могут пригодиться методы библиотеки NumPy для работы с массивами и numpy.random или scipy.stats — со случайными величинами. Для построения графиков используется Matplotlib." ] }, { "cell_type": "code", "execution_count": null, "id": "f19e76d1", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "8f1d7e6798cd47df9e761af092b080f6", "grade": false, "grade_id": "cell-80a1db70ab3d2c30", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def calc_likelihood(dist_name, theta_grid, sample):\n", " theta_grid = theta_grid.reshape((-1, 1))\n", " sample = sample.reshape((1, -1))\n", "\n", " if dist_name == \"normal\":\n", " return sps.norm(loc=theta_grid).pdf(sample).prod(axis=1)\n", " elif dist_name == \"expon\":\n", " ... # Ваше решение тут\n", " elif dist_name == \"uniform\":\n", " ... # Ваше решение тут\n", " elif dist_name == \"binomial\":\n", " ... # Ваше решение тут\n", " elif dist_name == \"poisson\":\n", " ... # Ваше решение тут\n", " elif dist_name == \"cauchy\":\n", " ... # Ваше решение тут\n", "\n", " assert False" ] }, { "cell_type": "code", "execution_count": null, "id": "3af95772", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "11e5939b607472ca700ba3774b00882c", "grade": true, "grade_id": "cell-9755b78c4d21f635", "locked": true, "points": 2, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "check_values = {\n", " 'normal': {'sample': [-1, 1], 'theta': [1.0, 2.0], 'likelihood': [0.0215, 0.0011]},\n", " 'expon': {'sample': [1, 2], 'theta': [6.0, 7.0], 'likelihood': [0.0168, 0.0133]},\n", " 'uniform': {'sample': [0.2, 0.8], 'theta': [1.8, 2.1], 'likelihood': [0.3086, 0.2268]},\n", " 'binomial': {'sample': [5, 5], 'theta': [0.6, 0.7], 'likelihood': [0.006, 0.0282]},\n", " 'poisson': {'sample': [5, 10], 'theta': [6.04, 7.03], 'likelihood': [0.0068, 0.0091]},\n", " 'cauchy': {'sample': [-0.5, 0.5], 'theta': [1.0, 2.0], 'likelihood': [0.0249, 0.0043]}\n", "}\n", "for dist_name, params in check_values.items():\n", " ans = calc_likelihood(dist_name, np.array(params[\"theta\"]), np.array(params[\"sample\"]))\n", " ref = np.array(params[\"likelihood\"])\n", " assert ans.shape == ref.shape\n", " assert np.allclose(ans, ref, atol=1e-4)" ] }, { "cell_type": "markdown", "id": "a870cef4", "metadata": {}, "source": [ "Посмотрим на графики функций правдоподобия:" ] }, { "cell_type": "code", "execution_count": null, "id": "006f6da9", "metadata": {}, "outputs": [], "source": [ "dist2samples = {\n", " \"normal\": [[-1, 1], [-5, 5], [-1, 5]],\n", " \"expon\": [[1, 2], [0.1, 1], [1, 10]],\n", " \"uniform\": [[0.2, 0.8], [0.5, 1], [0.5, 1.3]],\n", " \"binomial\": [[0, 1], [5, 5], [0, 5]],\n", " \"poisson\": [[0, 1], [0, 10], [5, 10]],\n", " \"cauchy\": [[-0.5, 0.5], [-2, 2], [-4, 0, 4]],\n", "}\n", "dist2grid = {\n", " \"normal\": np.linspace(-5, 5, 1000),\n", " \"expon\": np.linspace(0.01, 10, 1000),\n", " \"uniform\": np.linspace(0.01, 3, 1000),\n", " \"binomial\": np.linspace(0, 1, 1000),\n", " \"poisson\": np.linspace(0.1, 10, 1000),\n", " \"cauchy\": np.linspace(-5, 5, 1000),\n", "}\n", "dist2label = {\n", " \"normal\": r\"$\\mathcal{N}(\\theta, 1)$\",\n", " \"expon\": r\"$Exp(\\theta)$\",\n", " \"uniform\": r\"$U[0, \\theta]$\",\n", " \"binomial\": r\"$Bin(5, \\theta)$\",\n", " \"poisson\": r\"$Pois(\\theta)$\",\n", " \"cauchy\": r\"$Сauchy(\\theta)$\",\n", "}\n", "\n", "for dist_name in dist2samples.keys():\n", " label = dist2label[dist_name]\n", "\n", " plt.figure(figsize=(18, 5))\n", " grid = dist2grid[dist_name]\n", " for i, sample in enumerate(dist2samples[dist_name]):\n", " sample = np.array(sample).reshape((1, -1))\n", " likelihood = calc_likelihood(dist_name, grid, sample)\n", "\n", " plt.subplot(1, 3, i+1)\n", " plt.plot(grid, likelihood)\n", " plt.xlabel('$\\\\theta$', fontsize=16)\n", " plt.grid(ls=':')\n", " plt.title(label + ', sample=' + str(sample), fontsize=16)\n", " plt.show()" ] }, { "cell_type": "markdown", "id": "19d2ffff", "metadata": {}, "source": [ "Сделайте вывод о том, как функция правдоподобия для каждой модели зависит от выборки." ] }, { "cell_type": "markdown", "id": "d9fdd006", "metadata": {}, "source": [ "**Вывод:**\n", "\n", "<...>" ] }, { "cell_type": "markdown", "id": "2295b378", "metadata": {}, "source": [ "Является ли функция правдоподобия плотностью? Имеет ли она единственный максимум? Дайте ответы на эти вопросы в переменных следующей ячейке, записав в соответствующую переменную либо название распределения, для которого это свойство не выполняется (если таких несколько, можете вписать любое), либо `None`, если свойство верно всегда." ] }, { "cell_type": "code", "execution_count": null, "id": "7cc6d14a", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "fd0eaf749450363157a1a9957b47e9ff", "grade": false, "grade_id": "cell-8eb8f1f748b0130e", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "non_density = ...\n", "has_single_maximum = ...\n", "... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "cbd421e0", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "4759279a1b8dad4661c00e9d95d93176", "grade": true, "grade_id": "cell-f64f3adb82856fab", "locked": true, "points": 2, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "assert non_density in [None, 'normal', 'expon', 'uniform', 'binomial', 'poisson', 'cauchy']\n", "assert has_single_maximum in [None, 'normal', 'expon', 'uniform', 'binomial', 'poisson', 'cauchy']\n", "# А тут скрытые assert'ы :)" ] }, { "cell_type": "markdown", "id": "a698006e", "metadata": {}, "source": [ "### Задача 2.2\n", "\n", "Дана функция, которая по выборке $(X_1, \\ldots, X_n)$ и двум числам $\\mu_0, \\mu_1$ определяет, какое из двух распределений — $\\mathcal{N}(\\mu_0, 1)$ или $\\mathcal{N}(\\mu_1, 1)$ — более точно описывает выборку, путём сравнения функций правдоподобия:" ] }, { "cell_type": "code", "execution_count": null, "id": "cdb19e34", "metadata": {}, "outputs": [], "source": [ "def select(x, u0, u1):\n", " prob0 = sps.norm(loc=u0).pdf(x).prod()\n", " prob1 = sps.norm(loc=u1).pdf(x).prod()\n", " if prob0 > prob1:\n", " return 0\n", " else:\n", " return 1" ] }, { "cell_type": "markdown", "id": "934359f2", "metadata": {}, "source": [ "Пример работы для выборки размера 30 из $\\mathcal{N}(0.1, 1)$:" ] }, { "cell_type": "code", "execution_count": null, "id": "2462cbda", "metadata": {}, "outputs": [], "source": [ "np.random.seed(1)" ] }, { "cell_type": "code", "execution_count": null, "id": "6984314e", "metadata": {}, "outputs": [], "source": [ "select(sps.norm(loc=0.1).rvs(30), u0=0, u1=1)" ] }, { "cell_type": "code", "execution_count": null, "id": "4bd65921", "metadata": {}, "outputs": [], "source": [ "select(sps.norm(loc=0.1).rvs(30), u0=1, u1=0)" ] }, { "cell_type": "markdown", "id": "4c9a8dab", "metadata": {}, "source": [ "Однако она некорректно работает для выборок большого размера:" ] }, { "cell_type": "code", "execution_count": null, "id": "4c53d3ec", "metadata": {}, "outputs": [], "source": [ "select(sps.norm(loc=0.1).rvs(1000), u0=0, u1=1) # returns 1" ] }, { "cell_type": "code", "execution_count": null, "id": "3e39bc47", "metadata": {}, "outputs": [], "source": [ "select(sps.norm(loc=0.1).rvs(1000), u0=1, u1=0) # returns 0" ] }, { "cell_type": "markdown", "id": "dc88e6ab", "metadata": {}, "source": [ "Почему такое происходит?" ] }, { "cell_type": "markdown", "id": "b60cf87d", "metadata": {}, "source": [ "**Ответ:**\n", "\n", "<...>" ] }, { "cell_type": "markdown", "id": "23aeab53", "metadata": {}, "source": [ "Напишите исправленную версию функции, которая также выбирает подходящий параметр на основе значения правдободобия, но работает и для выборок большого размера.\n", "\n", "_Подсказка: обратите внимание на значения функций правдоподобия при маленькой и большой выборке. Нужно использовать некоторый метод класса `sps.norm` модуля scipy.stats._\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "4cd9ca17", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "ceae925b8392e9a8b60959eb6222550c", "grade": false, "grade_id": "cell-fb7b90ef4ab1ba0d", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def select_fixed(x, u0, u1):\n", " ... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "ff85b3a7", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "96a5864e1e007abf265c40f91449cea1", "grade": true, "grade_id": "cell-4205c4dacf84d4c0", "locked": true, "points": 1, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "x = sps.norm(loc=0.1).rvs(1000)\n", "assert select_fixed(x, u0=0, u1=1) == 0\n", "assert select_fixed(x, u0=1, u1=0) == 1" ] }, { "cell_type": "markdown", "id": "eb6a5a22", "metadata": {}, "source": [ "### Задача 3\n", "В этой задаче нужно визуально проверить *свойство состоятельности*.\n", "\n", "Пусть $X_1, ..., X_n$ — выборка из распределения $U(0, \\theta)$.\n", "\n", "Рассмотрим 5 оценок $\\theta$:\n", "- $\\widehat{\\theta}_a = 2\\overline{X}$\n", "- $\\widehat{\\theta}_b = \\max_i X_i$\n", "- $\\widehat{\\theta}_c = 2\\sqrt{\\overline{X^2}}$\n", "- $\\widehat{\\theta}_d = \\sqrt{3\\overline{X^2}}$\n", "- $\\widehat{\\theta}_e = (n + 1) \\min_i X_i$\n", "\n", "Дано множество выборок $X^1, \\dots, X^{300}$ из распределения $U[0, 1]$: $\\; X^j = (X^j_1, \\dots, X^j_{500}), 1 \\leqslant j \\leqslant 300$.\n", "
\n", "По каждой из них посчитайте оценки\n", "$\\widehat{\\theta}_{a,jn} = 2\\frac{X^j_1 + \\dots + X^j_n}{n}$,\n", "$\\widehat{\\theta}_{b,jn} = \\max(X^j_1, \\dots, X^j_n)$,\n", "$\\widehat{\\theta}_{c,jn} = 2 \\cdot \\sqrt{\\frac{\\sum_{i=1}^n {X_{ji}^2}}{n}}$ и т.д.,\n", "для $1 \\leqslant n \\leqslant 500$, то есть оценки параметра $\\theta$ по первым $n$ наблюдениям $j$-й выборки. При написании кода могут помочь функции `numpy.cumsum(axis=...)` и `np.maximum.accumulate(axis=...)` (см. пример в ноутбуке с семинара)." ] }, { "cell_type": "code", "execution_count": null, "id": "f672bc15", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "dafb90f3abb858172286faae7e5fd00f", "grade": false, "grade_id": "cell-e6345d3160da1130", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "x = sps.uniform().rvs((300, 500))\n", "\n", "estimations = [\n", " ..., # theta_a\n", " ..., # theta_b\n", " ..., # theta_c\n", " ..., # theta_d\n", " ..., # theta_e\n", "]\n", "\n", "... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "ef0b1cc6", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "e0fa21fce02737e5c242d8611c6f8fcc", "grade": true, "grade_id": "cell-d49a10ab75e010b3", "locked": true, "points": 2.5, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "assert np.allclose(estimations[0][42, 1], 2 * (x[42, 0] + x[42, 1]) / 2)\n", "assert np.allclose(estimations[1][42, 1], max(x[42, 0], x[42, 1]))\n", "assert np.allclose(estimations[2][42, 1], 2 * ((x[42, 0]**2 + x[42, 1]**2) / 2)**0.5)\n", "assert np.allclose(estimations[3][42, 1], (3 * (x[42, 0]**2 + x[42, 1]**2) / 2)**0.5)\n", "assert np.allclose(estimations[4][42, 1], 3 * min(x[42, 0], x[42, 1]))\n", "assert len(estimations) == 5\n", "assert all(estimations[i].shape == (300, 500) for i in range(5))" ] }, { "cell_type": "markdown", "id": "3fd1d919", "metadata": {}, "source": [ "Для каждой оценки $\\theta^*, \\widehat{\\theta}$ нарисуйте следующий график. Для каждого $j$ нанесите на один график зависимости $\\theta^*_{jn}$ или $\\widehat{\\theta}_{jn}$ от $n$ с помощью `plt.plot`. Каждая кривая должна быть нарисована *одним цветом* с прозрачностью `alpha=0.05`. Поскольку при малых $n$ значения средних могут быть большими по модулю, ограничьте область графика по оси *y* с помощью функции `plt.ylim((min, max))`. Примеры работы с графиками с помощью Matplotlib можно найти в этом ноутбуке." ] }, { "cell_type": "code", "execution_count": null, "id": "80359766", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "05df4a20e6d4f07a4e96192a7dded4bd", "grade": true, "grade_id": "cell-c8262f2f22b5ae50", "locked": false, "points": 0, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "... # Ваше решение тут" ] }, { "cell_type": "markdown", "id": "0421cb94", "metadata": {}, "source": [ "Укажите, для каких оценок, судя по графику, наблюдается свойство состоятельности:" ] }, { "cell_type": "code", "execution_count": null, "id": "6cd6e366", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "d53109dfcd17c059e618453a43a3da57", "grade": false, "grade_id": "cell-14575ad5a70123b1", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "consistent_estimators = {...}\n", "... # Ваше решение тут" ] }, { "cell_type": "code", "execution_count": null, "id": "bce64fde", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "a647f74c5e67d7713dd14c3cbc276b31", "grade": true, "grade_id": "cell-db93c0502a0559ea", "locked": true, "points": 2.5, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "assert isinstance(consistent_estimators, set)\n", "assert consistent_estimators <= {\"a\", \"b\", \"c\", \"d\", \"e\"}\n", "# А тут скрытые assert'ы :)" ] } ], "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.12.0" } }, "nbformat": 4, "nbformat_minor": 5 }