{ "cells": [ { "cell_type": "markdown", "id": "c57b2ce5", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Likelihood Function and Loss Function\n", "\n", "Feng Li\n", "\n", "School of Statistics and Mathematics\n", "\n", "Central University of Finance and Economics\n", "\n", "[feng.li@cufe.edu.cn](mailto:feng.li@cufe.edu.cn)\n", "\n", "[https://feng.li/statcomp](https://feng.li/statcomp)\n" ] }, { "cell_type": "markdown", "id": "055df123", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Likelihood Function\n", "\n", "- Given that $x_i\\sim N(\\mu,\\sigma)$ for $i=1,...,n$, the **likelihood\n", " function** is\n", "\n", " $$\\prod_{i=1}^n f(x_i,\\mu,\\sigma)$$\n", "\n", "- However the **log likelihood function** is more often used\n", " $$\\sum_{i=1}^n \\log f(x_i,\\mu,\\sigma)$$\n", "\n", "- Do you know why?" ] }, { "cell_type": "code", "execution_count": 6, "id": "42b9aabc", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "text/html": [ "0" ], "text/latex": [ "0" ], "text/markdown": [ "0" ], "text/plain": [ "[1] 0" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "-1926.43418681411" ], "text/latex": [ "-1926.43418681411" ], "text/markdown": [ "-1926.43418681411" ], "text/plain": [ "[1] -1926.434" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "NormLike <- function(mu, sigma, data)\n", " {\n", " out = prod(dnorm(x = data, mean = mu, sd = sigma))\n", " return(out)\n", " }\n", "\n", "logNormLike <- function(mu, sigma, data)\n", " {\n", " out = sum(dnorm(\n", " x = data, mean = mu, sd = sigma,\n", " log = TRUE))\n", " return(out)\n", " }\n", "\n", "set.seed(123)\n", "data = rnorm(1000, mean=2, sd=1)\n", "NormLike(mu = 1, sigma = 1, data=data) \n", "logNormLike(mu = 1, sigma = 1, data=data) " ] }, { "cell_type": "markdown", "id": "bd3ee8fd", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "\n", "### Walking APP example: How long do you walk every day?\n", "\n", "- Here is a list about my past six days walking statistics. Can you\n", " estimate how long do I walk everyday? and what is the variation?\n", "\n", "\n", "![walking](figures/walking.png)" ] }, { "cell_type": "code", "execution_count": 7, "id": "68acebc3", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAFoCAMAAAC46dgSAAAC/VBMVEUAAAABAQECAgIDAwME\nBAQFBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8QEBARERESEhITExMUFBQVFRUW\nFhYXFxcYGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIjIyMkJCQlJSUmJiYnJyco\nKCgpKSkqKiorKyssLCwtLS0uLi4vLy8wMDAxMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6\nOjo7Ozs8PDw9PT0+Pj4/Pz9AQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhJSUlKSkpLS0tM\nTExNTU1OTk5PT09QUFBRUVFSUlJTU1NUVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1e\nXl5fX19gYGBhYWFiYmJjY2NkZGRlZWVmZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29w\ncHBxcXFycnJzc3N0dHR1dXV2dnZ3d3d4eHh5eXl6enp8fHx9fX1+fn5/f3+AgICBgYGCgoKD\ng4OEhISFhYWGhoaHh4eIiIiJiYmKioqLi4uMjIyNjY2Ojo6Pj4+QkJCRkZGSkpKTk5OUlJSV\nlZWWlpaXl5eYmJiZmZmampqbm5ucnJydnZ2enp6fn5+goKChoaGioqKjo6OkpKSlpaWmpqan\np6eoqKipqamqqqqrq6usrKytra2urq6vr6+wsLCxsbGysrKzs7O0tLS1tbW2tra3t7e4uLi5\nubm6urq7u7u8vLy9vb2+vr6/v7/AwMDBwcHCwsLDw8PExMTFxcXGxsbHx8fIyMjJycnKysrL\ny8vMzMzNzc3Ozs7Pz8/Q0NDR0dHS0tLT09PU1NTV1dXW1tbX19fY2NjZ2dna2trb29vc3Nzd\n3d3e3t7f39/g4ODh4eHi4uLj4+Pk5OTl5eXm5ubn5+fo6Ojp6enq6urr6+vs7Ozt7e3u7u7v\n7+/w8PDx8fHy8vLz8/P09PT19fX29vb39/f4+Pj5+fn6+vr7+/v8/Pz9/f3+/v7////9SYPv\nAAAACXBIWXMAABJ0AAASdAHeZh94AAAag0lEQVR4nO2dd3wURfvAn7tQ0mmhhNCJBAOIUiQC\n0ozSAlhQAVGw0VVEfaUpKvqiwk/s/iy82EV9FQuigD0WFKUqXaVIJ/QSQjKfd7bc7cze5rLZ\nnbs9h+f7x+3ezszzzO73btvdzQFBpAa87gASWVCw5KBgyUHBkoOCJQcFSw4KlhwULDkoWHJQ\nsOSgYMlBwZKDgiUHBUsOCpYcFCw5KFhyULDkoGDJQcGSg4IlBwVLDgqWHBQsOShYclCw5KBg\nyUHBkoOCJQcFSw4KlhwULDkoWHJQsOSgYMlBwZKDgiUHBUsOCpacsILnAsCP2mwmQGdCngNI\njkavrCi6r0nFpFdLK30MkveXFaEgPz9/tcgunT4tMprt4OVK615wcV5e3jflyOiUx2ln4MVS\nCndXgX8p06VXZic2v2yhtvCnIS0Smt7+d7DS1TRAN4vGXKOip3KqJbeZVRxaEsIzAActMnER\nGOaDwUCL7tkKzpXwGHnZTO4Fn6aV3goXRRC5ANXH5pdSOBzid9HJPX5tpa5QXuNT4tT5Wkv1\nOl+AtWCu0YEO2pM+B8wlIWzICGxpLhMXgcUk2Nw9W8G5Eh4mr3PBJadDdw/REnw2wK2llW30\nwwg6+cYHkNS2Ou3QfYS8o2yVNgkAdfeqdYqyrQVzjUh/Oj2rKX0Yai7hObz6vhqgb2k+ExuB\ng93sw0O6Zy84W2KCyetcsAUb8mmlB348Ei5OGYTszCxpDjCltLJxWj8vBmixi5y4DKDiSZIO\nMKKI7G4PcI9aZyak9LISzDVaBRC3WDmiQ+WTfAnPy/rWU7c0l4mLwKG8OxQWxUH1bSHdsxWc\nK+Fh87KZnO2iS97r2SS+Se4rRYQM1FIq1U7N6t849YKbV+mtj93dKjn3d9okhz6ZBpBJ3ju3\nGa01p1vDyvU7Pq68JMbTqCs7+iq1nENOTu+a2vT6HWx6Np6eZppRWvCvHjXqdPu/U3T2UDJk\nKYvqALxOJ+toxWV/AFRQUnwOUFupsyMFHrvJSjDbiIwFuF1ZuG7dulN8SeiWCWxpPhMXwYpd\nNQE+JObu2Qs+t3TBFnnVTI4ElwzQE7U5wghedY42W/HeEqXF9iz10DHWEPyaDxqRkzl621ZH\nVMH1q6rPZnZXJxlM17l4IYI/r6stydpNyH9BPcU6Qp/+QqfH6fS/C+lOTan3N32ynk6HQOvT\nVoK5RiQNYIl1Cc+OJUuWTNe3NJ+JjWBJHsB4cyO7wdkSExZ51UxlCjYwBCsntM0uz/Epx5Lg\nMfh4EzrToH1l+jhHaawcE2olKQ11wbXp0aMRmUgXNO9aG9RdznilOClRi5+mPNwfTM7HW5vf\nAOC6/K2B0v300FixzXm08BL11fs+XXZ6xYoVJ4i6t4Klv9Fd1mGinVp9ScjX4PuBWAnmGu2l\nDysfOC+57V1HTSUWW+djfUtzmbgIVrwF0LjQ1MhucK6ExyKvlsmR4N4A1xL1WfWSoOApAP6X\n6MuvHVVVoJ0czC45dYshGCqNfuEN5VxpqhailyZ4avGxcYrKdWRVMsBlweSmeKZjMI2btlpd\nCfidtAZgrjaK6Usr83QRfW1dX0h2tqU13iRFreBmYimYa7SaVm6hrmzmFr4kjAMuk1UElqP0\nmPqauZHd4FwJT2hePZMjwTRbxrNbyfGvvvqqKCiYChiipqJv7PfUi85OyuZpbgimS0nJW2+9\nRU8HD3cBaKcKrlmsHmjgcaLsRpkzOVM8k+CGAJOU6SWZmW+SWgCGgZNX0ff2IkIeokGrtlB2\nALCAnnuk7Q8nWG+0VKld/WKaCvqbwpXugMtkEYHj3/TYVGxuZDc4X8IRmlfPZPcY3JwVPFUV\n3mz0u8dIcBddWAHgXbVqM4AHCaEXJQ8pz6YEBadokYry77+ytdJhTXBbumif3vtRjGBzPF7w\nCWo9eP+hOA5SgyXr2wAkKA0Le6q9pO9u+HVnKvyHhBEcaPSDsr8pUDeOfmQMhgsluKXZTKER\nOEroC3M+MTeyG5wv4QjJG8jkSHDh5Oramzr1haBg5V34g1q1B8CNpKQiwFzl2fNBwU3V0uXU\nvC/z6t4BwfRRFfwp4QWb4pkErwftBEilCKBaYP55ekDPXKltzCe7VUkfodTc9QTA+Tk5OfSs\nMjXnudD1DDZaSSt/QafF9AT6HVO4EIwtzWQKicBDD6XVi4i5kd3gphKWkLyBTI4E02369UTl\nNQW+Vew7WDvXzFJPdmmimcqzB5jLJMpJeuo0iPZ0UlmCzfF4wYdpC2O/WQXi1BN3cmwoXX7t\nYXYdlgEkFs9mjjQhF9NMowI6px7C2ql3RS3CMZi3tJrJFMHMdfSIyi1QG9kNblmiEZI3kMmJ\n4CNr1qyhi3ZQd/B08Bicpd+7+c2vvoa6que3hOTwgpVjxQainsKHF2yOZzoG19XO1Ui/li3f\nJWcBqPcFi+nlW4XntQpF06ZN20ynkwH6kHCC2UZK0g+IemFNT5y5klCC50FsJj6CmUJ6yvRR\naPdsB2dLTJjyBjM5Ebxef/vspkfST1TB9GxXeU/6XyZkZ3t6UrBHOxjMJcX3Ai94sbqYfOQr\nU7ApnknwcLpbpsekN2jLPwg9D/paWUj3xBDcA9Mr97ydJR8nAAQ/gtKPwfl0oxmbiGtEXwlN\nV5KCS2nOY3wJ34jb0lwmLoK50Y+0yVGL7tkNzpXwjbi8TCZHu2iaNK7LkL6p9EKX7r+qAjQc\n9wc51ojWbtqR9gTokZkco/toqFcFTIJ3KffvW7eifqFleMGmeCbBO1Loye35dH8EVxLy/wAz\nlIVZxvv0U3K/cgRJUfIE7xTpgmfQhX8FI3GNDtFrC6ju13JyJXwjbktzmbgI5kb0JZOuz3KN\n7AbnSvhGXF4mkyPBa2vo6x2vvHGuAq3acu06DCpOUQ+IX6onYpUu5wWT29U6TQbTvd/qsILN\n8Uz3oj9M00rb05fYnz64iC7axeyIP1WvXxUabgi2sRLMNyI/1dJmx5aYSsI44DMxEUIajTZW\nkGtkOzgpVTCXl8nk7CTr4OwuTRLS2oxXby/sHZae0Pw3OlP4SF6jlA43Bk46N1+fVaNv/rMm\nwcVPnJPUZsKhD0C5jxZWsCme+cOG3eMvrJbe4wX1zKM3+LdpuyXGFZl3UXpi9pQCo4mVYHOj\n3Xc3T2g6cFFISRgHpkxGhJBGFwEMCz5hGtkPXrpgNi+TKeJf2Zmq7kQjzGf6RbddHoVD5U8S\n0UYRCx4xwWMzMzsc1z6FfTBSOQzaQn3rz26sGVDDQY6INopY8IgJVk5BL/9sUS5AlT8jlcPg\ne1+p3+YJZdutDl5zEW0UweARE1x8pX74Sno/UilYhlt/IGDJjDqPlJQ7QUQbRTB4BI/BS/qf\nnVjngjt2Ry4Dw/758y2/iWbFofJvzAg3imBw/F605KBgyUHBkoOCJQcFSw4KlhwULDkoWHJQ\nsOSgYMlBwZKDgiUHBUsOCpYcFCw5KFhyULDkoGDJQcGSg4IlBwVLDgqWHBQsOShYclCw5KBg\nyUHBkoOCJQcFSw4KlhwULDkoWHJQsOSgYMlBwZLjXnCBvT9NQbzBueATT17/4Abyfl1IHlDa\nv3ch3uNY8AHlf6Zq/1I5tXtLqF1Qdn3EGxwLvhMmrFqcmdSAvnvfhDtEdgkRiWPB2coQowu0\nISK7nWsuPVxgwI6neKJAICec9v1MwrHghFH0YTu8rcyPTjQVbvIxI7X6ioyCDBBIXad9P5Nw\nLLiJMkLz8VErlPkr0syla5YFmQOFxvKku+cJY6L5ZYVY4Fjw1RU/DMxuSrD82wGd7zjBT68W\nxjMo2AaOBW9O9LVV//Vh9a1VfF+GqYiCPcX5dfDGy2s/pUyfg9pvh6uHgj3F1Z0s9R7Wpu/C\nD8WNgj0l8veiUbCnoGDJQcGSg4IlBwVLDgqWHBQsOShYclCw5KBgyUHBkoOCJQcFSw4reK6D\nv7AtGxTsKaxgiL/87ePCM6BgT2EFP93VD8lDPy7PPy3bAAV7Cn8M3vkUdVz95i9E/hoFBXtK\nyEnWzqe6+CH9th+FZUDBnhJ6Fr1iWmPlS8fN3hWUAQV7Ci+46IvbGgKkj1r0y4Rk389iMqBg\nT2EFv3ttNYCmd36v/nX4rzBRTAYU7CncZRK0nrYy8ORQ2qNiMqBgT2EFz9wciQwo2FP4Y/CG\nxfThuXVCM6BgT+EE3+brTB8r+CaUCMyAgj2FFTwHOi6gk4Xd4SWBGVCwp7CCu5+l3aUsym4n\nMAMK9hRWcNWR+syYFIEZULCnsIKb99Zn+jYTmAEFeworeETcfHW6MG64wAwo2FNYwfsbQe70\nF2f089XaKTADCvYU7jJpy7V+5XOGvmtFZkDBnmL6NGlP/utLttltu3udPn7O3u1haqFgT3H+\npbvl5wDUmavO9goXBQV7CqfmnUG5OmU33BTvz+0TD08r8yg4dmHVvAiQnKZRdsNBvk/oHj0z\nXrlxjYJjF1ZNi9R8+w0b91Qe1yf0I1aCd2wO8k7EBCdsFojgrxrGDIyakkq3lKNhyk3qZCp8\nYyF4Ezfk4EmjQKTgCSKHRYTyrPs/CUbNSd/t5WjYOVudHK3fotDiHbw1Cu/gMZUXiuNikfd2\nYglWTddGB+03nATj1HfmAhh0wptj8JgEcbFWDzgTBG9p1Wrepn0qZTc8cSGk5CkzUyGjJgqO\nWbhPk5KCRyQbLQ9MbK7tpedmha2Pgj2FVXOTQblilPy5JEwpCvaUf/DPR1GwHUyCj636QXQG\nFOwpnOC/Lq9ID6f3XBPus4Nyg4I9hRW8oz507A7kUcjYITADCvYUVvBYeIW8RhfMjRsjMAMK\n9hRWcMPuRBVM+p8lMAMK9hRWcNJIXfDoJIEZULCnsII7nK8LbtNWYAYU7Cms4OnwQLEieDpM\nEpgBBXsKK/h0F8i8AMa0hVYi/zQOBXsKdx1cOLsBANSYclhkBhTsKeZblUd+2y84Awr2FLwX\nrXFGCB5qIDADCvYUfowOnZRMgRlQsKewgk+q7FvSKWGBwAwo2FOsjsHHsmoI/BIpCvYUy5Os\nu2CruAwo2FMsBd9WWeBopCjYUywEl3xd5RyBGVCwp7CCkzUqA8wVmAEFeworOE/nuvkiM6Bg\nT8E7WRoo2DEo2FNYwfU4OgvKgII9hRU8KgN8ddvW80GjzpTLBGVAwZ7CCv7Wf8nvdLKuZ8Zf\nAjOgYE9hBfdrrP1r0vEmAwVmQMGewgquPUyfuaGewAwo2FPM34tWyU0XmAEFeworeJDvfXX6\ngb+/wAwo2FNYwX/V8F/10sI5V/kTVpZav/ygYE/hbnSs6KF+oaNluN9zlxsU7CmmO1lr3pn1\nyg8i/7kQBXsM/gBc48wQjD8Alw+nPwB/siqHuTgqA6HFrOC9IsdY3OuuL05/AL7x1sqQ0jKI\nqTQqQxnGruCaIsdYtDEwbDic/wD8U8grvTAag5HGruCk+8SNsXi/yzF1XfwAvFkYwQxn4DE4\nlgZNdvED8GvsfaCIgl0hUjD+AFwQsSoYfwAuiFgVjD8AF0SMCj763Pf4A3AhxKhgknSNu1jW\noGBXiBQ8pqaNgcDLDQp2hUjBRSNbzdt4+KiCu6AcKNgVIgXXqRNXjhHf7YKCXSFS8HADd0E5\nULArhF4mRQQU7ApRgse97Hq1SgEFu0KUYFBHTppTvn/jsAUKdoVYwcMjsMdGwa5AwWJAwY5B\nwa5AwWJAwY5Bwa5AwWKQX3DDQZTGMEjD9ToaoGBXCBPM43odDVCwK0QJXsbjeh0NULAr8F60\nGFCwY1CwK1CwGFCwY1CwK1CwGFCwY1CwK1CwGFCwY1CwK1CwGFCwY1CwK1CwGFCwY1CwK1Cw\nGFCwJYe27yhzWDwU7AoPBa++rg4AxGUMzg9bDQW7wjvB43yQ3qFPn5x6AGG/Lo+CXeGZ4Keh\n56/a3JqrYVaYiijYFZ4J7phVFJgtubCTuXSN8eWQOZzgu+cJ48rK4mLN69xvmTgSBK7lRK8E\npw4z5idXMRVu8jHf7/IVGQUZIgf5OzOo69SQhvN3cPPTwfnuIe/gwwUGh5jlJwoEsltgrH37\nBAYT2bECl0NauTgG916lza0fAo+46wQSOZyfRY8CqN+5/4AujQGGlwjsESIUF9fBywen0UNE\nXPrgr8R1BxGNuztZB7buLN8fPFT0+pTln0dFV4aicC+aI3G2uIuR2QniYi3rJ/QySeRaxvyH\nDRxJH4uL9bGtUa1tInRgoVhaSxSsgYLFEEurzoGCxRBLq86BgsUQS6vOgYLFEEurzoGCxRBL\nq86BgsUQS6vOgYLFEEurzoGCxRBLq86BgsVQbZG4WIuqiYtFRowQGCyW1jLKgv8U+OfTxX+K\ni0UKCgQGi6W1jLJgJNqgYMlBwZKDgiUHBUsOCpYcFCw5KFhyULDkoGDJQcGSg4IlBwVLDgqW\nHBQsOWe44CNzt3ndhQgTVcHPdKrS6Rk3AU5OvjC1yeBN5mDO4w6Hj8UE++ai1PSrBHVs/4Ts\nxOwJBWKCRVPwKMi6rhmMcx7g4IWQfdMlvoTlfDDncd8BTbDrYG9VqjtkQFyNLSJiFTSBbiO6\nQuZBEcGiKXg59CoiRZf4VjuOMAnG0scF/tZcMOdxt1dPVgW7DralQgeq4wUYJqJjk+Fp+jgb\n7hWyllEUPBhW0sdf4DrHEZqnnFQmubCbDeY4bkmPxpNVwa6DTYAflHiPPSsgFukLe+jj33Cp\niGDRFJxWT52k13EcITtPnfSBdWwwx3Ef9X87QxXsOljd+sFZ9x27D96gj6/AQyKCRVHwAdAG\nW+oAh90F2hNfu4gJ5jju8kqTiCrYdbAjcOGKfrXqD9woIBY90+hWcfC9gyvkHhYRLIqCt0J/\nddoHtruKsz4T/sMGcxr3ePa5hZpg18G2QdPkVjf08if+LKJj5KUKytAcrwroGImm4J0wQJ32\ngR0uohy9JyH+KS6Y07hj49cQTbDrYH8ATCwhZLHvPBEd+zf0X3lsRV+YJSJYFAUXx3VRpzlx\nLr4W/kkDyFvHB3MYdwk8RnTBroPtghrqsH+XwG73Hdsff/YpOik8K/GQ+2BRPclKb6JO6mc4\nD3EPtPg6JJizuDOD4xS96DpYcXw7dToKfnHfse9htDq9CX52HyzKl0nr6eMaGOw4wlwYFBi5\nlgnmLO7iUQodoPeofPfBeqWqI0p29R91H+tvfWesXC25DhZVwV/BUHqteDV86zRASVZGcGRO\nJpibuNplkutgn8FYuuecB3kiOtY6Tvnt2kJ/eyFrGc1blcOhx+QucKPj9n9CzV4ae7lgLuJq\ngt0HGw6tRlwM6dtExFqV4us5OtdXZa2IYFEVXPJwx9SOjzpv/3nwsLmdC+Yiri7YfbCZnVOy\nxxWIibXj5uzE7JG7xAQ7wz8ulB8ULDkoWHJQsOSgYMlBwZKDgiUHBUsOCpYcFCw5KFhyULDk\noGDJQcGSg4IlBwVLDgqWHBQsOShYclCw5KBgyUHBkoOCJQcFSw4KlhwULDkoWHJQsOSgYMlB\nwZKDgiUHBUsOCpacWBM84s6wxQLGd974pNsINsI5SHKw1nLXnbEgxgR/W2Uf4YaFPnBLy9Su\njwXLA+M7M7BjSBujNrNLeW6tqk3ZyMxoz0F2MeMsmbtRRjhjKUcgCxOZTfJwu9NWCVwSY4Lb\n30W4YaG31YXcEa3ger04ML4zAzuGtDFqM7uUZ1FlbduzkZnRng0Kumk0hI/M3SgjnLGUI5iF\nicwmOZbyit3NVA5iS/B3oAyHzAwLnQdvE1I8Bj5ViwPjO7MwlZlRm5mlHNdkAWjbnonMtAvl\nSKNLzd0oIxyzlCUkixqZnb2hXRmbxwmxIHj/yLOr9pijzA1tpTwaw0If9XdT5o6n9FQmwfGd\nWZgxpJlRm5mlHJfl5aWo256NzLQLZWStPYTvRlnhjKUcIVmUyNzsEvjJqgfuiAHBWxrF9RqR\nCbfRtU9TR6w3hoVepg/r17aScngKju/MwowhzYzazCw101Ld9mxkpl0Ii+E9vnLZ4YylHOYs\namRu9niF+0vviFNiQPC1yuqd6ujbQH6FV43FyrDQu6CXMns6DbYx4ztboFRmRm1mloZU1bY9\nEzm0ncGpzC585bLDMUtZzFm0yPxs266WK+cK7wXv9V+kTBZ0XkxehqXBxeqw0OQc/xd0fgrA\nWmZ851DUysyozWwIM/q2NyKHtGN4Qt2vst0oMxy7lMGcRY/MzQ6pbrVy7vBecD5MD8w+DIHX\ntz4sNFmaENdv5HnJTWAzM76zGb0yM2ozG8KMvu2NyKZ2LIfSBpgqlx2OXcpgyhKIzM3eAoXm\nZq7xXvDrynWmxiTYqc0EhoWm78Ir6tXss6or7GPGdzYRqMyM2syFMBHY9sHIfDuOx2CRqXLZ\n4bilBqYswcjs7GTYY27mGu8FL4EZgdlZ+jmRMSy0TsMa3PjOHMHK7KjNoSGC8NueRuba8Zzd\noNhUuexwFktDesdFNmbH+kJPGdziveBt0E+ZLKzwHH0z5yuzzLDQLz1Ld2tkKYzjxndmYSoz\nozYzS83o296IzLbj+QammrthIxyzlIXLEozMzV5dq5Q+u8B7waSvbyEhRT1868h6eJ7ww0IP\nhZfpCWjnuMCxT9tFn9p3IFCBrWyM2swuZSqr6NueicyM9sxXHg/5Id0oOxyzlK3MZGEic7Mt\neoffUk6IAcFra8Xljc2G2+ls/RsIPyz0H9X8nYc1qPhyoK4meAm0DCwwjSGtj9rMLmUqq+jb\nno1sjPbMVz47/iQxV7YRzljKVWbGlDYis7OH/DPtbzW7xIBgsmPYWcltXlB2cWMblZiGhd4w\nsE5yl8+DVUMEc5WDozazS0sxwkUOjvbMVd4GXQKzRmU74awFM2NKM5GZ2fkWt2VcEwuCDX63\nOWL9vrblCCq4cuRyD7y4HJXtEluCSS97Q9YvsvxcIDqVI5Z7b+VPyhHZLjEmeG3yHzZqLe1k\np1ZkKkcu9x159uvaJ8YEk5kRuN/+z+BgztZIhI01wYhgULDkoGDJQcGSg4IlBwVLDgqWHBQs\nOShYclCw5KBgyUHBkoOCJQcFSw4KlhwULDkoWHJQsOSgYMlBwZKDgiUHBUsOCpYcFCw5KFhy\nULDkoGDJ+R9mRYU1spn4owAAAABJRU5ErkJggg==", "text/plain": [ "Plot with title “Histogram of c(294, 262, 196, 79, 191, 677)”" ] }, "metadata": { "image/png": { "height": 180, "width": 240 } }, "output_type": "display_data" } ], "source": [ "hist(c(294,262,196,79,191,677))" ] }, { "cell_type": "markdown", "id": "b9bd295e", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### The likelihood function\n", "\n", "- We assume everyday’s walking steps ($x_i$) are independent, and\n", " $x_i$ follows standard normal distribution $\\sim N(\\mu,\\sigma)$, the corresponding likelihood function is\n", "\n", " $$\\prod_{i=1}^n f(x_i,\\mu,\\sigma)$$ which can be easily written in R\n", "\n", "- **The scope** Find a proper combination of $\\mu$ and $\\sigma$ that\n", " maximizes the loglikelihood function." ] }, { "cell_type": "markdown", "id": "a5944b92", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Conditional likelihood function\n", "\n", "- Fix other parameters\n", "\n", "![conddensity](figures/conddensity.png)\n", "\n", "Left: fix variance to allow $\\mu$ to change with likelihood function. Right: fix mean to allow $\\sigma$ to change with likelihood function." ] }, { "cell_type": "code", "execution_count": 8, "id": "f9071724", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "x <- c(294, 262, 196, 79, 191, 677)\n", "mu = 260:300\n", "sigma = 180:220\n", "\n", "parMat <- expand.grid(mu, sigma)\n", "muALL <- parMat[, 1]\n", "sigmaALL <- parMat[, 2]\n", "\n", "myLogLike <- matrix(NA, 1, length(sigma))\n", "for(i in 1:length(sigmaALL))\n", "{\n", " myLogLike[i] <- logNormLike(mu = muALL[i], sigma = sigmaALL[i], data = x)\n", "}" ] }, { "cell_type": "markdown", "id": "b1916648", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "```R\n", "persp(as.vector(mu), as.vector(sigma), \n", " matrix(myLogLike, length(mu),), \n", " theta = 90, phi = 30, expand = 0.5, \n", " col = \"lightblue\", xlab = \"mu\", \n", " ylab = \"sigma\", zlab = \"log likelihood\", \n", " ticktype = \"detailed\")\n", "\n", "filled.contour(as.vector(mu), as.vector(sigma), \n", " matrix(myLogLike, length(mu),), \n", " xlab = \"mu\", ylab = \"sigma\")\n", "```" ] }, { "cell_type": "markdown", "id": "876dcc3b", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "- 2D and 3D loglikelihood function\n", "\n", " ![3dlike](figures/3dlike.png)\n" ] }, { "cell_type": "markdown", "id": "dd211cfe", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "- Are $\\mu$ and $\\sigma$ we obtained the best combination?\n", "\n", " ![2dlike](figures/2dlike.png)\n" ] }, { "cell_type": "markdown", "id": "7f7f59d6", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Likelihood function for linear regression\n", "\n", "- Assume you want to make a regression model\n", " $$y_i = \\beta_0 + \\beta_1 x_i + \\epsilon_i$$\n", " where $\\epsilon_i \\sim N(0, \\sigma^2)$\n", "\n", "- What is the (log) likelihood function?\n", "\n", "- What are the unknown parameters?\n", "\n", "- How do we estimate the parameters?\n", "\n", " - Write down a likelihood function with respect to the unknown\n", " parameters.\n", "\n", " - Use an optimization algorithm to find the estimates of the\n", " unknown parameters.\n" ] }, { "cell_type": "code", "execution_count": 9, "id": "60704413", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAFoCAIAAAAAVb93AAAACXBIWXMAABJ0AAASdAHeZh94\nAAAgAElEQVR4nO3de1yOd+MH8M/dOToppVJCq+jA5lRakeNojMgYG8YzmrPs2Wab/Rg/h3iM\nYQ4PNqw5zWGOj/FgLOSwnEKFhFKOEZ3r/v3R/eskUe7u7/euz/u1P7q/99V1fUY+Lt/re1+X\nQqlUgoiI5KMjOgAREZWNBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkET\nEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGR\npFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmK\nBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQ\nRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1E\nJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESS\nYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkW\nNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkET\nEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESSYkETEUmKBU1EJCkWNBGR\npFjQRESSYkETEUmKBU1EJCkWNBGRpFjQRESS0hMdQDucO3cuNzdXdAoiqhJ6enrNmzcXnaIM\nLOiXO336dOvWrUWnIKIqdOrUqVatWolOURoL+uWys7MBZGVlGRgYiM5CRBWTlQVDw/I2yM7O\nNjQ0LPhjLhvOQRNRNZSejs8/h5UVjI3RpAk2bxYdqFJY0ERUDY0ahblz8fAhlErExOD997F7\nt+hMFccpDiLSVs+eYcMGxMejYUMMGAATE9X43btYs6b0xvPm4d13NRzwdbGgiUgrXb2KgAAk\nJqpeTp2KrVuxfz9iYqBUlrH9lSuaTKceLGgi0kqffFLUzgASE+Hnh5ycF27v7KyBUGrGOWgi\n0j6ZmYiIKD1YTjsDGDGi6uJUlZp+Bp2dnf3rr7+Wv8ImNjZWY3mI6FXk5SE//+WbKRRQKmFm\nhmnTMHhw1cdSt5pe0CkpKWFhYVlZWeVsk5aWBiA7O5vroIkkUbs23nwTZ868ZLOrV6FQoEED\n6OpqJJa61fSCdnR0vHTpUvnbLF++PCQkRDN5iOgVLV+OgAA8fap6qaeHUrdj8PBA48aaz6VO\nNb2giUhLtWyJmBgsX47r19G4MXr3Rr9+uHZN9W7duli7Vmg+dWBBE5G2srfHtGlFLy9exMaN\niIlBgwbo3x916ohLpiYsaCKSwpEj+OornDmDunUxcCC+/Ra1a1dsD0ZGGDKkasIJwoImIvHO\nnEGXLihYTnX7NsLCEB+PTZtExxKN66CJSLy5c1FqsevmzYiJEZRGGixoIhLv8OEyBl+2wKr6\nY0ETkWBHjyIlpYxxJyeNR5EMC5qI1C83F4mJr/RhPwD79pUxWL8+mjVTbyjtw4ImInXKzMSk\nSTA1hYMDzM0xderLa3rXrjIG33sPejV+EQMLmojU6YsvMH8+MjMB4OlTTJuGpk1x9eoLtz9x\nAufOlTHetWtVJdQiLGgiei3FT5CfPcPSpaU3iI1FYGDRZ7JLef6mdAAcHNCrl5ryaTMWNBFV\n0pYt8PKCgQHq18e0aUhNhbd32ff8jIvDnj1l76TMJ7p26waFQp1RtVSNn+MhokrZuRPBwaqv\nk5IwdSr27EF09Au3f9FdewMCyrjPUZcuakqp5XgGTUSVMXt26ZGTJ8vbPj0d8+dj0yY8e1Zi\n3NMTc+aUuB3o8OF4/301pdRyPIMmosq4fLkCGxsaYtYs1deOjtizB56eRe+GhuKdd7BnD7Kz\n0a4d/P3VmVOrsaCJqDJsbfHo0cs309GBrS2SkopGbt3CwIE4f77EZh4e8PBQc8JqgFMcRFRh\naWm4e7f0YJm393zvvTIu9124UOJ5r/QiPIMmoiJXrmDnTqSlwcICjo5wdy99Yrt5M1atQkwM\nHjwo/b1lnlCfOIGMjDLGX7TqjopjQRPVdDt3Yu5cXL0KY2MkJCAvr8S7wcEID4eBARITMW8e\nFiyo2M6zs+HtjT/+KDFoaQkXl9eNXRNod0Hn5+dfu3YtJyfH1dVVjx8LpZoqMxNxcbC2hq1t\nhb9340YMGFDeBr/9Bh0d3LyJEycqk83bG/PmwccHaWmqEV1dLFsGHU6vvgKt+UWaMmXK6tWr\nC1/m5uaGhYWZm5u7urp6eHiYmJiMHDny8ePHAhMSCTF/Pqyt0awZ7OzQtStu3Xr5t+TlIToa\nkZF4+hRff/3y7TdtqmQ7GxkhLAzu7oiORmgounfHJ58gMhL9+lVmbzWRUksACAgIKHw5duxY\nAHXq1AkODh45cqSPjw8Ad3f3zMxMtR962bJlANLS0tS+Z6LXtH69EijxX9u2ytzcog0yM5VZ\nWSW+5dQpZdOmqo3NzEp/eyX+UyjK2E+tWsoPP1RevqzhX4/KyMrKAhARESE6SBm05gy6uOjo\n6MWLF7dp0yYuLm7z5s3Lli07fvz4qlWrLl26NHPmTNHpiDTnp59Kjxw/jitXAOD8ebRvj1q1\nYGQECwsMG4YHDxAdjU6dipYwP3mihgwdOiAiAlZWRSN16+LkSaxbhyZN1LD/mkwr522PHTum\nVCoXLFhgVeyHYtiwYatXr967d++04o/5Jape8vMRHo79+5Gfj86dER9fxjbx8ahbF127Ft0F\n//Fj/PQTtmxBWhqUSjVH+uwzeHoiNhY//4z4eDRujCFDYGmp5qPUTFpZ0Ldv3wbg5eVVatzL\nyys8PFxEIiJNUCrRvz9++031Mjy8xCekCy1aBDu7Mp5RopbzZQB6ejAwQEYG3NwwfTq6dwcA\nS0uEhqpn/1RIKwvaxcUFQEJCgkfJJZrJyckNGzYUk4mo6u3aVdTOBUotiStQak2bGhkYwNcX\nYWFo3RpZWWXfiI7USJvmoGNjY2fMmPHbb7+5ublZW1vPmDGj+LunTp3avXt369atRcUjqmrH\njgk7tKkpFi1CVhYOHULBHzK2swZozRm0o6Pj7du3p0yZUjiyYcOGESNGdOjQAcDkyZMXLlxo\nZmY2depUYRGJ1Cc/H1ev4u5dGBoiKwuenrCwELB2uHdvzJqFx4/h7g5TU00fnbSmoG/evJmR\nkXH16tXY2Ni4uLi4uLjY2NjCD6ds3769bt2669atc3R0rOienz17lp2dXc4G6enplQxNVCkx\nMRgwAGfPFo0YGcHXF3/9pdEYlpZYubLE8gzSMIVS7dd0Rbh06VKTJk10Kn6Cce3aNVdX1/xX\nePjwkydPTHkKQVXv6VM0aPBKN4qrUq6u+OUX1IQpw+zsbENDw4iICF9fX9FZStOaM+jyubu7\nV+4bnZ2dz549W/4Z9NatW2fOnKngE3hII0aPFt/Ojo64fJmfxhavmhT063h+uV4pp0+f1kwS\nIgC7dgk4qJlZ0SI8ExOEh7OdpcCCJhJDqSzjRsk5OUhN1XQSf39s24Z16xATA0dHDBmC+vU1\nnYHKxIIm0qiMDHz3Hdaswf378PTEjBkIDFS9lZ+P/v3xChdEXpezMx4+hL4+3nwT/ftj8GDo\n6WHChCo/LlWUdhT04sWLiy+wK98j4RN4RC8WEoK1a1VfR0WhVy/s34+AANy8iW7dKvagv8rp\n2xebN5dx8k4S0o6C7tatW1xc3PLly7OyskxNTZ2cnEQnIqqMxMSidi6Qm4vx4+Hvj2XLyv5Y\noBrVqYMFC/Dhh2xnraEdBf3GG28sXLgwMDCwW7du7du337lzp+hERK/q/n188gn++APZ2WWv\nKT5/vvQTVKtIw4YYPFgTByJ10Y6CLvDOO++4urqKTkFUAbdvo1mzomVzz9/ASJPMzUUenSpB\nmwoaQOvWrfm5PpLZkyfYsQNJSUhKwt69iItT/+09K+3dd0UnoArSsoL+5ZdfREcgeqEzZ9Cj\nB5KTRecoS79+mDhRdAiqIC0raCJpKZX44APx7WxsjJYtoaeHpk3h4wNDQzx4gJYt4e0tOBhV\nAguaSA1yczF9OuLiBMeoUwenTsHZWXAMUhcWNFElnTiBlSuRmAgbG1y7hogIwXnmz8eYMdDX\nFxyD1IgFTVQxDx5g9mxs346rV0VHKWbSJE4xV0O8IQrRq8rNxa5dcHHBvHnC2tnSEj/+iKCg\nopFatTB3LubNE5OHqhTPoIleyfXr6N0bFy4IjrF0Kd5/H59+ivPn8fffsLJChw4wMRGciqoI\nC5rolQwaJL6dfXzw/vuqr5s1Q7NmQtNQ1WNBE5WWno4lS3D8OMzMEBSEXr2QkoITJwSnMjUV\n+dBYEoIFTVTCvXto3hx37qherlkDCwsYGwvNBADo3Jk3OapxeJGQqMj+/WjUqKidC6Smlh6p\nUkZG0HvuxMncHGFhmstAkmBBEwFAaio2bcJ77+HZM5ExbGywfj127ICnJxQK6OrC3h6fforo\naLzxhshgJASnOKgmevgQa9ciPh7Ozhg8GNu2YfRoZGQITvXddwgNRe3aANC9O9LTYWTEZwPW\naCxoqnEuXICvL54+Vb2cPBmZmZp40FT5GjbEpEmoVatopPjXVDPxb2eqWWJi0KZNUTsDSE/X\nUDvb2eHJE3TtWsZbCgV27GAjU2kvL+g1a9Y8KXwgO5E2e/AAvr7IzBRz9KAgmJpi1Sq4uZV+\n65tv4OUlIhPJ7eUFPXTo0Hr16vXt23fz5s0ZwmfpiF7Bgwfo3x+WlqhdGzY2aNMGPXvC2xv2\n9nj4UFiqHj0AwMEB0dFYuRKtW6NBA7Rvj59/xrRpwlKRzF4+B71kyZJNmzZt375969atJiYm\nvXv3HjBgQNeuXfV51yyST3IyfvgBc+ciN1c1kp6Oe/eEZgJ0dfHNN+jevejl8OEYPlxoJtIG\nLz+DHjVq1OHDhxMTExcvXtyyZctff/21R48etra2I0aMOHToUL7waytE/+/QIbi6YtasonaW\nwZtvIiYGU6eKzkFa6FUvEtra2o4ePbqwqT09PVetWtWxY0cHB4cJEyZERkZWaUqi8u3di4AA\ndOmCtDTRUYrR0UFICE6e5B30qZIqvIrD1tbWz8+vY8eOTk5OAO7cubNw4UIfHx83N7ctW7ZU\nQUKi0nJz8cMPeOst2NqiQwe0aYPAQPz5J/LyxOQxNlYtXi5kYoKvv8b161i6lHfQp8p71XXQ\nubm5R48e/f3337dv356QkADAzs4uJCSkT58+VlZW4eHhK1as6Nev38mTJ1u1alWVgdXs3r17\n48ePzy33n8TXr18HoJTn4cw1zP37OHgQjx/D21t1/7YvvsD8+ap3U1IERlPJyMCUKdi4EbGx\nANC2LVasgKen6Fik/V5e0Fu2bPn999937dr16NEjAM7Ozp999lmfPn18fHwU/3/vlhYtWnz4\n4YctWrTYsmWLdhW0oaFh48aNyy/oglWGCt6oRoRduzB4MB49Ur20tYWfHyT8p5qeHmJikJIC\nfX1YWopOQ9XFyws6ODgYQPPmzcePHx8UFNTsBfegdXZ2rlu3rpWVlZoDVjEzM7MZM2aUv83y\n5cv37dunmTxU3L17+OgjpKYWjSQn47ffxAV6sYIbZdSrJzoHVS8vL+h58+YFBQU1bty4/M3M\nzMzuCV/NRNVLRESJdpbEO+/g0CFkZxeNNGiAwEBxgaj6evlFwkmTJr20nYmqgmztXK8e5s7F\nnj3YsweurqrBtm2xezcsLIQmo2qKN0siST15gsREkQEUCqxejdhYODigTx+YmhYt1ejUiTPO\npAksaJLRihUYMwY5OSIzBAdj6NDyNuCMM1U1FjRJRKnEjh04cABLlkDsssZ+/bBihcgARGBB\nkzyystClC44eFRyjY0ds3Ii6dQXHIALvB03C7dsHX1+YmsLERHw729lhzx62M8mCZ9CkUWvX\nYuFC3L8PJyc4OODECcTHi8yzciWsrfHvfyM9Hf36ISREZBiiUljQpCEZGWjfHqdOqV7evCk0\nDQDAzAzDhkGhwHvviY5CVBYWNFWJnBz85z+4cQMuLujSBQoFOnUqamdJLF8OfoCfZMaCJvXI\nzsbatYiKgqUl2rXDuHG4ckX1lqUlJk/G8eOajlSnTtFNPAwN8e676N0b+/bh2jU0aIAxY+Dv\nr+lIRBXCgiY1ePoUb7+N8+dVLxWKEovkHj7EP/8pINW8eXB2xvnzsLVF9+4wMQGAjz4SkISo\ncljQpAaff17UzoDgJcwFmjfHxx9DoUD79qKjEFUWl9nR68rNxbp1gjOUuim+hQWOH+f8Mmk9\nFjS9rrNn8fSpyAC6uti3D/37o1492Nqif3/Ex8PYWGQkIrXgFAeVISMD+/bhzh14eKBduxJv\n5eTg2bMSN28Te4Nma2ssXIgOHdChg8gYRFWBZ9BU2oULcHdHUBBGjUL79ujcGenpAJCUhOBg\n1K6NOnXg5obduwFg1izMmSMs6po1SEzEBx8IC0BUpVjQVIJSiQ8+wI0bRSP//S++/hq5uejT\nB1u2qO4wFxuLoCAcOIBvvhEUFLC2Rp8+fCQrVWcsaCrhxg1ER5ce3L0bkZGIjCwxmJODjz9G\nfr7GokGn2E+rhQV++UW1co6ouuIcNJVQ5kNM7t1DeHgZ47dvq/noffrA2Rk3b8LfH/n5mDCh\n6C8AQ0McOoTMTNVnYd59F9bWaj46kWxY0FSChweMjZGRUWIwNRVLl1btcXV18fPP+PDDEoOd\nOmHxYty4gTfewLhxqgez8mIg1RwsaCohJwdmZqULuioYGSEzU/V1y5Y4dAimpqW3cXfHjz9W\neRIiaXEOukb7z3/w7rvw8kK/fjh9GgAOHUJKSpUft107JCVh0yZ8/z327cOpU2W0MxHxDLrG\nyc9XXW1bvRrDh6sGL17Etm2YOxf792sig5kZ6tRBv36aOBaR9tLuM+j8/Py4uLhLly7l5uaK\nziK71FSMGYO6dWFsDF9fHD6MCRNKbJCXh9BQ7NuniTAXLmjiKETaTmsKesqUKatXry58mZub\nGxYWZm5u7urq6uHhYWJiMnLkyMePHwtMKLOC1c1LluDBA2Rn4/hxdO6MtLQytlTvsjl9/bLX\nWnBCg+hVaE1Bz5gxY12xW/KEhoZ+8cUX+vr6wcHBI0eOfOutt1asWOHr65uVlSUwpLSiovCf\n/5QYycur8oMGBCA1FTt2lPFW//5VfnSiakBrCrq46OjoxYsXt2nTJi4ubvPmzcuWLTt+/Piq\nVasuXbo0c+ZM0emkk5SEr77S6BH19LB3Lw4dQq1a8PHBqFEl3m3VCl9+qdE8RFpKKy8SHjt2\nTKlULliwwMrKqnBw2LBhq1ev3rt377Rp0159V3l5ebt27crOzi5nmzNnzlQ+q2g3bqBFi6IH\ni1QpHR34+sLXF99+i9q1i8aXLEHfvvj9dzx9ioAADBwIXV1N5CHSdlpZ0Ldv3wbg5eVVatzL\nyyu8zE+8lburTz/9NLNwRW5ZCqZNlDLchb6CcnIwdqyG2tnYGKtWvfC+RR07omNHTcQgqk60\nsqBdXFwAJCQkeHh4FB9PTk5u2LBhhXbl5OSUlJRU/jbLly8PCQlRaMPt33NzEROD9HR4eODI\nEfzjH0hM1NChf/0VvXtr6FhENYQ2zUHHxsbOmDHjt99+c3Nzs7a2njFjRvF3T506tXv37tat\nW4uKJ9yJE/DwgKcn2rSBvT2CgqqkncucnejQAe+9p/5jEdV0Si3h6Oj4/DnswYMHC9798ssv\njY2Nraysbt68qfZDL1u2DEBaWpra96xG9+8r7eyUQNX+Z2mpXLdO+e23yoYNlYaGSiMjpaOj\n8rPPlI8fi/7/J6qsgjnMiIgI0UHKoDVTHDdv3szIyLh69WpsbGxcXFxcXFxsbKyenir/9u3b\n69atu27dOkdHR7E5NSwzE8OGYcsWlHuZs8L09VX3fS40YgSCg9GmDczNAaAiF2KJqJK0pqAB\nGBsbe3l5PX9tEMCWLVuaNGmio6NNMzaVEBWFK1fg6AhfX9XHtbt2xdGjaj6KuTlu3sTp05g1\nC5cvo0EDjB6NgQP5DFYiTdOmgi6Hu7u76AhV69kz9OuHvXtVL1u2xLZtOHhQ/e0MoE4dmJlx\n3QWReNWkoKu9zz8vamcAZ87ggw9w9myVHKtbtyrZLRFVVDWfE6g2tm0rPRIRgWfPKraT3r1h\nZvaSbdq2RVhYxXZLRFWEZ9BaIDsb9++/7k4UCsyahTNnSj+1pMDQofDygrs7unZFdZ/JJ9Ia\nLGjZPXyITp1Kr6moBKUS+voYNAh2dvjyS5w5U3Tjur59sWoVe5lIOixo2YWGqm2u+fRpODuj\nY0ecPIlbt7BnD1JT0aYNn/JHJCkWtFwSEnDrFlxdYWOjGvnjj8rsR6HA8/cOuXu36GtHR4wc\nWamIRKQp/GeteHfu4M8/ceYMevZEw4bw94edHUJCkJsLpRLp6S/fg4EBateGuTn8/LB2Lf78\nE0FBZWxWr57asxNRFeIZtEi5uRg9GitXln6OSX4+li+HjQ2++w4+Pi9/DFVODq5fR/36RSOP\nH2Pr1hLbWFoiIEA9sYlIM3gGLdL06Vix4oVPmfrlFwBYsODla+OUSkRFlRjp2RPTp8PAQPWy\nXj2sX180bUJEWoFn0GKcO4clS1D+zasTEtCyJRo1wk8/4dQpnDmDiIgXzng8X77ffIOhQ3Hy\nJExM0LYtHwNIpH1Y0ALs2YP33nv5UwHz8/H33/j7b2zZgg0bMGsWUlOxfDm2bMGpUyW2bNwY\nzZuXsQcHBzg4qC02EWkYpzgEGDWqws9snTABACws8MUXiIxESEjRsmVnZ2zeDENDNYckIuF4\nBq1pyclISHjJNs8vkktORmKi6jKgQoGlSxEainPnYGUFX1+2M1H1xILWNBMT6Oi88MIgAAsL\nAEhNLTGoq6saL+TiAheXKshHRNLgFIemmZjA3/+F71pYYN26Mh7u16VLiedkE1FNwDNoAX76\nCd26ITZW9fLNNzF3LuLiYGGBzp1hbQ1/f9y4gcOHVRu0bIlVqwRlJSJxWNACNGqECxewdy9u\n3ICLC955B7q66Ny5aANzcxw6hMhIxMbCyQl+fryTEVFNxIIWw8AAvXq9ZBtvb3h7ayQNEUmJ\nJ2ZERJLS1jPoJ0+epKWl6ejo1KtXr9o/K5aIaiYtq7aLFy8OGTLEzs7O3NzcwcHB3t7ewMDA\nwcFh4MCBERERAoPl5SE+vsLPoCIiKoc2FfTYsWObNWu2du1ahULh7e0dGBgYGBjYunVrhUKx\nfv16Pz+/Tz75REiwBQtgaYnGjWFqioED8fChkBREVN1ozRTHjz/+uHjx4nfeeWfWrFlvvfVW\nqXejo6OnT5++cuXKpk2bhoaGajLY2rWYOFH1tVKJ9euRloadOzUZgYiqJ605gw4PD3dzc9u1\na9fz7QzAw8Nj/fr1/v7+W0vdBbnqrVhRemTXLiQmajgFEVVDWlPQFy9e9PHx0dN74Sm/QqHw\n9/e/ePFihXZ7/fp1Y2NjRblCQkIK9l/mHq5de9VBIqIK0ZopDk9Pz8jIyLy8PF1d3Rdtc/z4\ncU9PzwrttlGjRvv27cvOzi5nm+jo6AkTJujr65f5bpMmSE4uY5CI6DVpTUEPGjRo9OjRPXv2\nnDNnjpeXV6l3Y2Njp02bdujQobCwsArtVqFQtGvXrvxtatWqVc67n31W9JnsAkOG8NklRKQG\nWlPQo0aNunDhwrJly/bu3evo6Ojk5GRpaalQKB49enTr1q34+HgAQ4cO/eyzzzQc7N13sWkT\npkxBbCzq1MHw4Zg6VcMRiKh6UihL3XhYbmfPng0LC9u/f//9+/cLRnR1dW1sbAICAkaOHNm+\nffuqOOixY8fefvvtrKwsg8Jn/JUlK4v3ZSbSPtnZ2YaGhhEREb6+vqKzlKY1Z9AF3nzzzV9/\n/RVAampqWlqavr6+jY1NVX+SsKCXDdm+RNVX+adfomjZGbQo586dy83NfZUtJ06caGxs/NFH\nH1V1pFeRmZn5ySeffPfdd40aNRKdBQAOHDhw4MCB2bNniw6iMm3atBYtWvTs2VN0EABITk7+\n5z//+cMPP9SpU0d0FgDYvHnz9evXv/jiC9FBVMaPHz927Ngq+s3S09NrXuZjPYVTklr16dNn\n3LhxolOopKWlATh16pToICpLlixxd3cXnaKIj4/PrFmzRKdQiYmJAZCYmCg6iMrXX3/dpUsX\n0SmKODo6rl27VnQKTdOaddBERDUNC5qISFIsaCIiSbGgiYgkxYImIpIUC5qISFIsaCIiSbGg\niYgkxYImIpIUC1rNDAwM5PlQv56eno6Ojjx5pPrFgWR5DAwMFAqFPHn09fXlCQPJfrM0hvfi\nULP79+8bGBiYmZmJDqJy/fr1xo0bi06hkp2dfffuXQcHB9FBVO7cuWNhYWFsbCw6iIpUv1lP\nnz5NT0+3kebW5jdv3rS3ty/nmUrVEguaiEhSnOIgIpIUC5qISFIsaCIiSbGgiYgkxYImIpIU\nC5qISFIsaCIiSbGgiYgkxYImIpIUC5qISFIsaCIiSbGgiYgkxYImIpIUC5qISFIs6Gru6tWr\nixcvFp2CXtXTp0/XrFlz+/Zt0UFICixotVm6dKmfn5+FhYWfn9/SpUtFx1FZtGjRlClTRKdA\nVlbW119/3a5dO3Nzc2dn54EDB167dk1UmPj4+IEDB7q4uNSuXdvLy+vzzz9//PixqDCljB07\ndujQoefOnROYwdHRUfEcgT9FR48e7dy5s7m5ub29ff/+/QX+5AigJHUICQkB4ObmNnjwYFdX\nVwBjxowRHUr5xx9/GBoaWlhYiI2Rmprq7+8PwN3d/R//+EfXrl0VCoWxsXFUVH1BzMUAAAmO\nSURBVJTmw8TFxdWuXVtPT69jx44hISHe3t4APDw8MjIyNB+mlM2bNxf8qdy1a5eoDOnp6QqF\nwt7ePqCkVatWCcmzYcMGAwMDe3v7gQMH9urVS1dX18rKKiEhQUgYzWNBq0FUVBSAbt265eTk\nKJXKnJycgg66cOGCqEiDBg1yc3Mr+NMuvKAnT54MYPTo0YUju3fv1tHRad68uebD9O3bV6FQ\n7Nixo3Bk4sSJABYtWqT5MMXdvn3b0tLSxMREbEGfP38ewPTp00UFKC4hIUFPT8/b2zs1NbVg\n5N///jeAIUOGCM2lOZziUIOwsDAAc+bMKXhgmp6e3qxZs5RK5dy5c0VFSk9Pd3Fx6dGjh6mp\nqagMhbZt22Zqavqvf/2rcCQwMLBjx47nzp27e/euhsP89ddfLVq06NmzZ+HIxx9/DODvv//W\ncJLilErl4MGDzc3Nx40bJzAGgNjYWABNmjQRG6PAwoULc3NzFyxYYG5uXjAyfPjw77//3sfH\nR2wwjalZT2CsIvv373dwcGjWrFnhSIsWLezs7P744w9RkbZu3VrwhZeXl/ArTjo6Ou3btzc0\nNCw+WPCE5kePHmnysaT5+flTpkxxcnIqPpiSkgLgjTfe0FiM5/3rX/86fPjwn3/+GRERITAG\ngLi4OABOTk7h4eFxcXEODg6+vr7u7u5CwmzYsMHR0bF4HSsUigkTJggJI4boU3it9+jRIwBv\nv/12qfGCyc0nT54ISVXI09NT+BTH8+7evWtkZFSvXr2CSSEh0tPTExMT9+zZ4+LiUq9evdjY\nWFFJoqKiDAwMJk+erFQqZ8+eDaFTHMOGDQNgbW1dWBE6Ojpjx47V/O9UWloaAH9//7Nnz/bs\n2dPGxsbR0TE4ODguLk7DSQTiFMfrKvgxsrKyKjVeMPLkyRMBmeQWGxvr6+ubmZk5e/bsgkkh\nIUJDQ+vXrx8YGJiUlFRQ00JiZGRkDBo0yN3dferUqUIClFIwxdGpU6fz58+npaX99ddfLVu2\nXLRo0fz58zWcJDU1FUBSUpKfn9+NGzd69Ojh4eGxdevW5s2bnz59WsNhhBH9N4TWu3PnDoBe\nvXqVGg8MDASQlJQkJFUhqc6gnz59+u233xobGxsZGS1evFhsmLNnz27cuPF///d/GzRoYGho\nuH37diExRo8ebWRkdPHixYKXws+gjxw5cvDgweIj9+7dq1OnjomJSV5eniaTXL9+vaCjvvzy\ny/z8/ILB/fv3KxSKt956S5NJBGJBv668vDxdXd127dqVGvfx8dHV1dXwz/Tz5CnoPXv2NGjQ\nAECPHj2uXLkiOk6RxMREU1PT+vXra/7QBw4cAPD9998Xjggv6DIFBwcD0PAsUHJyMgArK6vc\n3Nzi4127dgWQkpKiyTCicIrjdeno6NjY2Dx/IS4xMdHW1lZHh7/CAPA///M/gYGBpqamf/75\n586dOwuXAGrYtWvXli9ffvHixeKD9vb2rVq1SkxMLLicoElnz54FMHHixMLPg3z55ZcAevTo\noVAoVq1apeE8L1IwX5eTk6PJg1pbWxsZGTVq1EhXV7f4eOPGjQEIv/StGVzFoQYBAQHr16+P\njY0t+IgKgOjo6Fu3bn3wwQdig0lizZo133333YABA9asWVOweEOUlJSUkJCQcePGLVy4sPj4\nvXv3TExMCtdyaUzz5s0LPuJUKCoqKjIysnv37k5OTppf63bp0qW+ffsGBQXNnDmz+Pi5c+cM\nDQ0Lf7w1Q0dHJyAg4NixY5mZmUZGRoXjly9f1tHREfV3vKaJPoWvDg4fPgzgww8/LHiZn5/f\nv39/AEePHhUbTCnBFEd+fr6bm1v9+vVl+Khedna2jY2Nubn5tWvXCgc3bNiAsq4iCCF2iiMv\nL8/R0dHY2PjkyZOFgwUn8iNGjNB8nn379gEYPXp04VThxo0bAfTo0UPzYYTgGbQatG/ffujQ\noT///HNSUpKPj89ff/115MiR4cOH+/n5iY4mXkJCQkxMjLW1dVBQ0PPvrlu3rm7duhoLo6+v\nv2jRogEDBnh5eQUGBtrY2Fy+fPnQoUP16tVbsmSJxmJIS0dHZ926dX369Hn77bd79Ohha2t7\n/vz5iIiIpk2bzpkzR/N5unbtOnTo0CVLlhw5cqRt27bx8fH79++3s7OT5143VU703xDVRH5+\n/pw5c3x9fc3MzHx9fefOnSs6kYrwM+j//ve/5fz43b59W/ORDh482K1bNysrq1q1ajVv3jw0\nNPThw4eaj1EmGS4SJiQkfPzxx56eniYmJq1atZoyZYrYf/3MmzfPz8/P1NTU3d19zJgx8vxm\naYBCqVRq6u8CIiKqAK4xICKSFAuaiEhSLGgiIkmxoImIJMWCJiKSFAuaiEhSLGgiIkmxoImI\nJMWCJiKSFAuaiEhSLGgiIkmxoImIJMWCJiKSFAuaiEhSLGgiIkmxoImIJMWCJiKSFAuaiEhS\nLGgiIkmxoImIJMWCJiKSFAuaiEhSLGgiIkmxoImIJMWCJiKSFAuaiEhSLGgiIkmxoImIJMWC\nJiKSFAuaiEhSLGgiIkmxoImIJMWCJiKSFAuaiEhSLGgiIkmxoImIJMWCJgKA6OhoQ0PDDh06\nFI7k5OR4eXlZWVklJycLDEY1GQuaCAA8PDy++uqrw4cP//TTTwUjc+fOvXjx4g8//GBrays2\nG9VYCqVSKToDkRSys7NbtmyZlJR05cqVx48fe3l5denSZceOHaJzUc3FgiYqEhkZ6evrO2DA\ngOTk5KioqOjoaDs7O9GhqObSEx2ASCLe3t7jx4///vvvAaxdu5btTGLxDJqohKtXr7q4uNSu\nXTspKcnMzEx0HKrReJGQqIRJkyYZGBg8e/Zs8uTJorNQTceCJioSHh6+Y8eO2bNnBwcHL126\n9NixY6ITUY3GKQ4ilZSUFA8Pj4YNG0ZGRqakpDRt2tTBwSEqKsrAwEB0NKqheAZNpDJq1KjU\n1NQVK1bo6ura29vPnDnz0qVLs2fPFp2Lai6eQRMBwMaNGwcMGDBp0qR58+YVjOTn57dt2/bc\nuXNRUVFNmzYVG49qJhY0EZGkOMVBRCQpFjQRkaRY0EREkmJBExFJigVNRCQpFjQRkaRY0ERE\nkmJBExFJigVNRCQpFjQRkaRY0EREkmJBExFJigVNRCQpFjQRkaRY0EREkmJBExFJigVNRCQp\nFjQRkaRY0EREkmJBExFJigVNRCQpFjQRkaRY0EREkmJBExFJigVNRCQpFjQRkaRY0EREkmJB\nExFJigVNRCQpFjQRkaT+DyPjRYuycX2IAAAAAElFTkSuQmCC", "text/plain": [ "plot without title" ] }, "metadata": { "image/png": { "height": 180, "width": 240 } }, "output_type": "display_data" } ], "source": [ "beta0 <- 1\n", "beta1 <- 3\n", "sigma <- 0.5\n", "\n", "n <- 1000\n", "x <- rnorm(n, 3, 1)\n", "\n", "y <- beta0 +x*beta1 + rnorm(n, mean = 0, sd = sigma)\n", "plot(x, y, col = \"blue\", pch = 20)" ] }, { "cell_type": "code", "execution_count": 10, "id": "706f5887", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Warning message in dnorm(x = y, mean = mean, sd = sigma, log = TRUE):\n", "“NaNs produced”\n", "Warning message in dnorm(x = y, mean = mean, sd = sigma, log = TRUE):\n", "“NaNs produced”\n", "Warning message in dnorm(x = y, mean = mean, sd = sigma, log = TRUE):\n", "“NaNs produced”\n" ] } ], "source": [ "logNormLikelihood <- function(par, y, x)\n", " {\n", " beta0 <- par[1]\n", " beta1 <- par[2]\n", " sigma <- par[3]\n", " mean <- beta0 + x*beta1\n", " logDens <- dnorm(x = y, mean = mean, sd = sigma, log = TRUE)\n", " loglikelihood <- sum(logDens)\n", " return(loglikelihood)\n", " }\n", "\n", "optimOut <- optim(c(0.2, 0.3, 0.5), logNormLikelihood,\n", " control = list(fnscale = -1),\n", " x = x, y = y)\n", "\n", "beta0Hat <- optimOut$par[1]\n", "beta1Hat <- optimOut$par[2]\n", "sigmaHat <- optimOut$par[3]\n", "yHat <- beta0Hat + beta1Hat*x" ] }, { "cell_type": "code", "execution_count": 11, "id": "ea97b1ad", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAFoCAIAAAAAVb93AAAACXBIWXMAABJ0AAASdAHeZh94\nAAAgAElEQVR4nO3deVhV1f7H8fdhRkHECUVxyHkgrTSVnNM0tcyc54lDltlg/Rpuw81rt0yb\nvGUlB5zI1Mw0c0yvVqZW1lUTS8FZQVFywoHx7N8fIJOIoMjZ6Of1+Dxx1l57n+9G/LDbZ+21\nLIZhICIi5uPk6AJERCRvCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSIiEkpoEVETEoB\nLSJiUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JA\ni4iYlAJaRMSkFNAiIialgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQ\nIiImpYAWETEpBbSIiEkpoEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0\niIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAJaRMSkFNAiIialgBYRMSkFtIiISSmgRURMSgEt\nImJSCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSIiEkpoEVETEoBLSJiUgpoERGTUkCL\niJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAJaRMSkFNAi\nIialgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSI\niEkpoEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0i\nYlIKaBERk1JAi4iYlAJaRMSkFNAiIialgBYRMSkFtIiISSmgRURMysXRBZQMO3bsSE1NdXQV\nInJTuLi4NG3a1NFV5EEBfW2//fZbixYtHF2FiNxEW7dubd68uaOryE0BfW3JyclAUlKSm5ub\no2sRkcJJSsLdPb8OycnJ7u7u6f/MzUb3oEXkFnTxIi+8QPnyeHrSoAGLFjm6oOuigBaRW9AT\nTzB1KqdOYRjs2UP//qxY4eiaCk8BLSIl1YULhIfz6quEhXH+fFb7iRPMmQPgxfmaHHTCDrz7\nroOqvAG6By0iJdLevXToQExMxss33uDrr1m7lj17wG7vylortodZ5kpKDFVHM/OP3Q84tN7r\noYAWkRLJas1KZyAmhjZtqJgSO5qZYwivyUGASlCfqj/FLDb6DA3YA/4OKvY66RaHiJQ8iYls\n2pT10pm0Hqz4KqXXIWpM4rWaTgd5ABbBEfgRnsaL8y82KXk3oW/3K+jk5OQvvvgi/xE2UVFR\nxVaPiBREWhp2O0B1Do9m5mhmBnAEoAqMgmCola33WYDWbUte3JW8iotWXFzclClTkpKS8umT\nkJAAJCcnaxy0iEmUdk99pubyTvtsXVnjTBrO0BWs0DNbql2AL2EG/AK+vjz0kCMrvi63e0AH\nBAT8+eef+feZMWPG2LFji6ceEbmGAwcIC2PWrHePHctoGQ//BwHZ+mwDG3yRce1Mo0aEh1Oh\nQrHXeqNu94AWkZIhOZlly7DZWLcu4+5GOn/4z+WvE2ABhMJvAIabu2Vgb6xWOnbEYin+km+c\nAlpEzC06mrAwZs/mxIk8th6HCKgJETAfzgOcqdzA/Umr52PDS+JVc3YKaBExhR9/5B//4Pff\nqVCBwYN5/aXk0uu+ITSU//4Xw7jqbnYYfvlrd3f6PUxISNn77y+hl8y5KKBFxPF+/50uXUgf\nTuV99C+/KWH2aXNJigeoA2NgIGyCoVfZv0kTrFaGDcPXt9hqLgYKaBFxvKlTcUpOHMYiK7a2\nbAQwYABYoROkXw37gTOkZdutVCkGDMBqpXVrR1R90ymgRcTRIiM7f2v7lAhfTgPUh2AYARWz\n9fkR/p0tnZs1IySEwYPx8Sn2couPAlpEHCQxkW+/Tb/LHGwYuFy+ZG53+ZIZiIe5YIPdAHh6\n0rMnISF07uywsouRAlpEil5qKnFxVKmCU57TSWzfTmgoX3zB2bNZjdPgictfG7ABbLAE0h8j\na9ECq5WBA/H2vsm1m4gCWkSKUmIir7zCJ5+QmIiXF889x+uvX47phAQWLMBmY+vWPPZMX/Xz\nBMyGMIgGwMeHMUOwWmnWrJhOwEwU0CJSlF58kf9cfnLk/HkmTmT+fNa9vTVgtY0FC0hIuOqe\nEyAcdkMywKn6rcu9aGXAAEqVKo66TUkBLSI3xG7Puo9x4QKffpq1yYezQ5hnjbIF9NkO4AWj\nwR/eg0tXHCgN/uA0vhEMW+ZnXftXE26Focw3RNONish1WryYwEDc3KhalYkTOXOGli1JSQFo\nzZaZjI7FfzrjmrGd5jADYiEcJsGQPI62kbbDmVuVmKeZVuuhJrfEgyY3SlfQInI9vv2Wvn0z\nvo6N5Y03WLmS2F2nnyLCiq0JkQA+MBiscFe2PQ/Cj1mv0nwrfHR2eKg9+C8aZjZ26XLzT6Ak\nUECLyPWYPDnHy3b8aP3V1pevPEgEaAkhMABKX+6RBishFFZBGlgsdOhASIhz79726e5RL2SN\ncR4zhv79i+9EzEwBLSLX46+/ACoQP5y5VmwNMgYqQ3lYBS2ydT0E4TATYgDinf0YPaLCS8HU\nrZu+fcIEunZl5UqSk2nXjrZti/M8TE0BLSKFZxiPeK9/4HRYb5a4k3O9ixaX0zkFloEN1oId\nO07r6GzDuizt4fpb3f6om2Onxo1p3LjYqi8xFNAiUhjHjzNnjj00bObhvXl3WAuvQjJEwHGA\nOGf/MEaFM+bA5XWodu4kJoaqVYur5hJLAS0iWXbv5ttvSUigbFkCAmjU6PKFrd3O2rXYbPal\ny5zSUpwsUA8OQMoVh0iDf6f/13kVD37la11u7/732dxRc/78TT+XW4ACWuR29+23TJ3K3r14\nenLoEGlpObZae8R+2mKm8+xwDh4EnCrCSAiGevAV9MvjgEcImMnocMYcIaCcheYt+e67HB3K\nlcu8/yz5KdkBbbfb9+3bl5KSUq9ePReXkn0uItctMZHoaCpWpHLlQu+7cCEDB+bR7kxaN1Zb\nsfVYscJ5RSpO0Bms8AhkLp6cc1aMVFxW0MOGdTXd0nBOb2zZknffpVWrrEcInZ357LOrzNEh\nuRglxKuvvhoeHp75MiUl5Z133vHy8ko/C3d395CQkDNnztyMt/7ss8+AhISEm3FwkRv03nuG\nl5cBBhhduhiHD197l9RUIzLS+PlnIyHBqF07Y9/MPwEcfoN/HiYg43UVjH9g7MMwsv3ZhfEM\nhnfGPgeo+SqT/InJdSgPD2PnTsMwjMOHjQkTjAcfNKxW47ffbva3pHCSkpKATZs2ObqQPJSY\ngAY6dOiQ+XL8+PGAr69v3759H3vssVatWgGNGjVKTEws8rdWQItpzZ+fO15btzZSU7M6JCYa\nSUk5dtm61WjYMKNzmTJZO7qQ0ouly+mRinNG010YizGSs+XyBYzZGPdl7JOM6yL6drWsKVsm\nLVcZpUoZQ4caf/1VzN+P66GALgLZAzoyMtJisdx7773x8fGZHcLDw4HXX3+9yN9aAS2m9cAD\nuQMajMhIwzCMHTuMdu0MJyfDYjF8fIxRo4z4eCMyMkcop/+pxf43eSUG/9wborJF8w6MJzHK\nZmyKou6LTPbjOBidOhk7dxrly2ftV6FCRg0lgpkDukTet928ebNhGB9++GH58uUzG0ePHj1z\n5sxVq1ZNnDjRgbWJ3FR2O/PmsXYtdjudO3PgQB59DhygQgUeeIC4uIyWs2eZNYvFi0lIyFp/\n1ZWUXnxjxdaZdU7Y8zjQKigHS8AGvwIk4b6E3jasG+hoXJ7K6PnnadKEqChmz+bAAe64gxEj\nKFeuqM/8tlQiA/ro0aNAYGBgrvbAwMB58+Y5oiKR4mAYDBjAV19lvJw3D2fnPLp99BFVqmSl\nc6Zz5zK+qMPeYMJGMtuPKzpl9zQ8nfHlbhrYsM5leDwVXFzwdOPSJerXZ9IkHnwQoFw5Jky4\n3hOTqyiRAV23bl3g0KFDjXM+e3T8+PGaNWs6piaRm2/58qx0TpdrSFy6XGPaMrmT1JslVmwd\n2WDBwAUegpGQBCPymv8TEvH4ir42rD/SDnBzo0MQU6bQogVJSbi73+AJyTWUpICOiop68803\nGzRoUL9+/YoVK7755pvz58/P3Lp169YVK1YMGzbMgRWK3FSbN1/njg3YbcU2nLkViAeoDWNg\nFGQOy3sffs6xy04CbVg/Z+hpfAFvb956iyefzOqgdC4GJSagAwICjh49+tprr2W2LFiwICQk\npGPHjsDLL788bdq0MmXKvPHGGw4rUaTo2O3s3cuJE7i7k5REkyaULVvoscMeJPZjkRVbWzYC\nuEEvCIH7s63KeibrFjNwkVJf0j+UkC20Bh55hLff5uxZGjW6rdYCNIsSE9CHDx++dOnS3r17\no6KioqOjo6Ojo6KiMh9OWbp0aYUKFSIiIgICAgp75AsXLiQnJ+fT4eLFi9dZtMh12bOHgQPZ\nvj2rxcODoCB++qmgR2hCpBXbMCJ8OQ1QD4JhJFTM1mkThMKijJsb27grjOB5DDmLT/r2cuUI\nCyPbJ/FS3CxG5me6Jdmff/7ZoEEDp8I/nLRv37569erZ7Xl9hJ3TuXPnvHUJITff+fNUr87p\n09ezbykuDmChFVtrtmS19oBvIPPjxL8hAmzwJ0AC3vMZFEbw1hwzhFKvHp9/Toscbbem5ORk\nd3f3TZs2BQUFObqW3ErMFXT+GjVqdH071q5de/v27flfQX/99ddvvfWWRSvwSLEYN+560vku\ntlmxDeYLH87m3lYVnMGAH8AGX5M+pf5WWtiwLmBgArmvPAIC+OsvPY3teLdIQN+IK4fr5fLb\nb78VTyUiwPLlhejsTcJAFoQQ2pyr/5SGwWHYD1EAZ/GZxxAb1u00y+xSpkzWIDwvL+bNUzqb\nggJaxDEMgyv/rywlhTNnCrR7C7ZasQ1kgTcJAIEQBAvhyt3tsBpgC61DCfmS/hcplX1727Ys\nWUJEBHv2EBDAiBGaqdksFNAixerSJf71L+bMIT6eJk148026d8/YZLczYAD5fyDiw9khzLNi\na8Z2gNIwEKzQEoBm8HjuXU7jO5fhYQRH0iS9pXZtTp3C1ZVmzRgwgOHDcXHhmWeK8jSlSJSM\ngP7444+zD7DL3+nr+3hFpFiMHcvcuRlfb9tGr16sXUuHDhw+TLduGQv95SmIzVZs/fmyFBcB\n7oYQGARlLveww9Ycu/xIOxvWr+ibiEdmY58+LFqUx8W7mFDJCOhu3bpFR0fPmDEjKSnJ29u7\nRo0ajq5I5HrExGSlc7rUVJ5+mrZt+eyzvB8L9OX0MCKs2JoQCVAGBoEV7sl+XJgJ4XAIIJ4K\ncxgRRvBuGuQ4lC8ffsjQoUrnksPRszUVwurVq4GePXsW8/tqNju5ESdPGo88YpQqZbi4GH5+\neUw+d7U/7fghgqGX8Mhq+jdGQrYZ5lIxlmE8hOGMAXYs/6XTQOa7k5jnAe+6y9HfC1PSbHZF\no2vXrvXq1XN0FSKFcPQod96ZNWzuygmMrlSB+BHMCSasAbtzbPCBf1z++jCEw0w4ChCH32xG\nhhG8lzr5HNnHp9D1i2OVpIAGWrRooef6xMzOnWPZMmJjiY1l1Sqioyngo2AWjI5ssGLrzRJ3\nkvLocRZehiYwD9aAHTtOa+liw7qMh1NwveZb9OhRyJMRRythAf355587ugSRq/r9d3r25Pjx\nwu3lR9wI5lix1WEvgDukkOf8zEzO+G8s/rMYFc6YA9Qq4Lv068ezzxauMHG4EhbQIqZlGAwa\nVIh0dsLemXUhhD7MMldSsEA7sEIf+Aua55HRaTivppsN6wp6pF7lH6+nJ/fcg4sLDRvSqhXu\n7vz9N/fcQ8uWN3Bu4iAKaJEikJrKpElERxeosz+xo5k5mpm1OABQEUZAMNS/3KMeuOeYoPkI\nAeGMmcnoI+Q3HZivL1u3Urv29Z2EmI4CWuQ6/fwzYWHExFCpEvv2sWnTNfo7k9aN1VZsPVjh\nQioW6ARW6A1ulzslwtfwfkY6p+Kygh6hhKyhaxp5rZ6Szfvv8+STuF77XrSUGApokcL5+28m\nT2bpUvbuLeguARwZQ/hoZgZwBMALnoRgyH6p+xfYYC78DXCAWjascxgRi39B3uK553SL+Rak\ngBYpqNRUVq9m+PCCzjbnQmoPVoQQ2pU1zmR7CmUqjL389SX4CkLhJ4AUXL+hlw3rOjrbyWO+\nonLlePNN1q5lyZKMllKlmDiR55+//vMS01JAixTI/v088gg7dxaocy0OjCF8FLP8ic1j8x4A\ndoINIjKmN4qmbhjBcxgRh18+R/70U/r35/HH+eMP/vc/ypenY0e8vAp7NlIyKKBFCmTIkGun\nsyspvfjGiq0z65zyHigHwIcQTvokdEm4f82jYQRvoKPBNR7BbtWK/v0zvr7zTu68sxD1S0mk\ngBbJ7eJFpk9nyxbKlKF3b3r1Ii6On3/Ob5c67A0mbCSz/YgDqA6joRa8CHkOvEtgNw1sWOcy\nPJ4KBanK2/v6F42VEkoBLZLDyZM0bcqxYxkv58yhbFk8PfPu7E5Sb5ZYsXVkgwUDF+gJVuh6\neYmpSJiaY5dEPBbRz4Z1I20LVVjnzprk6LajgBbJsnYtvXtz4UKOxjNn8phEvwG7rdiGM7cC\n8QC1IBhGQZVsnbbAoqxXkTSxYY1g2Gl886nBw4PUVFJTczT6+DBlSuHPR0o4BbQIwJkzfPcd\nI0aQmJhfNw8S+7EomLB2/AjgBg+DFTqTNebiLMyDUNgBcIHSCxkQRvAWWl+zjEqVmDEDd3de\neIFdu3Byws+PXr145RWtcnI7UkDL7ejUKebO5cABatdm+HCWLGHcOC5dym+XJkRasQ0jwpfL\ng+wqwq9QM1unzWCDL0mfUn87zUIJ+YLBZynQPHL/+hcTJlC6NMCDD3LxIh4eWhvwtqaAltvO\nzp0EBXH+fMbLl18mMfGqC02V5kJ/vrRia82W3NuqX07n0zAXbLALIAHv+QwKI3grLQpeVc2a\nPPccpbItFliq1NV7y+1BAS23lz17uPfeHPcxrjZ/7V1ss2IbzBc+nM27x+/QF1zgG0gE2EoL\nG9b5DDpPHiOTq1Rhzx769uW773JvslhYtkyJLLldO6DnzJnTu3fvMmXKXLOniMn9/TdBQde4\ny+xNwkAWhBDanN8APKEd/A9O5tV7McBZfOYxxIZ1O83yOXLv3nh7Ex5O587s2ZNj06uvEhhY\nuHOR28G172+NHDnSz8+vT58+ixYtupT/XToRc/j7bwYMoFw5SpemUiXuvZeHHqJlS/z9OXXq\nqnu1YGsoITFUDSWkOb8RCNMgBlbDirx32UxQ+uOC45iefzoDPXsCVKvGrl2EhdGiBdWr0749\ns2czceJ1n6vc0q65KNb06dPbt2/v5OQEeHl5DR06dPny5cnJycWwHpdJaE3CEuTYMePllw0X\nl0Is/efDmSeYvo1mGa9LY4zC2Jxt6T8D49sc+5zCdxpPNSaygG/h7Gz885+O/tbIVZh5TcKC\nLhp77Nixjz/+ODOpy5UrZ7Va169fn5aWdlPrMwMFdEmxfr3h7V2IaA5i0yxGXqDU5UVVMT7B\nOJMtl9MwVmP0wXDNWJX1B9oNJcKDSwV/l2bNjL17Hf2tkau7FQI6U3pSt2vXLj2pq1Sp8vTT\nT//88883oziTUECb38qVRvv2hrNzgRLTl1NPMW0nTbKaOmL8mvOSOQZjEkbNjA7HqPweExrw\nV8FzGQwnJ2PsWON2+r/NEsnMAV3oURyVK1du06ZNfHz8kSNHDhw4cOzYsWnTpk2bNq1evXpv\nvfVWnz59ivD2i0ieUlP55BNmzeLYMRo25MIFtm699l4WjLZstGLry1ce5Pyg8ANoCkAarIZQ\nWAmpGFg20NGGdQm9k3C/2pE9PXFyyvH8oZcXTz+N1UqNGtd3iiJQ8GF2qampGzdu/Oabb5Yu\nXXro0CGgSpUqY8eOffTRR8uXLz9v3rzQ0NB+/fr9+uuvzZs3v5kFF7GTJ08+/fTTqbmeq81p\n//79gFHAxZmlqMXHs349Z8/SsmXG/G0vvsj772dsjYu79hEqED+cuVZsDdidd49PIAS+hZmk\nT6kfh99sRtqw7uPa60ddusRrr7FwIVFRAK1bExpKkyYFOTmR/Fw7oBcvXvzNN98sX7789OnT\nQO3atZ9//vlHH320VatWlstzt9x9991Dhw69++67Fy9eXLIC2t3d/Y477sg/oM+dOwdYNFGN\nIyxfnmOC/MqVadOGxYsLtK8FoyMbrNh6s8SdpPy6hkIogB2ntXSxYV3GwykUYvEoFxf27CEu\nDldXypUr+H4i+brmTZD0bk2bNn3jjTd27NhxtW5nz56tUKHC1KlTi/IGjDnoHrSjnDhhlC1b\niHu+mX/8OP4ik6Opk3nX2XgKYztGJEblvPeJwX8Sr9bkwHW8HRjz5jn6myXXq2Tfg3733Xd7\n9+59xx135N+tTJkyJ0/mOZRf5Dpt2pTHNHL5cMLemXUhhD7MMldSANqCFfqBx+VOjXNM0JyG\n82q62bCuoEdqwe74de3Khg0kJ2e1VK9O9+6FqFOkgK79E/ncc88VQx0iVyp4OvsTO4pZwYTV\n5CBABRgOVmiQrVM0fALrM14dISCcMTMZfYSAAr6Lnx/PP8+ECWzYwBNP5LjjXLZsQUsVKTjN\nxSEmde4cMTHX6ONM2oOsCiasBytcSAVoB49Db7LGXCTBErDBBjBIxWU5PcMIXk23tIxJ9fNm\nsTBzJlFRVKvGo4/i7Z0xzxxw//264yzFQQEtZhQaypNPkpJy1Q4BHBlD+GhmBqSPukj3MHyT\nrdNusMFc0qfUP0CtGTwWwbBY/AtSQ9++jByZXwe//BZ3FSkCCmgxEcNg2TLWrWP6dPIc1uhC\nak+WBxPWjdXOpOXenAgGJMFXEAobAVJw/YZeNqzr6GwvwOQz6fr1IzT0Rk5FpAgooMUskpLo\n0oWNG/PeWosDYwhPn5noqof4DurA36TPDxpN3TCC5zAijkJc63bqxMKFVCjQOq4iN5cCWhxs\nzRomTmTnThITcy/EB7iS0otvrNg6s84JO87QFTpBGHk/dLKfJNyX0NuGdQMdDQo3er1KFVau\nxP2qzwyKFCsFtBSruXOZNo34eGrUoFo1fv6ZAwfy7lmHvVZsI5ldiRMAATAaxpAx5qIu9Mq9\ny24a2LDOZXg8Bb0ADgujYkVsNi5epF8/xo69vtMSuSkU0FJMLl2iffusSTMOH867mztJvVkS\nQmgHvrdg4ALdIQS6kTXm4jwszNolEY9F9LNh3UjbQpVUpgyjR2Ox8PDDhT0bkeKggJabIiWF\n1as5eJC6denSBYuF+++/xpRGDdhtxTacuRXSR13UhDEwmhxjLrZCGMyHBIBImtiwRjDsNL7X\nUeSMGegBfjEzBbQUjeRk5s5l2zbKlaNdO556it2X7xGXK8fLL7PlijVX03mQ2I9FVmxtyfb5\n4AwIzrbgz1n4AmywDeACpb+kfyghP9Mqn5J8fbMm8XB3p0cPHnmENWvYt4/q1XnySdoW7oJb\npLgpoKUInD/Pfffxxx8ZLy2WHIPkTp3i//4vj72aEGnFNowIX07n2OAEIy6n888QCgvhIsA2\n7rJh/YLBZ/G5ZlXvvkvt2vzxB5Ur8+CDeHkBDBt2Hecn4hgKaCkCL7yQlc6Q9xDmTKW50J8v\nrdhac5WLajv0h5YwHyIBEvCez6AwgrfSooAlNW3KqFFYLLRvX8A9RExHAS03KjWViIgC9byL\nbVZsg/nCJ32gcnW4mPGYX27LYBnAVlqEErKQAQl453NkV9ccjx2WLcuWLbq/LCWeAlpu1Pbt\nnD+fXwdvEgYxP5iwFmwF8IC+YIW2cArugHO5dzmLzzyG2LBec6lswNmZNWuYMYPvv8+4ZP7s\nMzw9r/uERMxCAS15uHSJNWs4dozGjWnXLsemlBQuXMgxedtXX131OC3YasU2iPlenAdoDFYY\nBpkTDLmRa8KizQTZsH5J/4uUKkipFSsybRodO9KxY0G6i5QkCmjJbedOHn6YgwczXt5/P8uW\nUaoUsbE89RTLlpGSQr16vP8+PXrw9tu8807uI/hwdgjzrNiasR2gFPQDK9yXrVMczIVPSf+A\n8DS+EQyzYY2kEEtFzZnDoEG4FmLlE5ESxdErBpQAt9WKKna70bhx7uVCnnnGSEkxWrbM0ejq\naqxdazg55WhszeZZjLxAqYzXFTA+wjiTbbXsNIzvMPphuGKAHcsPtBtKhAeXCruIScWKxu3x\ndyI3V8leUUVuKwcPsmtX7sYVK+jbl19+ydGYksKoUdjtAL6cHkaEFVuT9FEXmf4BT17++hjM\ngjA4AHCSiumXzLtzTKqfHyenjLcDypbl888zRs6J3KoU0JJDnouYnDzJvHl5tB89Sjt+tGLr\ny1ceJObRYx0Mgv+BDZZDKgaWDXQMJWQpjySRe1KiRx+ldm0OH6ZtW+x2nnkmK5Hd3dmwgcTE\njGdhevSgYsUbOlMR81NASw6NG+PpyaVLORrPnOHTT3O0VCB+OHOt2BrkPafcZSuhSsaXcfjN\nZmQYwXupc2VHZ2dmz2bo0ByN99/Pxx9z8CB16vDUU9SpA+jDQLmNKKAlh5QUypTJHdCZLBgd\n2WDF1psl7iQBtIQQaAyPZzyHnYsdp7V0sWFdxsMpZH2c5+FB4uVr7nvuYcMGvK8Y6NyoEZ98\ncsOnJFJiKaBva6tX89FHHD5Mgwa8+CLNm7NhA3FxefSszPGRzB5DeB32ApSFoWCFOy/36JM7\noGPxn8WoMIIPUjPX0dq1Y+lS1q0jJoZGjTJmUxKRXBTQtx27HScngJkzGTMmozEykiVLmDqV\ntWtzdHbC3oW1VmwPs8yVFIA2YIV+kPkkSAp8A5cvddNwXk03G9YV9Ei9yg9YmTL4+tKvX1Gf\nm8itpWQHtN1u37dvX0pKSr169VxcSva53GxnzvDqqyxYQEIC99zDW2/xzDM5Ov5LkeYAABj9\nSURBVKSlMWFCRnYD/sSOYlYwYTU5CFAOhkMINMy2zz4Ig9lwHOAIAeGMCWfMUarlX8zOnUV0\nViK3tBITaq+99lqtWrVGjx6d/jI1NfX999+fNGnS+fPnAXd39xEjRkyZMsXH59qTnN2GDINB\ng1i9OuPlli107kzaFWuuAhZ7Wg9WW7H1YIULl1egKg+7yVqlJBmWQiisB4NUXFbQI5SQNXRN\ny/lcoKsrZcty8mTud7nydrOI5MHRA7ELCujQoUPmy/HjxwO+vr59+/Z97LHHWrVqBTRq1Cgx\nMbHI3/oWeFDl99+v/dxHAIff4J9HqJbnNiMZw8DYjfE8RsWM9v3UeoU3/YnJ84AdOhgXLhhb\ntuSxadIkR39HRC4z84MqJTKgIyMjLRbLvffeGx8fn9khPDwceP3114v8rUt6QMfEGF27XjWX\nXUjpxdLl9EjFOb/8DsRojWHBgCTcFtH3AdY4kZb3MV2MVauyCnjiiRxbmzc3UlIc9+0QycnM\nAV1ibnFkt3nzZsMwPvzww/Lly2c2jh49eubMmatWrZo4cWLBD5WWlrZ8+fLk5OR8+vz+++/X\nX6ujHTzI3XdnLSySXS0OjCF8NDOrcAygCjwCW0ifPyO3nQDR1A0jeA4j4vC7souTE0FBBAXx\n+uuULp3VPn06ffrwzTecP0+HDgwejLPzlXuLSG4lMqCPHj0KBAYG5moPDAycl+cTb/ke6vHH\nH09MzOspuMvSf8Ea+c9Cb0opKYwfnzud3Uh+mGVWbJ1Z54QdJ+gKVugJrnAs5xqAACThvoTe\noYR8TweDvAfEeXoSHs6gQXlX0qkTnTrd8PmI3GZKZEDXrVsXOHToUOPGjbO3Hz9+vGbNmoU6\nVI0aNWJjY/PvM2PGjLFjx1pKwkjd1FT27OHiRRo35scfCQ4mJiZra12igwkbwRw/4gCqwWgY\nA9WzHeKHHAf8i4bpq7LGZ31EmLcvvuCRR4rqPEQESlZAR0VFvfnmmw0aNKhfv37FihXffPPN\n+fPnZ27dunXrihUrht3GS879/DMjRhAVBeDjQ1JSxqN67iQ9ytdWbB343oKBM3QHK3TPNhfz\nefgSZsCvAJfw/Iq+Nqwbyb2uqrNzHsM/Onbk4Ydv4qmJ3KYcfRO8oAICAq68hl2/fn361pde\nesnT07N8+fKHDx8u8rcuER8SxscbVark/rCuIX++x4STVMhq6oVxNNvknwbG7xhjMcpkdPiD\nwPH8x5dTeX76V66cERFhvP66UbOm4e5ueHgYAQHG888bZ886+vxFrpc+JCwChw8fvnTp0t69\ne6OioqKjo6Ojo6OiojIfTlm6dGmFChUiIiICAgIcW2cxS0xk9GgWLyb7x5yeXOrLV1ZsbdmY\ne4cJUBWAc/AF2OB/ABco/SX9bVi30BpwdYWUHPuFhNC3L/feS/pA88J8ECsi16nEBDTg6ekZ\nGBh45WeDwOLFixs0aOCU+RjcLWrbNnbvJiCAoKCMR/4eeICN2UK4CZEhhA7lc1/yGrcBvAGP\nwXewEC4AbOMuG9Z5DDlHmfQuPj4cPsxvv/H22/z1F9WrM24cgwdrugyR4laSAjofjRo1cnQJ\nN9eFC/Trx6pVGS/vuYclS1i/PiOdS3OhP19asbVmC4AzuJLn/MxsgA0ACXgvYGAoIb/RPFcX\nX1/KlNG4CxHHu0UC+pb3wgtZ6Qz8/juDBrF9O83YHkLoYL7w4SxAPQiGEVAKgjIGL+eylRY2\nrPMZdJ681yPp1u1mnIGIFJoCumRYsiTHS28SGm1a8AG2FmwF8IBHwQrtyRqmHJAjoM/is7nW\nkH+fsG660CyfN2rdmilTirZ2EblOCugSIDmZ+PiMr+/lVyu2gSzw4jxAQwiB4VDucm8DfoBP\nYGVGw2aCbFgXWfr/trLU47+zaWju4wMjRxIYSKNGPPAAt/qdfJESQwFtdqdOcf/9lEo5O5gv\nQghtlv4gtgUGwePQJlvXkzAHbBAFcBrfCIaFErKLxgAGrq4MGUKVKrz0Er//nrXcX58+hIcr\nl0VMRwFtdp8M3fz0dlt/vizFxazWnpD5TLsd1oMNlkIywEbahhLyFX0T8ch+qN9+o3ZtOnXi\n1185coSVKzlzhnvv1Sp/IialgDaXQ4c4coR69ajkepqICGy2VyMj8+h3EBLgAsyGMNgHEE+F\nOYwII3g3DSwWrpw75MSJrK8DAnjssZt0EiJSNBTQjnfsGFFReHnxxhusWG60ZWOIxdbf6SvX\ntKtP4bQTKkIK2DGwfO/UaY5r8CqP3vUC3f8RQo0aTJvG11/n3skvjxnoRMS8FNCOlJrKuHGE\nhWG3U4H44cx9l9D67MGAUjAIusKnsC6vnZOIw282I8MI3mfUObKP2VWzNp49mzugy5WjQ4eb\neTIiUtQU0I40aRK2UKMT64MJ680Sd5IAWkAIDCRjmHKF3AFtx2ktXWxYl/FwCq4ABtu2UTVb\nQD/0EJMmMWlSxiPgfn7MnUulSsVzWiJSNBTQjrFjBxFTj5f6cnYU4XXYC1AWhoAVmmbrdwiy\njUo+ZvEPN0aHM+YgNXMd8MrwffVVRo7k11/x8qJ1ay0DKFLyKKCLnd2+9a21R18Pfdv41jV9\nRqIgsEJ/KHW5Typ8CzZYA3bScI67q5v/P62e9/VwDXepuJiDW3Mc8o47aNqUK1WrRrVrrK8t\nIualgC5GMTHMmkV4eIuDB1tkNi6GR7P1OQBhMIv0VaiOEDCT0eGMSTkWcKwXZeHFF3nhBZ54\ngtDQjIHMtWvz5Ze4uxfvuYjIzaeAvvnS0li1CpuNlStJTc29tTUAyfAN2OC/YCcVlxX0sGFd\nTbe09En1jxMTk3GX2WLh00+ZMIEdOyhfnqAgpbPIrUkBfTMdPszMmYSHc/ToVfs8CC1gGZwA\n2M8d8zzGzHMbtedcley9nJ0pWzbHfnXrUrfuzShaRMxCAX0TpKayfDk2G2vWkJaGBdrAOfgj\nr847YAfJuC3jYRvW38p2nhPhtH8xe2bn6NWlS451skXkdqCALlIHDhAWxqxZHDsG4AcjIRjq\nQArUhiO594im7kr/4Ls+HHEy3m9kWT7vTMWKtG3LwYN8/31Gn3vuITy8GM9CRMxBAV0UDINl\ny/jkE9atw27HCbqCFR4mfZgywOmMFUwyuLvH3Nv75ztDPB/s8GQ3i7Mz7bJt9PFhwwZ++YWo\nKGrUoE0bzWQkcjtSQBeF4GBmzgTwh1EwBmpl27oTbBABZwBo2JDgYIYPr1qhQp98j9qyJS1b\n3rSaRcT0FNA37KefmDmTavAR9Mz2Hb0AC8EGPwPg4cGwflittG3rsFJFpEQpqQF97ty5hIQE\nJycnPz8/B68Vu20bwDh45HLL/8AGX8A5AJo0wWpl2DB8fR1UooiUSCXs1mZkZOSIESOqVKni\n4+NTrVo1f39/Nze3atWqDR48eNOmTY6pqU4dgIWwCWZAC7gHPoO00owaxebN7NzJU08pnUWk\nsErSFfT48eOnT59uGEaVKlVatmxZvnx54NSpU0ePHp0/f/78+fODg4NtNltxl9W16+H6Xapv\nX5u5uElaYDPnx0MYPBgfn+IuRkRuISUmoD/55JOPP/64a9eub7/99l133ZVr665duyZNmhQW\nFtawYcMJEyYUZ2FzP3cK2fPtWD5rzZZY/OczyK9Gi28fL84SROTWZDGuXHjDlO67776///47\nMjLSxSXvXyqGYbRv395ut//0009F+9YzZswYO3ZsQkKCl5fXlVvbtOHKmytHj+aY/FNETCs5\nOdnd3X3Tpk1BQUGOriW3EnMPOjIyslWrVldLZ8BisbRt2zYyzwWirm7//v2enp6WfI0dOzb9\n+HkeYd++gjaKiBRKibnF0aRJk19++SUtLc3Z2flqfbZs2dKkSZNCHbZWrVpr1qxJTp/W/ip2\n7dr1zDPPuLq65rm1QQOOH8+jUUTkBpWYgB4yZMi4ceMeeuihd955JzAwMNfWqKioiRMnbtiw\nYcqUKXnufjUWi6Vdu3b59ylVqlQ+W59/PuuZ7HQjRmjtEhEpAiUmoJ944omdO3d+9tlnq1at\nCggIqFGjRrly5SwWy+nTp48cOXLgwAFg5MiRzz//fDEX1qMHX37Ja68RFYWvL2PG8MYbxVyC\niNyaSsyHhOm2b98+ZcqUtWvXxsfHp7c4OztXqlSpQ4cOjz32WPv27W/Gm27evPm+++5LSkpy\nc3PLp1tSkuZlFil5zPwhYYm5gk7XrFmzL774Ajhz5kxCQoKrq2ulSpVu9pOE6bnsrvQVuXXl\nf/nlKCXsCtpRduzYkXrlYih5efbZZz09PYcNG3azSyqIxMREq9X6r3/9q1atWtfuffOtW7du\n3bp1kydPdnQhGSZOnHj33Xc/9NBDji4E4Pjx4//3f//3n//8x9ccD50uWrRo//79L774oqML\nyfD000+PHz/+Jv1lubi4NM1zWU+HM6RIPfroo0899ZSjq8iQkJAAbN261dGFZJg+fXqjRo0c\nXUWWVq1avf32246uIsOePXuAmJgYRxeS4ZVXXunSpYujq8gSEBAwd+5cR1dR3ErMOGgRkduN\nAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAK6iLm5uZnnoX4X\nFxcnJyfz1GOqbw4mq8fNzc1isZinHldXV/MUg8n+soqN5uIoYvHx8W5ubmXKlHF0IRn2799/\nxx13OLqKDMnJySdOnKhWrZqjC8lw7NixsmXLenp6OrqQDKb6yzp//vzFixcrmWZq88OHD/v7\n++ezptItSQEtImJSusUhImJSCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSIiEkpoEVE\nTEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQX0LW7v3r0ff/yxo6uQgjp//vyc\nOXOOHj3q6ELEFBTQRebTTz9t06ZN2bJl27Rp8+mnnzq6nAwfffTRa6+95ugqSEpKeuWVV9q1\na+fj41O7du3Bgwfv27fPUcUcOHBg8ODBdevWLV26dGBg4AsvvHD27FlHFZPL+PHjR44cuWPH\nDgfWEBAQYLmCA3+KNm7c2LlzZx8fH39//wEDBjjwJ8cBDCkKY8eOBerXrz98+PB69eoBTz75\npKOLMr777jt3d/eyZcs6towzZ860bdsWaNSoUXBw8AMPPGCxWDw9Pbdt21b8xURHR5cuXdrF\nxaVTp05jx45t2bIl0Lhx40uXLhV/MbksWrQo/V/l8uXLHVXDxYsXLRaLv79/h5zCw8MdUs+C\nBQvc3Nz8/f0HDx7cq1cvZ2fn8uXLHzp0yCHFFD8FdBHYtm0b0K1bt5SUFMMwUlJS0jNo586d\njippyJAh9evXT//X7vCAfvnll4Fx48ZltqxYscLJyalp06bFX0yfPn0sFsuyZcsyW5599lng\no48+Kv5isjt69Gi5cuW8vLwcG9B//PEHMGnSJEcVkN2hQ4dcXFxatmx55syZ9BabzQaMGDHC\noXUVH93iKAJTpkwB3nnnnfQF01xcXN5++23DMKZOneqoki5evFi3bt2ePXt6e3s7qoZMS5Ys\n8fb2fu+99zJbunfv3qlTpx07dpw4caKYi/npp5/uvvvuhx56KLNl1KhRwP/+979iriQ7wzCG\nDx/u4+Pz1FNPObAMICoqCmjQoIFjy0g3bdq01NTUDz/80MfHJ71lzJgxH3zwQatWrRxbWLG5\nvVZgvEnWrl1brVq1O++8M7Pl7rvvrlKlynfffeeokr7++uv0LwIDAx3+iZOTk1P79u3d3d2z\nN6av0Hz69OniXJbUbre/9tprNWrUyN4YFxcH1KlTp9jKuNJ77733/fff//DDD5s2bXJgGUB0\ndDRQo0aNefPmRUdHV6tWLSgoqFGjRg4pZsGCBQEBAdnj2GKxPPPMMw4pxjEcfQlf4p0+fRq4\n7777crWn39w8d+6cQ6rK1KRJE4ff4rjSiRMnPDw8/Pz80m8KOcTFixdjYmJWrlxZt25dPz+/\nqKgoR1Wybds2Nze3l19+2TCMyZMn49BbHKNHjwYqVqyYGRFOTk7jx48v/r+phIQEoG3bttu3\nb3/ooYcqVaoUEBDQt2/f6OjoYq7EgXSL40al/xiVL18+V3t6y7lz5xxQk7lFRUUFBQUlJiZO\nnjw5/aaQQ0yYMKFq1ardu3ePjY1Nj2mHlHHp0qUhQ4Y0atTojTfecEgBuaTf4rj//vv/+OOP\nhISEn3766Z577vnoo4/ef//9Yq7kzJkzQGxsbJs2bQ4ePNizZ8/GjRt//fXXTZs2/e2334q5\nGIdx9G+IEu/YsWNAr169crV3794diI2NdUhVmUx1BX3+/PnXX3/d09PTw8Pj448/dmwx27dv\nX7hw4b///e/q1au7u7svXbrUIWWMGzfOw8MjMjIy/aXDr6B//PHH9evXZ285efKkr6+vl5dX\nWlpacVayf//+9Ix66aWX7HZ7euPatWstFstdd91VnJU4kAL6RqWlpTk7O7dr1y5Xe6tWrZyd\nnYv5Z/pK5gnolStXVq9eHejZs+fu3bsdXU6WmJgYb2/vqlWrFv9br1u3Dvjggw8yWxwe0Hnq\n27cvUMx3gY4fPw6UL18+NTU1e/sDDzwAxMXFFWcxjqJbHDfKycmpUqVKV34QFxMTU7lyZScn\nfYcB/vnPf3bv3t3b2/uHH3749ttvM4cAFrN9+/bNmDEjMjIye6O/v3/z5s1jYmLSP04oTtu3\nbweeffbZzOdBXnrpJaBnz54WiyU8PLyY67ma9Pt1KSkpxfmmFStW9PDwqFWrlrOzc/b2O+64\nA3D4R9/FQ6M4ikCHDh3mz58fFRWV/ogKsGvXriNHjgwaNMixhZnEnDlz/vWvfw0cOHDOnDnp\ngzccJS4ubuzYsU899dS0adOyt588edLLyytzLFexadq0afojTpm2bdv2yy+/PPjggzVq1Cj+\nsW5//vlnnz59evfu/dZbb2Vv37Fjh7u7e+aPd/FwcnLq0KHD5s2bExMTPTw8Mtv/+usvJycn\nR/2OL26OvoS/FXz//ffA0KFD01/a7fYBAwYAGzdudGxhhglucdjt9vr161etWtUMj+olJydX\nqlTJx8dn3759mY0LFiwgr08RHMKxtzjS0tICAgI8PT1//fXXzMb0C/mQkJDir2fNmjXAuHHj\nMm8VLly4EOjZs2fxF+MQuoIuAu3btx85cuTs2bNjY2NbtWr1008//fjjj2PGjGnTpo2jS3O8\nQ4cO7dmzp2LFir17975ya0RERIUKFYqtGFdX148++mjgwIGBgYHdu3evVKnSX3/9tWHDBj8/\nv+nTpxdbGabl5OQUERHx6KOP3nfffT179qxcufIff/yxadOmhg0bvvPOO8VfzwMPPDBy5Mjp\n06f/+OOPrVu3PnDgwNq1a6tUqWKeuW5uOkf/hrhF2O32d955JygoqEyZMkFBQVOnTnV0RRkc\nfgX93//+N58fv6NHjxZ/SevXr+/WrVv58uVLlSrVtGnTCRMmnDp1qvjLyJMZPiQ8dOjQqFGj\nmjRp4uXl1bx589dee82x//fz7rvvtmnTxtvbu1GjRk8++aR5/rKKgcUwjOL6XSAiIoWgMQYi\nIialgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSI\niEkpoEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0i\nYlIKaBERk1JAi4iYlAJaRMSkFNAiIialgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuI\nmJQCWkTEpBTQIiImpYAWETEpBbSIiEkpoEVETEoBLQKwa9cud3f3jh07ZrakpKQEBgaWL1/+\n+PHjDixMbmcKaBGAxo0b/+Mf//j+++9nzZqV3jJ16tTIyMj//Oc/lStXdmxtctuyGIbh6BpE\nTCE5Ofmee+6JjY3dvXv32bNnAwMDu3TpsmzZMkfXJbcvBbRIll9++SUoKGjgwIHHjx/ftm3b\nrl27qlSp4uii5Pbl4ugCREykZcuWTz/99AcffADMnTtX6SyOpStokRz27t1bt27d0qVLx8bG\nlilTxtHlyG1NHxKK5PDcc8+5ublduHDh5ZdfdnQtcrtTQItkmTdv3rJlyyZPnty3b99PP/10\n8+bNjq5Ibmu6xSGSIS4urnHjxjVr1vzll1/i4uIaNmxYrVq1bdu2ubm5Obo0uU3pClokwxNP\nPHHmzJnQ0FBnZ2d/f/+33nrrzz//nDx5sqPrktuXrqBFABYuXDhw4MDnnnvu3XffTW+x2+2t\nW7fesWPHtm3bGjZs6Njy5PakgBYRMSnd4hARMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuI\nmJQCWkTEpBTQIiImpYAWETEpBbSIiEkpoEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIi\nJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAJaRMSkFNAiIialgBYRMSkFtIiI\nSSmgRURMSgEtImJSCmgREZNSQIuImNT/A7gJwljddLTOAAAAAElFTkSuQmCC", "text/plain": [ "plot without title" ] }, "metadata": { "image/png": { "height": 180, "width": 240 } }, "output_type": "display_data" } ], "source": [ "myLM <- lm(y~x)\n", "myLMCoef <- myLM$coefficients\n", "yHatOLS <- myLMCoef[1] + myLMCoef[2]*x\n", "\n", "plot(x, y, pch = 20, col = \"blue\")\n", "points(sort(x), yHat[order(x)], type = \"l\", col = \"red\", lwd = 5)\n", "points(sort(x), yHatOLS[order(x)], type = \"l\", lty = \"dashed\",\n", " col = \"yellow\", lwd = 2, pch = 20)" ] }, { "cell_type": "markdown", "id": "c51a5cb6", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Loss function\n", "\n", "- Every statistical learning algorithm is trained to learn a prediction. \n", "\n", "- These predictions should be as close as possible to the label value (ground-truth value). \n", "\n", "- The loss function measures how near or far are these predicted values compared to original label values.\n", "\n", "- Loss functions are also referred to as error functions as it gives an idea about prediction error.\n", "\n", "- Loss Functions, in simple terms, are nothing but an equation that gives the error between the actual value and the predicted value.\n", "\n" ] }, { "cell_type": "markdown", "id": "ffc72b09", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Mean squared error loss\n", "\n", "- Mean squared error (MSE) can be computed by taking the actual value and predicted value as the inputs and returning the error via the below equation (mean squared error equation).\n", "$$Loss(\\beta_0, \\beta_1; y, x) = \\frac{1}{N}\\sum_{i=1}^{N}(y_i-(\\beta_0+\\beta_1 x_i))^2$$\n", "\n", "- The equation is very simple and thus can be implemented easily as a computer program. Besides that, it is powerful enough to solve complex problems.\n", "\n", "- The equation we have is differentiable hence the optimization becomes easy. This is one of the reasons for adopting MSE widely" ] }, { "cell_type": "code", "execution_count": 12, "id": "3f6a8d35", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "loss <- function(par, y, x)\n", " {\n", " beta0 <- par[1]\n", " beta1 <- par[2]\n", " yHat <- beta0 + x*beta1\n", " out <- mean((y- yHat)^2)\n", " return(out)\n", " }\n", "\n", "optimLossOut <- optim(c(0.2, 0.3), loss, x = x, y = y)\n", "\n", "beta0Hat <- optimLossOut$par[1]\n", "beta1Hat <- optimLossOut$par[2]\n", "yHatLoss <- beta0Hat + beta1Hat*x" ] }, { "cell_type": "code", "execution_count": 13, "id": "9b1191e7", "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAFoCAIAAAAAVb93AAAACXBIWXMAABJ0AAASdAHeZh94\nAAAgAElEQVR4nO3dd3hUVf7H8fekBxJCCC1AKNJLBBUEIl0QpClNes9EVkQU/VnWsrK4iuC6\nsopKJkGKCIgI0hEWXBBQ0QUkCIReEgggLZTUub8/EtIYQoAkcwOf18PzOHPuuXe+N8CH651z\nz7EYhoGIiJiPi7MLEBERxxTQIiImpYAWETEpBbSIiEkpoEVETEoBLSJiUgpoERGTUkCLiJiU\nAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAJaRMSkFNAiIial\ngBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSIiEkp\noEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIK\naBERk1JAi4iYlAJaRMSkFNAiIialgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJQC\nWkTEpBTQIiImpYAWETEpBbSIiEkpoEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIiJqWA\nFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAJaRMSkFNAiIialgBYRMSkFtIiISSmg\nRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSIiEkpoEVETEoBLSJiUgpo\nERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAJa\nRMSkFNAiIialgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJSbswsoGnbs2JGSkuLs\nKkSkQLi5uTVs2NDZVTiggL65X3/9tUmTJs6uQkQK0NatWxs3buzsKnJSQN9cUlISkJiY6OHh\n4exaROTWJCbi6Zlbh6SkJE9Pz7S/5maje9Aiche6coWXXyYgAG9v6tRhwQJnF3RbFNAichd6\n5hkmT+bsWQyDvXt56imWL3d2TbdOAS0iRdXly0RG8sYbRERw6VJm+6lTzJwJ4MOlqhx2wQ58\n8IGTqrwDugctIkXS/v20aUNMTPrbt9/m229Zs4a9e8Fu78gaK7buLHEnOYaKI5j++57HnFrv\n7VBAi0iRZLVmpjMQE0OLFpRJjh3B9JFEVuUwQFmoTcUfYxYavQYF7YUKTir2NukWh4gUPQkJ\nbNqU+daV1C4s/yb5iSNUmcCbVV0O8xgsgGOwAcbiw6VXGhS9m9D3+hV0UlLSV199lfsIm+jo\n6EKrR0TyIjUVux2gMkdHMH0E04M4BhAIwyEUqmXpfQGgecuiF3dFr+L8FRcXN2nSpMTExFz6\nxMfHA0lJSRoHLWISxT1Tnq+6rN0BW0dWu5KKK3QEK3TNkmqX4WuYBj+Dvz/dujmz4ttyrwd0\nUFDQH3/8kXufadOmjRo1qnDqEZGbOHSIiAi++OKDEyfSW8bA/0FQlj7bwAZfpV87U68ekZGU\nLl3otd6pez2gRaRoSEpiyRJsNtauTb+7kaYC/Pva63iYB+HwK4Dh4Wnp1wOrlbZtsVgKv+Q7\np4AWEXPbt4+ICGbM4NQpB1tPwmyoCrNhLlwCOF++juezVu+nhxTFq+asFNAiYgobNvDXv/Lb\nb5QuzYABvPVqUvG13xEezn/+g2HccDc7DLn22tOTPt0JCyv56KNF9JI5BwW0iDjfb7/RoQNp\nw6l8j+8uNynCPmUWiWcAasBI6AebYNAN9m/QAKuVwYPx9y+0mguBAlpEnG/yZFySEgazwIqt\nJRsBDOgLVmgHaVfD5cAVUrPsVqwYfftitdK8uTOqLnAKaBFxtp07Oyy1fcaX/pwDqA2hMBTK\nZOmzAf6RJZ0bNSIsjAED8PMr9HILjwJaRJwkIYGlSwkPZ+3akYDbtUvmVtcumYEzMAtssAcA\nb2+6diUsjPbtnVV1YVJAi0j+S0khLo7AQFwcTiexbRs2G199xYULmY1T4Jlrrw1YDzZYBGmP\nkTVpgtVKv374+hZw7SaigBaR/JSQwOuv8+mnJCTg48OLL/LWW9diOj6eefMID+fXXx3smbbq\n5ymYARGwDwA/P0YOxGqlUaNCOgEzUUCLSH565RX+fe3JkUuXGD+euXNZ+97WoFU25s0jPv6G\ne46DSNgDSQBnazcv9YqVvn0pVqww6jYlBbSI3BG7PfM+xuXLfPZZ5iY/LgxkjjXaFtRrO4AP\njIAK8E+4et2BUuF3zuE/m8FLylnX7G7A3TCU+Y5oulERuU0LFxIcjIcHFSsyfjznz9O0KcnJ\nAM3ZMp0RsVSYyuhGbKcxTINYiIQJMNDB0TbQajCzKxA7linVujW4Kx40uVO6ghaR27F0Kb17\np7+OjeXtt1mxgthd555jthVbA6IA/GAAWOGBLHsehg2Z71L9S398Ycg0u3UPdTIaO3Qo+BMo\nChTQInI7Jk7MfG3BaMlG6y+23nzjRQJAUwiDvlD8WqdUWAHhsBJSwWKhbVusVtcePexTPfe9\nnDnGeeRInnqqUM/FtBTQInI7du8GKMPpIcwKJaJO+kBlCICV0CRL1yMQCdMhBuCMazlGDC39\naig1a6ZtHzeOjh1ZsYKkJFq1omXLwjwPU1NAi8itM4wnfdc9di6iB4s8yb7eRZNr6ZwMS8AG\na8COHZe1tLdhXZLavfZWj99rZtupfn3q1y+06osMBbSI3IqTJ5k50x4eMf3ofscd1sAbkASz\n4SRAnGuFCIZHMvLQtXWodu4kJoaKFQur5iJLAS0imfbsYelS4uMpWZKgIOrVu3Zha7ezZg02\nm33xEpfUZBcL1IJDkHzdIVLhH2n/dV3J49/4W5fZO/95IWfUXLpU4OdyF1BAi9zrli5l8mT2\n78fbmyNHSE3NttXaJfazJtNdZ0Ry+DDgUgaGQSjUgm+gj4MDHiNoOiMiGXmMoFIWGjfl+++z\ndShVKuP+s+SmaAe03W4/cOBAcnJyrVq13NyK9rmI3LaEBPbto0wZype/5X3nz6dfPwftrqR2\nYpUVW5fly12Xp+AC7cEKT0LG4snZZ8VIwW05XWxYV9EpFde0xqZN+eADmjXLfITQ1ZXPP7/B\nHB2Sg1FEvPHGG5GRkRlvk5OT33//fR8fn7Sz8PT0DAsLO3/+fEF89Oeffw7Ex8cXxMFF7tA/\n/2n4+BhggNGhg3H06M13SUkxoqKMn34y4uON6tXT9834FcTRt/nbUYLS3wdi/BXjAIaR5dcu\njOcxfNP3OUi113mnAjE5DuXlZezcaRiGcfSoMW6c8fjjhtVq/PprQf9Ibk1iYiKwadMmZxfi\nQJEJaKBNmzYZb8eMGQP4+/v37t376aefbtasGVCvXr2EhIR8/2gFtJjW3Lk547V5cyMlJbND\nQoKRmJhtl61bjbp10zuXKJG5oxvJT7B4GV1ScE1vegBjIUZSlly+jDED45H0fZJwX0DvjpbV\nJUuk5iijWDFj0CBj9+5C/nncDgV0Psga0FFRURaL5eGHHz5z5kxGh8jISOCtt97K949WQItp\nPfZYzoAGIyrKMAxjxw6jVSvDxcWwWAw/P2P4cOPMGSMqKlsop/2qxsF3eD2GCjk3RGeJ5h0Y\nz2KUTN8UTc1XmFiOk2C0a2fs3GkEBGTuV7p0eg1FgpkDukjet928ebNhGB999FFAQEBG44gR\nI6ZPn75y5crx48c7sTaRAmW3M2cOa9Zgt9O+PYcOOehz6BClS/PYY8TFpbdcuMAXX7BwIfHx\nmeuvupP8BN9ZsbVnrQt2BwdaCaVgEdjgF4BEPBfRw4Z1PW2Na1MZvfQSDRoQHc2MGRw6xH33\nMXQopUrl95nfk4pkQB8/fhwIDg7O0R4cHDxnzhxnVCRSGAyDvn355pv0t3Pm4OrqoNvHHxMY\nmJnOGS5eTH9Rg/2hRAxjRjmu65TVWBib/nI3dSMIncWQM5R2c8Pbg6tXqV2bCRN4/HGAUqUY\nN+52T0xuoEgGdM2aNYEjR47Uz/7s0cmTJ6tWreqcmkQK3rJlmemcJseQuDQ5xrRl8CSxB4us\n2Nqy3oKBG3SDYZAIQx3N/wkJeC2gjw3rRloCHh60CWHSJJo0ITERT887PCG5iaIU0NHR0e+8\n806dOnVq165dpkyZd955Z+7cuRlbt27dunz58sGDBzuxQpECtXnzbe5Yhz1WbEOYVZozANVh\nJAyHjGF5H8JP2XbZSbAN65cMOoc/4OvLu+/y7LOZHZTOhaDIBHRQUNDx48fffPPNjJZ58+aF\nhYW1bdsWeO2116ZMmVKiRIm3337baSWK5B+7nf37OXUKT08SE2nQgJIlb3nssBcJfVhgxdaS\njQAe8ASEwaNZVmU9n3mLGbhCsa95KpywLTQHnnyS997jwgXq1bun1gI0iyIT0EePHr169er+\n/fujo6P37du3b9++6OjojIdTFi9eXLp06dmzZwcFBd3qkS9fvpyUlJRLhytXrtxm0SK3Ze9e\n+vVj+/bMFi8vQkL48ce8HqEBUVZsg5ntzzmAWhAKw6BMlk6bIBwWpN/c2MYDEYTOYeAF/NK2\nlypFRARZvomXwmYxMr7TLcr++OOPOnXquNz6w0kHDhyoVauW3e7oK+zsLl686KtLCCl4ly5R\nuTLnzt3OvsW40pf5VmzN2ZLZ2gW+g4yvE/+E2WCDPwDi8Z1L/whCt2abIZRatfjyS5pka7s7\nJSUleXp6btq0KSQkxNm15FRkrqBzV69evdvbsXr16tu3b8/9Cvrbb7999913LVqBRwrF6NG3\nk86N2B5G+AC+8uNCzm0VwRUM+C/Y4FvSptTfShMb1nn0iyfnlUdQELt362ls57tLAvpOXD9c\nL4dfHS4RL1Iwli27hc4+XOrP3DDCG3PjP6URcBQOQjTABfzmMNCGdTuNMrqUKJE5CM/Hhzlz\nlM6moIAWcQ7D4Pr/K0tO5vz5PO3ehK1WbP2Y50s8QDCEwHy4fnc7rALYQvNwwr7mqSsUy7q9\nZUsWLWL2bPbuJSiIoUM1U7NZKKBFCtXVq/z978ycyZkzNGjAO+/QuXP6Jrudvn3J/QsRPy4M\nZI4VWyO2AxSHfmCFpgA0gr/k3OUc/rMYEkFoFA3SWqpX5+xZ3N1p1Ii+fRkyBDc3nn8+P09T\n8kXRCOhPPvkk6wC73J27va9XRArFqFHMmpX+ets2nniCNWto04ajR+nUKX2hP4eas8WK7Sm+\nLs5lgAchDPpDiWs97LA12y4baRlO2Df0TsAro7FXLxYscHDxLiZUNAK6U6dO+/btmzZtWmJi\noq+vb5UqVZxdkcjtiInJTOc0KSmMHUvLlnz+uePHAv05N5jZVmwNiAIoAf3BCg9lPS5Mh0g4\nAnCG0rMYYsO6hzrZDuXPRx8xaJDSuehw9mxNt2DVqlVA165dC/lzNZud3InTp40nnzSKFTPc\n3Ixy5RxMPnejXy3ZMJtBV/HKbPoHRnyWGeZSMJZgdMNwxQA7lv/Qrh9zPUlweMAHHnD2z8KU\nNJtd/ujYsWOtWrWcXYXILTh+nPvvzxw2d/0ERtcrzZkhzLJiq8OebBv84K/XXh+FSJgOxwHi\nKDeDYRGE7qdGLkf287vl+sW5ilJAA02aNNFzfWJmFy+yZAmxscTGsnIl+/aRx0fBLBhtWW/F\n1oNFniQ66HEBXoMGMAdWgx07LmvoYMO6hO7JuN/0I7p0ucWTEWcrYgH95ZdfOrsEkRv67Te6\nduXkyVvbqxxxQ5lpxVaD/QCekIzD+ZmZmP7fWCp8wfBIRh6iWh4/pU8fXnjh1goTpytiAS1i\nWoZB//63kM4u2NuzNozw7ixxJxkLtAIr9ILd0NhBRqfiuopONqzL6ZJyg7+83t489BBubtSt\nS7NmeHry55889BBNm97BuYmTKKBF8kFKChMmsG9fnjpXIHYYM0KJqMYhgDIwFEKh9rUetcAz\n2wTNxwiazohIRh4jt+nA/P3ZupXq1W/vJMR0FNAit+mnn4iIICaGsmU5cIBNm27S35XUTqyy\nYuvCcjdSsEA7sEIP8LjWKQG+hQ/T0zkFt+V0sWFdRadUHK2eksWHH/Lss7jf/F60FBkKaJFb\n8+efTJzI4sXs35/XXYI4NoLpI4kM4hiADzwLoZD1Unc32GAW/AlwiGo2rDMZGkuFvHzEiy/q\nFvNdSAEtklcpKaxaxZAheZ1tzo2ULiwPI7wjq13J8hTKZBh17fVV+AbC4UeAZNy/44kIQtfQ\nwY6D+YpKleKdd1izhkWL0luKFWP8eF566fbPS0xLAS2SJwcP8uST7NyZp87VODSSyOF8UYFY\nB5v3ArATbDA7fXqj/dRIu2SOo1wuR/7sM556ir/8hd9/53//IyCAtm3x8bnVs5GiQQEtkicD\nB948nd1JfoLvrNjas9bF8UA5AD6CSNImoUvEcxE9bFjX09bgJo9gN2vGU0+lv77/fu6//xbq\nl6JIAS2S05UrTJ3Kli2UKEGPHjzxBHFx/PRTbrvUYH8oEcOYUY44gMowAqrBK+Bw4F08e6gT\nQehMhp6hdF6q8vW9/UVjpYhSQItkc/o0DRty4kT625kzKVkSb2/HnT1JfJLFYYS3Zb0FAzfo\nClboeG2JqSiYnG2XBLwW0MeGdSMtb6mw9u01ydE9RwEtkmnNGnr04PLlbI3nzzuYRL8Oe0KJ\nGMrM0pwBqAahMBwCs3TaAgsy30XRwIZ1NoPP4Z9LDV5epKSQkpKt0c+PSZNu/XykiFNAiwCc\nP8/33zN0KAkJuXXzIqE331ixtWSjBQMP6A5WaE/mmIsLMAfCYQfAZYp/zVM2rFtoftMyypZl\n2jQ8PXn5ZXbtwsWFcuV44glef12rnNyLFNByLzp7llmzOHSI6tUZMoRFixg9mqtXc9ulAVFW\nbIOZ7c+1QXZl4BeomqXTZrDB13AFYDuNbFjnMPACeZpH7u9/Z9w4ihcHePxxrlzBy0trA97T\nFNByz9m5k5AQLl1Kf/vaayQk3HChqWJc6cv8UCJCuO4busrX0vkczAIb7AKIx3ce/WxYt9Ik\n71VVrcqLL1Isy2KBxYrduLfcGxTQcm/Zu5eHH852H+NG89c2YrsV20Dm+HHBcY/foDe4wXeQ\nALCVJjas8+gXj+/13QMD2buX3r35/vucmywWlixRIktONw/omTNn9ujRo0SJEjftKWJyf/5J\nSMhN7jL7cKk/c63YmqQt8OcNreB/cNpR74UAF/Cbw0Ab1u00yuXIPXrg60tkJO3bs3dvtk1v\nvEFw8K2di9wLbn5/a9iwYeXKlevVq9eCBQuu5n6XTsQc/vyTvn0pVYrixSlblocfpls3mjal\nQgXOnr3hXo35dRpPx1IhnLAmbCUYpkAMrILljnfZQvMRTK9IzGim5p7OQNeuAJUqsWsXERE0\naULlyrRuzYwZjB9/2+cqd7WbLoo1derU1q1bu7i4AD4+PoMGDVq2bFlSUlIhrMdlElqTsAg5\nccJ47TXDze0Wlv7z4/xf+PR/PJD+vjjGcIzNWZb+MzCWZtvnLP5TeK4eu/L4Ea6uxt/+5uwf\njdyAmdckzOuisSdOnPjkk08ykrpUqVJWq3XdunWpqakFWp8ZKKCLinXrDF/fW4jm5myezvDL\nFLu2qCrGpxjns+RyKsYqjF4Y7umrsm6g5WBmeXE175/SqJGxf7+zfzRyY3dDQGdIS+pWrVql\nJXVgYODYsWN/+umngijOJBTQ5rdihdG6teHqmqfE9OfsGP69kwaZTW0xfsl+yRyDMQGjanqH\n05T+J+PqsDvvuQyGi4sxapRxL/3fZpFk5oC+5VEc5cuXb9GixZkzZ44dO3bo0KETJ05MmTJl\nypQptWrVevfdd3v16pWPt19EHEpJ4dNP+eILTpygbl0uX2br1jzt2JKNVmx9WOBF9i8K/wUN\nAUiFVRAOKyAFA8t62tqwLqJHIp43Oqy3Ny4u2Z4/9PFh7FisVqpUuY3zE0mX14BOSUnZuHHj\nd999t3jx4iNHjgCBgYGjRo3q2bNnQEDAnDlzwsPD+/Tp88svvzRu3LggC85np0+fHjt2bEqO\n52qzO3jwIGDkcXFmyW9nzrBuHRcu0LRp+vxtr7zChx+mb42Lu/kRSnNmCLOs2Oqwx3GPTyEM\nlsJ00qbUj6PcTIZGELqPmjc9/tWrvPkm8+cTHQ3QvDnh4TRokJeTE8nNzQN64cKF33333bJl\ny86dOwdUr179pZde6tmzZ7NmzSzX5m558MEHBw0a9OCDDy5cuLBoBbSnp+d9992Xe0BfvHgR\nsGiiGmdYtizbBPnly9OiBQsX5mlfC0YbfggjvAeLPEnMrWs4hAPYcVlLexvWJXRPylyH6ubc\n3Ni7l7g43N0pVSrv+4nk6qY3QdK6NWzY8O23396xY8eNul24cKF06dKTJ0/Ozxsw5qB70M5y\n6pRRsuQt3PPN+FWOky/zfjQ1M+46G89hbMeIwijveJ8YKkzgjWocvI2PA2POHGf/sOR2Fe17\n0B988EGPHj3uu+++3LuVKFHi9GmHQ/lFbtOmTQ6mkcuFC/b2rLVi684SD5IAWoIV+oDXtU71\ns03QnIrrKjpFELqMril5u+PXsSPr15OUlNlSuTKdO99CnSJ5dPM/kS+++GIh1CFyvbyncwVi\nhzEjlIhqHAIoDUPACnWydNoHn8K69HfHCJrOiEhGHiMoj59SrhwvvcS4caxfzzPPZLvjXLJk\nXksVyTvNxSEmdfEiMTE36eNKaidWWbF1YbkbKQCt4C/Qg8wxF4mwCGywHgxScFtOFxvWVXRK\nTZ9U3zGLhenTiY6mUiV69sTXN32eOeDRR3XHWQqDAlrMKDycZ58lOfmGHSpxfCSRI4kMSht1\nkaY7fJel0x6wwSzSptQ/RLVIRn7B8Fgq5KWG3r0ZNiy3DuVyW9xVJB8ooMVEDIMlS1i7lqlT\ncTis0Y2ULiwPJeJxVrqSmnNzAhiQCN9AOGwESMZ9Cd1tWNfQwZ6HyWfS9OlDePidnIpIPlBA\ni1kkJtKhAxs3Ot5alcMjiRzB9ArE3vAQ30MN+JO0+UH3UyOC0BkMi+MWrnXbtWP+fErnaR1X\nkYKlgBYnW72a8ePZuZOEhJwL8QHuJHdniRVbB9a4YMcVOkI7iMDxQycHScRzMU+GE7aetga3\nNno9MJAVK/C84TODIoVKAS2FatYspkzhzBmqVKFSJX76iUOHHPeswf5QIoYxoxxxAEEwAkaS\nPuaiJjyRc5c91LFhncWQM+T1AjgigjJlsNm4coU+fRg16vZOS6RAKKClkFy9SuvWmZNmHD3q\nuJsniU+yOIzwtqy3YOAGnSEMOpE55uISzM/cJQGvb+htw7qRlrd0yVyiBCNGYLHQvfttnJBI\ngVNAS4FITmbVKg4fpmZNOnTAYuHRR28ypVEd9oQSMYRZZdIWL6kKI2EE2cZcbIUImAvxAFE0\niCB0FkPO4X8bRU6bhh7gFzNTQEv+SEpi1iy2baNUKVq14rnn2HPtHnGpUrz2Glu2ON7Ri4Te\nfGPF1ooNma3TIDTLgj8X4CuwwTaAKxT7mqfCCdtC81xK8vfPnMTD05MuXXjySVav5sABKlfm\n2Wdp2fIOTlik4CmgJR9cusQjj/D77+lvLZZsg+TOnuX//s/BXvXZFUb4YGb7cy7bBhcYei2d\nf4JwmA9XALbTyIZ1DgMv4HfTqj74gOrV+f13ypfn8cfx8QEYPPg2zk/EORTQkg9efjkzncHx\nEOYMxbjyFF9bsYWw2XEPOzwFTWEuRAFcwudLBkUy8lfyOldiw4YMH47FQuvWedxDxHQU0HKn\nUlKYPTtPPRux3YptIHP80gYqV4Yr6Y/55bQElgD8SuNwwubRLx7fXI7s7p7tscOSJdmyRfeX\npchTQMud2r6dS5dy6+DDpf7MtWJrwlYAL+gNVmgJZ+E+uJhzlwv4fcWAcMJuulQ24OrK6tVM\nm8YPP6RfMn/+Od7et31CImahgBYHrl5l9WpOnKB+fVq1yrYpOZnLl7NN3vbNNzc8TmN+tWLr\nz1zftFEX9cEKgyFjgiEPckxYtIXmNqzz6XuFYnkptUwZpkyhbVvats1Ld5GiRAEtOe3cSffu\nHD6c/vbRR1myhGLFiI3luedYsoTkZGrV4sMP6dKF997j/fdzHsGPCwP4yortgbRRF8WgD1jh\nkSyd4mAWfEbaF4Tn8J/NYBvWKG5hqaiZM+nfH3f32z5XEXNz9ooBRcA9taKK3W7Ur59zuZDn\nnzeSk42mTbM1ursba9YYLi7ZGpuzeTrDL1E8/X1pjI8xzmdZLTsV43uMPhju6ftsoOVgZnlx\n9VYXMSlTxrg3fk+kYBXtFVXknnL4MLt25Wxcvpzevfn552yNyckMH47dDuDPuUF8acUWzM5s\nnf4Kz157fQK+gAjSptQ/Q+m0S+bd1M1jbS4u6R8HlCzJl1+mj5wTuVspoCUbh4uYnD7NnDkO\n2o8fpyUbrdj6sMCLBAc91kJ/+B/YYBmkYGBZT1sb1kX0SCTnpEQ9e1K9OkeP0rIldjvPP5+Z\nyJ6erF9PQkL6szBdulCmzB2dqYj5KaAlm/r18fbm6tVsjefP89ln2VoC+HMIs6zY6rI7t8Ot\ngMD0l6coO4NhNqz7qXF9R1dXZsxg0KBsjY8+yiefcPgwNWrw3HPUqAHoy0C5hyigJZvkZEqU\nyBnQGSwYbfjBiq0n33qSCNAUwqA+/CX9Oewc7Lispb0N6xK6J+GR0e7lRcK1a+6HHmL9enyv\nG+hcrx6ffnrHpyRSZCmg72mrVvHxxxw9Sp06vPIKjRuzfj1xcQ56liNuKDNDiajJPoCSMAis\ncP+1Hr1yBvQJAr9geAShh6iW42itWrF4MWvXEhNDvXrpsymJSA4K6HuO3Y6LC8D06Ywcmd4Y\nFcWiRUyezJo12Tq7YG/PWiu27izxIAmgBVihD2Q8CZIM38G1S91UXFfT0YZ1GV1TbvAHrEQJ\n/P3p0ye/z03k7lK0A9putx84cCA5OblWrVpubkX7XAra+fO88Qbz5hEfz0MP8SaHNpEAABkM\nSURBVO67PP98tg6pqYwbl57dQCAnhvNFKBHV0kZdlIIhEEa2MRcHIAJmwEmA41SKZGQkI4+l\nT6p/Qzt35r5dRKAIBfSbb75ZrVq1ESNGpL1NSUn58MMPJ0yYcOnSJcDT03Po0KGTJk3y87v5\nJGf3IMOgf39WrUp/u2UL7duTet2aq4DFntqZ1VZsXVnmxrUVqAJgD5mrlCTBYgiHdWCQgtsK\nOtuwruTx1OzPBbq7U7Ikp0/n/JTrbzeLiAPOHoidV0CbNm0y3o4ZMwbw9/fv3bv3008/3axZ\nM6BevXoJCQn5/tF3wYMqv/128+c+KnHsb7x9lCAH24IwkjAMjD0YL2GUSW8/SLU3mFCR4w4P\n2KaNcfmysWWLg00TJjj7JyJyjZkfVCmSAR0VFWWxWB5++OEzZ85kdIiMjATeeuutfP/ooh7Q\nMTFGx443zGVXUrrz3VK6puCaW34HYzTHsGBAEu7f0Ksjq1xIddjXzc1YuTKzgGeeyba1cWMj\nOdl5Pw6R7Mwc0EXmFkdWmzdvNgzjo48+CggIyGgcMWLE9OnTV65cOX78+LwfKjU1ddmyZUlJ\nSbn0+e23326/Vmc7fJgHH8xcWCSrqhweSeQIplcgFiAQnoQtsN3RgXYC7KdGJCO/YHgc5a7v\n4uJCSAghIbz1FsWLZ7ZPnUqvXnz3HZcu0aYNAwbg6nr93iKSU5EM6OPHjwPBwcE52oODg+c4\nfOIt10P95S9/SUhw9BTcNWn/wBq5z0JvSsnJjBmTM53dSe7G0jDCO7DGBTsu0BGs0BXc4UT2\nNQABSMRzET0iCF1HuxutyurtTWQk/fs7rqRdO9q1u/MTErm3FMmArlmzJnDkyJH69etnbT95\n8mTVqlVv6VBVqlSJjY3Nvc+0adNGjRplKQojdVNS2LuXK1eoX58NGwgNJSYmc2sN9ocSMZSZ\n5dNGXVSCETASKmc5xH+zHXAPdSIIncnQM5lfETr21Vc8+WR+nYeIQNEK6Ojo6HfeeadOnTq1\na9cuU6bMO++8M3fu3IytW7duXb58+eB7eMm5n35i6FCiowH8/EhMTH9Uz5PEHiyyYmvLegsG\nrtAZrNA5y1zMl+BrmAa/ACTg9Q29bVg30CrHp7i6Ohj+0bYt3bsX4KmJ3KOcfRM8r4KCgq6/\nhl23bl3a1ldffdXb2zsgIODo0aP5/tFF4kvCM2eMwMCcX9bVZs8HvHia0plNT2AczzL5p4Hx\nG8YojBLpHaKoP5aP/Dnr8Nu/UqWM2bONt94yqlY1PD0NLy8jKMh46SXjwgVnn7/I7dKXhPng\n6NGjV69e3b9/f3R09L59+/bt2xcdHZ3xcMrixYtLly49e/bsoKCbPCJxl0lIYMQIFi4k69ec\nXiT05ptQIlqxwUL2W+fjoCIAF+ErsMH/AK5Q7GuesmHdTAjg7g7J2fYLC6N3bx5+mLSB5rfy\nRayI3KYiE9CAt7d3cHDw9d8NAgsXLqxTp45LxmNwd6lt29izh6AgQkLSH/l77DE2bszsUJ9d\nYYQPZrY/jsZtAG/D0/A9zIfLANtpZMM6h4EXSH/Gx8+Po0f59Vfee4/du6lcmdGjGTBA02WI\nFLaiFNC5qFevnrNLKFiXL9OnDytXpr996CEWLWLduvR0LsaVp/jaii2EzQCu4I7D+ZlZD+sB\nLuEzj37hhG2lSY4u/v6UKKFxFyLOd5cE9F3v5Zcz0xn47Tf692f7dhqx3YptIHP8uABQC0Jh\nKBSDEHA05cWvNLZhnUv/eBw/cN2pU0GcgYjcMgV00bBoUba3vsTX3TT/Q2wPp4268IKeYIXW\nZA5TDsoW0BcpsbnawH+csv54+YFcPqh5cyZNyt/aReQ2KaCLgKQkzpxJf92ErVZs/ZjnSzxA\nXQiDIVDqWm8D/gufwor0hp9oFk7YAstTW1cUH/UbPw7KeXxg2DCCg6lXj8ce426/ky9SZCig\nze7sWR59lGLJFwbwlRXbA2mz4lugP/wFWmTpehpmgg2iAc7h/yWDbFh3Egxg4O7OwIEEBvLq\nq/z2W+Zyf716ERmpXBYxHQW02X06eMtz221P8XXxtFEXabpCxjPtdlgHNlhM2pT6P9IinLBv\n6H01c1J9gF9/pXp12rXjl184dowVKzh/nocf1ip/IialgDaXI0c4doxatSjrfo4vvyQ8/I2o\nKAf9DkM8XIYZEAEHAP4kYCZDIwjdTV2LhevnDjl1KvN1UBBPP11AJyEi+UMB7XwnThAdjY8P\nb7/NsmW0ZGOYxdbX9Rv3lBss3QrshDKQDHYMLP91aTvTPXSFV89awZ6vhVGlClOm8O23OXcq\n52AGOhExLwW0M6WkMHo0ERHY7QTw5xBm/YGtLrsxwBv6Q0f4DNY62jmRU5SdwbAIQvcbNY8d\n4IuKmRsvXMgZ0KVK0aZNQZ6MiOQ3BbQzTZiALdxoww9WbD351pNEgCYQBv3AB4DSOQPajsta\n2tuwLqF7Eh4ABtu2UTFLQHfrxoQJTJiQ/gh4uXLMmkXZsoVzWiKSPxTQzrFjB7M+OOU9f8Ze\nImqyD6AkDAQrNMzS7whkGZV80hI43RgeQeghquU44PXh+8YbDBvGL7/g40Pz5loGUKToUUAX\nOrv9l3fXHnvL9p6xxCNt1EUIWOEpKHatTwosBRusBjt2XOIe6BT4ltWrRVfXSLfSCzm0Ndsh\n77uPhg25XqVKVKpUsGcjIgVHAV2ITpzgiy+IiHj40KGHMxoXQs8sfQ5BBHwBJwCOU2k6IyIZ\nmXSi8oknKQmvvMLLL/PMM4SHpw9krl6dr7/G07Nwz0VECp4CuuClprJ6NTYby5aRkpJza3MA\nkuA7sMF/wE4qrsvpYsO6ksdT0ybVP0lMTPpdZouFzz5j3Dh27CAggJAQpbPI3UkBXZCOHWP6\ndKZP5+jRG/Z5HJrAEjgFcJiqc7xGzvEYvvtixay9XF0pWTLbfjVrUrNmQRQtImahgC4AKSms\nWEF4OKtWkZqKBVrARfjdUecdsINk3JfSLZywrSU7zJztsn8hu2dk69WhQ7Z1skXkXqCAzleH\nDxMZyfTppC1EWw6GQSjUgGSoDsdy7nGA6ssDQxt9NCzuz/JDSzK7PWXK0LIlhw/zww/pfR56\niMjIQjwLETEHBXR+MAyWLuXTT1mzBrsdF+gIVugO7tf6nCPrXBp4eMQ83OPn+61enduN7mRx\ndc22OKufH+vX8/PPREdTpQotWmgmI5F7kQI6P1it6Ze4FWA4jCTbMOWdYIPZcB6A2rUJDWXo\n0IplyvR0dLAMTZvStGmB1SwipqeAvmM//khkJJXgY+ia5Sd6GeaDDX4CwMuLgb2wWmnVSqv7\niUheFNWAvnjxYnx8vIuLS7ly5Zy8Vuy2bQCj4clrLf8DG3wFFwGoXx+rlcGDKVXqBocQEXGg\niN3ajIqKGjp0aGBgoJ+fX6VKlSpUqODh4VGpUqUBAwZs2rTJOTXVqAEwHzbBNGgCD8HnkFKM\noUPZtImoKMaOVTqLyK0qSlfQY8aMmTp1qmEYgYGBTZs2DQgIAM6ePXv8+PG5c+fOnTs3NDTU\nZrMVdlkdOx6t1b7y9rUZi5ukNmjoOsrKoEH4+RV2MSJyFykyAf3pp59+8sknHTt2fO+99x54\nIOeyp7t27ZowYUJERETdunXHjRtXmIXN+tLl6eilo/i8GT+dIHAu/ctWfXjp6MIsQUTuThbj\n+oU3TOmRRx75888/o6Ki3Nwc/6NiGEbr1q3tdvuPP/6Yvx89bdq0UaNGxcfH+/j4XL+1RQuu\nv7ly/Hi2yT9FxLSSkpI8PT03bdoUEhLi7FpyKjL3oKOiopo1a3ajdAYsFkvLli2jHC4QdWMH\nDx709va25GrUqFFpx3d4hAMH8tooInJLiswtjgYNGvz888+pqamurq436rNly5YGDRrc0mGr\nVau2evXqpLRp7W9g165dzz//vLu7u8Otdepw8qSDRhGRO1RkAnrgwIGjR4/u1q3b+++/Hxwc\nnGNrdHT0+PHj169fP2nSJIe734jFYmnVqlXufYoVK5bL1pdeynwmO83QoVq7RETyQZEJ6Gee\neWbnzp2ff/75ypUrg4KCqlSpUqpUKYvFcu7cuWPHjh06dAgYNmzYSy+9VMiFdenC11/z5ptE\nR+Pvz8iRvP12IZcgInenIvMlYZrt27dPmjRpzZo1Z86cSWtxdXUtW7ZsmzZtnn766datWxfE\nh27evPmRRx5JTEz08PDIpVtiouZlFil6zPwlYZG5gk7TqFGjr776Cjh//nx8fLy7u3vZsmUL\n+knCtFz2VPqK3L1yv/xyliJ2Be0sO3bsSLl+MRRHXnjhBW9v78GDBxd0SXmRkJBgtVr//ve/\nV6uWc5FZp1i7du3atWsnTpzo7ELSjR8//sEHH+zWrZuzCwE4efLk//3f//373//29/d3di0A\nCxYsOHjw4CuvvOLsQtKNHTt2zJgxBfSb5ebm1tDhsp5OZ0i+6tmz53PPPefsKtLFx8cDW7du\ndXYh6aZOnVqvXj1nV5GpWbNm7733nrOrSLd3714gJibG2YWke/311zt06ODsKjIFBQXNmjXL\n2VUUtiIzDlpE5F6jgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQIiIm\npYDOZx4eHuZ5qN/Nzc3FxcU89Zjqh4PJ6vHw8LBYLOapx93d3TzFYLLfrEKjuTjy2ZkzZzw8\nPEqUKOHsQtIdPHjwvvvuc3YV6ZKSkk6dOlWpUiVnF5LuxIkTJUuW9Pb2dnYh6Uz1m3Xp0qUr\nV66UNc3U5kePHq1QoUIuayrdlRTQIiImpVscIiImpYAWETEpBbSIiEkpoEVETEoBLSJiUgpo\nERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JA3+X279//\nySefOLsKyatLly7NnDnz+PHjzi5ETEEBnW8+++yzFi1alCxZskWLFp999pmzy0n38ccfv/nm\nm86ugsTExNdff71Vq1Z+fn7Vq1cfMGDAgQMHnFXMoUOHBgwYULNmzeLFiwcHB7/88ssXLlxw\nVjE5jBkzZtiwYTt27HBiDUFBQZbrOPFP0caNG9u3b+/n51ehQoW+ffs68U+OExiSH0aNGgXU\nrl17yJAhtWrVAp599llnF2V8//33np6eJUuWdG4Z58+fb9myJVCvXr3Q0NDHHnvMYrF4e3tv\n27at8IvZt29f8eLF3dzc2rVrN2rUqKZNmwL169e/evVq4ReTw4IFC9L+Vi5btsxZNVy5csVi\nsVSoUKFNdpGRkU6pZ968eR4eHhUqVBgwYMATTzzh6uoaEBBw5MgRpxRT+BTQ+WDbtm1Ap06d\nkpOTDcNITk5Oy6CdO3c6q6SBAwfWrl077W+70wP6tddeA0aPHp3Rsnz5chcXl4YNGxZ+Mb16\n9bJYLEuWLMloeeGFF4CPP/648IvJ6vjx46VKlfLx8XFuQP/+++/AhAkTnFVAVkeOHHFzc2va\ntOn58+fTWmw2GzB06FCn1lV4dIsjH0yaNAl4//330xZMc3Nze++99wzDmDx5srNKunLlSs2a\nNbt27err6+usGjIsWrTI19f3n//8Z0ZL586d27Vrt2PHjlOnThVyMT/++OODDz7YrVu3jJbh\nw4cD//vf/wq5kqwMwxgyZIifn99zzz3nxDKA6OhooE6dOs4tI82UKVNSUlI++ugjPz+/tJaR\nI0f+61//atasmXMLKzT31gqMBWTNmjWVKlW6//77M1oefPDBwMDA77//3lklffvtt2kvgoOD\nnf6Nk4uLS+vWrT09PbM2pq3QfO7cucJcltRut7/55ptVqlTJ2hgXFwfUqFGj0Mq43j//+c8f\nfvjhv//976ZNm5xYBrBv3z6gSpUqc+bM2bdvX6VKlUJCQurVq+eUYubNmxcUFJQ1ji0Wy/PP\nP++UYpzD2ZfwRd65c+eARx55JEd72s3NixcvOqWqDA0aNHD6LY7rnTp1ysvLq1y5cmk3hZzi\nypUrMTExK1asqFmzZrly5aKjo51VybZt2zw8PF577TXDMCZOnIhTb3GMGDECKFOmTEZEuLi4\njBkzpvB/p+Lj44GWLVtu3769W7duZcuWDQoK6t279759+wq5EifSLY47lfbHKCAgIEd7WsvF\nixedUJO5RUdHh4SEJCQkTJw4Me2mkFOMGzeuYsWKnTt3jo2NTYtpp5Rx9erVgQMH1qtX7+23\n33ZKATmk3eJ49NFHf//99/j4+B9//PGhhx76+OOPP/zww0Ku5Pz580BsbGyLFi0OHz7ctWvX\n+vXrf/vttw0bNvz1118LuRincfa/EEXeiRMngCeeeCJHe+fOnYHY2FinVJXBVFfQly5deuut\nt7y9vb28vD755BPnFrN9+/b58+f/4x//qFy5sqen5+LFi51SxujRo728vKKiotLeOv0KesOG\nDevWrcvacvr0aX9/fx8fn9TU1MKs5ODBg2kZ9eqrr9rt9rTGNWvWWCyWBx54oDArcSIF9J1K\nTU11dXVt1apVjvZmzZq5uroW8p/p65knoFesWFG5cmWga9eue/bscXY5mWJiYnx9fStWrFj4\nH7127VrgX//6V0aL0wPaod69ewOFfBfo5MmTQEBAQEpKStb2xx57DIiLiyvMYpxFtzjulIuL\nS9myZa//Ii4mJqZ8+fIuLvoJA/ztb3/r3Lmzr6/vf//736VLl2YMASxkBw4cmDZtWlRUVNbG\nChUqNG7cOCYmJu3rhMK0fft24IUXXsh4HuTVV18FunbtarFYIiMjC7meG0m7X5ecnFyYH1qm\nTBkvL69q1aq5urpmbb/vvvsAp3/1XTg0iiMftGnTZu7cudHR0WmPqAC7du06duxY//79nVuY\nScycOfPvf/97v379Zs6cmTZ4w1ni4uJGjRr13HPPTZkyJWv76dOnfXx8MsZyFZqGDRumPeKU\nYdu2bT///PPjjz9epUqVwh/r9scff/Tq1atHjx7vvvtu1vYdO3Z4enpm/PEuHC4uLm3atNm8\neXNCQoKXl1dG++7du11cXJz1b3xhc/Yl/N3ghx9+AAYNGpT21m639+3bF9i4caNzCzNMcIvD\nbrfXrl27YsWKZnhULykpqWzZsn5+fgcOHMhonDdvHo6+RXAK597iSE1NDQoK8vb2/uWXXzIa\n0y7kw8LCCr+e1atXA6NHj864VTh//nyga9euhV+MU+gKOh+0bt162LBhM2bMiI2Nbdas2Y8/\n/rhhw4aRI0e2aNHC2aU535EjR/bu3VumTJkePXpcv3X27NmlS5cutGLc3d0//vjjfv36BQcH\nd+7cuWzZsrt3716/fn25cuWmTp1aaGWYlouLy+zZs3v27PnII4907dq1fPnyv//++6ZNm+rW\nrfv+++8Xfj2PPfbYsGHDpk6dumHDhubNmx86dGjNmjWBgYHmmeumwDn7X4i7hN1uf//990NC\nQkqUKBESEjJ58mRnV5TO6VfQ//nPf3L543f8+PHCL2ndunWdOnUKCAgoVqxYw4YNx40bd/bs\n2cIvwyEzfEl45MiR4cOHN2jQwMfHp3Hjxm+++aZz/+/ngw8+aNGiha+vb7169Z599lnz/GYV\nAothGIX1b4GIiNwCjTEQETEpBbSIiEkpoEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKQU0CIi\nJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iYlAJaRMSkFNAiIialgBYRMSkFtIiI\nSSmgRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQIiImpYAWETEpBbSIiEkpoEVETEoBLSJi\nUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBGAXbt2eXp6\ntm3bNqMlOTk5ODg4ICDg5MmTTixM7mUKaBGA+vXr//Wvf/3hhx+++OKLtJbJkydHRUX9+9//\nLl++vHNrk3uWxTAMZ9cgYgpJSUkPPfRQbGzsnj17Lly4EBwc3KFDhyVLlji7Lrl3KaBFMv38\n888hISH9+vU7efLktm3bdu3aFRgY6Oyi5N7l5uwCREykadOmY8eO/de//gXMmjVL6SzOpSto\nkWz2799fs2bN4sWLx8bGlihRwtnlyD1NXxKKZPPiiy96eHhcvnz5tddec3Ytcq9TQItkmjNn\nzpIlSyZOnNi7d+/PPvts8+bNzq5I7mm6xSGSLi4urn79+lWrVv3555/j4uLq1q1bqVKlbdu2\neXh4OLs0uUfpClok3TPPPHP+/Pnw8HBXV9cKFSq8++67f/zxx8SJE51dl9y7dAUtAjB//vx+\n/fq9+OKLH3zwQVqL3W5v3rz5jh07tm3bVrduXeeWJ/cmBbSIiEnpFoeIiEkpoEVETEoBLSJi\nUgpoERGTUkCLiJiUAlpExKQU0CIiJqWAFhExKQW0iIhJKaBFRExKAS0iYlIKaBERk1JAi4iY\nlAJaRMSkFNAiIialgBYRMSkFtIiISSmgRURMSgEtImJSCmgREZNSQIuImJQCWkTEpBTQIiIm\npYAWETEpBbSIiEkpoEVETEoBLSJiUgpoERGTUkCLiJiUAlpExKT+H1CAwo807wVHAAAAAElF\nTkSuQmCC", "text/plain": [ "plot without title" ] }, "metadata": { "image/png": { "height": 180, "width": 240 } }, "output_type": "display_data" } ], "source": [ "myLM <- lm(y~x)\n", "myLMCoef <- myLM$coefficients\n", "yHatOLS <- myLMCoef[1] + myLMCoef[2]*x\n", "\n", "plot(x, y, pch = 20, col = \"blue\")\n", "points(sort(x), yHatLoss[order(x)], type = \"l\", col = \"red\", lwd = 5)\n", "points(sort(x), yHatOLS[order(x)], type = \"l\", lty = \"dashed\",\n", " col = \"yellow\", lwd = 2, pch = 20)" ] }, { "cell_type": "markdown", "id": "bde63ca2", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Discussion\n", "\n", "- What are the advantages and disadvantages of using likelihood function and loss function?" ] } ], "metadata": { "celltoolbar": "Slideshow", "kernelspec": { "display_name": "R", "language": "R", "name": "ir" }, "language_info": { "codemirror_mode": "r", "file_extension": ".r", "mimetype": "text/x-r-source", "name": "R", "pygments_lexer": "r", "version": "4.3.2" } }, "nbformat": 4, "nbformat_minor": 5 }