{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Математическая статистика (ФБМФ, ФМХФ)\n",
"## Домашнее задание 2 — часть 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 баллов**."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Данная часть задания проверяется автоматически. Для выполнения задания используйте этот ноутбук в качестве основы, ничего не удаляя из него. Можно добавлять необходимое количество ячеек. Ячейки с assert'ами удалять и изменять нельзя, в противном случае соответствующее задание *не будет оценено*.\n",
"\n",
"\n",
"> *Примечание.* Рекомендуется работать с данным ноутбуком **локально** в Jupyter Notebook (например, используя Anaconda или альтернативные варианты). Если вы используете софт по типу Google Colaboratory, то перед отправкой боту данного ноутбука необходимо проверить, что **в ячейках с assert'ами и \"# Ваше решение тут\"** в метаданных присутствуют поля *metadata* с `nbgrader`. Можно открыть ноутбук с решением в текстовом редакторе (MS Word, Блокнот) и выполнить поиск по документу слова `nbgrader`. Если поиск показал ровно **17 совпадений** — можете отправлять файл боту. Если совпадений меньше, решение может быть не оценено. В таком случае попробуйте скачать файл в форме `ipynb` еще раз и перенесите решения в новый файл. **Внимание! Бот не проверяет решение и не проверяет наличие метаданных.** "
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import scipy.stats as sps\n",
"import matplotlib.pyplot as plt\n",
"%matplotlib inline"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Задачи типа A (достаточно на \"хор\")\n",
"\n",
"---"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Задача 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",
"* Асимптотически нормальную оценку $\\theta$;\n",
"* Оценку асимптотической дисперсии данной оценки.\n",
"* Асимптотический доверительный интервал на $\\theta$."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "66eba1b5fce307606ecfe5053baef4a7",
"grade": false,
"grade_id": "cell-7039e8c24d3f54cd",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def poiss_asymp_norm_est(x):\n",
" ... # Ваше решение тут\n",
" return est"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "0e684f9892ae2fec1f3cc980e641010c",
"grade": true,
"grade_id": "cell-b2fc43b38ef17db6",
"locked": true,
"points": 1,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"x = np.array([1, 5, 5, 3, 2, 4, 2, 8, 4, 6])\n",
"assert np.allclose(poiss_asymp_norm_est(x), 4, atol=1e-2)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "5b2c1c25940df78b5b9e53aa13b0961e",
"grade": false,
"grade_id": "cell-a3e2906346fdb6f1",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def poiss_asymp_var_est(x):\n",
" ... # Ваше решение тут\n",
" return est"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "2a73ef69388749fd7c1203b41dfda489",
"grade": true,
"grade_id": "cell-9512a59aebda5e17",
"locked": true,
"points": 1,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"x = np.array([1, 5, 5, 3, 2, 4, 2, 8, 4, 6])\n",
"assert np.allclose(poiss_asymp_var_est(x), 4, atol=1e-2)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "f824ca37d221b51dd3fad102e10e8a25",
"grade": false,
"grade_id": "cell-8a45074da2440c35",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def poiss_asymp_confint(x, confidence_level=0.95):\n",
" ... # Ваше решение тут\n",
" return l, r"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "137c147809b4619b43c327d1b4952717",
"grade": true,
"grade_id": "cell-3061e1b0ae69a32e",
"locked": true,
"points": 1,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"x = np.array([1, 5, 5, 3, 2, 4, 2, 8, 4, 6])\n",
"assert np.allclose(poiss_asymp_confint(x), [2.76, 5.24], atol=1e-2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Проверим, что при большом размере выборки получившийся интервал действительно имеет уровень доверия $\\approx 0.95$:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "8d4e470b4ca94fd1994919ff5012593c",
"grade": true,
"grade_id": "cell-e72a5cbc1565ecdf",
"locked": true,
"points": 2,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"np.random.seed(0)\n",
"dist = sps.poisson(5)\n",
"l, r = np.stack([\n",
" poiss_asymp_confint(dist.rvs(size=(500,)))\n",
" for _ in range(10_000)\n",
"], axis=1)\n",
"estimated_conf_level = np.mean((l <= 5) & (5 <= r))\n",
"assert 0.94 <= estimated_conf_level <= 0.96, f\"Неправильный уровень доверия {estimated_conf_level:.4f}\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"-----\n",
"\n",
"### Задача 2\n",
"Дана выборка $X_1, ..., X_n$ из распределения Лапласа с плотностью $p(x) = \\frac{1}{2} e^{-|x-\\theta|}$. Реализуйте функцию, вычисляющую асимптотический доверительный интервал для параметра $\\theta$."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "4f7b362bce1202a19d0ff79d284c4ab9",
"grade": false,
"grade_id": "cell-7f73d5b085973f0e",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def laplace_asymp_confint(x, confidence_level=0.95):\n",
" ... # Ваше решение тут\n",
" return l, r"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "39405e35c3726e5f5b314194454d03f8",
"grade": true,
"grade_id": "cell-e9e1de07590b727c",
"locked": true,
"points": 5,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"np.random.seed(0)\n",
"dist = sps.laplace(loc=2)\n",
"l, r = np.stack([\n",
" laplace_asymp_confint(dist.rvs(size=(500,)))\n",
" for _ in range(10_000)\n",
"], axis=1)\n",
"estimated_conf_level = np.mean((l <= 2) & (2 <= r))\n",
"assert 0.94 <= estimated_conf_level <= 0.96, f\"Неправильный уровень доверия {estimated_conf_level:.4f}\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"---\n",
"\n",
"### Задача 3\n",
"В этой задаче нужно визуализировать *свойство асимптотической нормальности*. Посмотрите также на этот ноутбук.\n",
"\n",
"*a).* Пусть $X_1, ..., X_n$— выборка из распределения $U(0, 1)$. Согласно центральной предельной теореме оценка $\\widehat{\\theta} = 2\\overline{X}$ является асимптотически нормальной оценкой параметра $\\theta$. Вам нужно убедиться в этом, сгенерировав множество наборов случайных величин и посчитав по каждому из наборов величину $Z_n = \\sqrt{n} \\left( 2\\overline{X} - \\theta \\right)$ в зависимости от размера набора. \n",
"\n",
"Сгенерируйте множество выборок $X^1, \\dots, X^{300}$ из распределения $U[0, 1]$: $\\; X^j = (X^j_1, \\dots, X^j_{500}), 1 \\leq j \\leq 300$.
\n",
"По каждой из них посчитайте оценки $\\widehat{\\theta}_{jn} = 2\\frac{X^j_1 + \\dots + X^j_n}{n}$ для $1 \\leq n \\leq 500$, то есть оценку параметра $\\theta$ по первым $n$ наблюдениям $j$-й выборки.\n",
"Для этих оценок посчитайте статистики $Z_{jn} = \\sqrt{n} \\left( \\widehat{\\theta}_{jn} - \\theta \\right)$, где $\\theta = 1$."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import scipy.stats as sps\n",
"import matplotlib.pyplot as plt"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "a4855b52c0139ea51a1b6740e1cde940",
"grade": false,
"grade_id": "cell-47adc93fea048529",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"np.random.seed(0)\n",
"n = 500\n",
"num_samples = 300\n",
"x = ... # shape [num_samples, n]\n",
"estimation = ...\n",
"z = ...\n",
"\n",
"... # Ваше решение тут"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "ca425b3d1d8f37ce1d5ac6e07f68b151",
"grade": true,
"grade_id": "cell-71f1ecd8f619bf50",
"locked": true,
"points": 1,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"assert x.shape == estimation.shape == z.shape == (300, 500)\n",
"assert ((0 <= x) & (x <= 1)).all()\n",
"assert np.allclose(((1/3 <= x) & (x <= 2/3)).mean(), 1/3, atol=1e-2)\n",
"assert np.allclose(estimation[1, 2], 2 * (x[1, 0] + x[1, 1] + x[1, 2]) / 3)\n",
"assert np.allclose(z[1, 2], 3**0.5 * estimation[1, 2] - 3**0.5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Для каждого $j$ нанесите на один график зависимость $Z_{jn}$ от $n$ с помощью `plt.plot`. Каждая кривая должна быть нарисована *одним цветом* с прозрачностью `alpha=0.05`. Сходятся ли значения $Z_{jn}$ к какой-либо константе?"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "07759d0e7861221eced2a0a70a516cd5",
"grade": false,
"grade_id": "cell-9166b9178e7bad98",
"locked": false,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def plot_trajectories(z):\n",
" # plt.plot и прочее. Не вызывайте plt.show()!\n",
"\n",
" ... # Ваше решение тут"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "b3553b3e731e4aa5a27c8578b248cddd",
"grade": true,
"grade_id": "cell-a696a84d5079684a",
"locked": true,
"points": 2,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"plot_trajectories(z)\n",
"lines = plt.gca().lines\n",
"assert len(lines) == z.shape[0], \"Неверное к-во линий. Не вызывайте plt.show()!\"\n",
"assert np.allclose(lines[42].get_ydata(), z[42]), \"Неверные координаты линий\"\n",
"assert lines[42].get_alpha() == 0.05, \"Неправильная прозрачность линий\"\n",
"assert len(set(line.get_color() for line in lines)) == 1, \"Цвета должны быть одинаковыми\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Для $n=500$ по выборке $Z_{1,500}, ..., Z_{300,500}$ постройте гистограмму и график плотности распределения $\\mathcal{N}(0, \n",
"\\sigma^2)$, где $\\sigma^2$ — асимптотическая дисперсия. Не забудьте сделать легенду."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "ad37665fa199cea55f28ff965755c0fd",
"grade": true,
"grade_id": "cell-1df4043dcd0714d8",
"locked": false,
"points": 2,
"schema_version": 3,
"solution": true,
"task": false
}
},
"outputs": [],
"source": [
"def plot_distribution(z):\n",
" # plt.plot и прочее. Не вызывайте plt.show()!\n",
"\n",
" ... # Ваше решение тут"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"deletable": false,
"editable": false,
"nbgrader": {
"cell_type": "code",
"checksum": "080c2ebe399333bca65a49882ea3ecc8",
"grade": true,
"grade_id": "cell-68f4511bd1ec1511",
"locked": true,
"points": 0,
"schema_version": 3,
"solution": false,
"task": false
}
},
"outputs": [],
"source": [
"plot_distribution(z)\n",
"axes = plt.gca()\n",
"total_hist_area = sum(hist_rect.get_width() * hist_rect.get_height()\n",
" for hist_rect in axes.patches)\n",
"assert np.allclose(total_hist_area, 1.0, atol=1e-1), \"Гистограмма должна быть нормированна\"\n",
"assert len(axes.lines) == 1, \"Нет графика плотности\"\n",
"pdf = axes.lines[0]\n",
"pdf_fragment = pdf.get_ydata()[np.abs(pdf.get_xdata()) <= 0.3]\n",
"assert len(pdf_fragment) > 0, \"Некорректный график плотности\"\n",
"assert ((0.6 <= pdf_fragment) & (pdf_fragment <= 0.7)).all(), \"Некорректный график плотности\"\n",
"assert axes.get_legend() is not None, \"Нет легенды\"\n",
"assert len(axes.get_legend_handles_labels()[1]) == 2, \"Не все графики подписаны\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Сделайте вывод о смысле свойства асимптотической нормальности. Подтверждают ли сделанные эксперименты теоретические свойства?"
]
},
{
"cell_type": "markdown",
"metadata": {},
"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.12.0"
}
},
"nbformat": 4,
"nbformat_minor": 4
}