{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Load Packages" ] }, { "cell_type": "code", "execution_count": 90, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "TaskLocalRNG()" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "using Printf, Statistics, StatsBase, Random, Distributions\n", "include(\"jlFiles/printmat.jl\")\n", "Random.seed!(678) #set the random number generator to this starting point" ] }, { "cell_type": "code", "execution_count": 91, "metadata": {}, "outputs": [], "source": [ "using Plots\n", "\n", "gr(size=(480,320))\n", "default(fmt = :svg)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Introduction\n", "\n", "This exam explores how autocorrelation ought to change how we test statistical hypotheses." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Task 1\n", "\n", "Code a function for simulating $T$ observations from an AR(1) series\n", "\n", "$\n", "y_t = (1-\\rho)\\mu + \\rho y_{t-1} + \\varepsilon_t \\sigma\n", "$\n", "where $\\varepsilon_t$ is N(0,1).\n", "\n", "That is, generate $y_1,y_2,...,y_T$ from this formula.\n", "\n", "To make also the starting value ($y_0$) random, simulate $T+100$ data points, but then discard the first 100 values of $y_t$.\n", "\n", "Generate a single \"sample\" using `(T,ρ,σ,μ) = (500,0,3,2)`. Calculate and report the average (mean) and the first 5 autocorrelations (hint: `autocor()`) of this sample. Redo a 2nd time, but with `ρ=0.75`." ] }, { "cell_type": "code", "execution_count": 92, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "SimAR1 (generic function with 1 method)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "function SimAR1(T,ρ,σ,μ)\n", " y = fill(NaN, T + 100)\n", " e = rand(Normal(0, 1), T + 100)\n", " y[1] = 1\n", " for i = 2:T + 100\n", " y[i] = (1-ρ)μ + ρ*y[i-1] + e[i]*σ\n", " end\n", " return y[101:end]\n", "end" ] }, { "cell_type": "code", "execution_count": 93, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "average from one sample with ρ=0 1.986\n", "\n", "autocorrelations with ρ=0\n", "\n", " 1 -0.016\n", " 2 0.027\n", " 3 0.004\n", " 4 0.004\n", " 5 0.008\n", "\n" ] } ], "source": [ "y1 = SimAR1(500,0,3,2)\n", "\n", "printmat(\"average from one sample with ρ=0\", mean(y1))\n", "\n", "printmat(\"autocorrelations with ρ=0\")\n", "printmat(1:5, autocor(y1)[2:6])" ] }, { "cell_type": "code", "execution_count": 94, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "average from one sample with ρ=0.75 3.187\n", "\n", "autocorrelations with ρ=0.75\n", "\n", " 1 0.766\n", " 2 0.580\n", " 3 0.495\n", " 4 0.441\n", " 5 0.372\n", "\n" ] } ], "source": [ "y2 = SimAR1(500,0.75,3,2)\n", "\n", "printmat(\"average from one sample with ρ=0.75\", mean(y2))\n", "\n", "printmat(\"autocorrelations with ρ=0.75\")\n", "printmat(1:5, autocor(y2)[2:6])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Task 2\n", "\n", "Do a Monte Carlo simulation. Use the parameters `(T,ρ,σ,μ) = (500,0,3,2)`.\n", "\n", "1. Generate a sample with $T$ observations and calculate the average. Repeat $M=10,000$ times and store the estimated averages in a vector of length $M$. (The rest of the question uses the symbol $\\mu_i$ to denote the average from sample $i$.)\n", "\n", "2. What is average $\\mu_i$ across the $M$ estimates? (That is, what is $\\frac{1}{M}\\sum\\nolimits_{i=1}^{M}\\mu_i$?) _Report_ the result.\n", "\n", "3. What is the standard deviation of $\\mu_i$ across the $M$ estimates? Compare with the theoretical standard deviation (see below). _Report_ the result.\n", "\n", "4. Does the distribution of $\\mu_i$ look normal? _Plot_ a histogram and compare with the theoretical pdf (see below).\n", "\n", "\n", "## ...basic stats (the theoretical results)\n", "\n", "says that the sample average of an iid (\"independently and identically distributed\") data series is normally distributed with a mean equal to the true (population) mean $\\mu$ and a standard deviation equal to $s=\\sigma_y/\\sqrt{T}$ where $\\sigma_y$ is the standard deviation of $y$.\n", "\n", "To compare with our simulation results, you could estimate $\\sigma_y$ from a single simulation with very many observations (say 10'000)." ] }, { "cell_type": "code", "execution_count": 95, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Average across the simulations: 2.001\n", "\n", "Std across the samples (with ρ=0) and in theory:\n", "simulations theory\n", " 0.134 0.133\n", "\n" ] } ], "source": [ "M = 10_000\n", "\n", "(T,ρ1,σ,μ) = (500,0,3,2)\n", "\n", "μi1 = fill(NaN, M)\n", "σi1 = fill(NaN, M)\n", "\n", "# Monte Carlo simulation\n", "for i = 1:M\n", " y = SimAR1(T, ρ1, σ, μ)\n", " μi1[i] = mean(y)\n", " σi1[i] = std(y)\n", "end\n", "\n", "# Theoretical results\n", "y1 = SimAR1(10_000,ρ1,σ,μ)\n", "σy1 = std(y1)\n", "s1 = σy1/sqrt(T)\n", "\n", "printmat(\"Average across the simulations:\", mean(μi1))\n", "println(\"Std across the samples (with ρ=0) and in theory:\")\n", "printmat([\"simulations\", std(μi1)], [\"theory\", s1])" ] }, { "cell_type": "code", "execution_count": 96, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAFACAIAAADrqjgsAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3deXwTdf4/8Pckk6P33bSlBQpIoYicAgrlBhUUcXVlv3iA64rrsYKs13qAP/ECUdFFBYFFRAUVEAQWBArlKDeUo1ylhbb0TNqmR+7M8fvjg7Mxadq0TTqT9P188AeZTGbemU5e+eQzn5mheJ4HhBBC0iMTuwCEEEKNk7/99tti19Ay+/bte/fdd4OCgrp37y52LQHixo0ba9as+eGHHzZv3lxbW9uvXz+vLLa4uPjixYsmkykmJoaiKHezVVVVnT17tq6uLioqSi6Xu5utoaEhNze3oqIiLCxMqVS6m81ms128ePHatWshISFBQUFtfQ/IS/bv379gwQKVStWjR4+m51y/fv2///3vrl27ajSa9qmtRXiez8vLu3r1Ks/zERER7bE+0Wk0GoVCUVxc7PrUp59+qlAonnvuOWHKF198AQCLFy/2fPnff//98uXLvVBoIDp37lx4eDgAqNXqqKio2bNnu5tz8+bNr7766vjx4yMjIwFg7NixTSxz2LBhwj7Wo0eP3377zXW28vLyBx54QCa7+TMuLi7uiy++cJ3NZDK98MILarWazBYcHPzSSy9ZrVan2TiO++STT6Kjo8lscrl82rRpOp3O4y2BfGjZsmUAsHDhQmHK7t27ly9frtVqneZ8/vnnAWDr1q3tW6BHdu7c6fgFM2zYsHPnzvl0jZII6KioKAAoKipyfWrx4sUA8PTTTwtTNm/ePHr06J9++snz5Xfp0oWiKC8UGoiefPJJAGg0GZ3ExcWR/ZJkpbuAzs/Pj46Opijqscce++qrr+bMmaNSqRQKRWZmpuNsBoOhb9++ADBu3LilS5e+/fbbZPlLlixxWuCDDz4IAH369Pnkk08WL16clpYGAI8//rjTbO+99x4AJCQkLFiw4N///ndGRgYADB482Gw2t2R7IJ/YunXr6NGj161bJ0z585//DAAnTpxwmlOyAb1nzx6aplUq1Zw5c7766qtHH30UAGJiYgoKCny3Uv8L6FbAgG7CiBEjAKCqqqrZOd99992ffvrp2rVrO3bsaCKgH3jgAQCYN2+eMGXDhg0A0KtXL4ZhhIkLFiwAgClTprAsS6ZcuHAhKCgoJCSkvLxcmG379u0A0LNnz4aGBjJFr9d37doVALKysoTZCgsLVSpVaGgo+e3J8zzDMBMmTACATz75xNNtgdqRfwU0wzCkZfDzzz8LE+fNmwcAf/rTn3y3XtoLvSTtS6fTFRcXp6SkxMfHCxP1ev2hQ4eKioo4jouNjR04cGCvXr0AoKGhIS8vz2azAcCpU6fIzGq1uk+fPo4L3LNnT2lpaXh4+LBhw2677TbXlbIsu3fv3gsXLgQHB48aNSotLe3GjRtarbZHjx5CP9TZs2d5nu/fv7/ZbN65c2dhYeEtt9xy7733kmfPnTtHulCDgoL69+8/bNgwp27ZgoKC2tra9PR0lUpF1hUeHj5p0iShJy4vL+/AgQMGg2HIkCF33nmnh5urpKRk37595eXlsbGxGRkZt9xyi/BUUVFRVVVVWVkZAFy+fJm0i/v27euuh/eNN94g/7ly5Yq71Wm12l9//TU8PPyVV14RJj744IMDBw48ffr0wYMHR48eDQA8z69cuRIAFixYIHRxpKenT58+fdWqVT/88MPcuXPJxBUrVgDAa6+9FhoaSqZERka+/PLLzz333IoVK0aNGkUmfvvtt1arddasWcIvULlc/s477+zevXvFihUvvvhi01vJbDYfOXKkqKhIp9MlJiaOHDmyS5cuwrNkL4qIiHDtP62vr7969arTU1ardf/+/VeuXOE4Li0tbezYsU6b9NSpU0qlsm/fvgaDYefOncXFxbfddtv48eMBgOO4M2fOXLx4sbKyMjQ0dNCgQYMHD2605suXL2dlZVmt1j59+owZM8ZsNl+5ciUmJoZ8ewk4jjt69OiZM2csFkvXrl0nTJgQFhbmtKiKioqjR48WFxfLZLL4+PghQ4Y4LcSRVqu9ceOG0wcwJyeH4zjHj4PNZjt//nxYWFjPnj0BoKqqqqioKDk5WaPRcByXk5Oj1+vJuxA+CAMGDBB2BuLEiRNHjx4FgGHDht1+++3uSnKVk5Nz4sQJvV6fkJBwzz33OJbaOgcOHLhy5cqAAQMeeughYeKrr766ZMmSLVu2aLXatq+icb7Lfs+1qAXt2ge9Zs2akJAQp/f10ksv8TyflZXl+pbT09OF1y5atMjpUNLkyZNramocaygrKxs0aJDjPHPnzp09ezYA/Prrr8Js0dHRERER2dnZiYmJZLZJkybxPH/58uXk5GSnGgYPHlxYWOi4lilTpgDAb7/95rgjBgcH79ixg2XZuXPnOu67Tz75ZLNblWXZl19+mab/9x1MUdSMGTOEn/ykc8PJjRs3ml1yEy3o9evXA8C9997rNJ2E+xtvvEEe5uXlAUBSUpLTbKStfc8995CHHMeR/nHHNjXP89euXQMAjUYjTBk5ciQAbN++3XE2juNIl3RJSUkTb2fp0qXBwcGOG0Eulz///PNCe7+2tjYoKCguLs5mszm9lhxjf++99xw3TkpKiuPSUlNTndqJFEV17dp1165dMTExZJ5HHnmE5/mTJ0+6HhnLyMioqKhwel+zZ892/IIfPHjwxo0bAeCxxx5znPPMmTOkH0kQGxu7efNmx3kWLVqkUCicVura0STYvXs3/PEjWVRURF7l+GPlt99+A4BnnnmGPHTsgzaZTK47HgAYDAb+9xb0Tz/9NHXqVMdnZ8yYwXGcu6oEZrN51qxZji8MDQ39/vvvm31h015//XUA+Ne//uU0/b777gOA9evXt3H57vh9QF+/fl2hUMTGxn777bdXrlwpLCw8cODA22+//dFHH/E8X1NTs3v3bo1GQ1HU7t8dPnyYvHbp0qXkc7527dqrV68eOHBg7Nix5CMhfDjtdjtJzL/85S+nT58uLi7+/vvv4+LikpKSwCWg1Wp1fHz8I488snHjxuzsbJIXx48fz8jIWLFiRXZ29tWrV/ft2zdt2jTyoXLc4UhAd+nSZfTo0du2bTtx4sRbb70lk8ni4uLmz58fHx+/YsWKkydP/vjjj+QLwOlj5uq1114DgO7du2/atKmgoGDXrl0DBgwAgGnTppEZLly4sHv3btJU3LZtG9k4Foul2b9XEwE9f/58APjHP/7hNH3VqlXg8GNwy5YtADB8+HCn2XJycgCga9eu5OGNGzcAICgoyGk2lmXJF091dTWZQtovly5dcpqTND93797dxNuZN2/ejBkzfvnllzNnzly4cGHdunXp6ekA8OGHHwrzkD/Zli1bHF/IcVz37t1lMplwfHvv3r00TUdERHz00UcnT57MyclZsGCBUqmMiYlx/JKgKCoiIiIiIuKpp57asmXLoUOHyEHUXbt2jR8//ptvvjly5EheXt7u3bvJ53/cuHGO6/3ggw8AoGfPntu3by8uLs7Ozh41ahTZIR0DOj8/PzIyUqlUvvLKK0eOHMnNzV22bFlkZKRCoTh27BiZ58iRI+QrZMOGDfn5+fn5+Xv37n355Zf/85//uNtcJpNJrVZ3795dmPKf//wHAGQyGWmREOQnlNAh4BjQLMvu3r2bfKd+8cUXwgeTfOhIQKempvbu3Xv9+vWnT59eu3ZtQkICAHz77bdN/B2J//u//xM+rUVFRV999VVISIhcLj906JAwT1lZ2W8eyM3NFV5COu5WrFjhtLo5c+YAwPz585strHUkFNADBw4c5oLERxMB/e2338IfmzCuGu2DNhgMZDRCdna2MNFisZBOAGHH+umnnwBg6NChQlcpz/O7du0iX85OAe0Yf00jcXzw4EGnKYMGDXLsqCWNCJqmL1y4IExct24dADz66KNNLL+kpEShUCgUCscjGDU1NaTJdvToUWEi6e0xGo2elE00EdDPPvssALz99ttO00kiZ2RkkIek4+L+++93mo0kckhICHlIeqVSUlJcV0SOKF6+fJnneY7jyM8L1yEBd911FwD88MMPnr87nufLy8sjIyMdG/j//e9/waW38cCBAwAwceJE8pBl2bS0NLlc7vhn5Xn+3//+NwA4Do8hjd+///3vzVbCsizpxhHCoqGhISwsjKZpobed53mj0ega0CTcV65c6bjAPXv2AMCECRPIw/fff9/D4HM0ZswYALh+/Tp5OH36dLVaPWXKlNDQUOFHxqBBg2QymTCKxnUUR9N90ElJSfX19cLEX375BQCmTJnSdGHHjx8HgDvvvNPx00oaB6NGjRKmrF27ttEmvJOZM2cKLyGHajZt2uS0RnIo5dlnn226sFaT0IkqeXl5F1xUVFQ0/SoS7qdOnbLb7S1a3b59+2pra0eNGuXYn6tSqUjv56ZNm8iUrVu3AsDs2bMdexgmTJhw6623NrrYV1991ZO133///QBA9idHs2fPdhwLTD6cd911F2nTOU68fv16E8vfunWr3W6fNm1at27dhIlRUVFPP/2047vzOrPZDACu40PJFKPR6Dgb6b5wnc1kMvE838RsAEC+XA0GAwBYLBaO4xqdk8wmrNdDCQkJQ4cOLSsrI18YAHDXXXclJydv27atqqpKmG3NmjUAMGPGDPLw+PHjV65cGT16NPkwC2bNmqVQKMi3miNPdhWZTEZyVthVsrKyGhoa7r33Xsde7+DgYPKXFej1+u3btycmJjp1ZI0bNy4tLS0rK4tsW7J9jh07Rjagh8aNGwcAmZmZAMDz/L59+zIyMiZPnmwwGEidtbW1Z86c6d+/f2xsrOeLdfT3v//dsa984sSJ0Nw+D7//RV588UXHT+tjjz2m0WgOHDgg/O1uv/32zz1AxmkQTe+xLd3BPCehg4QXLlzo3Lmz08SPP/74pZdeauJVY8aMSU1N3bRpU7du3SZPnjx69OiJEycKg2GbcOnSJQAYOHCg03TS3XzhwgXykPSWOuYj0adPn9zcXKeJFEX17t3bdV0nTpxYvHjxuXPniouLHTvgHD/tBDmiIiDtRMcje8LEyspK928OLl68CACkT8OR07vzOpVKBb/npqOGhgYAELr7m55NrVaTNqa72QCgvr5eWKBKpSK/kIxGI3mJ02zCGOpGcRy3evXq77///urVq5WVlY7f9NXV1aRDWSaTPfLIIwsXLly/fj1p4pnN5g0bNoSHhwtdpaR/pq6ujnQuOQoKCiosLHScEhIS0uiBuAMHDixZsiQ3N7ekpISEAiHsKk3skI4Pz549y3GcWq12LcZms9nt9tLS0h49ekydOvWNN9744osvdu7cOXny5DFjxowfP144HuvOuHHj3nzzzczMzCeffPLixYvl5eWzZ88mqb1nz57hw4fv3buXZVkypXXIkAlBcHBwWFhYs821/Px8AHA600qhUPTp06eysvLixYukXyUtLc1p+c0i+5VrEDvt2F4noYBunZCQkOzs7Hnz5m3atGn58uXLly+nafqhhx5asmRJ02cikY+967FXMoVsd/j9m5M0NBy5TgGA8PBw1yzYvn371KlTaZoeP3781KlTyQtzc3O/++47lmWdZnb6S5OccjqERVoHfJMXUXH37sg2Ed6d15EfNDU1NU7TyRTyLACQb9BWzwYAZAwAmUEmk0VERNTW1lZXVzt9N5PXNv2FPWvWrFWrViUkJEyePLlTp05ka//www/nzp1jGEaY7a9//evChQvXrFlDAnrTpk11dXWzZs0S/jqkpAsXLhQUFDitQi6Xh4SEMAwjHLNttGm5bt26Rx99VK1W33333Q899BBpnR0/fnzTpk1CJe52SGGjEbW1tQBQWlr69ddfu64oKirKarUCQGJi4okTJ958883t27eTZqNKpZo5c+aiRYsa/eFCDB48OCIiggxsJ+3ocePGde/ePTU1NTMzc/78+cJEd0toltM+DwAymazpfR4AysvL4fcWjCOnD7XBYGg26wEgPDxc+ASRzVtdXe00j9Me63V+H9AAkJiYuGLFiq+++ur06dOZmZmrV69ev379jRs3Dh061MSryA8orVbrNJ20TIW9k2z6srIyx3FXZIqH5b322mssy+7fv9+xL+XLL7/87rvvPFxCK5B359rKJvtlE5+9NiINk5KSEqfppK9AaLaQHwruZhN+RnTp0kWlUjU0NNTV1Tl2m1RWVtpstoiICHLsiLzk+PHjJSUlTr82yCqcfpc4unjx4qpVq9LT048cOeK4WUins6OePXsOHTr02LFj58+f79u3r1P/Bvy+zZ966qnPPvvM3eoErue+8zz/yiuv0DR9+PBhxzbge++959glRXZIkkSOSktLHR+SYkaNGiUcL3Gne/fu69ats1qtJ06c2L1796pVq5YvX15bW0sG5DSKpumRI0du3bo1Nzc3MzMzMjKS/FYbO3bs2rVrDQbDnj17lEqlU1dPOyC/frRardMXmNOHevPmzY899lizS5s5c+bq1avJ/9PS0rZv3+60keH3Hayl7XHPBUJAEzRNDxkyZMiQIXPmzOnZs2d2dnZ5eTkZ8KBQKHieZ1nWsXuX/CQ8efKk03JOnDgBAMLgpIEDB+7fv3///v133HGHMI/ZbCbDM5vFMAzpunEauSwMyvYR0kXuuhby7tx1oLcdOX/v4MGDjq1FANi7d6/wLACkp6fHxMRcu3atuLjYsV+LzEZ+hwKAXC6/88479+3bl5WVRXrtG10a+f/x48ezsrLI8SviypUrpaWlnTp1auKyLWfPngWASZMmOaaz1WpttBdoxowZx44dW7t27ezZs/fu3XvLLbc47hUkpMi4iFaorq4uKSnp16+f0y/006dPOz4knXL79+93ernTiNJ+/frJZLKTJ0/abLYmLl0iUKlUI0aMGDFixD/+8Y+uXbtu3rzZ6S/oZNy4cVu3bt21axcZ+0Q+WePGjVu1atWPP/6Yl5c3cuRI18GvjsjYPsffKN6Sk5Pj+JVstVpzc3MpihJ6gdLS0p577rlml+P4x83IyPjkk0/27dv35ptvChNJwwv+uCt6mY8OPrZIW4bZCSeYOSLxKhxlJh94p3HH5CI+8MehFCaTiXyYhQFVp0+fpigqISHBcSguOXQLjY2Ddi0mOjo6NDTUsc68vDzymXnllVeEiWQUBznbRfD9998DwOuvv+44kXSMOI5zclVRUaFUKhUKheOx/qqqKuGYqjDRu6M4eJ4nX0XffPONMOXkyZMymUyj0ZCjfwTpKHjhhReEKTqdLi4uTi6XO46WI4fgMzIyhOPyDMMMGTIE/jj4lPyZEhMT9Xq9MPFvf/sbALz66qtNvBcyqvfhhx92nLhw4ULy93UaY0AGRGs0mnfeeQdcxg4JJ5tt2LDBdUWOOwAZB+00A0lSjUbjONLx1KlTpEfr/fffJ1PIyS8A8OOPPwqznTt3juxRjqM4yB71zjvvNFGM68eHZdmYmBi5XN70gMvz588DABngL1wnoKKigqIo0mvvNJLHdRQHObbkeGIe4e5MwoiIiNjY2CZK4nmeHAEaMmSI41Co5cuXg8tQxZYyGo3x8fEymcxxl/jmm2+gsdGiXuT3Ab1gwYKBAwd+9tlnmZmZeXl52dnZTz31FADcfvvtwkv++c9/AsAdd9zx0UcfLV++XLiOBxnsFRsbu2LFitzc3J07dw4fPpz8LR2H6ZChjomJif/6178+/fTT+++/X6VSkdB33I3cBTS5lMSECROysrIuX768Zs2a5ORkMrjCdwHN/z4kuUuXLuvWrcvNzd2yZQvJYqerWHge0KtXr541a9asWbPuueceAEhKSpr1O8cRUYcOHVIoFEFBQYsXLz5y5Mjq1atJR8SqVascl1ZWVhYfH09R1Jw5c7Kzs3/55RfScnS8MBbP83a7nbQZp06dmpmZuXv37kmTJpFPhePfiOf5mTNnAsDgwYN//fXXgwcPkgF/nTp1avos9pqamrCwMIqi5s2bd+bMmZycnNdee02pVJIjeK6DwMiA6ODgYMfhz47vXaVS0TT9wgsvbNu27dy5c7t27Vq6dOnw4cMd31ejAc3zPBmGf//992dnZ1+6dGn58uXx8fFkVxECmuf5rKwspVJJ0/TMmTM//fTTOXPmhIWFkfPaHQO6sLCQdKFOnz5948aNZ8+ezcrKWrVq1eTJk8ePH0/m+fvf/56RkbFs2bL9+/fn5eXt27eP7K7Nnr7McZzQv3TlyhVhuvDT03HcMd9YQP/4448A0K1btwULFpCjR3a7nfdGQHfu3Hnq1KnZ2dkXL178+OOPg4KCHId+txo59zUhIWH16tVHjhxZvHgxWbLjOF2v8/uAXrVqlesvqZEjRzqeEVdTU/PAAw8Ip0s5nkm4dOlSpzNfH374Yce44XmeZdn3339fOOnrtttu27dv3xNPPAEAe/fuFWZzF9BlZWVOP1qnT59OktenAc1x3Pz58x1HNcjl8meeecbpOnCeBzRJwEZVVlY6zrlp0ybH43Iqlerjjz92XeCpU6ccu4wpipo1a5br2XplZWVCpwcxYcIE1yHPZrP58ccfd5wtPT39/Pnzzb6vLVu2OO4DwcHBq1evJkHsGtBC37Qw/NnJ0aNHXa/XmpiY6Pj95C6gCwoKnLrLZ82aRaLNMaB5nt+3b59wTYLY2Nj3339/27Zt4HDmnrBA1yN1kZGRb731Fpnhvffecxr3AgD333+/08m0jSKnhHTq1MlxImnNOA6IJlwDmmXZF1980XHLO55J2JaAPn78uOMImYiIiI0bNzb7djyxePFix80VHR39yy+/eGXJ7lC8BO6ocu3aNZZlU1NTXfu8amtrdTpdZGSkcGS2rq5Oq9XGxcUJxwHsdvvp06evX79eV1en0WjS0tIaHevG83xFRYXZbFYqlY7nXtfV1R04cKC0tDQ0NPSOO+5w11/JMExFRUVwcDCJnoyMjEOHDhUUFAgDjck4qkbHTrEsS84NUygUQ4YMSUtLMxqNlZWVkZGRQpBVVFSYTKbk5GTHHkODwUCOeDiNQ7h27ZpCoXA6pbhRVVVV+/fv1+l0ERERGRkZriedl5SU2Gy21NTUJi7ZTOh0OnfDP7p06eJ0KWeTyZSZmVlWVhYVFTVu3Djh680Jy7IHDhzIz89Xq9UZGRlNXALi5MmT58+fpyiqX79+rsMHBQUFBYcPH7ZYLD179szIyHC6toM7NTU1R48evXHjRnx8/JgxYyIjI7VarcFg6NSpk1N+cRxH/tCufxQBz/Pnz58/c+aMyWTSaDRdu3YlPcLCDNeuXaNp2nVQKQDY7fbDhw/n5+erVKrhw4enpqY2NDTodLro6GjXkRtVVVU2m02j0cjl8s8//3z27NkffPCB67i6goKC48ePkytud+7cefDgwY7ndpvN5lOnThUVFRkMhqSkpFtvvTU1NdWTjabX6/V6fVBQkHBtAwCor6+vqqpy+oiR6ZWVlbGxsa4DHqqrq+vq6gCA7IRVVVX19fUJCQlOAzkKCwspinI6Vu8kPT390qVLJpNJLpfv27evsLAwOjp6woQJjY65ap3q6urMzEy9Xp+UlDRu3DjX0SbeJYmA9jvnz58fMGBA586dCwoKms01hHzNbrcPHTo0Jyfn+PHjLbqoUIARAjpg7tUQOKM4fGf58uVXrlyZMmVK165d9Xr98ePH58+fz7LsG2+8gemM2l9hYeHs2bMff/zx3r17K5VK0tOak5Nz1113deR0DkgY0M2zWq2ffvrpp59+KkwJCQlZtGhRo1eDQ8jX5HL5jh07fv31V8eJU6dOJRctQoEEuzg8cuPGjRMnTlRUVLAsm5KSMmrUKN+dO4RQswwGw+HDh0tLS+vq6qKjo++44w6nM3Q6ppUrV1ZXV//zn/9sYgS3f8GARgghiZLQ1ewQQgg5woBGCCGJwoBGCCGJwoBGCCGJwoBGCCGJwoBGCCGJ8nJAkyt8tn05OPjPE+RyKmJX4QdadMO9Dgu3kifa+RPn5YBOS0vzyhW4TSYT7i7NYlmW3LgINc139/QMJLiVPGGxWHxxkwF3sIsDIYQkCgMaIYQkytOArqure+CBB9LS0rp06TJlyhRywxuEEEK+42lAy2Sy6dOn79y5Mzs7u3fv3pMmTcLDUwgh5FOtuVhSeXl5UlJSdXW16x0llEql0Wh0vF9D6xiNxqCgIA9vh9FhMQzDMIxarRa7EKlraGhwurEZcoVbyRNms1mhULTb1fJatprz589XVVWtXLly2rRp7u73g5DUcPpKwOhBfqhlAb1w4cKrV6+WlJR8+eWXjc7AcRy5kTDx888/Dxs2rBVlmUwmlmWxBd005ndiFyJhdqvx58/tl44z9/5NNfQusauRNKPRiHcIapYXW9Dk3vBNz9OaLo6cnJzhw4dfvHjR9RafSqWypKRE6OIICwtr3TvBLg5PYBdHswz7N1vycuixfzF/807s0wsUyT3Erki6sIvDE+3cxdGaBBwwYEBUVNTVq1cbfTbKQcDc1wD5JZ43HPw1fML/yeOTQ0f/yXBkh9gFIdQyngb0pUuXcnNzOY6z2+3Lli0zGAwDBgzwaWUItZHl8ilZcJiyay8ACOo/wnLuMODpqciveBrQWq32oYceCgkJiYyM/Oabb3799dfY2FifVoZQG1kunwq67U7yfzomURYRY72WK25JCLWIp10Qo0aNunz5Mjlwh0cSkF+w5p+NevgF4WFwvxHmc9mqHreJWBJCLdKyPmi5XI7pjPwCZ6xnqyuVKf+717Wqx222wksiloRQS+EwCRSYrPnnlN36gEwuTFEk97BXFPN2m4hVIdQiGNAoMNmuX1R16+M4hVIo6fhke9k1sUpCqKUwoFFgspVecx31rOrSy1Z0RZR6EGoFDGgUmOxl1xSdujlNVHTuaSvGgEZ+AwMaBSC2rhoomTwsymm6snNPW3GeKCUh1AoY0CgA2UsLlMndXafT8cmsXssz9vYvCaFWwIBGAche2kj/BgBQcloencBUlbV/SQi1AgY0CkPQETEAACAASURBVED28uuKxNRGn1JoUpjK4nauB6HWwYBGAciuLaXjkxt9itak2CtL2rkehFoHAxoFIEZXSsclNfqUIj6FqbzRzvUg1DoY0CjQsPU1lEIlCwpt9Flak8JoMaCRf8CARoGG0ZbScZ3cPavQpNi1JYC3PEb+AAMaBRqmqlThpn8DAChVkCwohK3VtWdJCLUOBjQKNIzO7RFCgo5JZKor2q0ehFoNAxoFGkbr9gghIY9JwIBGfgHvGYgCwStvv5e5N4v8//Pesk9/Onjd/L9nb03rsWbFV8JDOlrD1mBAIz+AAY0CwacfLGCe+RnkCgCIYb/YdevfGkB18zl9WenuDxxnpmMSLFdOt3+RCLUUBjQKFL3GAK2KYurZq6qG3pP+N12bD3v/MKM8JoGprmzn6hBqBeyDRgEl2a69oYxveh46OgG7OJBfwIBGASXFri1R/DGgzfVVNXrKgTI63lSjU9Ny8vDxJ58WqViEmoFdHCigdLJrSxSaP0wy6Tm5Elb871aEHEBZ3tMpn+Xnq5IhZ0tp0Q/tXSVCnsEWNAooyTZtqTKu2dluqBJS7NjLgaQOAxoFlBRb5Q2nFnRjihSaZJu2HepBqC0woFFASbZrS5o7SAgA5YrYJHtVO9SDUFtgQKOAkmzTOR8kbEwZHZdkx8txIKnDgEaBI5QzU8DVyRu/0KijckVMor26HUpCqC0woFHgSLRVlStiPZmzTIktaOQHMKBR4EhidBWeBXQ5HZPA1FB4VWgkbRjQKHAk2qvLPAtoq0xpkAXFcHW+LgmhtsCARoEj0a7zsIsDsBsa+QMMaBQ4Euw1FXS0hzOX07GJNhxphyQNAxoFjkR7tect6DIFHidEUocBjQIHdnGgAIMBjQJHAtOCLo4yZVwSgy1oJGkY0ChABHNWFW/T0+Eezl8pj9bYa3xaEkJthAGNAkQCU1NBx3g+f6UiRoNdHEjaMKBRgNAwNZUKT/s3AKBCEa1h9L6rB6G2w4BGAULD1FR63AENAEZZEA8QStl9VxJCbYQBjQJEPKOvVLSgiwMAtHSUBkw+qgehtsOARgFCw9RW0lEtekmlIkYDRh/Vg1DbYUCjAKFharQt6eIAgEpFtIbCFjSSLgxoFCDiGX2LDhICQKUiJh67OJCEYUCjAKGxt7wFjX3QSNowoFGA0LAtb0HT0RoZBjSSLgxoFAhCFHIAMMiCWvQqLR2tAbNvKkLICzCgUSCID1FVyls2hAPIQUIcxYEkDAMaBYL4YKVO0eKA1tJR8TiKA0kYBjQKBJoQVUsHQcPNkwkpNYV3JkQShQGNAkFckELb8i4OAKgCdTjt9XIQ8g4MaBQIYoOVVXREK16o5UMiac7r9SDkFRjQKBDEBqu0Le/iAIAqXh0u93o5CHkHBjQKBPHBCp28NS3oKgiKkGMfNJIoDGgUCOKCla1rQet4DGgkXRjQKBDEBau0dGQrXqiD4HAMaCRVGNDI/3FchFKul3t6N0JHOlBH0BjQSKIwoJHfYw21eivDUq3ZmXU8tqCRdGFAI7/H1et1JmvrXqvjgyKxBY2kCgMa+T22oUZnauWtBbEPGkkZBjTye2y9vspka91rLSDneOAseEUOJEUY0MjvcYbaKnMrAxoA6lgZ16D3Yj0IeQsGNPJ7bL1e19oWNADUMcBiQCNJwoBGfo8z1FabW9kHDdiCRhKGAY38HtvQ+j5oAKhngK3HgEZShAGN/B7XUKtrUx80hS1oJE0Y0MjvsQ366rYFNNtQ68V6EPIWDGjk5ziONxtrLUyrF1DHUtjFgaQJAxr5N9ZYJwsOY/nWn2yCXRxIsjCgkX/jGmplYa25jp2gnqFYAwY0kiIMaOTf2Aa9PKw1V4IW1LIUV6+HNrTBEfIRDGjk39regmZ4oBQqzmL0VkkIeQsGNPJvrKFWHtqmgAYAWVgkhwM5kPR4esd5m8327bffZmZm6vX6Xr16vfTSS8nJyT6tDCFPcA16Wdu6OABAHh7F1uvpeNylkbR42oLW6XQbN26cNGnSyy+/XFNTM2rUKLPZ7NPKEPIE21Arb1sXBwDIwqJwIAeSIE9b0J06ddqxYwf5/6hRoyIiIs6dOzd06FCfFYaQRzhDrazNXRzysCi8XhKSoNb0QZeUlFitVuziQFLA1uvl4W3t4pCFRmJAIwnytAUtsNvtTzzxxAsvvNCpUyfXZzmOGzJkCEVR5OHnn3/ev3//VpRlMplYlpXJ8BhmU5jfiV2ImNj6GrNM2aYlsCyjUDOVJQaDwVtV+SOj0Sh8cpE7ZrNZoVDQdIuT01VwcHCzEdey1bAs+8gjj4SFhS1cuLDRGWQy2bJly4Tqe/bsGRoa2qJVEBRFBQUFYUA3jaSzWq0WuxDx8HydqSFMk9SWZcjl8qDYBGP+2dbtqwGD5/kOvgU8IZfLvRXQnmjBaliWnTlzZn19/ZYtWxQKhbvZBg4c2MSzCHkRZzZSSjUlb+unRR4WyRlwmB2SHE/3bJ7nn3322aKioh07dqhUKp/WhJCHuAZ9G89SIWShGNBIijztQ7h06dLXX3995syZlJSU6Ojo6Ojobdu2+bQyhJrllTF2ACAPw4OESIo8bUGnp6fzeLECJDFeOUsFAChVEADwVjP5D0ISgUfhkB/zynnehDw0ksVeDiQxGNDIj3mrBQ03TybEgEbSggGN/Ji3+qDh5rkqGNBIWjCgkR/zynnehDwsAgdyIKnBgEZ+jG3Qy0MjvLIoWVgUtqCR1GBAIz/GNdR6qw9ajkOhkfRgQCM/5s0+aLxmP5IeDGjkr3irGXjOWyOX5WE4zA5JDgY08lesobaNt4t1JAuNxGv2I6nBgEb+yosd0HDzbG9sQSNpwYBG/opt8NoYOwCQhUTwFhPPduiLayOpwYBG/opr0HvrCCEAAEXJgsM4Y73XFohQm2FAI3/FGuq8cq1RgSw0AgdyIEnBgEb+ivPelZIIGd46FkkMBjTyV1xDrXdb0PKwSM5Q58UFItRGGNDIX7He7YMm10vCodBISjCgkb/y4pWSCDmeTIgkBgMa+Su2wZsnqsDNs72xDxpJCAY08k8cy5mNsuAwLy4Sb6qCpAYDGvkl1lAnDwkHmTd3YFlYJNeABwmRhGBAI7/ENei9O4QDAOShOMwOSQsGNPJLrKFO5qVL9QtkYZGcEVvQSEJosQtAqDW4hlp5WLRXFsUwjF7/e8NZRtdUlFKqYPKIoqjISC+30xHyHLagkV9ivdXFUXPjcPbBpC7dyb9r2poRt98uPIxLSFr00UdeWAtCrYItaOSXOEOtd+5GWFsOqUMsL2wjj3TXXgt7eYkluM/NZ7e+azAYvLAWhFoFW9DIL3n9PG9CR0fGMTjSDkkFtqCRf1j4+bKvV6wUHi7qH/FLycYjVW+QhyzLemUtVXRkLAY0kgwMaOQfVq5ceW3oPyCxF3kYbFqem3LvNXmnm09fHO6VtWALGkkKBjTyH4m9oMtA8t/YK9aqlGGgiPPuGqroyL6WAu8uE6FWwz5o5JdimPpq2svjoAFAR0fFMXiuCpIKDGjkf8I4o42iLZTS60vGPmgkKRjQyP/EMXVa2pvXsRNgHzSSFAxo5H/i7DU6hXdOI3SiU0RhCxpJBwY08j9xTJ3OBx3QAGChlBzIQjizLxaOUEthQCP/E8fU6GiftKABezmQlGBAI//juxY04EAOJCUY0Mj/xDE1Vb45SAgAVXRELIMXHUWSgAGN/E8sU6ejfXUV0Co6KhZb0EgaMKCR/4ln9DqftaB1dEQctqCRNGBAI/8T68uAxhY0kg4MaORnKJ6PZeqq8CAh6gAwoJGfieAMZpnKRil8tHw8SIikAwMa+ZlYptZH53kT2IJG0oEBjfxMvF1f5bMhHIAtaCQlGNDIz8Sytb4bYwcARlkQBXwwZ/HdKhDyEAY08jPxjF7rs/O8iSo6MgYb0UgCMKCRn/HpEA5CR0fGYzc0kgAMaORn4nw5CJrAy/YjicCARn4mzl6rU/iwDxrwgnZIMjCgkZ+JZWur5L4NaDyZEEkEBjTyM/FMja8PEmqxDxpJAwY08icUz0cz9dV0uE/XUknHaJgan64CIU9gQCN/EsU1NMiCGYr26Vq0iiiNHQMaiQ8DGvmTOLtep/DtEA4A0NLR8diCRhKAAY38SRzj29MICS0dGcPUy4D39YoQahoGNPInsUytTy/EQTAUXS8LiWbxZEIkMgxo5E809upKRUw7rEirjIq34UAOJDIMaORPNIy+0senERLYDY2kAAMa+RMNU13p40HQhJaOwpF2SHQY0MifaBh9+wR0JR2DLWgkOgxo5E/arQ+6UhGlsWMfNBIZBjTyJ/GM3qf3uxJgHzSSAgxo5DfCeTND0SaZuh3WhX3QSAowoJHf0PCGSt+fRkhgHzSSAgxo5DcS+IYKeXt0QAOAVhEZb6+l2mdlCLmBAY38hoar1yrbqQVtoxQmmTpSxrTP6hBqFAY08hsavqGyvVrQAFCpiI6nbe22OoRcYUAjvxHPNWjbqw8aALR0tEaOAY3EhAGN/IaGN1TS7diCpqM0cnu7rQ4hVxjQyG9o+Pp2G8UBAJWKmHgMaCQqDGjkNxL4hop2bEFr6Sjsg0biwoBG/oECiOMN7XMaIaGlozVyHMWBxNSCe7uVlpaePHkyLy9v/PjxAwYM8F1NCLkKV8hMoLTKlO22xkpFVLzcVtpu60PIRQsC+uGHH6ZpOi8vLzQ0FAMatbNYlaySCm3PNVbS0XiQEImrBQGdnZ0NAGPHjvVZMQi5FaOUVcrC23ONWkV0PG1tzzUi5AT7oJF/iFbJtO3bgrZQSisvU/HYDY1E04IWtCc4jnv44Ydlspu5/8Ybb/Tq1asVyzGZTDzPC8tBjWIYhmEYjuPELqQ9xKioCso3LWj3N+8uZ5Qqm8lkMvlkvRJjNpvlcrnYVUid2WxWKBQ07YXkVKvVzUaclwOaoqg///nPwp85KSlJpVK1YjkMw6hUKgzopsnlcrlc3rot7HfiVLITsjCfLNr9JZHKGVUExXWQLWyz2TrIO20LjuO8FdAU1fzFuHwS0AqFoo3LIdGDAd00nud5nu8grZ5opayS8k1Au1fOKlPA1kG2MPnQiV2F1Ml/1z6rwwRE/kGjlpe170FCAChnlCE8nquCRNOCgH799dcHDx588uTJDz/8cPDgwfv27fNdWQg5iVPJyqiIdl5pOasM5TCgkWha0MXx1FNPPfjgg8LD7t27+6AehBrBM/ZQmqqhQtp5veWMIgQDGomnBQGdmpqamprqu1IQcoetraqy8lwTh/N8o5xRBWNAI/FgHzTyA2ytTmtl23+95awqFPugkXgwoJEfYPU6rUWEgDZwMp4Hzmxo/1UjBBjQyC+wtTqtRZzzcYwyJVtbJcqqEcKARn6Ara3SWd2f8OdLGNBIRBjQyA8weq0ofdAAYKCUrF4ryqoRwoBGfkDELo4GuZrBgEYiwYBGfkCsURxAWtA1GNBIHBjQSOp4q5ln7A12cfqgDTIlq68UZdUIYUAjqWNqKulojVhrN8jUjF4n1tpRB4cBjaSOramUixfQRpmCa9DzLF62H4kAAxpJnbgtaA4oWVgUV1cjVgGoI8OARlLH6rXyKNECGgDoqHgGu6GRGDCgkdQx1WK2oAFAHhWPAzmQKDCgkdSxNZXyGFEDOjoeh0IjUWBAI6lj9JXyqHgRC6CjNWx1hYgFoA4LAxpJGm+z8DaLPDRSxBromESmBgMaiQADGkkaU11BR2vAg/sf+448JoGpKhexANRhYUAjSRN9CAcA0FHxnKGWZ+ziloE6IAxoJGmMroyOSRC5CJlMHhGD17RD7Q8DGkkaU11OxyaKXQXQsUnYy4HaHwY0kjSmqoyOSxK7CpDHJDDVGNCovWFAI0ljdGXyWPEDmo5JxJF2qP1hQCMJ4zi2VkdHi90HDUDjQA4kBgxoJF2MXisLjaQUSrELATo2kakqE7sK1OFgQCPpYnSltAT6NwCAjuvEVJcDL85NA1CHhQGNpIupLpfCEUIAoJRqWVAoW4tX7kftCgMaSRdbVU7HiD/GjqDjk+3aErGrQB0LBjSSLkZXSsd1EruKm+j4ZAYDGrUvDGgkXfbKGwpNithV3ETHdWJ0pWJXgToWDGgkUTzLsLU6uQROIyQU2IJG7Y4WuwCEGsfoSuXRCZRc1F2UY7b8N7OwpAIAYmX2J4K1c598WnhSqVT866UXu3fvLl59KMBhQCOJYiqLFZpkkYu4evicMvqcvT8A0BT/bPC2n+x9rSAnT6qzVk4cMxIDGvkOBjSSKHtlCR0vgQ7oLgNh5N8AgAG4kZ/T7Za7Lqm7kmeUBfvELAx1ABjQSCr+89265597Tni4eGT3Q6V1mx/5J3lo5cS8Zj9xWdWlp7VYCGiEfA0DGknFjp27zPf9PxgyjTxMLZm3fORMs6obeSh761bxSrspT9U5zVIMEWLXgToMDGgkJcpgCI4CABnw3ZjKaxHpIAu6+ZSod70irqi7/lm/R+wqUAeCw+yQFHWxldfQ4QYhnaUhT9W5p7VY7CpQB4IBjaQo3XJdgl29hcrEOKY2mLOIXQjqKDCgkRT1thReVHcTuwpnLCW7rkzqYb0hdiGoo8CARlLUW5ItaAC4ou7cy1IkdhWoo8CARlKUbrl2UZ0qdhWNuKDudqvlmthVoI4CAxpJTjhriGIMxUrx73Tl6nxQj77mfLGrQB0FBjSSnF7WoivqzhyIP67OVW5Qt96WQjnPiV0I6hAwoJHk3GbOP6/uIXYVjWuQhejoyG5WvO4oag8Y0EhyBpiu5ASliV2FW+eDuve1FIhdBeoQMKCR5PQ3Xz0bfIvYVbiVG9T9Vgxo1C4woJG0xDB1EayhQCmVO125OqPuOcB0RewqUIeAAY2kZYD5ytngnrwErrzhzpngnumW6yrOJnYhKPBhQCNp6W++mhPUU+wqmmKSqfNVKbdZcLAd8jkMaCQtg0yXJB7QAHAyuNdg4yWxq0CBDwMaSYiC4vqb8k4E9xa7kGacDE4fbMaARj6HAY0kpL+spkCVXC8PFbuQZpwISR9svCSTbj85ChAY0EhC7pTrDofeJnYVzaugo2vloT2DWLELQQEOAxpJyB209miw+Le28sTB0AF3hmNAI9/CgEZSoQS+n1x/PDhd7EI8sj9s4J0RjNhVoACHAY2k4hal7SwbZZAHi12IR46E9L0thJFx2IhGPoQBjaSij9KSySaKXYWnDLKgiyZ5ZF2Z2IWgQIYBjaSB59OVlj12vwloANhRo4itwov3Ix/CgEaSYCu5auVl17kwsQtpgV01iih9CW/De8giX8GARpJgOrb7lDVI7CpaRs9Q9WHx5vNHxC4EBSwMaCQ+3m4z5ew/YfGzgAaAqrhuptNZYleBAhYGNBKfKWe/smu6npOLXUiL1UR3sV2/wBnrxS4EBSYMaCQ+45EdoXfeLXYVrcHKaHWvweYzB8UuBAUmDGgkMnt5IavXqnsPEbuQVgoePNZ47Dexq0CBCQMaicx4ZEfIsLtB5q+7orr37TzHmnOPil0ICkC02AWgDsRisTz/4stW2//uRRJKsXNCypYaE2t/Pnji1Cm4Y4SI5bUSRUVMmlG39T9BfYaChG8Eg/wRBjRqP0uXLv3m4GV24J+EKUtCT6+xdF9qSwcAmW6neKW1Bm+3btmy5fr16wAwobo2652XCtXxwrNBQUFPP/20SqUSr0Dk9zCgUbuSpdzGjvwb+f8Q04VhNw6Pu+VjkAUBgGzfF5yotbWUuezqj/buG+qqAWBrUMyH8fl/LZYx/M1GNHVo9fjx49PT/ePaT0iaMKCROGie+aD0y7cTZxll/jf8+SaOZYdOZwc/BADZANcL/98znW75OH46eTI8d7uoxaFA4K9HZpC/m6P7sUQZvzN8mNiFeM3clNkP1e65ux6PFiKvwYBGIniqevN9tQde7vQPsQvxpip55KyU1z8o/aKntVjsWlCAaEFAMwyzadOmJUuWnDp1yncFoYD3RPXWGdXb/5L6vpaOFrsWLzsf1GNB4pPfFr7d35wndi0oELSgD/qBBx6oqqoaOXLk5MmTFy1a9Pjjj/uuLBSQVJz9/YSyO6p/nZb6frkiRuxyfGJT5GijTL26cMHyblK/9S2SPk8D+siRI8eOHbt+/XpISMjIkSOff/75Rx99VOa3JxegdsbWaM0Xjk6qPLqOC7m7+2f+ctuU1vktfNhlddevq2eEbviofuiE4AEj6fhksYtCfsnTgP7tt9/GjRsXEhICABMnTqyoqLh8+TIOIUJO6uvrZz77YsmNYgCIU0KfUOrWMOrWUEop4883wKrzZVnx49mATmeiSJnw571lP2b8qVthvnr/ZlYVYu3a19q5DxMSCXIaAEJCQnr16iV2mUjqPA3o0tLSpKQk8n+FQhEXF1dWVuYa0DzPf/jhh0LL+i9/+UtyciNtB9P+X4BvasxrbW3tuXPn3D1rMpmCgoIoN2dtmUym4GC3EWA2m4OC3I7rauK1PM+bTCbyFeXKarXK5XKabnx7+qhgnueNRmNoaOM/pVmWtdls7l5rsVhomnZXcLW+luN4cHNanJznE6MjHKdQwCs4RskzjLHhCb46tGdoBA0yCo7U8HsquQWX+HwjAADcqKAMJ2DnR40uljPVwfkdUF/Z+JtlWTi8Fi5lNl4TAPz2Cck+Z9p8sBrdrRSuH+Pqyt0+W1MM+ZTbgq1GOPMrVF1v9NmG+tp7/vYPAJBR1OCkyCm3HBvTNS4hVGVhOK3JWmxhrscmmkEGABwPVv5/P0Z54G02O69Q8W62v1Iui4iIaPSppvdSm81GUZRCoWi84IaG0NDQ1u2lrX6W4ziLxeLu2aY/Vkaj0d07bbakpnOgiSVzHNe3b9+YmKY66OSxSapb72hiBkKpVLrb2gJPA5qiKJ7nhYeO/3ei1+vl8pvXjWSYxm97zFuMwDUV0JzJwJuN7p6V2awcz7rrYKGsFp5yWx5ltfLgdtVNvJbnecpq4d106lB2OwfAu9n1KauFB87decAym7V1BfM8T1nNvNzN35jjwGZz91rKZuPtMt7Nrs8Z6mJiYt0FtNVmq7Pa/1AJUDaQW0BZD4pyThYdnmIDmYFSQDDEJcOU32erqe5st9s0CfpGF1s6YURoWFhEROPPXhufkZwiVyobf/by3Xf16t7Q6FNMckRe6PD0ro2/0BCbXlOj6dy58Wd1IRk8z8XHN/5syfgRkVERoaGNP5s/amjX1FSavrlL1ABsBAAAlYoNCWa0Rfm91bFqngUAGYAS/nfzWZZjzSZzUkQ4QON7hc1qdffp4HkeLGZ3eynY7TwAzzS+l8psFt5MudtLm/5YNbcPu/9YcRxls7l7lrLbeYpyt5fKbBaO4t1lXKs/VgBAWd1uQ57jwGLizWp3rwUA8N5NdjwN6MTExPz8fPJ/hmGqq6sTExu5fRxFUQsXLnT3FS1QTXmymRmMxtQZQdjH3TSGYRiGUaub3FcQQENDQ1iYP91MSxS4lTxhNpsVCoW7Rr3XeZqAEyZMyMzMNJvNALB3797Y2Fif9qAdPXq0tLTUd8sPDAUFBWfPnhW7CqkzGAy7d+8Wuwo/sGvXLqPR7c9WROTk5JCrr7QPTwN6xIgR/fr1u/vuu+fNm/fEE0/MmzdP6Mfwhc8+++zkyZO+W35g2LVr1zfffCN2FVJXXl7+yiuviF2FH5g7d65OpxO7CqlbuXLl3r172211LWiob9u27eeffy4tLd2wYcMddzTfBY4QQqgtWhDQCoVi+vTpvisFIYSQI+93dev1+mYPEjbLbrcbDAa9vvHj44gwmUxWqxW3UtPq6uo4jsOt1CyO4+rq6nBDNc1msxmNRq9spfDw8GY7iqkmBsy1wrRp07xyQMZut8vlchzF0TSWZXmeb7cDyn6K53m73a5UKsUuROpsNptCoWh2ZG4HxzAMRVFeOQK3f//+vn37Nj2PlwMaIYSQt2ATFSGEJAoDGiGEJAoDGiGEJAoDGiGEJEpyAwB4nl+xYsWWLVuioqJeeuml/v37i12RFJWUlGRmZp45cyYlJWXu3LlilyNR1dXVa9euPXz4sM1mGz58+HPPPdfE5c06skWLFp06daq+vv6WW26ZPXt29+7dxa5I0ubPn69QKN588812WJfkWtDLli1btGjRCy+8MGjQoLFjx1ZWNn7xyQ5u69atmzZtunjx4rZt28SuRboOHTp06tSpBx988K9//euGDRtmzJghdkUSVV9fP23atDlz5jAMM2LEiIaGxq8LiABgzZo1K1eu3LhxYzutj5eYnj17bty4kfz/vvvu++CDD8StR8q+/PLLMWPGiF2Ffzh58qRCoWAYRuxCpC4sLOzIkSNiVyFR5eXl6enpX3zxRf/+/dtnjdLq4mhoaMjLyxs+fDh5OHz48BMnTohbEgoM165dS0pK8ukVvvyawWAwm83//e9/IyMj+/TpI3Y5EvXcc8+98847Vqu13dYorS4O0qERFRVFHsbExFRUVIhaEQoEOp1u7ty577//vtiFSNeMGTO6dev2zDPPLFmyBK8K3aj169fbbLYHH3ywPVcqrRY0uc2MxWIhJ+aaTCbcV1Ab6fX6iRMnzpgxAy/11QTSqZqdnT1p0qTU1NQBAwaIXZG0VFdXv/XWW1lZWe28XmkFdHx8vFqtvn79er9+/QDg+vXrKSkpYheF/FhdXd1dd901bty4d999V+xa/MDw4cMHDRqUnZ2NAe3k/PnzRUVF5NIZNpvNYrFER0eXlJT4elyQtLo45HL5gw8+uHLlSgDQ6/UbN26cNm2a2EUhf2U0Gu+7775hw4YtXrxY7Fqkq7a2tqqqivz/6tWrOTk5zV7BpwMaPXq0zWarqampqan5+uuv+/btW1NT0w6jNqXVggaABQsWTJw4ceDAgeXl5VOmTBk7dqzYFUnRtm3bHn/8cavVarPZoqOj//SnP5FvNeToxx9/PHjwxqbMiQAAAMRJREFU4Pnz57/77jsy5dKlSxqNRtyqpKa4uHj06NFJSUkymayoqGju3LmjRo0Suyh0kxSvZsey7MWLF6OiopKTk8WuRaLI9bKFh0qlson7z3dYVqvVZDI5TomIiMBr2LqyWq0FBQU8z6empuK5PM0iXRzh4eHtsC4pBjRCCCGQWh80QgghAQY0QghJFAY0QghJFAY0QghJFAY0QghJFAY0QghJFAY0QghJFAY0QghJFAY0QghJFAY0QghJFAY0QghJ1P8Hf5UA1iVZg9QAAAAASUVORK5CYII=", "image/svg+xml": "\n\n\n \n \n \n\n\n\n \n \n \n\n\n\n \n \n \n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "histogram(μi1,bins=0:0.1:4,normalize=true,legend=false,title=\"Histogram of 10000 averages with ρ=0\")\n", "plot!(μi1->pdf(Normal(μ, s1), μi1))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Task 3\n", "\n", "Redo task 2, but now use `ρ=0.75` (the other parameters are unchanged)." ] }, { "cell_type": "code", "execution_count": 97, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Average across the simulations: 1.989\n", "\n", "Std across the samples (with ρ=0.75) and in theory:\n", "simulations theory\n", " 0.534 0.206\n", "\n" ] } ], "source": [ "M = 10_000\n", "\n", "(T,ρ2,σ,μ) = (500,0.75,3,2)\n", "\n", "μi2 = fill(NaN, M)\n", "σi2 = fill(NaN, M)\n", "\n", "# Monte Carlo simulation\n", "for i = 1:M\n", " y = SimAR1(T, ρ2, σ, μ)\n", " μi2[i] = mean(y)\n", " σi2[i] = std(y)\n", "end\n", "\n", "# Theoretical results\n", "y2 = SimAR1(10_000,ρ2,σ,μ)\n", "σy2 = std(y2)\n", "s2 = σy2/sqrt(T)\n", "\n", "printmat(\"Average across the simulations:\", mean(μi2))\n", "println(\"Std across the samples (with ρ=0.75) and in theory:\")\n", "printmat([\"simulations\", std(μi2)], [\"theory\", s2])" ] }, { "cell_type": "code", "execution_count": 98, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAFACAIAAADrqjgsAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3dd3wT9f8H8E922nTvlrS0hdIySkEoyKrsLwKijJ8KoggqylC/LhREwAEOUFFUZCjiRGXIUGS0yCpllFVauke6m3SmadYl9/vjo/G+yaWkbZK7pO/ngz/I3efu3kmTVy6f+9wdhyRJBAAAgH24TBcAAACAHp/pAszp9fpnnnmmR48eb731FtO1uAm9Xr9v375bt27J5XIul7t161a7rFapVObl5XG53Pj4eIlEYq2ZVqvNz89XqVR9+vQJCAiw1sxoNObn5zc2NkZGRkql0na2W1paWl1dHRgY2KdPny49AWA/BoNh8eLFERERb7/9dvstZTLZ+vXr77rrrqeffto5tXWUQqEoKCjw9vaOj48XCAQMV0M63ooVKwQCwdtvv205q6GhQSAQeHl5maa0tbUhhPr372/7+q9cubJt27a8vDw71Op2DAbDxIkTEUIcDsff3z84ONhay5ycnI0bNz700EO9e/fG743S0lLalm1tbcuXLxeLxbiZRCJZsWKFVqs1a2Y0Gjdt2uTv74+b8fn8uXPnKhQKyxUeOHCgZ8+epvdkSkpKbm6uZbPMzMwhQ4aYmsXHx6empnbkxQCOotVqEUIJCQmmKTKZbNu2befOnTNreeXKFYTQ7NmznVugTWpqambNmsXj8fAbLCgo6PPPP7/jUr/99ls7AfvFF1+YWi5fvpy2zcaNG62t3Bl70ARB6PV6g8FgOYskSb1ez+FwTFN4PN7YsWNjYmJsX/+RI0fWrVu3a9cu2KWylJ6efvLkybFjxx4+fNjLy6udlt98880HH3yAEBIIBDwej/bvhc2bN++3334bMGDAokWLDAbDjh07PvjgA7lc/vXXX1ObvfPOO2vWrAkPD3/77bf9/Px+/vnnn376qaio6MyZMyKRyNTs8OHDs2fPFovFL7/8cmxs7F9//fXLL7+MGzfu8uXLPXr0MDXLzc0dP358S0vLwoULhw8fnp2d/eWXX957771paWmjRo3q/AsE7IHL5Y4dO5b6LZuVlfX0008vW7bMVf46KpVq0qRJWVlZ48ePnzVrlkKh+Pzzz5ctW6bT6f773/+2s2BwcDDeB6LS6XRnzpxBCFk+/UGDBgUFBVGnUF83c3b55mnfiy++iBBat26d5az6+nqEkFAo7Mr6161bhxDatWtXV1birr766ivU7le0ycGDB7du3Xrp0iWNRhMeHo6s7EEfOnQIIZSQkNDa2oqnNDQ0REVFIYTOnj1ralZcXCwUCr29vYuKivAUvV4/YcIEhNAnn3xiaqbRaKKiojgcztGjR00TX3jhBYTQggULqNudMmUKQmjDhg2mKd999x1CKCkpyWg02vJSAGf6/fffEULLli0zm87aPeh33nkHITR9+nSDwYCnZGdne3h4eHp6VlVVdXRte/bsQQglJydTJ+I96MOHD9u+Htb1QZMkefXqVQ8Pj379+pkmGo3GCxcuFBYWyuXywMDAmJiYESNG4L2wW7duVVdXI4RKS0szMzNx+4SEBFOvqF6vP3PmTE5OjtFojIuLGz9+vOm3OVVhYWFqaqpare7bt+/EiRN1Ol1OTo6fn1+vXr1wg4qKitra2tjYWH9//8zMzIsXL2q12iVLluC11dXVXbx4USaT6XS66Ojo8ePH+/r6Utff1NRUVFQUGhoqlUqLiorS0tLUavXw4cOHDx+OG6hUqj/++EMmk0VGRk6bNq2dXl0qjUaTmppaWFjI5XL79et3zz338Pl//01bW1vz8vKysrIQQnK5HL84EREROHwtzZgxw5Yt7tixAyG0atUqU4X+/v4vv/zyc889t2PHjtGjR+OJu3fv1ul0S5cujY2NxVP4fP6bb76Zmpq6Y8eO5557Dk88duyYTCZLSUnB+Yu98cYbW7du3bNnz5YtW7y9vRFC5eXlx44dCwwMxF/22Pz58zds2HDjxo2LFy/efffd7Zd969Yt/FYRCoVJSUkjR47kcv89Qp6VlaXT6QYOHGjZ54hnJSUlmV5YhNDt27fT09Pr6+sjIiImTJhg9pKWlZUpFIo+ffp4eXllZGRkZmYaDIZnn30Wb7Gqqury5csymcxgMPTq1Wv8+PG0f+uWlpajR49WVFSEhYXde++9AQEB2dnZGo3mrrvuov7ixC/OqVOnamtrAwIC7rnnHlMPlYlerz937lxxcXFzc3NgYGDv3r2HDx9OfTpURqPx2rVrZh9AmUwml8v9/f1Nf02EUH5+vlKpTExMxDtY1I9tSUlJYWEhQqiurs70qYyMjAwJCaFuq66u7ujRo3K5PDY29t577/Xw8KAtyZJCoTh27FhdXR2fzx8xYsTQoUNtXLAdO3fuRAi9/fbbpjdGv379HnnkkZ07d/74448vvfRSh9aGf00uWrSoq2V19JuhEzq0B23ZB11SUpKYmGhWdlBQEJ5rClCqjIwMPDczM9Os30MqlZp1XBqNxpUrV1I/rgMHDjx8+DBCaNasWaZm+C+0e/fuadOmmVrW1taSJDlnzhzq4gghX1/fH374gbqVAwcOIIReeOGF1atXUz9gCxcuNBgMaWlp1F89MTEx1vp/qU6ePEntBEAIxcfHX7t2Dc89e/as5Svz1ltv3XG1JEla24M2GAw4Terq6qjTCwoKEEIRERGmKfiX3Z9//mm2uJ+fH0KouroaT3n22WfR/+4XY2PHjkUI/f777/jhrl27EEL/93//Z9YMv7VoD2+YlJWVWf6ETEpKKigoMLV5+OGHEUL79+83W7a0tJTL5cbFxZl20uvq6qZOnUpdlUgkWr9+PXUp/LHcu3fvuHHjTM10Oh1JkhMnTjSL18DAwIMHD5pt98SJE8HBwaY2Xl5e+/fvT0hIMK0HU6vVTz75JPW9x+VyFy9eTG1z/fr16Ohos6c/aNCgdl6xyMhIiURCPaiAn0hiYiK1mVQq9fb2xtsy64PGr6eZLVu2kJQ96F27dlF3laKjo6l/kXbs27fPdGADe/DBBzUajS3LWoPfwOHh4ZbbQgj95z//6dDaysvLeTyeh4dHY2MjdXon9qBdIKAnTZqEEHr22WevXr0qk8kyMzN37959//3347nnzp177LHHEEKvvPLKiX80NTWRJCmTyfCwgf/+9783b97Mzs5eu3Ytj8cTi8U3btwwrf/TTz9FCMXExBw8eFAmk2VkZEyePDkiIgLRBXRUVFR8fPz27dvT09P37NnT0tJCkuTkyZNXr159/Pjx27dvZ2Zmbty40dfXl8/nZ2ZmmhbHAd2zZ8+AgIAvv/zyypUre/fuxcMV1q9f7+vru3z58r/++uvs2bN4T9b0BK25ceOGWCzGu6XZ2dk3btzAYRcQEFBeXk6SZGNj44kTJxYvXowQWr58OX5liouLbfmTWQvo0tJSnBdm0wmCwDGBX3mSJAMDAxFC+fn5Zi0HDx6MEEpLS8MPcefdnj17zJo98cQTCKFNmzbhhytWrEAIrVixwqzZZ599hhCaN29eO8/l1q1bI0aM+PLLL8+dO1dYWHj69Gn8hunfvz9BELjN0aNHaV9zPJTonXfewQ/b2toGDhyIEHrooYdSU1Nzc3P379+P9wCox4JwQEdFRSUlJX399dfp6enfffedXq8nSXLkyJFvvfXWyZMnc3NzL1++/M4770gkErFYTH2hCgsLJRIJn89/9913CwsL8/Pz165dK5FIcCqZwtdoNOK3SkpKyu+//56bm3vs2DH8vbh8+XLT2vr27cvlcteuXXvz5k2ZTHbp0qVt27aZdR+ZWbBgAULozJkzpmctFou5XC6Hw6mpqcETb9++jRCaNm0afmgW0Ddv3ly/fj1CaMaMGaZPpUwmI/8J6KioKIlE8vbbb6enp6elpeHvvHHjxrVTlekvxePxYmJiDh06JJPJTp8+jXefn376aWqzEydOHLOBqT3uuBs5cqTZ5q5fv46rvWNhVLjT9dFHHzWbjgN63LhxQ4cOTUxMnDp16pYtW0xdhbScF9BSqfRuC/jFbSegjUajQCCIi4trZ/3W+qBxNj311FPUiWvXrkUITZ06FT/UaDRBQUEcDicrK8vURqvV4p9ylgHt7e1teo+2Y//+/QihRYsWmabggOZyuVeuXDFNPHjwIN4FeO2110wTVSpVQEAAj8dTq9XtbAK/p832iHGuUd+sGzduNO282M5aQF+6dAnv7Fgugr8LCwsLSZIkCAI/r/r6erNmOJF//vln/DApKQkhdOLECbNmOJFfffVV/HDhwoUIoffee8+sGe7pmzRpUoeeHUmSc+fORZQdfIPBEBkZyefzzf648fHxXC63rKwMP8ShY/aOqqqq8vHxCQ4ONu3E4YAOCQlpbm6+YyX4t/BLL71kmoIXX7NmDbUZ7penBjR+j6WkpJj6TEmS1Gg0cXFxPB4P/+1w79/48eNte1X+9u233yKE1q5dix+eOHECIYR3in/88Uc8EX81fvTRR/ih5SiO9vugEUJ79+41TVSr1eHh4RwOh3aQjwlBENHR0Vwu99atW6aJDQ0Nfn5+HA4Hv/cwax04Zkw/jHD/xowZM8y2WFFRgRDy8PC486v2D6PRiNPj1KlTZrNwQPN4vPDwcNNx8ri4OGrlZpx3oopcLs+2kJeX1/5SHA7H19e3tra2uLi4o1vEmfjaa69RJz7//PMeHh7Hjx9XqVQIofT0dIVCMXHixAEDBpjaCIXCpUuX0q5z4cKFoaGhd9z01KlTBQIBjjOqe+65hzpQLCUlBf/H9PFDCHl6eg4bNsxgMMhkMmvrV6lUx48f9/T0NHXmYvjJ4o+uI6jVaoSQj4+P5Szcd9Ha2mpqRtsSN8MvfjsrxD34Zs3MuvUtm9nu/vvvRwiZ/kBcLveRRx4hCAInPnbu3Lm8vLwJEybgQ6AIoe+//x4htHr1auqqwsPDZ86cKZfLTdGDLVmyhPaFar8ShNDhw4e5XK7ZX5ba+Y7hY6SrVq2idnGIRKInnnjCYDAcP34cIeTj48PlcvHBmztWYoKP5aampuKH+D/r1q0TiURmE3HLToiNjZ09e7bpoVgsTklJIUmypKSknaWuXLlSWlo6ZcqU/v37myb6+/vjPSHqcLdPPvnkUxuYupvafx+q1Wqj0WjjU0tNTS0uLo6JibnnnnvMZk2ePPnUqVMajaaqqkqlUp09e3bw4MEFBQUzZ860NmjKeQcJV65cifdeqRoaGvBv4XYsWrTogw8+6Nev36RJkyZMmDBp0iTqn8cauVwul8v9/PyohzUQQvhAB/5uuOuuu/Lz8xFC1OMhmLVNWLZECDU1NW3atOno0aOVlZW1tbWm6bgDh8qsQ9zPz08gEHh6epodPMH9j3V1ddYGDubl5REE0bdvX7PM6t27t6+vL37u1E5Me8Ff+ziFzbS0tCCE8HEe096BSqUyqxA3M3U+4paWCatUKk1ra2e7Zs2suX79+saNG69duyaTyajbUigUpv8vXLjw/fff37179/PPP4+n7N69GyGEf+8jhNra2nJzc0Ui0RdffGG2fnxArLS0lDqmivatIpfLP/jggxMnTlRVVVFD0/RWaWhokMvlERERZp8LqVTq6+vb3NxsmnLt2jWE0KFDh06dOkVtiTsfcGeUp6fnww8//OOPP8bGxk6ZMmX8+PGTJ0+mPWxDFRERkZCQkJGR0dLS4uPjk5qaGh8fHx8fP2zYsJMnTyKEDAbD6dOng4KCqLs1HRIfH282Be/3UD8+lnB8419dVIMGDUIIZWdnm6ZY28Gypv03mEgkMjvI1A48buqJJ54wO9iAELrvvvtM/+fxeKNHjz516lTfvn2zsrKOHz9+7733Wq6NdaM4LL377rtRUVHbtm07cuTIkSNHEEIJCQmffPLJ5MmT21kKv9ZmwYeFhoZmZ2fjlx5/c+LdOiqzoxAmZgMYEUJKpXLEiBG5ubmJiYlz584NCAgQCoUIobVr15p+6Zt4enqaTeFyubQTEULtfGm3/+yam5tbWlocEdC4H6OhocFsOkmSTU1N6J/XTSAQeHt7K5XK+vp6s4DGy5pOKcTtLb/JcDPTXwH/x3K7Zs1o4S5OkiTHjx8/Y8YM3DgvL2/Xrl3U3ZY+ffoMHz48IyPj5s2bAwcO1Gg0e/fu9fHxmTlzJm6A+ysIgti+fbvlVvz9/c12gizfKnK5PDk5uaysbMiQIY8++qi/v79AIDAajatWrTK9Vay9IfEmqAGNX/CffvqJtqXp/7t27erfv//XX3+9d+/evXv3IoSSk5M///zz5ORk668ZmjBhQm5u7tmzZ0eNGnX16tVnnnkGT1y3bl1RUVFDQ0NDQ8NDDz1ke2yZ6cR7HiGEe2ws39j4g4A/0Rge3HnHMkwjXqy9sW15g1E1NTUdPHiQy+Xi4xx35OvrO2XKlF27dl26dMlVA5rL5S5btmzZsmXl5eVpaWn79+8/fPjwfffdl5mZ2c4XOB6eVVdXZzmrpqYG/fNzBr/0+A9PVVlZaWN5O3fuzM3NffLJJ/H4M0ylUq1cudLGNXQCfna0uxv42Vn2BthFdHS0UChsbm5WKpW4Bqy6upogCHymIp7Sp0+fzMzMiooKs18wuFPP9MsgPj7+zJkzlq+2ZTNE90cpLy9HdLtjVKtWrdJqtcePH8dHm7Hdu3fjkSFUCxYsyMjI+Pbbbzdt2nTgwIGmpqannnrKFCX4+Xp5eSkUis4F05YtW8rKyl555RV8QhBWWVm5atUq00PcnVpVVWW2rNFoNHuX+vj4NDU13bx5s/0z44VC4apVq1atWlVQUHDq1KlffvklNTV18uTJOTk51gZcIoQmTJjw+eefp6am4lPMcFfGxIkT161bl5qaimOr0/0bnabX6xHdhxp/EKgdFAkJCZa7R5aMRiPez8XvNPyuo7LlDUb13XffqdXqqVOnRkZG2rgI/ulgrZvOBQLaJDIycsGCBQsWLFixYsXGjRsPHTqEAxrvsVruv4SFhdXU1BQWFlJHhtbX1xcXFwsEAvyi33XXXQihM2fOGI1G6qfur7/+srGqGzduoH8OoZjgsW6depY2wVcJKCoqamxspH695+XltbS0hIWFWe6+2YVAILj77rvPnDlz+vTp6dOnm6anpaUhSpc6QmjMmDGZmZl//fUXdWJOTk5NTU1UVJRp4NeYMWN27Nhx6tQpU8cCQkij0Zw/f57L5ZpGVY8ZMwYhdPr0abM/E94unmvNjRs3goKCqOmMEDINzqWaO3fuiy+++P3337/77rtm/RsIIS8vr7i4uIKCguzsbMtxn7agfatcvXqV+lAikSQkJNy+ffv69ev4lzt24cIFfCDOZPDgwTKZLD09/cEHH7Rl63FxcXFxcYsXL547d+6ePXtOnjz56KOPWms8btw4Ho+HA5rH4+Hu1OHDh+MeD1sCGo8otyUlOwr37VDh15C6u7Z06dJ2ToU1MfVC9OvXLygoqLi42GxQpuUbu334W79Dw5/xu8Lqt6ztRyc7rSvD7HQ6neVgho8++ghRDnN/8803iHLQ2WTZsmUIoSeeeII68fXXX0f/O6AKD/z6+uuvTVPy8vLwfpPlKA7qoWfqdOripj0O6oUvTOOgzRYXiUSWoy9xLpw+fZq0Dh9cMnvWeEHqKCv7juIgSRL/wB87dqzpCLher8dHPn/99VdTsytXrnA4nB49epgG3pH/DMZYtWqVaUpjY6Ovr69AIMjOzjZN3LJlC0JoypQp1O3iPxN1QN65c+c4HI5UKrW8DAiVVCoVCoXUEallZWX472s5xgCn57Zt23g8HnX4M7ZhwwaE0IwZM0zj80yUSqXp//jzaXmdEDz9wIEDpil6vR4PZOrTp49p4rvvvosQmjp1qmkrOp3OdMTJNIoDv6MSExOpm8Y0Gg1u1tbWhof3UeEj0jt27LB8raiSk5PxX5B6Oty0adMCAwM9PT3NRvJYjuLAPcKm4VIm1s4kxKdTHzp0qJ2S3n//fYQQh8OhDpNVKBT4WGhJSUn7z6h9eJQq9bODD+TweLzbt2+bJmZnZ3/yySf79u2zXAP+1g8MDKQdlK1UKi1H1P3xxx9cLpfH41kOSMXYHtAVFRWBgYH//e9/f/vtt5s3b2ZlZW3fvt3f35/P51+/fh23uXXrFofDCQgIWLVq1datW7dt24bPH6mursa/uJcsWZKRkXH16tUVK1bgPt+cnBzTFjMyMjw8PLhc7vz58z/++OOXXnrJz88P73DZEtD4lJaQkJDdu3fn5uaeOnVq+vTpPXr0EIvFDg3onJwcT09PHo+3cuXKq1evZmRkPPXUU/hbwXQaCNmRgM7Ly1v8D5xf8+bNww+p55vodDq8Zzdr1qy0tLQTJ07gkwBTUlLM4gzvoA0bNuzQoUNnz57F/ZhSqdRs7N3HH3+MEIqMjPzuu+8uXLiwYcMGoVDo4eFhOuMGS0tL4/F4Eolk8+bNFy5c2LlzJ/7jmp0QZAmXkZKSkpqampeX98MPP0RHR+OOF8uAxgOi8dM3DX82UavVw4YNQwiNGjXqu+++u3LlyoULF3788cfHH3+c+re2FtB43IVUKt2zZ09eXt7x48fHjx+PK6EGtFqtxq/w0KFD33vvvfXr1yclJSUkJODBJKbANRqN+OskISHhyy+/vHDhAh5c/8ILLwQGBuIv14yMjIiIiFWrVh05ciQ7O/v69esffvihWCz28vKqrKxs/3V79dVX8VcCdQzohx9+iCea7fdYBrRGowkICODz+cuXL9+yZcu2bdvw2LiuB3RUVJRUKv3pp5/y8/N///13fMzw2Wefbf/p3FFVVVVISAiHw3n++efPnz9/4MABvOalS5dSm+EBeffcc4/lGvAe4Ysvvki7/oyMDB8fnyeeeALn+7Zt2+bOnYt34amDLM2wPaDlcrlpkJNJcHAwdWeNJMlPPvmEesTMdCZhVlYWPrPApFevXpZX2Dp//rxp9Ju/v//atWvxT5vHHnvM1MZaQJMk+fLLL1OP2MbGxl67dg2PjTW1sXtAkyR59uxZsyPygwYNou6Kkh0J6HZ6dcwGIFdWVpr1KkyePFkul5utsK2tbf78+dRm/fv3NysPe/PNN6nnWIeEhJidgojt2bOHevTMw8Pjs88+u+PzqqurMzsgNmvWLHyGmGVA4wHRCCHq8GeqlpaWRYsWmQ2zFYlE1N9k1gLaaDTisfkm/fr1w3ua1IAmSVKhUDzyyCN4aAGPx5sxY0ZFRQUePEttptfrV69ebXamOI/HGzNmDB5QnJOTYzkqtGfPnqYThdqBB+qh/x2ljn+PI8qAaMwyoEmSPHr0aFxcnGm7ZmcSmm3O9oD+5JNP5s2bR31Gjz/+OPXkyU67evUqtWAOh/PUU0+ZrdlaQKvVatzZSN27p8rKyqIetsH8/Pzef/996kh2MxzS8XdUUSgUjY2NgYGBlpcDNhqNRUVFXC7XFDQkSRYWFopEImouFxcXZ2dn19TUiMXimJiY5ORk6uXQTFpaWvDAqR49epgaGI3GS5cu5eTkGAyG+Pj4kSNHWhvEXl9fr9FoQkND+Xz+V1999eSTT77xxhumy1LX19c3NzeHhobSXjmhqKgoMzOzubm5V69eKSkpfD6/tLSUw+GY+rNUKlVtba2vr6/Z8KmSkhIul2t2LrJcLlcqlREREbSXDaEiCOL8+fP5+fk8Hq9///7Jyclmx6+am5vr6+uDgoLuOCYXD8+knRUQEGA5ruDy5ctZWVlcLnfQoEHU3lIzhYWFuP+0T58+o0ePtnZ4ra6uLi0trbm5OSoqaty4cdaeuEqlSk1NxdeDnjBhgo2H141G48WLF3Nzc3k83tChQ/v169fW1lZTU+Pj42PZWV9TU9PW1sbn8y33DKjVpqen19TUSCQSqVQ6ZMgQ6suL/3zh4eG04/9yc3OvX7/e2trap0+fUaNG8Xg8fFDE8rCSWq2uq6sLCQnx8PBQKpW+vr5xcXGWpw60tLSkp6eXlZWJxeLw8PBBgwaZDe/JyckpKCjAzxd/fExX1GwHQRB4JH5UVJTpI0P+M1SZ+hHDCgoKhEKh5Vn1Go2mpqbGaDTiN6FWq62srJRIJGbfHPjzFRYWZjnAw+SDDz549dVXP//886VLl2ZlZWVmZhqNxuHDh9sy7tZGBoPhzJkzhYWFYrF4zJgxlmfJt7S01NbWenp6ml1lQafTVVRUcLlcy0WoK79x40ZRUVF9fb1QKIyNjR02bFg7zxchp/RBuxyCIEaOHInozgUCgBG4b+GZZ55huhAm4T1oW67R7DZcaRSHg8jl8gULFjz++OP9+/f38PC4ffv25s2b09PTx4wZY3kuEABOsHjx4r59+44aNSosLEwmkx05cuTDDz+USCSW5xMC9wYBjTgcTlpaGj46ZHLvvffu3r3b8lwgAJygvLycOqweIdSzZ89du3ZRe0hBd+CMPmj2a2trS09Pr6ioaGpq8vf3Hz58OL60IwCMMBqNN2/ezMnJUSgUQqGwX79+d999Nx7v351dvHjxr7/+mjx5Mh5z2R1AQAMAAEs572p2AAAAOgQCGgAAWAoCGgAAWAoCGgAAWAoCGgAAWAoCGgAAWMrZAV1bW9uJO8hZsv0WYd0TvD7twCfRMl0Fe8Gbp33OfH2cHdDLly83O2evc+yS8m4MXp926PV6fG8OQAtfURJY48wPF3RxAAAAS0FAAwAAS0FAAwAAS0FAAwAAS0FAAwAAS0FAg27H2KY0NsqZrgKAO4ML9oPuhCQbf/lUfe0M4vE8Ekf6P/gs4t757nwAMAX2oEE30vTbdqK+OnD1NwGvbjc0yZt//4bpigBoDwQ06C70FUXqa6cDF67mCEVcsWfAo6+2XT2tK8lmui4ArIKABt1F08EdPlMXcD288EOuxMfnP4+0HPuR2aoAaAcENOgW9NWlRF25Z/IE6kTP5An62nKdLJ+pqgBoHwQ06BZazx6SjJzG4f3PUXEOj+81Zobq/O9MVQVA+yCggfsj9Tr1tTOSEVMsZ3kOGafOSuTTp0EAACAASURBVCf1OudXBcAdQUAD96fJvSKQ9ub5BFjO4vkGCqS9NDkXnV8VAHcEAQ3cn/rGOY+k0dbmeg4Z33b1tDPrAcBGENDAzZEGQpNzyWPgSGsNPPoP0+ZfJw2EM6sCwBYQ0MDN6csLeAFhtP0bGNfLjx8coSu57cyqALAFBDRwc9qCG+LeA9tvI+6brLl92Tn1AGA7CGjg5rRFWaI4GwI694pz6gHAdhDQwJ2RBkJXelsYO6D9ZsKoPoaGOqOqxTlVAWAjCGjgzvSyPH6w1HR6t1VcrjA6QVuS45SiALAVBDRwZ5qCm3fs38BEsQN0xbccXQ8AHQIBDdyZtvCG6E5HCDFhrwHawixH1wNAh0BAA7dFGghdWZ4wtr8tjYVR8fpaGalVO7oqAGwHAQ3clq4sTxAayRVLbGnM4QsEEbG68kJHVwWA7SCggdvSyfKEPRNsby+M6qOT5TmuHgA6CgIauC19RaFA2sv29hDQgG0goIHb0pUXCqVxtrcX9oyHgAasAgEN3BOp0xoa6/hhUbYvwg8MJ7Uag7LRcVUB0CEQ0MA96SuLBKFRZrdQuQMORxgZpy8vcFhRAHQMBDRwT7ryAkFUB/o3MEGPXrrKIkfUA0AnQEAD96SvKBT26MARQkwQEU1UlTqgHAA6AwIauCddRaEgsuN70BGxuqoSR9QDQCdAQAM3ROp1hKJKENazowvyQ6SGhhq4hyxgCQho4Ib0VcX8EClHIOzoghwenx/cQ19T5oiqAOgoCGjghnQVhUJp784tK4iI0VeX2rMaADoLAhq4IaKqVBAR27llBeExeuiGBuwAAQ3ckL6mTBAe3bllBREQ0IAtIKCBG9LXlHXoHEIqQUSMvqrYvvUA0DkQ0MDdGFubEEI8b//OLc7zDUQIGVrghG/APAho4G701WWdGGBHJQiP1ldDLwdgHgQ0cDf6mi4HNHRDA3aAgAbuhqiRdboDGhOEw0g7wAoQ0MDd6GvLBF0M6IhoOE4I2AACGrgbfY2M37UuDn5IJFFXgUjSXiUB0DkQ0MCtGFubEEl2eggHxhV7csUSQ5PcXlUB0DkQ0MCt6GtkgtAu9W9g/NBIfV1F19cDQFdAQAO3oq8p44d3qX8D44dIidryrq8HgK6AgAZuhagps8setCBESsAeNGAaBDRwK/oaWRcHQWP8kEh9HexBA4ZBQAO3QtSV80OkXV8PH/agAQtAQAP3QWrVRo0aX0yji/gBocY2JalVd31VAHQaBDRwH/q6cn5wD8Th2GFdHA4/KIKQV9phVQB0FgQ0cB9EXSU/pIe91sYPkephIAdgFAQ0cB+EvFIQbLeAFoRGQjc0YBYENHAfhLyKb7+AhuOEgHEQ0MB9EPIKfrAdhnBgMNIOMI5ve9PMzMyLFy/KZLKHH3540KBBlg0OHDhw8eJF/H8Oh/Puu+/ap0YAbEPIq+zYBy0IkRLySkSS9jnqCEDHdSCgV6xYERoampaWNnDgQNqAPnbsWF1d3ZQpUxBCHHhPA+cyKBsRl8v19LbXCjkiD66Hl6FRzgsIsdc6AeiQDgR0amoqQmjIkCHttElOTl68eHFXiwKg4+zbAY3xQyP1deUQ0IApHQhoW5w4cUImk0VHRz/55JOBgXY4XwAAGxF1FQJ7nENIxQ+KIBTV9l0nALazZ0APGjSoV69ePj4+x44d27x58/Xr10NDQ83aqFSqdevWffHFF/jhAw88sGjRok5sS6VSQS9KO7rh66OuLOH6Bre2tt6xpU6nQwgJhcI7tjT6BKmry5AN63QnKpWK6RJYzV4fLk9PTy73DsM07BnQzzzzDP7P008/PXbs2J07d77++utmbUQi0axZs1JSUvDDnj17enl5dWJbJEl2bsFuohu+PtrmOs9e4zxseNa2BzSvR7QqI7+7vZIIoW74lG3nzA+Xnbs4TPr27VtTU0OzPT5/4MCBEydOdNB2QbdF1FXavw86KIJQVNl3nQDYrqvjoAsLC/fv34//X1n594ULqqqqjhw5MmzYsC6uHABbkSRRX+2IgDY01MDNCQFTOhDQixYt6tWrV3Z29gsvvNCrV6/09HSE0Llz51atWoUbJCUlJSYmjh49Oj4+furUqY888ohDSgbAgqFRzpX4cIQi+66WIxByJT6GRrg5IWBGB7o4Pv30U71eb3ro7e2NEHr44YdnzJiBp1RWVubl5alUqri4uKCgIPsWCkA7iPoqflC4I9aMezlgpB1gRAcCmrZfXCwWi8Vi/H+RSDRw4ED71AVARxDyKn5QhCPWjANa1IfmzCwAHA2uxQHcAVFf7bCADofjhIApENDAHRCKan6gY7o4gmEgB2AMBDRwB4Simue4Pmg5BDRgBgQ0cAeEwpEHCeurYaQdYAQENHB5BmUjhy/gejjk5C6OyIMr9jS0NDhi5QC0DwIauDyDwlFHCDHo5QBMgYAGLo9QVDuofwODE74BUyCggctz3Bg7jAcBDRgCAQ1cnuOOEGL8YOjiAMyAgAYuz+EBHRRO1ENAAwZAQAOXR8ireA49SBgYblDQXDsXAEeDgAauzahpIwk9z8vPcZvgenojDsfYpnTcJgCgBQENXJsB9284+P5evMAwuDkhcD4IaODaiHpHXYWDih8UTtRDQANng4AGrs3Rg6AxfmC4oR66oYGzQUAD10bU1/ACwxy9FT50cQAmQEAD12aor3FWFwfsQQNng4AGro1QVPMdvwfNCwwzQB80cDoIaODKjAZDs4Ln7/AbBvL9QwzKRtJAOHpDAFBBQAMXRjTWcX38OXyBw7fE5fF8Ag0NtQ7fEAAUENDAhTmnAxrjBYZBNzRwMgho4MKc0wGNwVBo4HwQ0MCFEQ21TtuD5geGwVBo4GQQ0MCFEfXVThgEjfGDwmEoNHAyCGjgwgyKGiecRojxA2EoNHA2CGjgwpxzIQ6MFxQO91UBTgYBDVyVUd1KGgiuxMc5m+OKJRwe39ja7JzNAYAgoIHrMtQ7r38Dg4EcwMkgoIGrIhTO69/AYCg0cDIIaOCqiHpnXGiUih8Ie9DAqSCggasi6mt4AU4aY4fBUGjgZBDQwFUZ6mucdhohxg8MI+rhchzAeSCggatyfhcHXHQUOBkENHBNRqOhud4JFxql4vuHGFoa4KKjwGkgoIFLIprkXC8/Z1xolIrL4/kGGhrrnLpR0I1BQAOXZKiv5geGOn+7MNIOOBMENHBJhBOvBE0Ft/cGzgQBDVyS88fYYXzYgwZOBAENXJKhvpofxEBA8wLD4FwV4DQQ0MAlEfXOu1Q/FZyrApwJAhq4JGdeqp+KHwiX7QfOAwENXA+pVZM6Dc/Lz/mb5kp8EGk0qludv2nQDUFAA9fz93XsOBxGts6DgRzAWSCggeshnH4VDio+HCcEzgIBDVwPUV/NY+IIIcYPDCMUsAcNnAECGrge51/HjooXGEY0QEADZ4CABq6HaGC2iyPcAAM5gFNAQAPXQ9TXMNnFERQOJxMC54CABq6GJA0NtfwABq6UhPECQg1NcmQ0MFUA6D4goIGLMbQ0cD28OEIRUwVweHyujz/RKGeqANB9QEADF0Moqhjs38D4AXDCN3AGCGjgYgiFs+90ZYkfBLf3Bs4AAQ1cjKG+msEhHBg/EI4TAmewNaBbWlpWr1597733Dh06tK6O/pY/VVVV06dPDwwMHDx48JkzZ+xXJAD/IhTV/KAIZmuAi44C57A1oNVqtUajmTNnTmZmpl6vp23z9NNP9+jRQyaTrVixYubMmSqVyn51AvA3or6GDV0cMBQaOIGtAR0aGrpp06b/+7//s9agqqrqzz//fOuttyQSydy5c2NiYvbt22enIgH4F6GoYuRCo1Rw0VHgHHbrgy4oKAgJCQkN/Xt0alJSUm5urr1WDgBm1LSRhJ6RC41ScSU+CCFjG1x0FDgW314ramho8PLyMj309fVVKBSWzZqbm6m74c8999z69es7sTmVSsVh6GqTLsFdXx9DVQnXP6S1a71nOp0OISQUCruyEo5/cEt5Ib9H766shJ2gc7J99vpweXp6crl32EW2W0AHBAQolUrTw+bmZtPeNJWvr++vv/46Z86cLm6OJEnq9wEw466vj1rdTIT06OJTs0tAa4N7iNpaPNzxRUYIueWbx16c+eGyWxdHr1695HK5XP736VW3bt2Ki4uz18oBwP6+VD8LwFBo4AQdCOiSkpLS0lKEkEwmKy4uxhN37NjxzTffIISkUumECRPeeecdvV5/8ODBvLy82bNn279e0L0R9dU8podwYPyAMBgKDRytA10c06dP12g0sbGx8+fPRwjl5+fzeLzS0lLTT8Xt27cvXLgwICBAKpX++uuvPj4+DikZdGOEotpj4Cimq0AIIV5QOHHjLNNVADfXgYDOzs62nEg9xBcVFZWammqHogCwwsCCQdAYXBUaOAGc6g1ch9FgaKnn+YcwXQdCCPEDQgwtDaSBYLoQ4M4goIHLIBrqeN7+HJ7dhh51CZfH8w0yNNQyXQdwZxDQwGUQiiqWHCHE+EFwPiFwLAho4DIM9TWMXyaJih8UQSiqmK4CuDMIaOAyCBZcaJSKB3vQwMEgoIHLYM9ZKhg/uAchr2S6CuDOIKCByyDklfxgVnVxwMmEwLEgoIGLIEmCbX3QgeGGhlpkNDJdCHBbENDANRhaGrhiD47Ig+lC/sURCLkSX6IJbu8NHAUCGrgGQl7JD+rBdBXm+MERBhjIARwGAhq4BkJRxaoOaIwfFEHIIaCBo0BAA9dAyKtY1QGNwXFC4FAQ0MA1EPJKfjD7ujhgDxo4EgQ0cA2EvJKVe9BwMiFwIAho4ArwGDsW9kEHRxD11YgkmS4EuCcIaOACWDjGDuMIxVyxxNBcz3QhwD1BQAMXwM4xdhj0cgDHgYAGLoCdY+wwfjAENHAUCGjgAtg5xg6Da9oBx4GABi6AULA3oAXBUqKugukqgHuCgAYugJ2DoDG46ChwHAhowHpsHWOH8UOkhKIKrmkHHAECGrCdoaWBKxKzcIwd9vc17RrrmC4EuCEIaMB2bD5CiAlCI6EbGjgCO+5gD4B1RF05PySyo0vV19c3NzfTztLr9WFhYUKhsMul/e3vbui+Q+21QgAwCGjAdkRdBT9E2qFFVCrVgORRSo2edq5e2fjaS8+9uW6dHYpDCOFu6Lpye60NABMIaMBezc3Ni5575TFe2Z8K46V128zm8jjos43rhw0bZrmgTqdrVNRpP6qlX++uJ95666233nyTdqZAKFLI63x8fGyvkx8iVd/KsL09ADaCgAbsdfv27WOnzr44e+AxvwdKOP5mcz2Ob8zPz6cN6DuoL0fTVqH719LO5K+QajSajgV0cA/ogwaOAAENWE3i7RuK1OX9ZiLEM5vFv/w9IyVZ4geEGlXNpFbN2qEmwEXBKA7AalGevHJRKGGRzuzC4cAlk4AjQEADVov14hYKO3aEkBH84B5EHZxPCOwMAhqwWk8Jt0jkCgEdItVDNzSwNwhowGoxEl6xkKVX4aDih0gJOQQ0sDM4SAhYLcaTu0Ps3D1okiwrK2ttbaWdGR4e7uFBcyRQECJVnTvs4MpAtwMBDVitnT1oktAdOXKkspKm51etVhsNhs5tUUtyx02fw+XTfDQIdev0qVN++X635Sx+iBT6oIHdQUAD9uKqlRojauZ50c7VyLJ+UcXuVdLdD7Ct2UAQnduoUa1Urc1E3iE08y7/2qo4RF+qpzfi8w0tjTwf8/HaAHQaBDRgL0GLokTVzo4wh0x+0HD3PJo5ihKU7uxR0oLQKKJWBgEN7AgOEgL24isVxa2d7KlwPkFYT31NGdNVALcCe9CAvQTK+jIVmy6E31J7+lx6cFQv2pnzevu9MHem15gZTi4KuDEIaMBe/GZ5CasCWl6iDYprm7+VdubtQwuIOpmTKwLuDQIasBdfqShpY1cXB0ckQcExtLOKtXxBs9zJ9QD3Bn3QgKVIrZqrVVepSaYLsVW9jkQIGVubmC4EuA8IaMBS+upSwifYSLpMQCOECJ9gfQ30cgC7gYAGLKWvLtX70g1GZjGdbzAM5AB2BAENWEpfU6b3c7GA1nsHERDQwH7gICFg2I5vf9q282vL6W/F8X4saVFrNM4vqdP0PsH6qptMVwHcBwQ0YNjGjz4u6D8fhfUxmx5h2H5SHKcnTjFSVefofIL1maVMVwHcBwQ0YIGYZBQ9lDohwNAiyN9ZGzyAi1wpoA1iCSJJY2sz18uX6VqAO4CABmzUV1OSK45muoqOMeq1P/zww2Jv7rGP3qsT+ZnNHTBgwLRp0xgpDLguCGjARvEaWZ4oiukqOkZVX7vzelNsAi+3su6H5v85v8bYXBv13c8Q0KCjIKABG/XRynLEMQjVMF1IRxgNxD1P3/Zv6K2t0EYs+Z9Zsuvk/qcZKgu4MAhowEbxmrIDvve4WEAjhBC6LY6Z0XzWfKqura6xZUjKJNpFuBy0ecO6UaNGObw44GogoAHrcEiyj7Y8XxSF0CWma+mwHHF0X00pF5FGxPl3akO5Sq25mvwS7SLi1M25ubkQ0MBSxwJapVLl5+dHR0f7+9NclbytrU2r1Zoe0rYB3VBZWdnEGXNaGhto5za2mN/9L4JQqLnCRr6P40uzvxaeVyPPu6euukQYQZ3OEXqQfSfQLsK7+otTSgOupwMBffz48fnz58fGxubn57///vtPPfWUWYMXX3zx22+/FYvFCCEul6tQKOxZKXBZp0+frhJJ25b8QDuX+94YsynxmrI8VxvCQZXtEdtPU2IW0AB0gq2nehuNxqVLl27evDkjI+PkyZMvvPBCY2OjZbM33nijoaGhoaEB0hlQcT28UXAM/T8Ox6xxgrY0T+hiQziocsSx/dUlTFcB3IGtAX3lyhWFQvHggw8ihO66665+/fodOkRz90yDwVBVVWU0suki68DVJKoLszzp71riEnLEMX01ENDADmwNaJlMFhUVxf/nXvSxsbFlZTQXhdm4ceOwYcN8fX3Xr19Pux6CIG7evHnyHwUFBZ2rG7ixRHVRltiFAzpbHNNfU8R0FcAd2NoHrVKpRCKR6aGHh0drq/mxnTfeeOOLL77gcrk3btwYP378gAED7r//frM2Wq12//79586dww/vv//+RYsWdaJuy60DKla9PhqNhrT5ss7eRlUg0Vwi7PH34/aWY+RS0XfeaIUgxNOoDTC0NPBsOs5JIqTRaJRKZZdrsw+VSmX736sbsteHy9PTk8fjtd/G1oAODQ2ldjrX19cnJiaatenR4+8PVVJS0uzZs9PS0iwDWiKRrFu3bs6cOTZutx3e3t5dX4kbY8/rIxaLORYdzdYkthXdFscYOP/8tmtvOVvXaVd33ijJ4eSKo/upS855Jdm4RrFYzJ6/F4fD8fLyYroKVnPaH8vWLo6BAweWl5dXV1cjhAwGw8WLF4cOHdpO+6qqKj8/88sRAHBHAzRFtzxcuH8DyxbH9NVCNzToKlsDOiIiYs6cOU888cTp06eXLl0qlUrHjBmDEPr111/HjRuH2yxbtmzfvn1Hjx594YUXTp8+/eijjzqqauC+BqiLs1w/oHPEMf1gIAfosg7cUWX79u1DhgzZsGGDSCT6448/8I/WqKio8ePH4wZSqfTHH3/87LPPjEbj1atXe/fu7ZCSgVtL1BTecuUjhFiOR2x/TTHTVQCX14ETVSQSydtvv202cfjw4cOHD8f/X7lypd3qAt2Sp1HTQy8vEEmZLqSr8kRR0bpqIanXcQRM1wJcGNyTELBIf01xnqgnwXH5S8ToOIIiUQ/YiQZdBAENWCRRXeQGHdDYNY/4wW35TFcBXJvL76oAllj/0ZZDB2lOLkUIKWqr9QHmgzJpDVAXZXom2LUuxlzzjB+jvI4Cma4DuDIIaGAf76xdrXn8GyT0oJl34lO+wUAz3cIATdGuwOl2rowh1zz6PFf3M9NVANcGAQ3sp08K8qA7d+7KPqRpvuPSYlLXU1eTJ+5p/8KYUCySBhiU/kTLHa+bSuo1v+w7kJNHf9mDkOCgV1952QEFAhcAAQ3YYlBb/m1xtNsMezAizg1x70Ga/FNe7Z3ShRDSybKOB8cd96TLcYNesGUtBHS3BQEN2CK5LSfTsy/TVdjTdc/4wW0FdwxohBBKnIJGL6SZrmtDf35g98KAq4BRHIAthrTlXnGvgL7m2WewOo/pKoALg4AGrMAhycFteZme8UwXYk/XPPsMasvnwJXhQGdBFwfoANrb6GBdTKFeuopWnmcdP6BLa2EZBc+vlSeO1ld38aoc7bzsPj4+d7xkJXBdENDAVn/8cfT+WbN5QjHtXK1W05WVJ7fdvixxq/4NDJ+u0vmANugJvji8J/3JOwadZtXK195cu6bTqwcsBwENbHXu/Hni3teIaVauuLKkS1fIHdKWm+nhJqeoUF31TBjSdns/6uwvAwNBqpXaL1X0c//cqFK1dLo2wH7QBw1YYWhbzhXPfkxXYX8XJQNGtGUxXQVwVRDQgHn+ZFuIvilf5MJ38rYmWxwbRDQHcdRMFwJcEgQ0YN5QY/lVz/h/b3PlRoyIc9mz7whuHdOFAJfkhh8J4HKGGCrcbIAd1QVJ4gheDdNVAJcEAQ2YN4YoTpcMZLoKR7kgSRzBgz1o0BkQ0IBh/iJerFFx1V2uMmrptjgmkKMJ9YDRyqDDIKABw0ZFeF/gxehd/y4q1hgR55IhZFgw3YVYAWgXBDRgWEqE91leDNNVONYFQ+iIEAho0GFuu9sCOufYsWPXr1+nnXUhIwP5pNh9i6MjfLby3fwG8OnGsEdDIaBBh0FAg3/V1NTc98Bsw7iltHPJ6/koxc4B3VtbgRAq4rr5jaFyjb4+Ao5UX1chCGG6FuBKIKDB/xB4eutnraedxSm6YPfNpbReO1vZgkLtvmJ2IREnraptcszFrwPvY7oW4EogoAGTxrReP1ClRHcxXYfjnahsW9Byye4B3dzcXFxcTDtLIBBERkbad3PAySCgAWP4JDFclb2ippXpQpzhbI36Y3WBj6G1hedlt5VW3d715+GfDh+nnamuk1VWlIeFhdltc8DpIKABY0aosvLFkfWa80wX4gwagzFD0v+e1muHfcfYbaVKOTnmSdXsd2lnSl7vpdfr7bYtwAQYZgcYc1/zud99RzFdhfOc9Bk2WXmR6SqAK4GABszgI8OklktHvbtRQB/3HjZWeZVPEkwXAlwGBDRgxqjWG6Wi8AphMNOFOI+C718ijBimymG6EOAyIKABM6Z3s/4N7IRP8sTWy0xXAVwGHCTsdmpqaz/86GPaWSqViiCc8QOcTxKTWy5+HDLPCdtilaM+I38sfWN96EK3vPg1sDsI6G5nweLlx2sFSDqAZl6zmuOU4/5jVNeLRNIqQZATtsUqhaLIWn7gmNZrf3kPYboW4AIgoLulwfejIbNopldkofO7nbD96c3nu2H/BvaL/8QHm05CQANbQEADZxOS+oktlzaFzGe6EGYc9B3zWu23fgZlE69L90G/IyOh37p1q6+vL+3c5OTk8ePHO7QA0HUQ0MDZZjSfveXZu1rg5hdIsqaF53XK664Hmk5/EzjdoRvSaDQfZDRyPGgOKpC1BfecyYCAZj8IaOBsj9cf+ThkLtNVMOln/4kra791dEAjo8EwZQUKkNLMunaILP3esVsH9gCHkoFTJbfl+BpaT3kPZboQJp33GhRANPdTlzBdCGA72IMGTrWw/sg3gdONiMN0IUwyIs4v/pMebfxjpccypmqQ11Rt376ddhZJko888oiXl/0u6gQ6CwIaOE+ovmFM6/VXezCWSuyxK2DamYJnPgt+sFLAxLmU5Tduy2pe/Jn+lBnDjT/i4uKgh5oNIKDd0969e+vq6mhnlZSWoAgnl/O3BQ1/7Pcbq+RKmNk8mzTyfX4OmPSU4rd14U8xsHlCywmLV839gnamr+I/JEk6uSJACwLaDZ0+fXrpijX6ZPoDcaSsEo10ckUIIeRlVM9tPD4nhv7amN3QtoAH0gqXfRY8p4npSgBrQUC7J1FYjMbKnav4134zOrkahBBCy+W/pnoPLRL1YGLjbCQX+B/wG7u4/uAHTFcCWAtGcQBnkArJuQ3HNnbXk1Os2Ro0++HGEwEiHtOFAJaCPWjgDCui1V8HPlQrCGC6EHapFgTu9Ru/MrHmBaYroSI06pVvbgj8bBvt3LuHDF67eqWTS+q2IKBdlUqlqq2tpZ1VVVXFqmM8d3nqhnobXgqayXQhbLQxZH5q4L4xZN1Zpisx0dRXXo5fgEL708xrqsza/RUEtNNAQLuqeU89m3ryBFcgtJylV7XoA6KcXxItPkm8Hdm4USZW3yViuhY2UnNFa240vnv3tYmkTsOh+WsyI2406ks3zK46F2V+5fRqui8IaPa6cePG3aNTNK0ttHM5nr7kK2lImkgzL/UzXgZbTuR9veabWj3vN7mA6ULY668a7XXC/8W6HzeEPs50LYBdIKDZ6+LFi9xhDyIrg1U5L0awqRuD3tSW9P8oM+4tCyCRkulaWG2dJumPpr9yRdH7/cYyXQtgEQho4CjRuur1lVsXRK9pvvgxvM/apzCKHo7ZsKdkFclBB3zHMl2OdYS2WdX2nwcepJ3J5aB1r700fPhwJxflxuCDw7CmpqaGhgbaWQqFwsnF2FEvbeWusrc+Cp170yOO6VpcQ7EwYn70Wz+VvKFHgiOsvZtBY6VK2Xo8jO5uDwgJ07++LzMTAtqOIKAZNmT0+LqGJkR38SBNUx39fU9Yb1zrlY8qPn0v9NGf/ScxXYsryRdFzY9+88fSNRxkPOw7hulyrBB6oqGzaefwik47uRa3BwHNMFlBDvGJAvHpRjh89bjL/Xk4JLmkft/j9UeeiHr9qmc80+W4ntvi6Ed6vvm9bK2WIzzu42K7ogZF+bJlXy5bRn8xLB//wOYGF/5RyAiXSwDAXh5G3YcVG6W6uvtiP4JzUjotxyNmQdSab8ve3jjuLQAADFxJREFUCiUafhLc5Yy7rNuJobEKPfIpGvsMzTyjoXUJXL+0wzpwqjdBEPv379+8eXNmZqa1NpcuXdq8efNvv/1mMBjsUR5wGb39PA9UrtNwhHNi34V07qIsj94PxbwzRZlxqvL1B/qEchH7B+wAh+jAHvTMmTMVCkVKSsq0adM++OCDxx57zKzB9u3b161b99hjj/3www8//PDDr7/+atdSWa20tLS+vp52lsFguHLlCp9P/1Ibja792ROTumlN5x9qOhE3c/Bm33G7Q+l7J0FH5YuiHun55t2K9BWJ+csKnt0Y+ugJ72Ekx2VvdECSpIePxNfqN/eq1159feWrzqzIJdga0BcuXLh48WJJSYlEIklJSVm+fPn8+fO53H93wAmCePPNN3fv3j1p0iSlUtmzZ89r164NHjzYMWUzQKvVVlZW0s5qaGi4d/Zcg9iHdq6ystgYGifsmUQ712h01Z8aA9UFDzeenN5y7qpH/NeBM05u+ozY9AvTRbmbDHH8rH2ZEzZ++HLd9y/XfX/Fs2+OOPa2ODpP3LOV68F0dR1BGsm25rbN1fRz/9z00Rc7vvmJfpfO11ty4shv/v7+DiyPrWwN6GPHjk2YMEEikSCEJk+eXFNTk5ub269fP1ODmzdvtrS04LsweHt7jxs37s8//2RhQLe2tlq7kj1CqKyszFrnzJavvj1+/ARf5Gk5y6DTaPR68sMc2gW5740xDp6p+c+L9Js8vevORTPKx9AapauL1NdGxggiA0oiy96K0tdEamurBEF7/SdM7r2lhh+AEEIu/lOAzVK9k9O8hg5V3x6gLhqoLpjbeLy3tlzO98tJRLncupyWjBxxdLkg1AX2rz2thGxtQVPYoIYpr9DO9Ng68/PPPw8JCaGdGxkZGR9Pfziaw+FER0dz2P+yWGdrQFdWVkZE/H0fDoFAEBwcXFVVRQ3oqqqq0NBQHu/vCydGRERUVVVZrker1e7fvz8vLw8/HDp06NixYy2btf21H7Xb73bh4uUKuvUjhDRqdXlllbW/iZYwEKpma6vl8AUkoaedJRWIF40bhbzpblDU1ohk19H55+jXKdWTmr/Q+VL6TQ6JRhkvIQ7dwQCvMg6HQ1pbbWIImfcxKjM/8CJApKdWhnohdGEp/YL9vcn6PejCSfyQxyG9OAaEkIBDehIqNKGX6PLjYo4BIeTPIyL5Gh4iZYRHhUEkQ6qSquozugAZEVJORLWRPITKEDJdmYFEJz9FArrhKOU3jIQG/bmRth7SQKAL36PcUzTzSq6QSoW1BY3qFnTzd9RE97NG1YQMemsLogYZKuRYnYsQSv0Cielu+1KWaWxrsFoPoUcZP6HCdJp55TdJVaPVBduaUdafqJVuhIOuDZFG9OdGEqHLCP1zi6pwHgqL5mv6Fd3sG6Z/KP/LvgKVL1ffYBBqSI4W8TQkR2Pkav7TV5u/UssVaUiuluRoSMpFTfklnIhW0to7JCmIrNyJmvYihHQkV01S3pxKOSfey+p7sicilcfQebo9FdKIhkZb+4wgv2qOQExWfEO/2viA6oNfW9n3RjcFoiN6Lf08vlAs9jCFkllBYg+JtEc47XJ6vV6v13t60uyNIYSEfP79902zUs7fPO+Zhe70xSAUCu/45WFrQHM4HOpdcGjviGPWgHbbJEm2trY2Nf19Ewm1Wk27OVKjQu1ekI1UNQd70l9ZRsMxBsRGSqzc8rKxocHPrx+HSxOIhF6vamvz9fWlXbClucXTk88XtNIU48tv9OgV4E8zCyHUJuzJ5XLFHvRzG/r2CQhoo52l9YjQE3ovL/oFG+NifH2NXJ75XCPiqL0Dmpo4Qb703zRNUZEenmKR8O+5JOJoERcvqCP9qtrCQr28lYiDEJIh3jkkUHP4SIAQQvIQPofLTQziJiIdQjqz1eZOnpgQp0KI5rm0BAxSKpU9ejTS1lM0YUxkFE8opJmrDZdWRI3oFU2/YOWkMd4+Pj4+NHNJkszjTU6wsmCtxyiSNIaF0c/NmzSuT7yGwzF/ggghZWBiU1NjZCT9giUTR4dHCMVimrn6HmFl4aN7W6mnatJoiZfE15d+7m3Df/paWbDeY1ClXk94hWUhJCINYoGBTxp5iBQgIx+RtWJlpLeXFxf5kiQfkXz07/tB5xHS1ibxs/IOqZdG+PoK+Xw9QoiHSD769wYPRl9Rs1eMv5W3uiq+J5/PF4mtvdXjrX1GNOJIg4GQWH2r9/Tz8+vUZ7bZ09OTL6C5DgxpJBubGgOsZIiqVc8TisUe9HOFAgGpVtHOsjtbAzo8PLywsBD/nyCI+vr68PBwswZyudxoNOKO6draWtr+DbFY/NBDD82ZM6f9zYnuf7L9BikTHvb29rax+G5IqVTC62ONTqdDCAmFrLl0HMu0trbCLb3b4cwPl63D7CZNmpSamop3eNPS0oKCghISEhBCdXV1crkcIZSYmCiRSM6cOYMQUqlUp06dmjTJgWeRHThwwHErd3VyufzsWfZcXph1bt++nZNDf8AAIITOnTvXznEacPDgQaPRSbeNszWgR48enZSUNGXKlDVr1ixcuHDNmjW4Z2fFihVvvPEGQkgoFK5evfrRRx9ds2bN5MmTR44cmZyc7Li6FyxY4LiVu7qcnJx334V7s1q1f//+vXv3Ml0Fe73//vtZWVlMV8FeixYtctp5Hh0YB33kyJFff/21srJy7969I0aMwBOXL19u6oNftmxZUlLShQsXnnvuudmzYTwsAAB0SQcCWiAQzJs3z2zi0KFDqQ9Hjx49evRoO9QFAADdnrOvxUGSpEqlamykPyrdIXZZiVtSKpUEQcDrY41arTYajfD6WEMQhFKphNfHGpIkGxsbBXSDQzrEx8fHyhDAf3FoB8w5zldfffXSSy9x6UbMdIhWqxWJ4B539IxGo8Fg6PobyF3hDsQ7fja6Lb1ez+Pxuv4hdVf2Cp/Tp08nJtLdso7C2QENAADARvAlCQAALAUBDQAALAUBDQAALAUBDQAALOV6t7xqbGxcv359VlbWwIEDX3/9dT8/P6YrYpHGxsYTJ05kZmYaDIZNmzYxXQ7rXLt2bc+ePdnZ2X5+fvPmzZs6dSrTFbHL119/nZaWJpfLe/bsuWTJEhZeLpgNiouL33///YceeghfXdmhXG8P+uGHHy4vL1+5cqVMJps7dy7T5bDLpUuXdu7cWVZW9u233zJdCxt99dVXYrF4yZIlY8aMmTt37v79+5muiF3q6+unT5/+yiuvREREpKSkFBUVMV0R65AkuXjx4sOHDzvnbHgXG2aXnZ2dnJwsl8slEolKpQoODr569Sq+bBMwOXPmzJw5c+B6N+175ZVXqqqqfvjhB6YLYanBgwe/9NJL8+fPZ7oQdtm6deutW7cKCgqmTZv2/PPPO3pzLrYHnZmZOWjQIHxjF4lEkpSU1M4dbAFoR1FRkVQqZboK1mlra2toaPjjjz+qq6tHjRrFdDnsUllZuWXLlvXr1zttiy7WB11bW0u9NVlgYGB1tbU7LQBg1b59+86ePfvll18yXQjrrF27dvv27SqV6sMPP4yJiWG6HHZ58sknN2zY4MzjXi62By2RSDQajelhW1sbXJYedFRaWtqSJUsOHz5s7TZ33dnGjRubm5tv3bq1ceNGuOo61e7du729vR944AFnbtTF9qCjoqJKSkpMD0tKSqKiohisB7ics2fPzp07d+/evXfffTfTtbBXQkLC1KlTT58+PXPmTKZrYYtTp0798ccfAQEBCKHW1tbz589fv3591y7H3vTZxfagJ06c2NLSkpqaihBKTU1tbW11wkgX4DYuXLgwZ86cn376KSUlhelaWEer1ZaXl+P/y+XyU6dO3fFSPt3KN99809ra2tDQ0NDQMHbs2A0bNjg6nZHL7UGLxeIvvvjiwQcf7NOnT0FBwZdffgnXtKMqLCwcNmwYQRCtra0BAQH9+vU7d+4c00WxyJo1axobG023xBw1atThw4eZLYk92traBg8eHBQUJBaLi4uL586d+/jjjzNdVHfnYsPsMKVSWVhY2Lt3b+iANmMwGFpaWkwPeTyej48Pg/WwTWtrq17/762s+Xw+vIWoDAZDUVGRVquNioqydqtsgBBqbW0VCARO2Dt0yYAGAIDuwMX6oAEAoPuAgAYAAJaCgAYAAJaCgAYAAJaCgAYAAJaCgAYAAJaCgAYAAJaCgAYAAJaCgAYAAJaCgAYAAJaCgAYAAJb6f/0GnIw7sSOrAAAAAElFTkSuQmCC", "image/svg+xml": "\n\n\n \n \n \n\n\n\n \n \n \n\n\n\n \n \n \n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "histogram(μi2,bins=0:0.1:4,normalize=true,legend=false,title=\"Histogram of 10000 averages with ρ=0.75\")\n", "plot!(μi2->pdf(Normal(μ, s2), μi2))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Task 4\n", "\n", "You decide to test the hypothesis that $\\mu=2$. Your decision rule is \n", "\n", "- reject the hypothesis if $|(\\mu_i-2)/s|>1.645$ with $s=\\sigma_y/\\sqrt{T}$\n", "\n", "With this decision rule, you are clearly assuming that the theoretical result (definition of $s$) is correct.\n", "\n", "Estimate both $\\mu_i$ and $\\sigma_y$ from each sample.\n", "\n", "In what fraction of the $M$ simulation do you reject your hypothesis when $\\rho=0$ and when $\\rho=0.75$? For the other parameters, use `(T,σ,μ) = (500,3,2)` (same as before)." ] }, { "cell_type": "code", "execution_count": 99, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Frequency of rejections:\n", " with ρ=0 with ρ=0.75 \n", " 0.098 0.536\n", "\n" ] } ], "source": [ "(T,σ,μ) = (500,3,2)\n", "\n", "rejection1 = 0\n", "rejection2 = 0\n", "\n", "# Count how many times we reject the hypothesis\n", "for i = 1:M\n", " # rejections for ρ = 0\n", " s1 = σi1[i]/sqrt(T)\n", " if abs((μi1[i] - 2)/s1) > 1.645\n", " rejection1 += 1\n", " end\n", "\n", " # rejections for ρ = 0.75\n", " s2 = σi2[i]/sqrt(T)\n", " if abs((μi2[i] - 2)/s2) > 1.645\n", " rejection2 += 1\n", " end\n", "end\n", "\n", "println(\"Frequency of rejections:\")\n", "printmat([\"with ρ=0 \", rejection1 / (M)], [\"with ρ=0.75 \", rejection2 / (M)])" ] } ], "metadata": { "@webio": { "lastCommId": null, "lastKernelId": null }, "kernelspec": { "display_name": "Julia 1.7.1", "language": "julia", "name": "julia-1.7" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.7.1" } }, "nbformat": 4, "nbformat_minor": 4 }