QuantitativeMethods_I.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Quantitative Methods Part I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "import math\n",
    "from scipy.stats.mstats import gmean \n",
    "import statistics\n",
    "\n",
    "import csv\n",
    "from urllib.request import urlopen\n",
    "from bs4 import BeautifulSoup\n",
    "\n",
    "from matplotlib import pyplot as plt\n",
    "%matplotlib inline\n",
    "\n",
    "import scipy\n",
    "import scipy.stats as stats\n",
    "from scipy.stats import norm, kurtosis, skew"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The Time Value of Money\n",
    "## Interest Rates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Differentiate:\n",
    "- real risk-free interest rate\n",
    "- inflation premium (The sum of the real risk-free interest rate and the inflation premium is the nominal risk-free interest rate.)\n",
    "- default risk premium\n",
    "- liquidity premium\n",
    "- maturity premium"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Time value of money equates cash flows that occur on different dates. Cash flows in the future must be discounted at appropriate interest rates to find the equivalent present value. \n",
    "\n",
    "Interest rates can be thought of as `required rates of return`, `discount rates`, or `opportunity costs`. They are determined by the interaction between the supply of funds from savers and the demand for funds from borrowers.\n",
    "\n",
    "Using a `building block approach`, we can calculate an interest rate, $r$, as the sum of the `risk-free rate` and any applicable premiums to compensate for various risks:\n",
    "\n",
    "$$r=r_f+M+I+L+D$$\n",
    "\n",
    "where\n",
    "\n",
    "$r_f$ is the `real risk-free rate` that savers require as compensation for lending today rather than consuming.\n",
    "\n",
    "$I$ is the `inflation premium` that compensates lenders for expected inflation. Adding this to the real risk-free rate gives the nominal risk-free interest rate.\n",
    "\n",
    "$D$ is the `default risk premium` that lenders charge for the possibility that borrowers will not meet their obligations on time and in full.\n",
    "\n",
    "$L$ is a `liquidity premium` that compensates the investor for the extra cost of converting the investment to cash.\n",
    "\n",
    "$M$ is the `maturity premium` earned by investors who are willing to lend over longer periods, which requires greater exposure to interest rate risk."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Future Value of a Single Cash Flow (Lump Sum) \n",
    "The future value of an investment has three components:\n",
    "- Principal is the amount of the original investment\n",
    "- Simple interest is the interest earned on the principal\n",
    "- Compounding is the return on reinvested cash flows (or \"interest on interest\")\n",
    "\n",
    "### Frequency of Compounding\n",
    "#### Annual Compounding\n",
    "$$FV_1 = PV(1 + r)$$\n",
    "<div style=\"text-align: right\"> Equation 1</div>\n",
    "\n",
    "$$FV_{N} = PV(1 + r)^N$$\n",
    "<div style=\"text-align: right\"> Equation 2 </div>\n",
    "#### Non-Annual Compounding (Future Value)\n",
    "The formula used in the previous example assumes annual compounding. However, investments are often quoted in terms of a stated annual rate with more frequent compounding. In these cases, we have to adjust the stated annual interest rate ($r_s$) to be consistent with the number of compounding periods per year ($m$):\n",
    "\n",
    "$$FV_{N}=PV(1+\\displaystyle\\frac{r_{s}}{m})^{mN}$$\n",
    "<div style=\"text-align: right\"> Equation 3 </div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Stated and Effective Rates\n",
    "Rather than quote the periodic monthly interest rate, financial institutions often quote an annual interest rate that we refer to as the **stated annual interest rate or quoted interest rate**. \n",
    "The expression for the future value of a lump sum in $N$ years with continuous compounding is<br><br>\n",
    "$$FV_{N}=PVe^{r_{s}N}$$\n",
    "<div style=\"text-align: right\"> Equation 4 </div>\n",
    "\n",
    "Distinction between the stated annual interest rate and the `effective annual rate (EAR)`:\n",
    "#### For discrete compounding\n",
    "$$EAR=(1+Periodic\\ interest\\ rate)^m-1$$  \n",
    "<div style=\"text-align: right\"> Equation 5 </div>\n",
    "i.e., \n",
    "$$EAR=(1+\\frac{r_s}{m})^m−1$$\n",
    "\n",
    "#### For continuous compounding\n",
    "$$EAR=e^{r_{s}}−1$$\n",
    "<div style=\"text-align: right\"> Equation 6 </div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Future Value of a Series of Cash Flows, Future Value Annuities\n",
    "- An `annuity` is a specified number of level cash flows. \n",
    "- `Ordinary annuities` start one period from now while,\n",
    "- an `annuity due` begins immediately. \n",
    "- A `perpetuity` is a stream of cash flows that never ends, with the first cash flow occurring in one period.\n",
    "\n",
    "### Equal Cash Flows—Ordinary Annuity\n",
    "The future value of an annuity is the sum of the future values of each payment. If N annuity payments of amount A are made and accumulated at interest rate r, the future value is:\n",
    "\n",
    "$$FV_{N} = A [(1+r)^{N-1} + (1+r)^{N-2} + ... + (1+r)^{0}]$$<br><br>\n",
    "$$FV_{N} = A \\Big[\\displaystyle\\frac{(1+r)^N-1}{r}\\Big]$$\n",
    "<div style=\"text-align: right\"> Equation 7 </div>\n",
    "\n",
    "### Unequal Cash Flows\n",
    "Shortcuts to take when the cash flows are close to even; these shortcuts will allow you to combine annuity and single-period calculations: (TBC)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## PV\n",
    "### PV of a Single Cash Flow (Lump Sum)\n",
    "#### Annual Compounding\n",
    "$$PV=FV_{N}(1+r)^{−N}$$\n",
    "<div style=\"text-align: right\"> Equation 8 </div>\n",
    "\n",
    "#### Non-Annual Compounding \n",
    "In general, with more than one compounding period in a year, we can express the formula for present value as<br><br>\n",
    "$$PV=FV_{N}\\Big(1+\\frac{r_{s}}{m}\\Big)^{−mN}$$\n",
    "<div style=\"text-align: right\"> Equation 9 </div>\n",
    "where\n",
    "\n",
    "$m$ = number of compounding periods per year\n",
    "\n",
    "$r_{s}$ = quoted annual interest rate\n",
    "\n",
    "$N$ = number of years\n",
    "### PV of a Series Cash Flows\n",
    "#### PV of a Series of Equal Cash Flows\n",
    "##### Ordinary Annuity\n",
    "The present value of an `ordinary annuity` with $N$ equal payments of $A$ is:\n",
    "$$PV=\\displaystyle\\frac{A}{(1+r)}+\\frac{A}{(1+r)^{2}}+\\frac{A}{(1+r)^3}+…+\\frac{A}{(1+r)^{N−1}}+\\frac{A}{(1+r)^N}$$\n",
    "<div style=\"text-align: right\"> Equation 10 </div>\n",
    "where\n",
    "\n",
    "$A$ = the annuity amount\n",
    "\n",
    "$r$ = the interest rate per period corresponding to the frequency of annuity payments (for example, annual, quarterly, or monthly)\n",
    "\n",
    "$N$ = the number of annuity payments\n",
    "\n",
    "The formula is shortened as<br>\n",
    "$$PV=A\\displaystyle\\Big[\\frac{1−\\frac{1}{(1+r)^N}}{r}\\Big]$$\n",
    "<div style=\"text-align: right\"> Equation 11 </div>\n",
    "\n",
    "##### Annuity Due\n",
    "An `annuity due` is like an ordinary annuity except that the first payment is received immediately (at $t=0$). It can be thought of as an ordinary annuity with $N−1$ payments plus one payment immediately. For example, the following two payments are equivalent:\n",
    "1. 5-year ordinary annuity plus one immediate payment now; and\n",
    "2. 6-year annuity due.\n",
    "\n",
    "##### Perpetuity (with Level Payments)\n",
    "A `perpetuity` is an ordinary annuity (one in which the first payment is in one period) with no ending. Examples of perpetuity: \n",
    "- value dividends from stocks\n",
    "- certain government bonds \n",
    "- preferred stocks \n",
    "\n",
    "$$PV = A\\sum_{t=1}^\\infty{\\frac{1}{(1+r)^t}}$$\n",
    "<div style=\"text-align: right\"> Equation 12</div>\n",
    "As long as interest rates are positive, the sum of present value factors converges and\n",
    "$$PV=\\frac{A}{r}$$\n",
    "<div style=\"text-align: right\"> Equation 13</div>\n",
    "\n",
    "*Note that Equation 13 is valid only for a perpetuity with level payments.*\n",
    "\n",
    "##### PV Indexed at Times Other than t = 0 \n",
    "One type of PV that are index at the time other than t=0 is a `deferred ordinary annuity` is like an ordinary annuity except that the first payment is not made until more than one year in the future. The present value at time zero can be calculated in two steps:\n",
    "1. Calculate the present value at the time one period prior to the first payment.\n",
    "2. Discount this present value back to time zero.\n",
    "\n",
    "Perpertuity can also be indexed at a later date. Consider a level perpetuity of £100 per year with its first payment beginning at t = 5. What is its present value today (at t = 0), given a 5 percent discount rate? \n",
    "#### PV of a Series of Unequal Cash Flows\n",
    "When we have unequal cash flows, we must first find the present value of each individual cash flow and then sum the respective present values. For a series with many cash flows, we usually use a spreadsheet."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1645.4049495837637"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Step 1: calculate t at the previous 1 point before payment starts\n",
    "pv_4 = 100/0.05\n",
    "pv_4\n",
    "# Step 2: calculate pv at t=0\n",
    "pv_0 = pv_4/(1+0.05)**4\n",
    "pv_0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "User case: constructing an annuity as the difference between two perpetuities with equal, level payments but differing starting dates. <br>\n",
    "\n",
    "Given a 5 percent discount rate, find the present value of a four-year ordinary annuity of £100 per year starting in Year 1 as the difference between the following two level perpetuities:<br>\n",
    "\n",
    "Perpetuity 1\t£100 per year starting in Year 1 (first payment at t = 1)<br>\n",
    "Perpetuity 2\t£100 per year starting in Year 5 (first payment at t = 5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "354.5950504162363 354.5950504162363\n"
     ]
    }
   ],
   "source": [
    "# Solution 1: Substracting the PVs of the 2 perpetuities\n",
    "pv1_0 = 100/0.05 # PV of Perpetuity 1 at t=0\n",
    "pv2_4 = 100/0.05 # PV of Perpetuity 2 at t=4\n",
    "pv2_0 = pv2_4/(1+0.05)**4 # PV of Perpetuity 2 at t=0\n",
    "pv_s1 = pv1_0 - pv2_0\n",
    "\n",
    "# Solution 2: PV of an ordinary annuity of n=4 using equation 11\n",
    "r = 0.05\n",
    "a = 100\n",
    "n = 5-1\n",
    "rate = 1-(1+r)**(-n)\n",
    "pv_s2 = a*rate/r\n",
    "\n",
    "print(pv_s1, pv_s2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Solving for Rates, Number of Periods, or Size of Annuity Payments\n",
    "### (Compound) Interest Rates or Growth Rates\n",
    "Depending on the context, the interest rate, $r$, could be considered a growth rate $g$.\n",
    "$$g = (\\frac{FV_{N}}{PV})^{1/N} – 1$$\n",
    "<div style=\"text-align: right\"> Equation 14</div>\n",
    "\n",
    "### Number of Periods\n",
    "Based on the following formula (as for rates) and use a calculator:\n",
    "$FV_N=PV(1+r)^N$\n",
    "\n",
    "### Size of Annuity Payments (Combining FV and PV of Annuities)\n",
    "Given the present value, number of periods, and interest rate, the annuity payment $A$ can be solved for using the following formula:\n",
    "$$PV=A\\displaystyle\\Big[\\frac{1−\\frac{1}{(1+r)^N}}{r}\\Big]$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Present Value and Future Value Equivalence\n",
    "A lump sum can be seen as equivalent to an annuity, and an annuity can be seen as equivalent to its future value. Thus present values, future values, and a series of cash flows can all be considered equivalent as long as they are indexed at the same point in time.\n",
    "$$PV\\ of\\ Annuity=(FV\\ of\\ Annuity)(1+r)^{−N}$$\n",
    "\n",
    "### The Cash Flow Additivity Principle\n",
    "Amounts of money indexed at the same point in time are additive. \n",
    "\n",
    "Suppose cash flows are \\\\$4 at the end of the first year and \\\\$24 (actually separate payments of \\\\$4 and \\\\$20) at the end of the second year. Rather than finding present values of the first year’s \\\\$4 and the second year’s $24, we can treat this situation as a \\\\$4 annuity for two years and a second-year \\\\$20 lump sum. If the discount rate were 6 percent, the \\\\$4 annuity would have a present value of \\\\$7.33 and the \\\\$20 lump sum a present value of \\\\$17.80, for a total of \\\\$25.13."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example 1:\n",
    "\n",
    "Toyota Motor Corporation, one of the largest automakers in the world, had consolidated vehicle sales of 8.96 million units in 2018 (fiscal year ending 31 March 2018). This is substantially more than consolidated vehicle sales of 7.35 million units six years earlier in 2012. What was the growth rate in number of vehicles sold by Toyota from 2012 to 2018?\n",
    "\n",
    "Example 2:\n",
    "\n",
    "You are interested in determining how long it will take an investment of €10,000,000 to double in value. The current interest rate is 7 percent compounded annually. How many years will it take €10,000,000 to double to €20,000,000?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.03356258254589006"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# example 1\n",
    "sales1 = 7.35\n",
    "sales2 = 8.96\n",
    "n = 2018 - 2012\n",
    "\n",
    "# Applying equation 14\n",
    "g = (sales2/sales1)**(1/n) -1\n",
    "g"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "10.244768351058712"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# example 2\n",
    "# Using logarithm \n",
    "fv = 20 # unit: millions\n",
    "pv = 10\n",
    "r = 0.07\n",
    "#1.07**n = 2\n",
    "math.log(fv/pv, 1+r)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example 3. \n",
    "\n",
    "You are planning to purchase a \\\\$120,000 house by making a down payment of $20,000 and borrowing the remainder with a 30-year fixed-rate mortgage with monthly payments. The first payment is due at t = 1. Current mortgage interest rates are quoted at 8 percent with monthly compounding. What will your monthly mortgage payments be?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100000 360 0.006666666666666667\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "733.764573879378"
      ]
     },
     "execution_count": 53,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculating monthly mortgage payments\n",
    "monthly_rate = 0.08/12\n",
    "down_payment = 20000\n",
    "price = 120000\n",
    "n = 30*12\n",
    "pv = price - down_payment\n",
    "print(pv, n, monthly_rate)\n",
    "\n",
    "# Applying equation 11, not equation 7 because the bank is calculating annuity (monthly payment) that has a present value of 10,000\n",
    "numerator = 1 - 1 / (1+monthly_rate)**n\n",
    "denominator = monthly_rate\n",
    "m = numerator/denominator\n",
    "payment = pv/m\n",
    "payment"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example 4. \n",
    "\n",
    "An individual wants to retire with a specified retirement income.<br>\n",
    "Jill Grant is 22 years old (at t = 0) and is planning for her retirement at age 63 (at t = 41). She plans to save \\\\$2,000 per year for the next 15 years (t = 1 to t = 15). She wants to have retirement income of $100,000 per year for 20 years, with the first retirement payment starting at t = 41. How much must Grant save each year from t = 16 to t = 40 in order to achieve her retirement goal? Assume she plans to invest in a diversified stock-and-bond mutual fund that will earn 8 percent per year on average."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "54304.2278549568 401653.2501156755 981814.7407449294 580161.4906292539 8689.636530502932\n"
     ]
    }
   ],
   "source": [
    "yearly_rate = 0.08\n",
    "\n",
    "# Annuity one (a1): saving from t=1 to t=15, saving for 15 years and compounding until t=41\n",
    "a1_a = 2000\n",
    "a1_n = 15 # annuity period \n",
    "a1_m = 41-15 # lump sum period from t=15 to t=41 (or t=16???)\n",
    "a1_fv_1 = a1_a * (((1+yearly_rate)**a1_n -1) / yearly_rate) # FV of firstly batch of savings at t=15, equation 7\n",
    "a1_fv_2 = a1_fv_1 * (1+yearly_rate)**a1_m # ???\n",
    "\n",
    "# Annuity two (a2): pension from t=41 onwards\n",
    "a2_a = 100000\n",
    "a2_n = 20\n",
    "a2_pv = a2_a * ((1-1/(1+yearly_rate)**a2_n) / yearly_rate) # PV of desired pension fund at t=41\n",
    "\n",
    "# Annuity three (a3): saving from t=16 to t=40\n",
    "a3_fv = a2_pv - a1_fv_2\n",
    "a3_n = 40-16\n",
    "a3_a = a3_fv / (((1+yearly_rate)**a3_n-1) / yearly_rate) \n",
    "\n",
    "print(a1_fv_1, a1_fv_2, a2_pv, a3_fv, a3_a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Note: in the reading, the solution discounted the FV (at t=41) of a2 to PV (t=16), and then calculate the annuity of that PV; alternatively, I solved for the FV of a1, a2, and then use the FV (at t=41)to calculate the annuity for a2 at t=16, the result slightly differs (\\\\$8,342.87 in reading).* <br> \n",
    "Follow up: find out if there is any miscalculations."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Statistical Concepts and Market Returns\n",
    "## Introduction\n",
    "Four properties of return distributions are discussed in this reading:\n",
    "1. `Central Tendency` – where returns are centered\n",
    "2. `Dispersion` – how far returns are spread out from the center\n",
    "3. `Skewness` – whether the distribution is symmetric or lopsided\n",
    "4. `Kurtosis` – how likely are extreme outcomes\n",
    "\n",
    "#### The Nature of Statistics\n",
    "Statistics can refer to both data and method. Statistical methods include the following:\n",
    "- Descriptive Statistics: Turn data into information by summarizing\n",
    "- Statistical Inference (`inferencial`): Make forecasts based on probability theory\n",
    "\n",
    "#### Populations and Samples\n",
    "A `population` is all members of a specified group. `Parameters` are descriptive measures of the population:\n",
    "- size: $N$\n",
    "- mean: $\\mu$\n",
    "- std: $\\sigma$\n",
    "- variance: $\\sigma^2$\n",
    "- proportion: $p$\n",
    "\n",
    "A `sample` is a subset of the population. A process of sampling is used to try and make the sample representative of the population. A sample statistic is a descriptive measure of a sample.\n",
    "- size: $n$\n",
    "- mean: $\\bar{X}$\n",
    "- std: $S$\n",
    "- variance: $S^2$\n",
    "- proportion: $\\hat{p}$\n",
    "\n",
    "#### Measurement Scales\n",
    "- `Nominal scales` are used to categorize but not rank data. An example is the style of a mutual fund (small-cap vs. large-cap). The scale categorizes the funds based on their investment strategies but does not rank them.\n",
    "- With `ordinal scales`, data is sorted into categories and ordered based on some characteristic. An example includes Morningstar star ratings for mutual funds.\n",
    "- `Interval scales` provide ranking and ensure the difference between scale values are equal. An example is the Celsius temperature scale. There are no true zero points, so ratios are meaningless.\n",
    "- `Ratio scales` are like an interval scale, but with a true zero point. This allows ratios to be meaningful. An example is money measurement."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Frequency Distribution and Graphic Presentation of Data\n",
    "#### The Frequency Distribution\n",
    "A `frequency distribution summarizes` the data into a few intervals, which is helpful when working with large data sets.\n",
    "\n",
    "When constructing a frequency distribution\n",
    "1. first sort the data and calculate the `range` (maximum value – minimum value). \n",
    "2. Determine the desired `number of intervals`, $k$, which makes the width of each interval equal to $\\frac{range}{k}$. \n",
    "3. Then, construct a table that shows the number of observations in each interval, which is the absolute frequency. Each observation will only fall in one interval since the intervals will not overlap.\n",
    "\n",
    "The above approach can be modified to make the intervals stop and start at rounded numbers (e.g., whole numbers). The choice of $k$ is somewhat arbitrary but should be chosen to best summarize the data.\n",
    "\n",
    "The `relative frequency` is $\\frac{absolute\\ frequency}{total\\ number\\ of\\ observations}$ \n",
    "\n",
    "The `cumulative relative frequency` adds up the relative frequencies as you move from the first to the last interval. \n",
    "\n",
    "The `frequency distribution` illustrates where most of the data lie and shows the shape of the distribution.\n",
    "#### The Histogram\n",
    "A histogram uses the data from a frequency distribution to create a bar chart.  The height of each bar represents the frequency. If no observations were in a particular interval, the bar height would be zero.\n",
    "#### The Frequency Polygon and the Cumulative Frequency Distribution\n",
    "A `frequency polygon` plots each interval's midpoint on the x-axis and the absolute frequency on the y-axis. The points are connected with `line segments`.\n",
    "\n",
    "A line graph can also be used to illustrate the cumulative frequency distribution. This shows how many observations lie below a certain value. Most observations will lie on the steep slope."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Measures of Central Tendency\n",
    "A measure of central tendency indicates where the data is centered. These are the most common statistical measures. The most common central tendency measures are arithmetic mean, median, mode, weighted mean, and geometric mean.\n",
    "\n",
    "### The Arithmetic Mean\n",
    "The `arithmetic mean` is the most common measure of where the data are centered. It is the sum of the observations divided by the number of observations. Both population means and sample means can be computed.\n",
    "\n",
    "The `population mean` is the arithmetic mean for the entire population. All observations are included in calculating this parameter. For a population with $N$ members, the population mean is: \n",
    "$$μ=\\frac{\\sum^N_{i=1}{X_i}}{N}$$\n",
    "\n",
    "The `sample mean` or average is the arithmetic mean for a sample. For a sample with n observations, the sample mean is:\n",
    "$$\\bar{X}=\\frac{\\sum^n_{i=1}{X_i}}{n}$$\n",
    "\n",
    "`Cross-sectional data` are observations that capture characteristics of different units at a specific point in time. For example, the average dividend yield among FTSE 100 companies is a cross-sectional mean.\n",
    "\n",
    "`Time-series data` are observations of the same unit at different points in time. For example, the average month-end dividend yield for the FTSE 100 over the past 10 years is a time-series mean.\n",
    "\n",
    "The arithmetic mean is like the center of gravity for an object. The sum of the deviations from the mean is zero.\n",
    "\n",
    "$$\\sum_{i=1}^n{(X_i−\\bar{X})}=0$$\n",
    "\n",
    "The value of the mean can be skewed by extreme values, or `outliers`. For example, the mean of 3, 8, 9, 10, and 100 is 26. One method of addressing outliers is to use a trimmed mean, which is calculated by discarding a certain percentage of the highest and lowest values. For example, with a sample of 100 observations, a 2% trimmed mean would be the arithmetic mean without the highest value and the lowest value. With a Winsorized mean, no observations are excluded from the calculation, but any outliers' values are adjusted to either an upper or lower limit.\n",
    "\n",
    "The mean is easy to work with and uses all the information about the observations' size and magnitude. Since the mean is sensitive to extreme values, it can be misleading.\n",
    "### Weighted Mean\n",
    "A portfolio’s return is a weighted mean return computed from the returns on the individual assets, where the weight applied to each asset’s return is the fraction of the portfolio invested in that asset.\n",
    "\n",
    "The weighted mean formula is:\n",
    "\n",
    "$$\\bar{X_w}=\\sum_{i=1}^n{w_iX_i}$$\n",
    "where:\n",
    "\n",
    "$w_i$ is the weight\n",
    "\n",
    "$\\sum^n_{i=1}{w_i=1}$\n",
    "\n",
    "**Weights less than zero are permissible, which would represent a `short asset` in a portfolio.**\n",
    "\n",
    "Weighted averages are used to calculate past portfolio or index returns. They can also be used to calculate future expected returns by weighing various scenarios.\n",
    "### The Geometric Mean\n",
    "The geometric mean is used to average rates over time or compute growth rates. It is often used to average portfolio returns from different time periods.\n",
    "\n",
    "\\begin{equation*}\n",
    "G = \\sqrt[n]{x_{1}x_{2}...x_{n}} = \\sqrt[n]{∏_{i=1}^n{X_i}}\n",
    "\\end{equation*}with $X_{i}≥0$ for i=1,2,…,n.\n",
    "\n",
    "One must be added to each return when calculating the geometric mean. For T holding period returns, the geometric mean return is:\n",
    "\n",
    "$$R_G=\\sqrt[T]{∏_{t=1}^T{(1+R_t)}}−1$$\n",
    "\n",
    "The geometric mean is always less than or equal to the arithmetic mean. The two are equal when all observations have the same value. Greater dispersion in the values of observations will increase the difference between the geometric and arithmetic means.\n",
    "\n",
    "**The geometric mean return represents an investment's compound growth rate over multiple periods. By contrast, the arithmetic mean is used to measure performance for a single period.**\n",
    "\n",
    "### Harmonic Mean\n",
    "The harmonic mean is not as commonly used but is appropriate for some circumstances. It is the reciprocal of the average of the reciprocals, which is calculated as follows:\n",
    "$$\\bar{X_H}=\\frac{1}{\\Big(\\frac{1}{n}\\Big)\\Big(\\displaystyle\\sum_{i=1}^n{\\frac{1}{X_i}}\\Big)}$$\n",
    "where $X_i>0$.\n",
    "\n",
    "The observation's weight is inversely proportional to its magnitude. This would be appropriate if you want the average price paid for a security when investing the same dollar amount for several time periods (also known as the `cost averaging technique`).\n",
    "\n",
    "The harmonic mean is always less than the geometric mean (which is less than the arithmetic mean) unless all the observations are the same.\n",
    "### The Median\n",
    "The median is the middle item of a sorted list. Assume there are $n$ items:\n",
    "- If $n$ is odd, then the median is the $(\\frac{n+1}{2})$th term.\n",
    "- If $n$ is even, then the median is the mean of the $(\\frac{n}{2})$th and $(\\frac{n+2}{2})$th terms.\n",
    "\n",
    "The median is not sensitive to extreme values because it focuses only on the position, not the magnitude. The median is more complicated to calculate (less mathematically tractable) than the mean.\n",
    "### The Mode\n",
    "The mode is the most frequently occurring value in a distribution. Some distributions have more than one mode, while others have none. A distribution with just one mode is called unimodal.\n",
    "\n",
    "Data grouped in intervals have modal intervals. This is the highest bar in a histogram.\n",
    "\n",
    "The mode is the only measure of central tendency that can be used with nominal data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Other Measures of Location: Quantiles\n",
    "### Quartiles, Quintiles, Deciles, and Percentiles\n",
    "If we arrange the observations in ascending order, then the quantile is a value at or below which a stated fraction of the data is found.\n",
    "\n",
    "There are many common quantiles used in practice. Distributions are often divided into four quartiles, five quintiles, ten deciles, or one hundred percentiles.\n",
    "\n",
    "The $y_{th}$ percentile is the value at or below which $y$ percent of the observations lie. \n",
    "\n",
    "The values of quartiles, quintiles, and deciles are based on percentiles. For example, the 90th percentile is equal to the ninth decile, and the 50th percentile is equal to the second quartile. Note that the 50th percentile is also the median.\n",
    "\n",
    "The location of the $y_{th}$ percentile ($L_y$) in a list of n observations ranked in ascending order (lowest to highest) is calculated as follows:\n",
    "\n",
    "$$L_y=(n+1)\\frac{y}{100}$$\n",
    "If $L_y$ is not an integer, then linearly interpolate between the two values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 233,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Year</th>\n",
       "      <th>Selected American Shares (SLASX) %</th>\n",
       "      <th>T. Rowe Price Equity Income (PRFDX)%</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2013</td>\n",
       "      <td>34.90</td>\n",
       "      <td>31.69</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2014</td>\n",
       "      <td>6.13</td>\n",
       "      <td>7.75</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2015</td>\n",
       "      <td>2.69</td>\n",
       "      <td>-7.56</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2016</td>\n",
       "      <td>11.66</td>\n",
       "      <td>18.25</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2017</td>\n",
       "      <td>21.77</td>\n",
       "      <td>16.18</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Year  Selected American Shares (SLASX) %  \\\n",
       "0  2013                               34.90   \n",
       "1  2014                                6.13   \n",
       "2  2015                                2.69   \n",
       "3  2016                               11.66   \n",
       "4  2017                               21.77   \n",
       "\n",
       "   T. Rowe Price Equity Income (PRFDX)%  \n",
       "0                                 31.69  \n",
       "1                                  7.75  \n",
       "2                                 -7.56  \n",
       "3                                 18.25  \n",
       "4                                 16.18  "
      ]
     },
     "execution_count": 233,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# geometric mean\n",
    "# The most recent five years of total returns for two US large-cap value equity mutual funds in 2018\n",
    "d = {'Year': [2013, 2014, 2015, 2016, 2017], \n",
    "     'Selected American Shares (SLASX) %': [34.90, 6.13, 2.69, 11.66, 21.77], \n",
    "     'T. Rowe Price Equity Income (PRFDX)%':[31.69, 7.75, -7.56, 18.25, 16.18]}\n",
    "df = pd.DataFrame(data=d)\n",
    "df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 234,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Calculate the geometric mean return of SLASX\n",
    "slasx = df['Selected American Shares (SLASX) %'].tolist()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation*}\n",
    "R_{t} + 1 = \\sqrt[t]{(1+R_{1})(1+R_{2})...(1+R_{t})}\n",
    "\\end{equation*}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 235,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Redefine returns to make them positive by adding 1.0\n",
    "g_slasx = list(map(lambda i: i/100+1, slasx))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 236,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.14858526227425917"
      ]
     },
     "execution_count": 236,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gmean(g_slasx) -1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 237,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "15.429999999999998"
      ]
     },
     "execution_count": 237,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the arithmetic mean return of SLASX and contrast it to the fund’s geometric mean return\n",
    "df['Selected American Shares (SLASX) %'].mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 238,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.1249999599109004"
      ]
     },
     "execution_count": 238,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the geometric mean return of PRFDX\n",
    "g_prfdx = list(map(lambda i: i/100+1, df['T. Rowe Price Equity Income (PRFDX)%'].tolist()))\n",
    "gmean(g_prfdx) - 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 239,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "13.262"
      ]
     },
     "execution_count": 239,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the arithmetic mean return of PRFDX and contrast it to the fund’s geometric mean return\n",
    "df['T. Rowe Price Equity Income (PRFDX)%'].mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 242,
   "metadata": {},
   "outputs": [],
   "source": [
    "# harmonic mean practice\n",
    "price_per_share = [12.1, 13.5, 17, 11.98, 15]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 243,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Harmonic Mean is 13.67310629317005 \n"
     ]
    }
   ],
   "source": [
    "print(\"Harmonic Mean is % s \" % (statistics.harmonic_mean(price_per_share))) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Real time data application - Quantile**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 244,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Web scraping Euro Stoxx 50 data\n",
    "\n",
    "html = urlopen(\"https://www.dividendmax.com/market-index-constituents/euro-stoxx-50\")\n",
    "soup = BeautifulSoup(html, \"html.parser\")\n",
    "table = soup.findAll(\"table\")[0]\n",
    "rows = table.findAll(\"tr\")\n",
    "\n",
    "with open(\"editors.csv\", \"wt+\", newline=\"\") as f:\n",
    "    writer = csv.writer(f)\n",
    "    for row in rows:\n",
    "        csv_row = []\n",
    "        for cell in row.findAll([\"td\", \"th\"]):\n",
    "            csv_row.append(cell.get_text())\n",
    "        writer.writerow(csv_row)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 245,
   "metadata": {},
   "outputs": [],
   "source": [
    "euro_stoxx = pd.read_csv('editors.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 246,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Company</th>\n",
       "      <th>Exchange</th>\n",
       "      <th>Price</th>\n",
       "      <th>Currency</th>\n",
       "      <th>\\nMarket Cap\\n\\n</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>🇫🇷AIAir Liquide S.A</td>\n",
       "      <td>Euronext Paris (France)</td>\n",
       "      <td>€133.50</td>\n",
       "      <td>EUR</td>\n",
       "      <td>€41.7bn</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>🇩🇪ALVAllianz SE</td>\n",
       "      <td>Frankfurt Stock Exchange (Germany)</td>\n",
       "      <td>€166.52</td>\n",
       "      <td>EUR</td>\n",
       "      <td>€75.2bn</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>🇧🇪ABIAnheuser-Busch In Bev SA/NV</td>\n",
       "      <td>Euronext Brussels (Belgium)</td>\n",
       "      <td>€46.79</td>\n",
       "      <td>EUR</td>\n",
       "      <td>€75.2bn</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>🇱🇺MTLArcelorMittal</td>\n",
       "      <td>Luxembourg Stock Exchange (Luxembourg)</td>\n",
       "      <td>€12.94</td>\n",
       "      <td>EUR</td>\n",
       "      <td>€19.5bn</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>🇳🇱ASMLASML Holding NV</td>\n",
       "      <td>Euronext Amsterdam (Netherlands)</td>\n",
       "      <td>€320.25</td>\n",
       "      <td>EUR</td>\n",
       "      <td>€134.5bn</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                            Company                                Exchange  \\\n",
       "0               🇫🇷AIAir Liquide S.A                 Euronext Paris (France)   \n",
       "1                   🇩🇪ALVAllianz SE      Frankfurt Stock Exchange (Germany)   \n",
       "2  🇧🇪ABIAnheuser-Busch In Bev SA/NV             Euronext Brussels (Belgium)   \n",
       "3                🇱🇺MTLArcelorMittal  Luxembourg Stock Exchange (Luxembourg)   \n",
       "4             🇳🇱ASMLASML Holding NV        Euronext Amsterdam (Netherlands)   \n",
       "\n",
       "     Price Currency \\nMarket Cap\\n\\n  \n",
       "0  €133.50      EUR          €41.7bn  \n",
       "1  €166.52      EUR          €75.2bn  \n",
       "2   €46.79      EUR          €75.2bn  \n",
       "3   €12.94      EUR          €19.5bn  \n",
       "4  €320.25      EUR         €134.5bn  "
      ]
     },
     "execution_count": 246,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "euro_stoxx.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Quantiles in practice: determing or approximating the value of a percentile (when the location is not a whole number). "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\\begin{equation*}\n",
    "L_{y} = (n+1) * \\frac{y}{100}\n",
    "\\end{equation*}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 247,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Make a function for calculating Ly with a given column and desired percentile\n",
    "def loc_y(column, percentile):\n",
    "    n = len(column)\n",
    "    loc = (n+1)*percentile/100\n",
    "    return loc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Situation I\n",
    "When the location, $L_{y}$, is a whole number, the location corresponds to an actual observation. \n",
    "<br>$X_{i}$ is defined as the value of the observation in the ith ($i = L_{60}$) position of the data sorted in ascending order (i.e., $P_{30} =8.16$)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Situation II\n",
    "When $L_{y}$ is not a whole number or integer, $L_{y}$ lies between the two closest integer numbers (one above and one below), and we use **linear interpolation** between those two places to determine $P_{y}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Making a function to calculate $P_{y}$ given an array (observation) and the percentile of interest. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Practice based on CFA textbook example data:**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 248,
   "metadata": {},
   "outputs": [],
   "source": [
    "his_euro_stoxx = pd.read_csv('euro_stoxx_2018_nov.csv', delimiter=';')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 249,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>No.</th>\n",
       "      <th>Company</th>\n",
       "      <th>Market Cap (Euro Billion)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>RWE</td>\n",
       "      <td>9.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>Carrefour</td>\n",
       "      <td>12.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>E.ON</td>\n",
       "      <td>15.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>Inditex</td>\n",
       "      <td>15.5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>Unibail Rodamco</td>\n",
       "      <td>16.1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   No.          Company  Market Cap (Euro Billion)\n",
       "0    1              RWE                        9.9\n",
       "1    2        Carrefour                       12.2\n",
       "2    3             E.ON                       15.5\n",
       "3    4          Inditex                       15.5\n",
       "4    5  Unibail Rodamco                       16.1"
      ]
     },
     "execution_count": 249,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "his_euro_stoxx.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 250,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 1. Turn the column into a list and sort ascending\n",
    "# 2. Get location\n",
    "# 3. Calculate Py based on the location: if y is integer, direct indexing; if not an integer, use linear interpolation\n",
    "def p_given_percent(column, percentile):\n",
    "    list1 = column.tolist()\n",
    "    sorted_list = sorted(list1)\n",
    "    n = len(sorted_list)\n",
    "    loc = (n+1)*percentile/100\n",
    "    if loc.is_integer():\n",
    "        print('Ly is an integer,', percentile, '% of data lies at or below Py =', sorted_list[(int(loc)-1)])\n",
    "    else:\n",
    "        i = int(loc)\n",
    "        P_approx = float(sorted_list[i-1]) + (float(loc) - i) * (float(sorted_list[i]) - float(sorted_list[i-1]))\n",
    "        print('Ly is not an integer, Py is approximately', P_approx, sorted_list[i-1], 'and', sorted_list[i], 'bracket Py from below and above, respectively.')\n",
    "        \n",
    "# WIP        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 251,
   "metadata": {},
   "outputs": [],
   "source": [
    "#?isinstance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 252,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ly is not an integer, Py is approximately 16.16 16.1 and 16.7 bracket Py from below and above, respectively.\n"
     ]
    }
   ],
   "source": [
    "# 10th percentile     \n",
    "p_given_percent(his_euro_stoxx['Market Cap (Euro Billion)'], 10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 253,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ly is not an integer, Py is approximately 104.7 104.7 and 104.7 bracket Py from below and above, respectively.\n"
     ]
    }
   ],
   "source": [
    "# 90th percentile     \n",
    "p_given_percent(his_euro_stoxx['Market Cap (Euro Billion)'], 90)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 254,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ly is not an integer, Py is approximately 37.5 36.6 and 38.4 bracket Py from below and above, respectively.\n"
     ]
    }
   ],
   "source": [
    "# Second quartile\n",
    "p_given_percent(his_euro_stoxx['Market Cap (Euro Billion)'], 50)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 255,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ly is not an integer, Py is approximately 20.52 20.1 and 22.2 bracket Py from below and above, respectively.\n"
     ]
    }
   ],
   "source": [
    "# First quintile \n",
    "p_given_percent(his_euro_stoxx['Market Cap (Euro Billion)'], 20)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 256,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "37.5"
      ]
     },
     "execution_count": 256,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Value of the median\n",
    "np.median(his_euro_stoxx['Market Cap (Euro Billion)'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 257,
   "metadata": {},
   "outputs": [],
   "source": [
    "# How many quintiles are there, and to what percentiles do the quintiles correspond? \n",
    "# Answer: 4, P20, P40, P60 and P80"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Real time Euro Stoxx data application. Example 1 with $L_{y}$ as a whole number. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution 1.** Step by step"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 258,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Change Price data type as float, get rid of € symbol\n",
    "euro_stoxx['Price'] = euro_stoxx['Price'].str.replace('€', '')\n",
    "price = pd.to_numeric(euro_stoxx['Price']).tolist()\n",
    "price = sorted(price)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 259,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "30.0"
      ]
     },
     "execution_count": 259,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Obtaining the location (i) at 60th percentile\n",
    "loc_y(euro_stoxx['Company'], 60)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since $L_{y}$ is a whole number, $P_{60}$ = $X_{30}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 260,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "46.77"
      ]
     },
     "execution_count": 260,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Get value of X30 (= P60)\n",
    "price[29]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution 2.** Using `p_given_stoxx()` funtion"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 261,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ly is an integer, 60 % of data lies at or below Py = 36.70\n"
     ]
    }
   ],
   "source": [
    "# Using p_given_percent() funtion:\n",
    "p_given_percent(euro_stoxx['Price'], 60)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "60% of the values in this data set lies at or below €44.78 in price.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution 3.** Using Numpy's `percentile` function"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To verify the above calculation for future efficiency, I proceed with using Numpy's `percentile` function to obtain the value. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 262,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "45.603999999999985"
      ]
     },
     "execution_count": 262,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Obtain the value at 60th percentile\n",
    "np.percentile(price, 60)\n",
    "# But this is not in the observation, what went wrong?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example 2 with $L_{y}$ as a fraction.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution 1.** Step by step"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 263,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The location is at 37.5\n"
     ]
    }
   ],
   "source": [
    "# Obtaining the location (i) at 75th percentile\n",
    "print('The location is at', loc_y(euro_stoxx['Company'], 75))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The value of $P_{y}$ need to be approximated between $X_{37}$ and $X_{38}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 264,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Obtain values on both sides of Ly\n",
    "X37 = price[36]\n",
    "X38 = price[37]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 265,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "80.41"
      ]
     },
     "execution_count": 265,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Linear interpolation\n",
    "P75_approx = X37 + (37.5 - 37) * (X38 - X37)\n",
    "P75_approx"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution 2.** Using `p_given_percent()` function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 266,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ly is not an integer, Py is approximately 52.84 52.34 and 53.34 bracket Py from below and above, respectively.\n"
     ]
    }
   ],
   "source": [
    "p_given_percent(euro_stoxx['Price'], 75)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution 3.** Verify with Numpy `percentile`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 267,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "74.3"
      ]
     },
     "execution_count": 267,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Obtain the value at 75th percentile\n",
    "np.percentile(price, 75)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Quantiles in Investment Practice\n",
    "Quantiles are used in portfolio performance evaluation as well as in investment strategy development and research."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**User case 1. Ranking performance—for example, the performance of portfolios.**<br>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**User case 2. Investment research<br>**\n",
    "Analysts refer to a group defined by a particular quantile as that quantile. For example, analysts often refer to the set of companies with returns falling below the 10th percentile cutoff point as the bottom return decile. Dividing data into quantiles based on some characteristic allows analysts to evaluate the impact of that characteristic on a quantity of interest. For instance, empirical finance studies commonly rank companies based on the market value of their equity and then sort them into deciles. The 1st decile contains the portfolio of those companies with the smallest market values, and the 10th decile contains those companies with the largest market value. Ranking companies by decile allows analysts to compare the performance of small companies with large ones."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Investment style based on liquidity**<br>\n",
    "Ibbotson et al. (2018): Investment style based on liquidity—buying stocks of less liquid stocks and selling stocks of more liquid stocks.<br>\n",
    "Computed:\n",
    "- yearly measures of liquidity as the annual share turnover (the sum of the 12 monthly volumes divided by each month’s shares outstanding);\n",
    "- size as the year-end market capitalization; \n",
    "- value as the trailing earnings-to-price ratio as of the year end;\n",
    "- momentum as the annual return. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 1. assigned one-fourth of the total sample with the lowest liquidity in a year to Quartile 1 and the one-fourth with the highest liquidity in that year to Quartile 4."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 2. treating each quartile group as a portfolio composed of equally weighted stocks, they measured the returns on each liquidity quartile in the following year (so that the quartiles are constructed “before the fact”.) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Measure of Dispersion\n",
    "The most common measures of **(absolute) dispersion**: `range`, `mean absolute deviation`, `variance`, and `standard deviation`.\n",
    "### Range\n",
    "The range is the difference between the maximum and minimum values:\n",
    "\n",
    "$$Range=Maximum\\ Value−Minimum\\ Value$$\n",
    "\n",
    "- Pros: The range is easy to compute, \n",
    "- Cons: it is limited as a dispersion measure because it uses only the highest and lowest values. The interquartile range, which is the difference between the 25th and 75th percentiles, is less sensitive to extreme values.\n",
    "\n",
    "### Mean Absolute Deviation\n",
    "The mean of the deviations around the mean will always be zero. Therefore, it would not be a useful measure of dispersion. The `mean absolute deviation (MAD`) adjusts for this in the following formula:\n",
    "\n",
    "$$MAD=\\frac{\\displaystyle\\sum^n_{i=1}{|X_i−\\bar{X}|}}{n}$$\n",
    "\n",
    "- Pros: It uses all the observations in the sample, which makes it better than the range. \n",
    "- Cons: it does not work as well with portfolio optimization software as the variance and standard deviation measures.\n",
    "\n",
    "### Population Variance and Population Standard Deviation\n",
    "$$\\sigma^2 = \\frac{\\sum_{i=1}^{N} (X_{i} - \\mu)^2}{N}$$\n",
    "$$\\sigma = \\sqrt{\\frac{\\sum_{i=1}^{N} (X_{i} - \\mu)^2}{N}}$$\n",
    "\n",
    "*In practice, only used for well-defined group that can be considered to represent the population.*\n",
    "\n",
    "### Sample Variance and Sample Standard Deviation\n",
    "Often, only a sample of the population is available to analyze. Summary measures are then called statistics.\n",
    "\n",
    "The sample variance is:\n",
    "\n",
    "$$s^2 = \\displaystyle\\frac{\\sum_{i=1}^{n} (X_{i} - \\bar{X})^2}{n-1}$$\n",
    "\n",
    "This is very similar to the formula for the population variance, except X¯ is used instead of μ, and the divisor is n−1 rather than N. The quantity n−1 is called the degrees of freedom. Since the sample mean must be estimated first, there are only n−1 independent deviations from it.\n",
    "\n",
    "The sample standard deviation is:\n",
    "$$s = \\displaystyle\\sqrt{\\frac{\\sum_{i=1}^{n} (X_{i} - \\bar{X})^2}{n-1}}$$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 268,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Mean Absolute Deviation practice\n",
    "# Changing Price column as numeric data type\n",
    "euro_stoxx['Price'] = pd.to_numeric(euro_stoxx['Price'])\n",
    "# Obtaining the mean absolute deviation using Pandas mad function\n",
    "euro_stoxx['Price'].mad()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 270,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "7516.714429571011"
      ]
     },
     "execution_count": 270,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Obtaining the variance\n",
    "euro_stoxx['Price'].var(ddof=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 271,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "86.69898747719613"
      ]
     },
     "execution_count": 271,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Obtaining the standard deviation \n",
    "euro_stoxx['Price'].std(ddof=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Mean Variance and Standard Deviation\n",
    "Present the arithmetic mean and standard deviation together when summarizing data. When dealing with data that represent a time series of percent changes, presenting the geometric mean—representing the compound rate of growth—is also very helpful."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 272,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "7673.312646853741"
      ]
     },
     "execution_count": 272,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Obtaining the variance, with ddof = 1 (default value for Pandas std function)\n",
    "euro_stoxx['Price'].var(ddof=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 273,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "87.59744657724757"
      ]
     },
     "execution_count": 273,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Obtaining the standard deviation, with ddof = 1 (default value for Pandas std function)\n",
    "euro_stoxx['Price'].std(ddof=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Semivariance and Semideviation for Downside Risks\n",
    " Variance and standard deviation take into account returns above and below the mean. They are symmetric risk measures. Investors really only care about `downside risk`.\n",
    "\n",
    "Semivariance is the average squared deviation below the mean. Semideviation is the positive square root of semivariance.\n",
    "$$s_{semi}^2 = \\frac{\\displaystyle\\sum_{\\forall X_{i}\\leq\\bar{X}}^{n} (X_{i} - \\bar{X})^2}{n-1}$$\n",
    "\n",
    "$$s_{semi} = \\sqrt{\\frac{\\displaystyle\\sum_{\\forall X_{i}\\leq\\bar{X}}^{n} (X_{i} - \\bar{X})^2}{n-1}}$$\n",
    "\n",
    "Similar concept for `Target semivariance and target semideviation`. `Target semivariance` is the average squared deviation below a target value (e.g., 10%). Target semideviation is its positive square root.\n",
    "\n",
    "$$Target\\ Semivariance=\\frac{\\sum_{for\\ all\\ X_i≤B}{(X_i−B)^2}}{n−1}$$\n",
    "where $B$ is the target.\n",
    "\n",
    "Variance and semivariance can rank prospects differently if the distributions are `asymmetric`.\n",
    "### Chebyshev’s Inequality\n",
    "According to `Chebyshev’s inequality`, for any distribution with **finite variance**, the proportion of the observations within $k$ standard deviations ($\\sigma$) of the arithmetic mean ($\\bar{X}$) is at least $1 − \\frac{1}{k^2}$ for all $k > 1$.\n",
    "\n",
    "A two-standard-deviation interval around the mean must contain at least 75% of the observations, and a three-standard-deviation interval around the mean must contain at least 89% of the observations, no matter how the data are distributed.\n",
    "### Coefficient of Variation (CV)\n",
    "A scale-free measure of relative dispersion, by expressing the magnitude of variation among observations relative to their average size, the CV permits direct comparisons of dispersion across different data sets.\n",
    "\n",
    "$$CV = \\displaystyle\\frac{s}{\\bar{X}}$$\n",
    "Note that the coefficient of variation has no units of measurement."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 274,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compute the calculation with stock data\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 277,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Market</th>\n",
       "      <th>Arithmetic Mean Return (%)</th>\n",
       "      <th>Standard Deviation of Return (%)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Australia</td>\n",
       "      <td>5.3</td>\n",
       "      <td>11.9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Hong Kong SAR</td>\n",
       "      <td>5.6</td>\n",
       "      <td>15.8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Japan</td>\n",
       "      <td>15.7</td>\n",
       "      <td>16.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>South Korea</td>\n",
       "      <td>4.2</td>\n",
       "      <td>8.6</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          Market  Arithmetic Mean Return (%)  Standard Deviation of Return (%)\n",
       "0      Australia                         5.3                              11.9\n",
       "1  Hong Kong SAR                         5.6                              15.8\n",
       "2          Japan                        15.7                              16.3\n",
       "3    South Korea                         4.2                               8.6"
      ]
     },
     "execution_count": 277,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "index_returns = {'Market': ['Australia', 'Hong Kong SAR', 'Japan', 'South Korea'], \n",
    "     'Arithmetic Mean Return (%)': [5.3, 5.6, 15.7, 4.2], \n",
    "     'Standard Deviation of Return (%)':[11.9, 15.8, 16.3, 8.6]}\n",
    "index = pd.DataFrame(data=index_returns)\n",
    "index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 278,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Obtaining the coefficent of variance for each index\n",
    "index['CV'] = index['Standard Deviation of Return (%)'] / index['Arithmetic Mean Return (%)']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 279,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Market</th>\n",
       "      <th>Arithmetic Mean Return (%)</th>\n",
       "      <th>Standard Deviation of Return (%)</th>\n",
       "      <th>CV</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Japan</td>\n",
       "      <td>15.7</td>\n",
       "      <td>16.3</td>\n",
       "      <td>1.038217</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>South Korea</td>\n",
       "      <td>4.2</td>\n",
       "      <td>8.6</td>\n",
       "      <td>2.047619</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Australia</td>\n",
       "      <td>5.3</td>\n",
       "      <td>11.9</td>\n",
       "      <td>2.245283</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Hong Kong SAR</td>\n",
       "      <td>5.6</td>\n",
       "      <td>15.8</td>\n",
       "      <td>2.821429</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          Market  Arithmetic Mean Return (%)  \\\n",
       "2          Japan                        15.7   \n",
       "3    South Korea                         4.2   \n",
       "0      Australia                         5.3   \n",
       "1  Hong Kong SAR                         5.6   \n",
       "\n",
       "   Standard Deviation of Return (%)        CV  \n",
       "2                              16.3  1.038217  \n",
       "3                               8.6  2.047619  \n",
       "0                              11.9  2.245283  \n",
       "1                              15.8  2.821429  "
      ]
     },
     "execution_count": 279,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Ranking CV by ascending order\n",
    "index.sort_values(by=['CV'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Determine whether there is more difference between the absolute or the relative riskiness of the Australia and South Korea markets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 280,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.3000000000000007"
      ]
     },
     "execution_count": 280,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Std (consider using %)\n",
    "11.9-8.6"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 281,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.20000000000000018"
      ]
     },
     "execution_count": 281,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# CV (consider using %)\n",
    "2.25-2.05"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Symmetry, Skewness, and Kurtosis in Return Distributions\n",
    "### Symmetry and Skewness\n",
    "An asymmetrical distribution is skewed. \n",
    "- A `positively skewed` distribution (`long right tail`) has frequent small losses and a few extreme gains. The mode is less than the median, which is less than the mean.\n",
    "- A `negatively skewed` distribution (`long left tail`) has frequent small gains and a few extreme losses. The mean is less than the median, which is less than the mode. Investors should be concerned if returns have this distribution.\n",
    "\n",
    "**Investors should be attracted by a POSITIVE skew (distribution skewed to the right) because the mean return falls above the median. Relative to the mean return, positive skew amounts to a limited, though frequent, downside compared with a somewhat unlimited, but less frequent, upside.**\n",
    "\n",
    "Skewness is the average cubed deviation from the mean divided by the cubed standard deviation. The formula for sample skewness is given in the curriculum, but it is unlikely that you will be asked to use it on the exam. It is far more likely that you will be given a skewness measure and asked to interpret it.\n",
    "\n",
    "If the value of the skewness measure is greater than zero, the distribution is positively skewed. A negative skewness measure indicates a negatively skewed distribution. As a general rule, skewness measures are considered significant if their magnitude exceeds ±0.5."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 282,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Year</th>\n",
       "      <th>Return (%)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2008</td>\n",
       "      <td>-35.75</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2009</td>\n",
       "      <td>25.62</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2010</td>\n",
       "      <td>15.15</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2011</td>\n",
       "      <td>-0.72</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2012</td>\n",
       "      <td>17.25</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>2013</td>\n",
       "      <td>31.69</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>2014</td>\n",
       "      <td>7.75</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>2015</td>\n",
       "      <td>-7.76</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>2016</td>\n",
       "      <td>18.25</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>2017</td>\n",
       "      <td>16.18</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Year  Return (%)\n",
       "0  2008      -35.75\n",
       "1  2009       25.62\n",
       "2  2010       15.15\n",
       "3  2011       -0.72\n",
       "4  2012       17.25\n",
       "5  2013       31.69\n",
       "6  2014        7.75\n",
       "7  2015       -7.76\n",
       "8  2016       18.25\n",
       "9  2017       16.18"
      ]
     },
     "execution_count": 282,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Making a Pandas Dataframe of 10 years of annual returns on the T. Rowe Price Equity Income Fund (PRFDX).\n",
    "RoR = {'Year':[2008, 2009, 2010, 2011, 2012, 2013, 2014, 2015, 2016, 2017], \n",
    "     'Return (%)':[-35.75, 25.62, 15.15, -0.72, 17.25, 31.69, 7.75, -7.76, 18.25, 16.18]}\n",
    "rate_of_return = pd.DataFrame(data=RoR)\n",
    "rate_of_return"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 284,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAO+klEQVR4nO3df6zddX3H8efLtqAZTKa90a60XIxkGy4q7AYxLgvRuQESumWY1Cz+mksTI1ETk6VKApFliWwJJgqxaQYRDMHfc3WWIE4W9Q+qt12pLZV5dS4t6aSCFggOV/feH+eLO7ve23Nue+49px+fj+Sbfn987vm+2p688j3f8/1+b6oKSdLp7znjDiBJGg0LXZIaYaFLUiMsdElqhIUuSY1YPa4dr127tqanp8e1e0k6Le3evftHVTW10LaxFfr09DSzs7Pj2r0knZaS/Mdi2zzlIkmNsNAlqREWuiQ1wkKXpEZY6JLUCAtdkhoxsNCTPDfJN5M8mORAkg8uMObMJJ9KMpdkV5Lp5QgrSVrcMEfozwCvrapXAK8ELk9y6bwx7wB+XFUvBT4M3DTamJKkQQYWevU81S2u6ab5D1HfBNzRzX8WeF2SjCylJGmgoe4UTbIK2A28FLi1qnbNG7IeOARQVceTHANeCPxo3utsAbYAbNy48dSSSxqZ6a1fGtu+f/ChN4xt360Z6kvRqvp5Vb0SOBe4JMnvnszOqmp7Vc1U1czU1IKPIpAknaQlXeVSVT8B7gcun7fpEWADQJLVwPOBx0YRUJI0nGGucplKck43/zzg9cB35g3bAby1m78G+Gr5y0olaUUNcw59HXBHdx79OcCnq+qfktwIzFbVDuA24BNJ5oDHgc3LlliStKCBhV5V+4CLFlh/fd/8fwFvHG00SdJSeKeoJDXCQpekRljoktQIC12SGmGhS1IjLHRJaoSFLkmNsNAlqREWuiQ1wkKXpEZY6JLUCAtdkhphoUtSIyx0SWqEhS5JjbDQJakRFrokNcJCl6RGWOiS1AgLXZIaYaFLUiMsdElqhIUuSY2w0CWpERa6JDViYKEn2ZDk/iQPJTmQ5D0LjLksybEke7vp+uWJK0lazOohxhwH3ldVe5KcDexOcl9VPTRv3Ner6qrRR5QkDWPgEXpVHamqPd38k8BBYP1yB5MkLc2SzqEnmQYuAnYtsPnVSR5Mck+Sly3y81uSzCaZPXr06JLDSpIWN3ShJzkL+Bzw3qp6Yt7mPcB5VfUK4KPAFxZ6jaraXlUzVTUzNTV1spklSQsYqtCTrKFX5ndV1efnb6+qJ6rqqW5+J7AmydqRJpUkndAwV7kEuA04WFU3LzLmxd04klzSve5jowwqSTqxYa5yeQ3wZuDbSfZ26z4AbASoqm3ANcA7kxwHfgpsrqpahrySpEUMLPSq+gaQAWNuAW4ZVShJ0tJ5p6gkNcJCl6RGWOiS1AgLXZIaYaFLUiMsdElqhIUuSY2w0CWpERa6JDXCQpekRljoktQIC12SGmGhS1IjLHRJaoSFLkmNsNAlqREWuiQ1wkKXpEZY6JLUCAtdkhphoUtSIyx0SWqEhS5JjbDQJakRFrokNcJCl6RGDCz0JBuS3J/koSQHkrxngTFJ8pEkc0n2Jbl4eeJKkhazeogxx4H3VdWeJGcDu5PcV1UP9Y25Arigm14FfKz7U5K0QgYeoVfVkara080/CRwE1s8btgm4s3oeAM5Jsm7kaSVJixrmCP0XkkwDFwG75m1aDxzqWz7crTsy7+e3AFsANm7cuLSkkjRC01u/NLZ9/+BDb1iW1x36S9EkZwGfA95bVU+czM6qantVzVTVzNTU1Mm8hCRpEUMVepI19Mr8rqr6/AJDHgE29C2f262TJK2QYa5yCXAbcLCqbl5k2A7gLd3VLpcCx6rqyCJjJUnLYJhz6K8B3gx8O8nebt0HgI0AVbUN2AlcCcwBTwNvH31USdKJDCz0qvoGkAFjCnjXqEJJkpbOO0UlqREWuiQ1wkKXpEZY6JLUCAtdkhphoUtSIyx0SWqEhS5JjbDQJakRFrokNcJCl6RGWOiS1AgLXZIaYaFLUiMsdElqhIUuSY2w0CWpERa6JDXCQpekRljoktQIC12SGmGhS1IjLHRJaoSFLkmNsNAlqREDCz3J7UkeTbJ/ke2XJTmWZG83XT/6mJKkQVYPMebjwC3AnScY8/WqumokiSRJJ2XgEXpVfQ14fAWySJJOwajOob86yYNJ7knyssUGJdmSZDbJ7NGjR0e0a0kSjKbQ9wDnVdUrgI8CX1hsYFVtr6qZqpqZmpoawa4lSc865UKvqieq6qlufiewJsnaU04mSVqSUy70JC9Okm7+ku41HzvV15UkLc3Aq1yS3A1cBqxNchi4AVgDUFXbgGuAdyY5DvwU2FxVtWyJJUkLGljoVfWmAdtvoXdZoyRpjLxTVJIaYaFLUiMsdElqhIUuSY2w0CWpERa6JDXCQpekRljoktQIC12SGmGhS1IjLHRJaoSFLkmNsNAlqREWuiQ1wkKXpEZY6JLUCAtdkhphoUtSIyx0SWqEhS5JjbDQJakRFrokNcJCl6RGWOiS1AgLXZIaYaFLUiMGFnqS25M8mmT/ItuT5CNJ5pLsS3Lx6GNKkgYZ5gj948DlJ9h+BXBBN20BPnbqsSRJSzWw0Kvqa8DjJxiyCbizeh4AzkmyblQBJUnDWT2C11gPHOpbPtytOzJ/YJIt9I7i2bhx40nvcHrrl076Z0/VDz70hrHt+1fNOP+ftXL8fx6dFf1StKq2V9VMVc1MTU2t5K4lqXmjKPRHgA19y+d26yRJK2gUhb4DeEt3tculwLGq+qXTLZKk5TXwHHqSu4HLgLVJDgM3AGsAqmobsBO4EpgDngbevlxhJUmLG1joVfWmAdsLeNfIEkmSTop3ikpSIyx0SWqEhS5JjbDQJakRFrokNcJCl6RGWOiS1AgLXZIaYaFLUiMsdElqhIUuSY2w0CWpERa6JDXCQpekRljoktQIC12SGmGhS1IjLHRJaoSFLkmNsNAlqREWuiQ1wkKXpEZY6JLUCAtdkhphoUtSI4Yq9CSXJ3k4yVySrQtsf1uSo0n2dtNfjj6qJOlEVg8akGQVcCvweuAw8K0kO6rqoXlDP1VV1y5DRknSEIY5Qr8EmKuq71fVz4BPApuWN5YkaamGKfT1wKG+5cPduvn+LMm+JJ9NsmGhF0qyJclsktmjR4+eRFxJ0mJG9aXoF4Hpqno5cB9wx0KDqmp7Vc1U1czU1NSIdi1JguEK/RGg/4j73G7dL1TVY1X1TLf498DvjSaeJGlYwxT6t4ALkpyf5AxgM7Cjf0CSdX2LVwMHRxdRkjSMgVe5VNXxJNcC9wKrgNur6kCSG4HZqtoBvDvJ1cBx4HHgbcuYWZK0gIGFDlBVO4Gd89Zd3zf/fuD9o40mSVoK7xSVpEZY6JLUCAtdkhphoUtSIyx0SWqEhS5JjbDQJakRFrokNcJCl6RGWOiS1AgLXZIaYaFLUiMsdElqhIUuSY2w0CWpERa6JDXCQpekRljoktQIC12SGmGhS1IjLHRJaoSFLkmNsNAlqREWuiQ1wkKXpEZY6JLUiKEKPcnlSR5OMpdk6wLbz0zyqW77riTTow4qSTqxgYWeZBVwK3AFcCHwpiQXzhv2DuDHVfVS4MPATaMOKkk6sWGO0C8B5qrq+1X1M+CTwKZ5YzYBd3TznwVelySjiylJGmT1EGPWA4f6lg8Dr1psTFUdT3IMeCHwo/5BSbYAW7rFp5I8fDKhV9Ba5v8dJv+zxy9lPk2cjrnNvDKay3yKPXLeYhuGKfSRqartwPaV3OepSDJbVTPjzrEUp2NmOD1zm3llmHl4w5xyeQTY0Ld8brduwTFJVgPPBx4bRUBJ0nCGKfRvARckOT/JGcBmYMe8MTuAt3bz1wBfraoaXUxJ0iADT7l058SvBe4FVgG3V9WBJDcCs1W1A7gN+ESSOeBxeqXfgtPm9FCf0zEznJ65zbwyzDykeCAtSW3wTlFJaoSFLkmNsNAXkOSvk+xLsjfJl5P8Zrc+ST7SPeJgX5KLx531WUn+Lsl3ulz/kOScvm3v7zI/nOSPx5mzX5I3JjmQ5H+SzMzbNpGZYfCjMCZFktuTPJpkf9+6FyS5L8l3uz9/Y5wZ+yXZkOT+JA9174v3dOsnNjNAkucm+WaSB7vcH+zWn989CmWuezTKGcsepqqc5k3Ar/fNvxvY1s1fCdwDBLgU2DXurH05/whY3c3fBNzUzV8IPAicCZwPfA9YNe68XbbfAX4L+Bdgpm/9JGde1eV5CXBGl/PCcedaJOsfABcD+/vW/S2wtZvf+uz7ZBImYB1wcTd/NvBv3XthYjN3mQKc1c2vAXZ1/fBpYHO3fhvwzuXO4hH6Aqrqib7FXwOe/eZ4E3Bn9TwAnJNk3YoHXEBVfbmqjneLD9C7XwB6mT9ZVc9U1b8Dc/Qe5zB2VXWwqha6W3hiMzPcozAmQlV9jd5VZ/36H9NxB/AnKxrqBKrqSFXt6eafBA7Suwt9YjMDdH3wVLe4ppsKeC29R6HACuW20BeR5G+SHAL+HLi+W73QYxDWr3S2IfwFvU8ScPpk7jfJmSc52zBeVFVHuvn/BF40zjCL6Z7YehG9o92Jz5xkVZK9wKPAffQ+xf2k7yBrRd4nv7KFnuQrSfYvMG0CqKrrqmoDcBdw7XjT9gzK3I25DjhOL/fYDZNZ41G9cwETd91ykrOAzwHvnfdpeWIzV9XPq+qV9D4ZXwL89jhyrOizXCZJVf3hkEPvAnYCNzDcYxCWzaDMSd4GXAW8rnvjw4RnXsRYMw8wydmG8cMk66rqSHe68NFxB+qXZA29Mr+rqj7frZ7ozP2q6idJ7gdeTe+U7OruKH1F3ie/skfoJ5Lkgr7FTcB3uvkdwFu6q10uBY71fRQcqySXA38FXF1VT/dt2gFsTu+XkJwPXAB8cxwZl2CSMw/zKIxJ1v+YjrcC/zjGLP9PktC76/xgVd3ct2liMwMkmXr2qrIkzwNeT+/8//30HoUCK5V73N8QT+JE7whhP7AP+CKwvv7v2+xb6Z0f+zZ9V2aMe6L3xeEhYG83bevbdl2X+WHginFn7cv1p/TOLT4D/BC4d9Izd9mupHcFxveA68ad5wQ57waOAP/d/Tu/g95jrf8Z+C7wFeAF487Zl/f36Z1O2df3Pr5ykjN3uV8O/GuXez9wfbf+JfQOROaAzwBnLncWb/2XpEZ4ykWSGmGhS1IjLHRJaoSFLkmNsNAlqREWuiQ1wkKXpEb8L7l/QPEXgv3UAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the distribution\n",
    "plt.hist(rate_of_return['Return (%)']);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 285,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-1.4033657997639521"
      ]
     },
     "execution_count": 285,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the skewness Sk\n",
    "rate_of_return['Return (%)'].skew()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Kurtosis in Return Distribution\n",
    "    If a return distribution has positive excess kurtosis (leptokurtosis) and we use statistical models that do not account for the fatter tails, we will underestimate the likelihood of very bad or very good outcomes.\n",
    "\n",
    "`Kurtosis` measures if a return is more or less peaked than a normal distribution. A normal distribution has a kurtosis value of 3. \n",
    "- `mesokurtic`: a distribution that is peaked just like a normal distribution\n",
    "- `leptokurtic`: A distribution that is more peaked than the normal distribution has a $K_E > 3$ (`fat tails`) and `positive excess kurtosis`. *(Most equity returns are `leptokurtic`, meaning that extreme returns are more common than would be expected if they were normally distributed.*)\n",
    "- `platykurtic`: A less peaked distribution. $K_E < 3$ and `negative excess kurtosis`.\n",
    "\n",
    "Sample `excess kurtosis` formula:\n",
    "\n",
    "`Excess kurtosis` = kurtosis - 3. This measures kurtosis relative to the kurtosis of a normal distribution.\n",
    "\n",
    "The formula for calculating kurtosis is given in the curriculum, but it is more likely that you will be given kurtosis measures and asked to interpret them. Consider the following examples:\n",
    "- A distribution with kurtosis of 4 is leptokurtic and has excess kurtosis of 1.\n",
    "- A normal (mesokurtic) distribution has kurtosis of 3 and excess kurtosis of 0.\n",
    "- A distribution with kurtosis of 2.5 is platykurtic and has excess kurtosis of -0.5."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 287,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nOzdeXxU1d348c+ZSSb7vkBIAglkJQkEDEsVFEUREUFciiJVf2pd6tJqVaz2sWrr8zy2PrbWUq2oVesuAiJLraIoKPsmOwQSICGQfU9mMjPn98fNxCwzIcvMJJmc9+uVF5l779z7nYR858z3nHuOkFKiKIqiDHy6vg5AURRFcQ6V0BVFUTyESuiKoigeQiV0RVEUD6ESuqIoiofw6qsLR0ZGyoSEhL66vKIoyoC0Y8eOUilllL19fZbQExIS2L59e19dXlEUZUASQpxwtE+VXBRFUTyESuiKoigeQiV0RVEUD9FnNXRFUVyrqamJgoICGhsb+zoUpQd8fX2Ji4vD29u7y89RCV1RPFRBQQFBQUEkJCQghOjrcJRukFJSVlZGQUEBiYmJXX6eKrkoiodqbGwkIiJCJfMBSAhBREREtz9dqYSuKB5MJfOBqye/O1VyUZRmG4+WsjWvrOXxBUmRTBoZ0YcRKUr3qBa6ojR7fPle/vpVLi99nctfv8rlqc8O9HVIA1pgYGCPnztt2jS333g4Z84cMjMz7e5bv349ISEhZGdnk52dzTPPPNOy78UXXyQzM5OMjAz+8pe/tGz/+OOPycjIQKfTue21qBa6ogCNTRZOVdTz4KUp/PLSZJ5dfYC3Np3AYpXodaps4emWLVt2zjegqVOnsmrVqjbb9u3bx5IlS9i6dSsGg4GZM2cye/ZskpKSyMzMZNmyZdx1112uDL0N1UJXFOBYSS1SQvIQ7Y86OToIk9lKQUV9H0c28NXW1jJ9+nTGjx9PVlYWn376KQD5+fmkpaVx0003kZ6eznXXXUd9fcef9z333ENOTg4ZGRn87ne/a9m+bds2zj//fMaOHcvEiROpqanBYrHwyCOPMGHCBMaMGcM//vGPLsX3wgsv8Nvf/rbbr+3gwYNMmjQJf39/vLy8uOiii1i2bBkA6enppKamdvucvaFa6IoC5BbXApAUrSX0Uc3/Hj1by4iIgD6Ly1me/mw/B05XO/Wco4cF87urMs55nK+vL8uXLyc4OJjS0lImT57MnDlzADh8+DCvv/46F1xwAbfddht///vfefjhh9s8/9lnnyU8PByLxcL06dP54YcfSEtLY/78+Xz44YdMmDCB6upq/Pz8eP311wkJCWHbtm0YjUYuuOACZsyY0enQv//6r//i17/+Nf7+/p2+jk2bNjF27FiGDRvG888/T0ZGBpmZmTzxxBOUlZXh5+fHmjVryMnJ6cJPzzVUC11R0BK6XidIaE7etsSeW1Lbl2F5BCkljz/+OGPGjOHSSy+lsLCQs2fPAhAfH88FF1wAwMKFC9m4cWOH53/00UeMHz+ecePGsX//fg4cOMDhw4eJiYlhwoQJAAQHB+Pl5cV//vMf3n77bbKzs5k0aRJlZWUcPXrUYWy7d+/m2LFjzJs3r9PXMH78eE6cOMGePXu4//77ufrqqwGtFb5o0SJmzJjBzJkzyc7ORq/X9+jn5Ayqha4o2Fri/hi8tDZOiJ830UE+LS33ga4rLWlXeffddykpKWHHjh14e3uTkJDQMr66/dC89o/z8vJ4/vnn2bZtG2FhYdx6662djs2WUvLSSy9x+eWXdym2TZs2sX37dhISEjCbzRQXFzNt2jTWr1/f5rjg4OCW72fNmsUvfvELSktLiYyM5Pbbb+f2228H4PHHHycuLq5L13YF1UJXFLSWeFJU206xpOhAjnpIQu9LVVVVREdH4+3tzddff82JEz/O/nry5Ek2bdoEwHvvvceUKVPaPLe6upqAgABCQkI4e/Ysa9euBSA1NZWioiK2bdsGQE1NDWazmcsvv5yXX36ZpqYmAI4cOUJdXR0AaWlpHWK75557OH36NPn5+WzcuJGUlJQOyRzgzJkzSCkB2Lp1K1arlYgIbUhrcXFxy2tZtmwZCxYs6PHPqrdUC10Z9JosVvJL67g8Y0ib7cnRgXyysxAppbpBpxduuukmrrrqKrKyssjJyWmTWFNTU1m8eDG33XYbo0eP5p577mnz3LFjxzJu3DjS0tLalGcMBgMffvgh999/Pw0NDfj5+fHll19yxx13kJ+fz/jx45FSEhUVxYoVKygtLW1JyF31yiuvAHD33XezdOlSXn75Zby8vPDz8+ODDz5o+T9x7bXXUlZWhre3N4sXLyY0NBSA5cuXc//991NSUsKVV15JdnY2n3/+eY9/jl0huvsinSUnJ0eqBS6U/iC3uIZLX/iWP88fy7xxP35c/temfP7r0/1s+s0lxIT49V2APXTw4EHS09P7OgyH8vPzmT17Nvv27XP5tVatWsXx48d54IEHXH4tZ7L3OxRC7JBS2u15VS10ZdCz1cmTo4PabE9qfpxbXDsgE7ryo9mzZ/d1CG6haujKoHf0rJbQR0a1HZ7YMtJF1dFdIiEhwS2t88FEJXRl0MstqSUuzA9/Q9sPrJGBBkL9vVXHqDJgqISuDHpHz9a2tMZbE0KQFBWoWujKgKESujKoWaySY3aGLNokRauErgwcKqErg1phRQNGs7VlDpf2kqIDKa8zUV5ncnNkitJ9KqErg1puSQ2A3ZJL6+2qlT7w5OXlMWnSJJKSkpg/fz4mU8c35bKyMi6++GICAwO577772ux74okniI+P7zAL4yuvvEJWVhbZ2dlMmTKFAwf6zzTLKqErg5pthEtSVJDd/clDtO1Hi2vcFpMCZrO51+dYtGgRDz74ILm5uYSFhfH66693OMbX15ff//73PP/88x32XXXVVWzdurXD9gULFrB37152797No48+ykMPPdTrWJ1FJXRlUMstriUqyIcQf/srqw8L8cXfoFct9G7Kz88nPT2dn//852RkZDBjxgwaGhoAbUKsyZMnM2bMGObNm0dFRQWgLWrxq1/9ipycHF588UWmTZvGgw8+SE5ODunp6Wzbto1rrrmG5OTkc051K6Xkq6++4rrrrgPglltuYcWKFR2OCwgIYMqUKfj6+nbYN3nyZGJiYjpsbz2vS11dXb+6i1jdWKQMakeLa0l2UG6B5pEuntAxuvYxOLPXueccmgVX/K/D3UePHuX9999nyZIl/PSnP+WTTz5h4cKF3Hzzzbz00ktcdNFFPPnkkzz99NMtK/2YTKaW1X0+++wzDAYD27dv58UXX2Tu3Lns2LGD8PBwRo0axYMPPtgyn0p7ZWVlhIaG4uWlpbi4uDgKCwud9tIXL17MCy+8gMlk4quvvnLaeXtLtdCVQUtKSe45EjpodfQjZ1XJpbsSExPJzs4G4LzzziM/P5+qqioqKyu56KKLAK3l/O2337Y8Z/78+W3OYZs3PSsri4yMDGJiYvDx8WHkyJGcOnXKTa+ko3vvvZdjx47x3HPP8Yc//KHP4mhPtdCVQaugooFao5m0mOBOj0sfGsyynYVU1JkICzC4KTon66Ql7So+Pj4t3+v1+paSS2cCAtrerWs7h06na3M+nU7XaZ09IiKCyspKzGYzXl5eFBQUEBsb292XcE433HBDhwnF+pJqoSuD1qEzWqs7daj9DlGbtJigNscrPRcSEkJYWBgbNmwA4F//+ldLa72npk+f3qGcIoTg4osvZunSpQC89dZbzJ07t1fXsWm9YMbq1atJTk52ynmdQSV0ZdA6VKQtyZY6pPOEbkv4h844dwm3weqtt97ikUceYcyYMezevZsnn3yyx+eyWq3k5uYSHh7eYd9zzz3HCy+8QFJSEmVlZS2LUKxcubLNNRMSEnjooYd48803iYuLaxmG+OijjxIXF0d9fT1xcXE89dRTAPztb38jIyOD7OxsXnjhBd56660ex+9savpcZdC6992d7DtdxTePXNzpcVJKcv7wJZemD+G568a4Kbre6+/T5zrDvn37eOONN3jhhRf6OhSX6O70uaqFrgxaB89Uk3aOcgtoH9/TYoJUC70fyszM9Nhk3hMqoSuDUoPJQn5pHWlDO+8QtUkbGszhszVYrH3ziVZRukIldGVQOlpcg1XSpRY6aMc1Nlk5UVbn4sgUpedUQlcGpUNF2oiVcw1ZtElvPk6NdFH6sy4ldCHETCHEYSFErhDisU6Ou1YIIYUQdgv2itJfHDpTg5+3nuHh/l06Pik6EJ1QCV3p386Z0IUQemAxcAUwGrhRCDHaznFBwC+BLc4OUlGc7dCZalKGBqHXdW0eDl9vPYmRAS1DHRWlP+pKC30ikCulPC6lNAEfAPZG6P8eeA5odGJ8iuJ0UkoOFlWT3sX6uU1aTLBqobvQ9u3beeCBBwAwGo1ceumlZGdn8+GHH/ZJPJ1Nrdvaxx9/TEZGBjqdDntDsU+ePElgYGCbGR1ffPFFMjMzycjIaJnHxhm6cut/LNB60oQCYFLrA4QQ44F4KeVqIcQjTotOUVygpMZIRX1T2ztErRbY9wmUHvlx25AMSJ8LOq3dkz40iNU/FFFrNBPoo2bNcLacnBxycrRq7a5duwBtZsauslgs6PV6p8Vjm1p33759nS5mnZmZybJly7jrrrvs7n/ooYe44oorWh7v27ePJUuWsHXrVgwGAzNnzmT27NkkJSX1OuZed4oKIXTAC8Cvu3DsnUKI7UKI7SUlJb29tKL0yMHmVnba0GCQEo78B16ZAst+Dt8+Dxv+T/v341thyTQ4/s2PxwOHVSv9nPLz88nMzGx5/Pzzz7fcaTlt2jQWLVrExIkTSUlJaZkGYP369cyePZvi4mIWLlzItm3byM7O5tixY6xbt45x48aRlZXFbbfdhtFoBLS7PBctWsT48eP5+OOPSUhI4De/+Q3Z2dnk5OSwc+dOLr/8ckaNGsUrr7zSrdfQ2dS6raWnp5Oammp334oVK0hMTCQjI6Nl28GDB5k0aRL+/v54eXlx0UUXsWzZsm7F5khXmhmFQHyrx3HN22yCgExgffO8wEOBlUKIOVLKNp8/pJSvAq+CdqdoL+JWlB6z1cHTIr3h/RvhyFoIS4Tr34LRc0EIsFph78fw1e/h7TmQcQ1pF2sfjQ+dqea8EWF9+RK67bmtz3Go/JBTz5kWnsaiiYt69Fyz2czWrVtZs2YNTz/9NF9++WXLvujoaF577TWef/55Vq1aRWNjI9OmTWPdunWkpKRw88038/LLL/OrX/0K0Cbi2rlzJwCPPfYYw4cPZ/fu3Tz44IPceuutfPfddzQ2NpKZmcndd9/d+xfeRbW1tTz33HN88cUXbcotmZmZPPHEE5SVleHn58eaNWtaPpn0VlcS+jYgWQiRiJbIbwAW2HZKKauASNtjIcR64OH2yVxR+otDZ2qIDfImbM1dWjK/7BmYdA94tZpJUaeDsfO1BP/di7D+v4nV6Qn2ubZlyKPSc9dccw3w47S6nTl8+DCJiYmkpKQA2pS7ixcvbknonU25W1tbS1BQEEFBQfj4+FBZWUloaKiTX419Tz31FA8++GCHJezS09NZtGgRM2bMICAggOzsbKeVis6Z0KWUZiHEfcDngB54Q0q5XwjxDLBdSrnSKZEoipscOl3Jc4ZX4fAXcMWfYNKdjg/29oVpi0DvhVj3DM8FNPJGUf+ZLrWretqS7ikvLy+sVmvL48bGtmMlbFPh6vX6Xi8356wpd5cvX87TTz8NwGuvvdbrVvOWLVtYunQpjz76KJWVleh0Onx9fbnvvvu4/fbbWyYLe/zxx4mLi+vVtWy61LMjpVwDrGm3ze4UaVLKab0PS1Fcw2i2cF35P5ii/wKmPd55Mm9tykPQUMEV37/E8TpfrNbz0XVxyONgNGTIEIqLiykrKyMwMJBVq1Yxc+bMHp0rNTWV/Px8cnNzSUpKcsqUu2lpaRw61LYENW/ePObNm9er87Zm6xsArbXeerRMcXEx0dHRnDx5kmXLlrF582anXFPdKaoMKqe3r+F2/RryRy6Aix7t+hOFgMt+T17cHO5hKWf2f+O6ID2At7c3Tz75JBMnTuSyyy4jLS2tx+fy9fXln//8J9dffz1ZWVnodLpe1cJLS0vp6iyzjqbWveOOO1qGKC5fvpy4uDg2bdrElVdeyeWXX37O81577bWMHj2aq666isWLFzutDKSmz1UGD1MdtX+ewNk6K16/+I4RQ+2vR9mZA/mnCf7nVIKDQwn+5aa2dfd+ZjBMn9sTq1at4vjx4y1j3vszNX2uojiy/n8JbCjkD7q7GD6k44IIXZEUN5SnLbcTXJOrdZYqA87s2bMHRDLvCZXQlcGhaA9sWsznPpfTOOwnNA+x7TaDl46zQy9kk99F8O0fofTouZ+kKG6iErri+axWWPkA0j+CJ2qvJysupFenyxgWwm8abkJ6+8Fnv9JuTuqn+qqkqvReT353KqErnu/gp1C0m8KJj1Nq8SdjWNemzHUkKzaE/MZAKib/Bk5shGPrnBSoc/n6+lJWVqaS+gAkpaSsrOycd6m2pyakUDyb1Qrf/BEiU/jebxpwgKzY3rXQbc/fHDqLWSGLYf1zMGq6NhKmH4mLi6OgoAA1zcbA5Ovr2+3x6SqhK57t0CooPgDXLOGH47UE+niREBFw7ud1ImVoIN56wQ9FDcya8iCsfgiOr4dRnS827W7e3t4kJib2dRiKG6mSi+K5pNRa5+GjIOMa9hVWkzEsuNc3BPl46UkZEsT+01UwbiEEDYNvnuvXtXRlcFAJXfFch9fC2b1w4SOY0XGwqJrMXpZbbLJiQ9hbWIXUG2DKg3ByE+RvdMq5FaWnVEJXPJOUWqs5LAGyrudocS1Gs7XX9XObjNgQKuubKKhogPE3Q+BQ7XqK0odUQlc8U/5GKNqtzcGi92JfYRWAU1vogFZ28faFCx6A/A1wepdTzq8oPaESuuKZtr0GfmEw5qcA7CusIsCgZ2Rk7zpEbdKa1yPd2/xGwbiF4O0P2153yvkVpSdUQlc8T3WRNrpl3ELw9gNgb2EVo53QIWrj660nOTqQvYXNi0b7hmhvHnuXQkOFU66hKN2lErrieXa+DVYznPf/ADBbrBxwYoeoTVZsCPsKq368cSfndjA3wO73nXodRekqldAVz2Ixw443tRt9IkYBcORsLY1NVsb08pb/9sbEh1JeZ9I6RgFixkDcRNj+uhrCqPQJldAVz3JkLdSchgl3tGzadUorgYwf7tx1QMfFa3NY7zzZqsQy4Q4oy4U8NV+64n4qoSueZdtrEBIPKT8uMrDrZCURAQaGh/s79VJpQ4Pw89az62TljxtHzwX/CC0ORXEzldAVz1F2TLsF/7xbQPfjoru7TlYwbnhoj6fMdcRLr2NMXAi7TrVK6N6+WmfsoTVa56yiuJFK6Irn2PMBCB1kL2zZVFXfxLGSOsY5udxiM254GAdOV9HYZPlx4/hbQFpg78cuuaaiOKISuuIZrFb44QNIvAiCY1o27y7QWs+2erezjRseSpNFsv909Y8bI0ZBbI72BqMobqQSuuIZTm2GypMw9oY2m3edrEAIbUSKK9jeKHadbDf2fOwNULwfzux1yXUVxR6V0BXPsOcD8A6AtNltNu86WUnqkCACfVwzU3R0sC+xoX5t6+gAGdeAzlu10hW3UgldGfiaGmH/Cki/CnwCWzZbrbKlQ9SVxg0PZdeJdi30gAhInqHV0S1ml15fUWxUQlcGviP/BmMVjJ3fZvPx0jqqG82Mi3dNh6jNuOFhnK5q5ExVY9sdY+dD7VnIW+/S6yuKjUroysD3w4cQFKN1iLZiq2u7o4UOsPtUu1Z6ykxtjpc9H7r0+opioxK6MrDVlcHR/0DW9W3GngPsOlVJkK8Xo6ICHTzZOTKGBWPQ69reYATg5aPV0g+tAmONS2NQFFAJXRnoDn6qTcQ1Zn6HXbtOVpIdH+q0GRYd8fHSkxEb3DGhgxZXUz0c/rdLY1AUUAldGej2L4eIZBiS0WZzndHM4TPVLruhqL1x8WH8UFhJk8Xadkf8JG3N0f3L3RKHMriphK4MXLUl2spEGVdDu9v6d52sxCphvIvr5zbnjQijscna9gYjAJ0ORs+B3C+hsdr+kxXFSVRCVwaugytBWiFjXoddW/PK0AnISQh3SygTEsNarttBxjywGOHI526JRRm8VEJXBq4DK7RyS/ToDrs255WTGRvishuK2osO8mVkZABbjpd33Bk3UZVdFLdQCV0ZmFrKLfM6lFsamyzsPlXJRDe1zm0mJoazNb8ci7Xd4hY6nTatriq7KC7WpYQuhJgphDgshMgVQjxmZ//dQoi9QojdQoiNQoiOTSZFcaaWcsvVHXb9UFCFyWxl0sgIt4Y0aWQ4NY1mDp+xM0Qx42pVdlFc7pwJXQihBxYDVwCjgRvtJOz3pJRZUsps4I/AC06PVFFas41usVNusdWxJyS4Z4SLzcTEiDbXb0OVXRQ36EoLfSKQK6U8LqU0AR8Ac1sfIKVs/TkyAFALKiquU1sCJ76zW24B2JJXTtrQIEL9DW4NKzbUj9hQP7bm26mjq7KL4gZdSeixwKlWjwuat7UhhLhXCHEMrYX+gL0TCSHuFEJsF0JsLykp6Um8igKHPtPKLaPndthltljZcaKCiYnurZ/bTEoMZ2teOdLeItG2ssvR/7g/MGVQcFqnqJRysZRyFLAI+K2DY16VUuZIKXOioqKcdWllsDm4CsISO9xMBLD/dDX1JkufJfSJieGU1po4XlrXcWfcRAiIhoOfuT8wZVDoSkIvBOJbPY5r3ubIB0DHnipFcYbGKsj7FtJn2y23bM3Tyh3uHuFiY3sjscXRhk4HabO0sktTY8f9itJLXUno24BkIUSiEMIA3ACsbH2AECK51cMrgaPOC1FRWjn6BVibIO0qu7u35JWTGBlAdLCvmwPTJEYGEBnoYz+hgxa3qRbyvnFvYMqgcM6ELqU0A/cBnwMHgY+klPuFEM8IIeY0H3afEGK/EGI38BBwi8siVga3g59pZYu4CR12Wa2SbfnlfdY6BxBCtNTR7UqcCoYgVXZRXKJLt9FJKdcAa9pte7LV9790clyK0lFTo1auyLpOK1+0c/hsDVUNTUzoo/q5zcTEcFbvLeJUeT3x4f5td3r5QMoMOLwWrJYOU/4qSm+oO0WVgSPvG61c4aDcsvFoKQAXJLn3hqL2zh+lXX9jbqn9A9JmQ30pnNrixqiUwUAldGXgOPiZVq5InGp398bcUkZFBRAT4ufmwNpKig5kSLCP44SedCnoDdpoHUVxIpXQlYHBatHKFCkztLJFO0azhS15ZUxN7vvhsEIIpiRF8X1uKdb287oA+AbDyGnaSkb2xqsrSg+phK4MDKe2aGWKtCvt7t5xooLGJisXJEW6OTD7piZHUlHf1HF+dJu0K6HyBJzd597AFI+mEroyMBxaDTpvSLrM7u7vckvR6wSTR/Zth6jN+UnnqKOnzgIEHFpjf7+i9IBK6Er/JyUcXgOJF2rlCjs2Hi1lXHwoQb7ebg7OvuggX9KGBrEx18EUF4HNQy8Pq4SuOI9K6Er/V3oUyo9D6hV2d1fWm/ihsKrflFtspiRFsi2/gsYmi/0DUq+Aot1Q1dmN14rSdSqhK/2frRWbOsvu7k3HypBSq1v3J1OSIzGZrWyzN/si/NgfcOTf7gtK8WgqoSv93+E1EDMWQjpM8gnAhtxSAn28GBvvngWhu2piYjgGva5lfHwHkSkQPlKVXRSnUQld6d9qS+DUVoetc9Dq55NHhuOt71//nf0NXowfEcoGRwldCO115X0LRjurHClKN/WvvwBFae/o54B0WD8/WVbPyfJ6pvSz+rnN1OQoDhRVU1prtH9A6hVgMcGxr9wbmOKRVEJX+rfDayE4DoaOsbv7myPFAExN6fsbiuyx1fU3HHUw2iV+MviFaa9TUXpJJXSl/2pq0FquqVfYnfscYN2hYhIjAxgVFejm4Lomc1gIUUE+rDtYbP8AvRckX64tHm0xuzc4xeOohK70X3nfQlO9w3JLvcnM98fKuCQt2s2BdZ1OJ7gkNZpvjpTQZLHaPyj1Cmgoh4Kt7g1O8TjC7tqHbpCTkyO3b9/eJ9dWuq6ysZKHv32YOpOdJdVcfvET0FABMWOw1/aoamwiv7SOkVGBBPl0aSboPlHV0ER+WR2jogIJtBentEDRHu1mo+A4t8cX7BPMn6f9GX9v/3MfrPQ5IcQOKWWOvX39969A6ReOVh5lS9EWsiKzCPEJcd+FJVD/A/iFg6/92/nPVtagk4K44EhHFZl+IcggyS8Bo9GXuBAHpSGfUKivguhMt8ZW3ljO96e/50T1CdIj0t16bcX5VEJXOmWymAB4ZMIjjIse574LF+6E796Def8FY2/osFtKyaT/XseFCWH8/bLz3BdXD91csJWC0/W8vGCa/QO2LoE1D8PVn0Bksv1jXGBj4Ubu+fIejBYHo3CUAUXV0JVO2RK6QW9w74UPrwWhg+QZdnfvP11NcY2RS9KGuDeuHpqeFs3x0jqOl9TaPyDlcu1fN4928dFrUxE3WZvcel3FNVRCVzpltGotN4POzQn9yFptSJ+//XLLuoPFCAHTUvvncMX2bB23Xx1yMNoldDgMyXL7NADeOm0yM9VC9wwqoSudarJoLTdbS84tKk/Bmb2QOtPhIV8dOkt2fCiRgW6Mqxfiw/1JGRLoOKGD9npPboJ6B3O/uIDt92r7JKYMbCqhK52ytdzcWnKxtVId3O5fXNPInoIqpvfj4Yr2XJI2hK155VQ3OihvpF4B0gpHv3BbTCqhexaV0JVO9UkN/fBaCB/lsHPwq+abdAZK/dxmeno0Zqtk/WEHd43GjIPAIW6drMtbr5VcTFaV0D2BSuhKp2wJ3W0lF2MN5G9weDMRwJp9ZxgR4U96TJB7YnKS8cPDiAz04d/7iuwfoNNBykzIXQdm9yRY2+9V1dA9g0roSqdsLTe3dYoe+0qbrKqTxSy+zy3liswYRH8efG6HXieYmTmErw+VUG9ycJt/6hVgqoETG90Sk+33qkounkEldKVTRosRgcBL56ZbFg7/G3xDtREudnxx4Cxmq37ZPKoAACAASURBVGRW1lD3xONkszJjaGiy8I2jskviReDl57bhi7ZSmkronkEldKVTJosJg97gntawxax1iKZcrk1aZcfafWeIDfUjK9aNd6060cTEcMIDDKzZd8b+AQZ/GHWxltDdMC2HLaGrkotnUAld6ZQtobtFwVZtkioHo1uqG5vYcLSEWVlDB1y5xcZLr+PyjCF8dfBsJ2uNzoKqU3B2n+vj0XmhF3rVQvcQKqErnTJajO7rED20GvQGSJpud/eXB87SZJHMyopxTzwuMisrhjqThW+POCi7pFwOCDjkntEuBr1BJXQPoRK60qkma5N7OkSl1IbrJUwFH/ujV9bsPcOwEF+y+9naod01eWQEof7erHVUdgmMhrgJbhu+aNAb1LBFD6ESutIpo8XonpJL6REoPw5p9sstNY1NfHu0hJkDcHRLe956HTNGD+HLA2cxmh2UXdJmQdFuqCp0eTw+Oh/VQvcQKqErnXJbDd3WGk2xP1zxq0PFmMzWATu6pb0rsmKoMZrZcMTBAtK2foQjrh/t4q33VgndQ6iErnTKZDG5p4Z+aA3EjIWQWLu7l+8qZFiIL+OHh7k+Fje4YFQkYf7erNjtoAUemQLhI91SR/fR+6hRLh5CJXSlU0aLsWVGPpepLYaCbZB6pd3dJTVGNhwtZe64WHS6gV1usTF46Zg9ZhhfHDhrf24XIbRWet630Fjt2lhUp6jH6FJCF0LMFEIcFkLkCiEes7P/ISHEASHED0KIdUKIEc4PVekLJqsbWuhH/g1Ih3eHrtxzGotVcs04+633gWre+FiMZiv/3uugczR1Flib4Ng6l8ahOkU9xzkTuhBCDywGrgBGAzcKIUa3O2wXkCOlHAMsBf7o7ECVvuGWksuhNRASD0Oz7O5evquAzNhgkocMrLlbzmVcfCiJkQEs21Vg/4D4SdoSfC4uuxh0BlVy8RBdaaFPBHKllMellCbgA2Bu6wOklF9LKeubH24G3L/SreISJoupZUY+lzDWavO3pF2JvYVBj56tYV9hNfPGed5/KSEEV2fHsvl4OYWVDR0P0Htpn1qOfO7Sybp89D4t894rA1tXEnoscKrV44LmbY7cDtjtmhdC3CmE2C6E2F5S4uCmCqVfcfmNRblfgsUIabPt7l62qxC9TjBn7DDXxdCH5jWXkVbsctA5mjYbjFXaDJQuYtCrFrqncGqnqBBiIZAD/Mnefinlq1LKHCllTlTUwFg6bLBrsjS5dtjiodVaWWH4Tzrssloln+4qZGpyJFFBA2Nlou4aHuFPzogwlu8qRNqbu2XUxeAdoP2cXETV0D1HVxJ6IRDf6nFc87Y2hBCXAk8Ac6SU6u3eQxitRtfdKWo2aeWE1Fl2J+PanFfG6arGllasp5o3Ppbc4lr2FdoZzeLtp02FcGg1WK0uub6PXt1Y5Cm6ktC3AclCiEQhhAG4AVjZ+gAhxDjgH2jJvJNFE5WBxqWdovkbtHJCmv3hiku3FxDo48WM0Z5xM5Ejs7OGYfDS8dH2U/YPSJsNtWegcIdLru+t81YlFw9xzoQupTQD9wGfAweBj6SU+4UQzwgh5jQf9icgEPhYCLFbCLHSwemUAcald4oeWqWVE0Zd3GFXZb2JVXuLuHrcMPwMetdcv58I8ffmyqwYVuwqtL/wRcoM0HnBoc9ccn3VQvccXaqhSynXSClTpJSjpJTPNm97Ukq5svn7S6WUQ6SU2c1fczo/ozIQmK1mLNLimoRutWrD8ZKma2WFdj7ZWYjJbGXBxMFxS8ONE4dTYzSzao+d5en8wrRJyw6ucskc6erGIs+h7hRVHHLpAtGFO7QyQvpVHXZJKXlvywmy40MZPSzY+dfuhyYkhJEUHci7W0/aPyB9NpQfg5LDTr+2rVPUbqesMqCohK445NIFog99ppURkmd02LU1r5xjJXUsmDTc+dftp4QQLJg4nD2nKtl/uqrjAbZpEVxQdvHR+2CVVszSwTqnyoChErrikG0om9PncpESDqyExAvBr+Pc5u9vPUmQrxdXjfHMseeOXDs+Dh8vHe/ba6UHx0DcRO3n5mS2UUzq5qKBTyV0xSHbyAent9DP/AAVeTD66g67KupMrNl3hnnjYj2+M7S9EH9vrhwTw4pdp6kz2mktj56r/ezKjzv1umpdUc+hErrikMtKLvtXgNDbvTv0k50FWmfoICq3tHbTpOHUGs2s3HO6487RzTNu7F/h1GuqhO45VEJXHLIldKfO5SIlHFihlVsCItrsMlusvPl9PhMSwkgbOjg6Q9sbPzyM0THB/PO7vI6dlKHxEJsDBz516jVtb9iq5DLwqYSuOOSSksvZfVrJYPTcDrv+c+AsBRUN3D5lpPOuN8AIIbh9SiJHztay4aid1YxGz9WWpivPc9o1bW/YqoU+8KmErjjUZNVabE699d9WbrEzXPH1jXkMD/fnstFDnHe9AeiqscOIDvLhtY12krbtjdCJrXQfnfaGreZzGfhUQlccsrXYnDYO3VZuSZgCAZFtdu08WcGOExXcdkECeg9ZlainDF46bjk/gW+PlHD4TE3bnWEjYNh47eforOs1/37VzUUDn0roikNOv7Ho7H4oy7Vbbnl9Yx5Bvl5cnxNv54mDz4KJw/H11vGGo1b66V1Qke+Ua6mE7jlUQlcccvoolwOfgtBBetuZIQoq6lm7t4gFE4cT4NNx1sXBKCzAwLXj41i+u5DS2na17Yzm4Z5OGpNu+/2qGvrApxK64pCtpuqUFrqUsH85jLgAAtvOhf/P7/IRQnDL+Qm9v44HuW1KIiazlbc3nWi7IywBYrJh/zKnXKelha5q6AOeSuiKQy01dGd0ihbthrKjkHVdm82ltUbe23KSuWOHMSy04yRdg9moqEBmjB7Cm9/lUd3Ybkhh1nVa2aU0t9fXUSUXz6ESuuKQU0sue5eCzrtD/fy1DXk0mi3ce0lS76/hgR6Ynkx1o5m3v89vuyPzWkDAvqW9vobtDVuVXAY+ldAVh5zWKWq1aAk9eYY2FWyz8joTb2/K56oxwxgVFdi7a3iozNgQpqdF89rGPGpbTwcQPEwbLfTDR72eUtf2hq1a6AOfSuiKQ04btpi/UZsqt1255Y2NeTQ0WbhPtc47df/0ZCrrm3hnc7taetb12pS6p3f16vyq5OI5VEJXHDJZTOiEDi9dL0ee7P0YDIGQMrNlU1V9E29+n8+szBhShgT1MlLPlh0fyoUpUSz59njbFY1GzwG9Qfv00wuqU9RzqISuOOSU9USbGrXhdelXgcG/ZfPr32klBNU675pfTk+irM7Eu5tbTa3rF6aVsfZ9opW1ekjV0D2HSuiKQyarE9YTzf1CWwi6VbmltNbI6xuOMzNjKOkxg3MSru46b0Q4U5Mj+fv63LYjXrKu08pZ+Rt6fG69To+X8FKTc3kAldAVh0wWU++HLO79GPwjIXFay6aX1h2l0WzlkZmpvTv3ILNoZhoV9U28sv7YjxtTZmrlrB8+7tW5DXqDaqF7AJXQFYeMFmPvWuj15XB4rTbETq/V4fNK63h3y0lumBCvRrZ0U2ZsCHOzh/HGd3mcqWrUNnr7aXfeHvgUTHU9PrdK6J5BJXTFIZOllyWXvUvBYoJxC1s2Pf/5YQxeOn55abITIhx8Hp6RitUKf/7iyI8bxy0EU02vpgIw6A0ts2sqA5dK6IpDve4U3f0ODM2CmDHaw1OVrN5bxB1TRxId5OukKAeX+HB/Fk4ewcc7TnH0bPNMjCPOh7BE2P1uj89r0KkWuidQCV1xyGTtRQ296Aco2gPjfgaAlJL/Xn2QyEADd144eBewcIb7L0kiwMeLZ9cc1FY1EgLG3aR1jPZwvVEfvY8ah+4BVEJXHOpVDX33u9oY6azrAVixu5Ct+eX8ekYqgWpGxV4JCzDwy+nJrD9cwpcHi7WNYxcAAna/16NzGvQGldA9gEroikNNlqaelVzMRvjhQ0i7EvzDqW5s4tnVhxgbH8p8Nd+5U9xyfgIpQwJ5auV+GkwWCImFpOlaQu/BmHSV0D2DSuiKQ0aLsWcLRB9eCw0VkK11hv75iyOU1Rn5w9xMdIN8NSJn8dbreGZuJoWVDby8vnnGxeyboLoQjq/v9vl89D6qhu4BVEJXHDJajD1roe96B4KGwaiLOXC6mre+z+emScPJigtxfpCD2OSREVydPYxXvjlOXmmd9onILwx2/avb5/LWe6sWugdQCV1xqMna1P1O0YoTkPslZC/Ago7frthLqL+Bh2eom4hc4fFZ6Ri8dDz56T6k3gBj5sOh1VBb0q3z+Oh81FwuHkAldMWhHnWKbn9DG3WR8/948/t8dp6s5IlZ6YT6O2ldUqWN6GBfHrk8lQ1HS/l4ewHk3K6N/d/5VrfOo2ronkEldMWhbt9Y1NQIO9+G1FnkNYXxp88PcUlaNNeMj3VdkAo/mzyCSYnh/H7VAYoM8ZB4kfbGajGf+8nNVEL3DCqhKw51+8ai/cugoRzrhJ/z6NI9GPQ6/ueaLIRQHaGupNMJ/njdGMxWyWOf7EVO/LnWOXpkbZfPoTpFPUOXEroQYqYQ4rAQIlcI8Zid/RcKIXYKIcxCiOvsnUMZWKSU3Z9tcesSiEzlzdPD2ZZfwZNXZTAkWN0R6g4jIgJYNDOVb46UsLQmE4LjtN9HFxn0BlVD9wDnTOhCCD2wGLgCGA3cKIQY3e6wk8CtQM/ualD6HbM0Y5XWrneKFuyA0zs5m7aQ5z4/zCVp0VyrSi1udfNPEpiUGM7Tq49QMXoh5H0DJYe79FxVcvEMXWmhTwRypZTHpZQm4AOgzUq/Usp8KeUPgNUFMSp9wDY3dpdLLtuWIA2B3LknmSBfb567dowqtbiZTif4v5+ORSfggSNZ2qiXba916bkGnZbQZS/XJ1X6VlcSeixwqtXjguZt3SaEuFMIsV0Isb2kpHvDqhT3stVTu3RjUc1Z2LeMzYGX8UOplb/MzyYqqJcrHSk9Ehfmzx+vG8uG04J9IZdod442VJ7zeT56HyQSs7XrHalK/+PWTlEp5atSyhwpZU5UVJQ7L610ky2hd6mFvuVlpMXEb4qmcM9Fo5iSHOni6JTOzMwcys0/GcFjRReCqRa2v37O59j6SlTH6MDWlYReCLSegCOueZviwbpccmmswrL1NT6Xk4gYPpqHLktxQ3TKuTw+Kx05dAzfMRbLppehqaHT49VC0Z6hKwl9G5AshEgUQhiAG4Cez6SvDAhdLbk0bHoNvamGd7yu4aUbx+GlVyNh+wNfbz1/v2k8bzAXfX0Jxu3vdHq8rfNbdYwObOf865NSmoH7gM+Bg8BHUsr9QohnhBBzAIQQE4QQBcD1wD+EEPtdGbTieraWmo/OcQvdbKzHuOFvfGfN5Fc/u55hoX7uCk/pgoTIAG658Wfsto6iet3/YTU7XpGopYWuEvqA1qXmlJRyjZQyRUo5Skr5bPO2J6WUK5u/3yaljJNSBkgpI6SUGa4MWnE92x92Z+PQP3//RUKt5TT95AFyEsLdFZrSDRemRlM69h6izEX8++N/ODzOVlpTNfSBTX0+Vuw6V0J/5/tjpB9/k0K/VKbN/Kk7Q1O6afq82yg2DGf4wSV8sv2U3WNUDd0zqISu2GVrqdlL6Gv2FrFr9auM1J1h6JWPa5NxKf2W0OkJu/xRMnX5fLn8DdYfLu5wjCq5eAaV0BW7bH/Y7Ue5bD5exiMfbGORzzKsQ7PRZ8y193Sln/HOvhFLRDKP+Szl3ne2s/tU27Hptk5RVXIZ2FRCV+yyffRu3UL/oaCSn7+9nbuDNhBtLUZ36X+p1vlAofdCf8lvGWE9yY1+m7jtzW0cPlPTstv2xq1a6AObSuiKXS0ll+aW296CKha+toUhvhZ+oVsOI6bAqOl9GaLSXelzIGYsi3yW4yvMLFiymSNntaSuSi6eQSV0xa7WJZe9BVXc9Npmgny9WTpuD/r6Epj+pGqdDzQ6HUx/Eu+aU3z6k1z0OsGCJZs5erZGJXQPoRK6YpftD/tQUT0LX99CkK83H/0shdCdL0PKTBg+qY8jVHpk1HQYMYWonX/lg1szEUJw45LNnCzVxqirGvrAphK6Ypethn7nW3sI9ffmgzsnE7vzT9rcINN/18fRKT0mBFz2DNQVM3L/3/jgzsl463Xc9+4eQFtHVhm4VEJX7Np1SpsNc2RUCEvvPp/4+oOw4y2YdDcMaT8dvjKgxJ0H42+GzS8zSp7ik3vOJzooEIDdBWoW1IFMJXSlDatV8pcvj/D5gQKQej74+flEBXjB6ocgcAhM67BglTIQTX8KDIGw+mGGhfjy3h0XALB0Zz5Lvj2u5kUfoFRCV1rUGc3c+95O/vLlUZKH+uHv7UOQrzfseBOKdsPlz4JvcF+HqThDQARc+js4sRH2LiUqUGuhpwz15dk1B/n1R3tobLL0cZBKd6mErgBwsqyea1/+ns/3n+G3V6bzk1HB2siH2hJY9wwkTIXMa/s6TMWZxt8Cw8bDf55A11iNl86LS0dH8OvLUli2q5D5r26mqKrzaXeV/kUldIXVPxRx5V83cLqygX/+v4ncMXUkTdYmLaF/9oA2l/as59UwRU+j08PsF6CuFNY+io/eB5PVxP3Tk3ll4Xnknq1h1osbWHfwbF9HqnSRSuiDWGOThd+u2Mu97+1kVHQgqx+YykUp2kpSRosRQ5MRDq/RPppHp/VxtIpLDBsHFz0KP3yIj/xxuOrMzKF8dv8UYkL8uP2t7fxh1QFMZrVkcH+nEvogtfNkBVf+dQPvbD7JXReO5OO7f0J8uH/LflNDBT71ZZB4IUy6pw8jVVxu6sMQex7exhpMjT/O8TIyKpBlvzifm38ygtc25jHnbxvZf7qqDwNVzkUl9EGmscnC/6w9yHUvf0+DycK/bp/Ib2al4916pSGLGVPhDgwIuPpl7Q5DxXPpveCaJfhYJcYT34H1x5a4r7eeZ+Zm8trNOZTVmZj7t+/48xdHVGu9n1J/qYPI14eLufwv3/KPb44zf0I8nz94IVOT7SzWve4pTMZKDKEjICTO/YEq7hcxCkNgNE11xfDtnzrsvnT0EL548EJmj4nhxXVHmfXXDXx/rLQPAlU6oxL6IFBQUc9d/9rO//vnNvRC8M7tk/ifa8ZoQxLb2/0+fP8SxqBh+AQNc3+wSp8xBA7BGBwL6/8bDn7WYX+ov4G/3DCO12/JwWi2sGDJFh54fxdnqhr7IFrFHq++DkBxnar6Jhavz+XN7/LR6eCRy1O5Y2oiPl56+084tU0b1ZIwFVNIIH7nWCBa8SwGnQFjZBI0eMGyu+D2RBia2eG46elDuCApkr+vP8Yr64/xnwNn+PnUkdx54Uj7jQTFbVQL3QPVGc288s0xLvzT1yzZcJw52cP46tfTuPfiJMfJvKoAPrwJgmLgp29jsjZ1ukC04nl89D40SQvMf1e7gez9G7X7EOzw9dbz0GUpfPnQRVw2eigvfZXLtD+t5/WNeeqGpD6kEroHqTWaWfx1LlOe+4r/XXuIsfGhrL5/Ks9fP5ZhoX6On1hzBt6aA6Z6uPED8A/HZDF1ukC04nm89d7abIvBMXDDu1BXAm/Phfpyh88ZHuHPSzeO49N7LyBlSBC/X3WAKc99zWsbjtNgUond3VTJxQMUVTXw5nf5vLf1JDWNZqalRvHA9GTGDw8795NrS7RkXnMGfra8ZeItldAHH9uNRQDEngc3vg/vzYd/XQ03rwS/UIfPHRsfyvt3Tmbz8TL+uu4of1h9kMVf57Jw8gh+9pMRRAf5uulVDG4qoQ9QUkq25JXz7paTrN1bhFVKrsiK4c6pIxkb7/gPr426Mq0FVnkSFi5tM8e50WJUCX2QMegNbRe4GHWx1lL/YAG8c432hu8b0uk5Jo+MYPLICLbnl/Pqt8f529e5/OOb48weE8NNk0cwfngoQt1x7DIqoQ8wxTWNrNx9mve3nuRYSR3Bvl7ccn4Ct56f0ObGoHMqzYX3roeqQljwISRMabPbZDV1WCBa8WwGnaHjikXJl8H1b8FHP4M3Zmr/V0KHn/NcOQnh5CSEk1daxz+/y2PZzkKW7SokbWgQN0yIZ052LOEBqsHgbCqhDwC1RjPrDp5l+a5CNhwtxWKVjBseyp+uG8PsMcPwMzjo6HQkfyN8cJM2l8ctn9ldfUiVXAYfH72P/RWL0mbBwk/gw5thySVaP0tcTpfOmRgZwDNzM1k0M42Ve07z3paTPPXZAf6w+iDTUqOZNy6Wi9Oi8DeoVOQM6qfYT1XUmVh/pJg1e8/wzZESTGYrw0J8uevCkVwzPpak6KDun1RK2PYa/Ps3EJ4ICz7S/u1wmNRKLjqV0AeTDiWX1kZOgzu+gHevhzevhCtfgOwFXZ6wLcDHixsnDufGicM5WFTN8l2FrNhVyJcHz+LnrefitChmZsZwUUoUIX5q6GNPqYTeT1itkgNF1Ww4WspXh86y40QFVglDg31ZMHE4V46J4bzhYeh0Paw/1pyFlffB0f9o60pe9zr42e80NVvNAKrkMsh0mtABolLhjnXw8a3w6S/gyL/hqhfBP7xb10mPCSY9JphFM9PYklfG2r1nWLvvDGv2nkGvE0xICGN62hCmpkSSOiRI1dy7QSX0PiKl5FhJLVvyytl8vJzvc0spq9P+mEbHBHPvxUlckhbN2LjQnidx7UKw7xNY+yiY6uCKP8KEn3c6P4vtY7cquQwuBr0Bk9WElNJxEg2MgltWwvd/ha+ehVNb4crnIW12t6dX1usE54+K5PxRkTw1J4NdJytYd6iYrw4W8+yag7AGooJ8uGCU1tE6MTGcxMgAleA7oRK6m9Q0NrGvsJpdpyrYeaKSXScrWhJ4dJAPF6ZEMSUpkinJkQwJdtIQr4Id8Plv4NQWiMmGef/o0jS4tqFrKqEPLrZPZC1z4Tui08OUB2HUJbD8bvhwobYAysz/gaFZPbq2XidaOlIXzUzjdGUDG3NL2Xi0lI25pazYfRrQEvz44aGMHx7GuOFhZAwLJsBHpTEb9ZNwMiklxTVGDp2p4WBRNQdOV7OvsIq8sjpsyzSOjAxgWmo0ExPDmJQYwYgIf+e2Ogp3wHd/hQMrICAa5rwE2c2doF1g+9itSi6Di63PpMtDVmPGwl3faksUfv0svDIVsq6HCx7ocWK3GRbqx09z4vlpTjxSSo6X1rHleDnb8svZebKCz/dri24IAaOiAsmKDWF0TDBpMUGkDQ0mMtAwKFvyKqH3kNFs4VR5A3mldeSX1nGspJbc4lqOFtdS1dDUctywEF8yYkO4elwsWbEhZMeHEuaK4Vpmo1Yf3/yKtk6kT7A2z/WUX4FP9zpQbQndW6c6pwYTWxLvtI7ent4bJv4csq6DDf8H296AvR/ByIth0t2QNF07pheEEIyKCmRUVCALJmlDJktrjew+Wcm+01XsK6zi+2OlLN9V2PKcMH9vkqODGBUdyKioABIjta+4MH8MXp57g7xK6A7UGs2cqWrkTFUjpysbKKhsoLCigYKKek6V11NU3UjrhdHDAwwkRQdy5ZgYkqMDSY8JJm1oEKH+LixbWJrg5GatRr5/OTRWQnAszHgWxt/c4wWdbTV01UIfXGy/724ldBu/MJjxB5j6a9j+T9jyCrw/H/wjtLVoM66BuAna3OtOEBnow6Wjh3Dp6CEt28rrTBw6U83Bohpyi2vILa5l7b4iKut/bGDpBMSE+BEf7kdcmD+xoX7EhvkRG+rH0BBfhgb7DugSzsCNvJssVkl1QxMV9SYq6k2U1Wr/ltaaKKkxUlZnori6kZIaI8U1RmqN5jbP1wkYEuxLbKgfk0dFMDzcnxER/iRGBpIYEUCIvxtasxYzFB+Agq1w7Gs4/g2YasDLD9Jnw5gbtOFlvfyjsf1Bqxr64OLd3JK2Oxa9q/zCYOpD8JP74Ng62PMB7HgLtr6q3WU68mLt/2j8RIhK63IZsCvCAwwtnaytVdSZyCurI6+kjhPlWoPsZHk9G46WUFxjbNMwAwjy8SIq2IfoIB+ig3yJDPQhItBAVKAPYQEGwpu/wvy9Cfb17t2gBSfr0l++EGIm8CKgB16TUv5vu/0+wNvAeUAZMF9Kme/cUDVnqhrJK62j3mSmzmShzmimttFMrdFMTaOZmsYm7V9jE1UN2ld1g5nqxqYOvzibIB8v7RcW5EN6TDAXpvgwNMSXmOZ37GHN795tVvVxJUsTVBdCxQkoOQwlB6H4IBTtgaZ67ZjgOMi6FpIug5EXdbus0hnVKTo4tbTQrT1oobfnZYDUK7SvxiqtAZL7pfZ1YIV2jCFQ66yPTtOSe1QqhCVA0DCnteQBwgIMhAUY7M5tZDJbKapqoMj2abyqgeJqI8U1jRRXG9lTUElZralDA89GCAjx0xJ7iJ83wX5eBPt6E+TrRZCvN4E+XgT5ehHg0/xl0BPg48XIqACXzG9zzp+aEEIPLAYuAwqAbUKIlVLKA60Oux2okFImCSFuAJ4D5js9WmDF7kL+d+0hu/sCDHqCWn6YXkQF+pAUFUiInzch/to7aqi/N2H+BiICfAgPNBARYMDX23mtBKxWsJjAYtQSs7kRmhrB3KDNZmiqAWMtGGugoaL5q1ybJKv2rPZVfRpkq5nqfEMgKl0ro8RN0O7SCx3R7WFiXdUybFHdWDSo2H7fPSq5dMY3BDKu1r6khPLjULAdCrZB0W7Y86H2d2Ej9FrpMGgIBA6BgChtrLtfGPiGao0X25e3n/YJ1dsPvHxAb9D+1Xl3aelEg5eOEREBjIgI6PS4BpOFsjoj5XUmyupMVNSZqKxvorLeRGXDj43HqoYmzlYbWxqW9Q5mnPzD1ZksnDyiWz/GrujK2+BEIFdKeRxACPEBMBdondDnAk81f78U+JsQQkjpqE3ccwG1ixmT/CUIaJ3O2qc2K1AhoaIeqEf73OBQJ2HK9vtbfd/y8mTb77tFaB87dV5g0EOkP0Snax1Jem/Q+2j7gbm1DwAABLpJREFUAIyH4dhhOPZON6/RPfVm7VOAaqEPLrbf98PfPIyv3k2zI0YGQmQOWJuaG0JN2pe1Cazl0FACdWawWujR31ZLomjOEMJe1hAdE4jdc3WBLwT6QiAQ07ypQ/aQ4FszF/h9187ZDV1J6LHAqVaPC4D2k3+0HCOlNAshqoAIoM2ig0KIO4E7AYYPP/cEP3aDCYshqbzzGd9+vGB3d4h237b75bdsEz/+xxDNrQAhtO+FTtuv02ktDaHTErbQ/5i4bV96g1M/WjrT+cPOJyUspa/DUNwoMzKTq5Oupq6prq9D6UjSnPSbwGpu9WXRPs1aLSCtzd9btSdIa/OX5MdGl2zVSJP2G2yy/YUdBdT1ze0NCx3atQO7ya3ZREr5KvAqQE5OTo9a75ec/wiXnP+IU+NSFAWCDEH8/gLntxoV9+lKL18hEN/qcVzzNrvHCCG8gBDOUeRQFEVRnKsrCX0bkCyESBRCGIAbgJXtjlkJ3NL8/XXAV66onyuKoiiOnbPk0lwTvw/4HG3Y4htSyv1CiGeA7VLKlcDrwL+EELlAOVrSVxRFUdyoSzV0KeUaYE27bU+2+r4RuN65oSmKoijd4bmTGiiKogwyKqEriqJ4CJXQFUVRPIRK6IqiKB5C9NXoQiFECXCiTy7eO5G0uwN2EBhsr3mwvV5Qr3kgGSGljLK3o88S+kAlhNgupczp6zjcabC95sH2ekG9Zk+hSi6KoigeQiV0RVEUD6ESeve92tcB9IHB9poH2+sF9Zo9gqqhK4qieAjVQlcURfEQKqEriqJ4CJXQe0EI8WshhBRCRJ776IFLCPEnIcQhIcQPQojlQojQvo7JVYQQM4UQh4UQuUKIx/o6HlcTQsQLIb4WQhwQQuwXQvyyr2NyFyGEXgixSwixqq9jcRaV0HtICBEPzABO9nUsbvAFkCmlHAMcAX7Tx/G4RKsF0a8ARgM3CiFG921ULmcGfi2lHA1MBu4dBK/Z5pfAwb4OwplUQu+5PwOP0v2VawccKeV/pJTm5oeb0Vat8kQtC6JLKU2AbUF0jyWlLJJS7mz+vgYtwcX2bVSuJ4SIA64EXuvrWJxJJfQeEELMBQqllHv6OpY+cBuwtq+DcBF7C6J7fHKzEUIkAOOALX0biVv8Ba1BZu3rQJypfy453w8IIb4E7C3N/QTwOFq5xWN09nqllJ82H/ME2kf0d90Zm+J6QohA4BPgV1LK6r6Ox5WEELOBYinlDiHEtL6Ox5lUQndASnmpve1CiCwgEdgjhACt/LBTCDFRSnnGjSE6laPXayOEuBWYDUz34PViu7IguscRQnijJfN3pZTL+joeN7gAmCOEmAX4AsFCiHeklAv7OK5eUzcW9ZIQIh/IkVIOxFnbukQIMRN4AbhISlnS1/G4ihDCC63TdzpaIt8GLJBS7u/TwFxIaK2St4ByKeWv+joed2tuof//9u7YBmEYCMPofwMwER01K7AFo6RgBlZgiXSM4zSmTwGKOL03wUmWPp3swvcxxvXoWb7BHTp7LElOSV5VtVbV4+iBfmE+/H4+RH8neXaO+XROcktymWe7zs2VP2RDB2jChg7QhKADNCHoAE0IOkATgg7QhKADNCHoAE1sZRTbxZ52p+4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# An example illustration from SciPy documentation\n",
    "x = np.linspace(-5, 5, 100)\n",
    "ax = plt.subplot()\n",
    "distnames = ['laplace', 'norm', 'uniform']\n",
    "\n",
    "for distname in distnames:\n",
    "    if distname == 'uniform':\n",
    "        dist = getattr(stats, distname)(loc=-2, scale=4)\n",
    "    else:\n",
    "        dist = getattr(stats, distname)\n",
    "    data = dist.rvs(size=1000)\n",
    "    kur = kurtosis(data, fisher=True)\n",
    "    y = dist.pdf(x)\n",
    "    ax.plot(x, y, label=\"{}, {}\".format(distname, round(kur, 3)))\n",
    "    ax.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculating kurtosis and sample excess kurtosis using SciPy's `kurtosis` function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 288,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.7639701249639383"
      ]
     },
     "execution_count": 288,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the Kurtosis of T. Rowe Price Equity Income Fund (PRFDX) returns\n",
    "kurtosis(rate_of_return['Return (%)'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 289,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.7639701249639383"
      ]
     },
     "execution_count": 289,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the Sample Excess Kurtosis of T. Rowe Price Equity Income Fund (PRFDX) returns\n",
    "kurtosis(rate_of_return['Return (%)'], fisher=False) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The above result is different from that from using $K_{E}$ formula, which is 2.32. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using Geometric and Arithmetic Means\n",
    "Geometric means are preferred when reporting historical returns over **multiple periods**. This corresponds to what would have to be earned each year to match the cumulative returns.\n",
    "\n",
    "The arithmetic mean should be used to estimate the return over a **one-period horizon**. It is also better for estimating a forward-looking equity risk premium. The arithmetic mean is also better for predicting ending wealth, so it should be used for discounting.\n",
    "\n",
    "The geometric mean is approximately the arithmetic mean minus half the variance. We can express this relationship mathematically as follows:\n",
    "\n",
    "$$R_G≈R_A−\\frac{σ^2}{2}$$\n",
    "\n",
    "As implied by the above formula, the difference between the geometric and arithmetic means increases with higher volatility levels.\n",
    "\n",
    "A semilogarithmic scale plots an arithmetic scale for time on the horizontal axis and a logarithmic scale for value on the vertical axis. On a logarithmic scale, 10-fold increases are equally spaced. So, the increase along the vertical axis is the same for a \\\\$9 increase in value from $1 to \\\\$10 as it is for a $90 increase from \\\\$10 to \\\\$100. Graphically, a constant growth rate would be shown as a linear increase on a semilogarithmic scale rather than the exponentially upward sloping increase that would appear if both axes were plotted on arithmetic scales."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*$R_{G}$ for historic performance with semilogarithmic scales* "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 240,
   "metadata": {},
   "outputs": [],
   "source": [
    "# y-axis in log scale\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*$\\bar{X}$ for expected returns*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Probability Concepts\n",
    "## Probability\n",
    "The return on a risky asset is a random variable because the outcomes (possible values) are uncertain. An `event` is a specified set of outcomes.\n",
    "\n",
    "Three approaches are used to estimate probabilities:\n",
    "- `Subjective probabilities` are based on personal judgment.\n",
    "- `Empirical probabilities` are derived from relative frequencies from historical data. This only works if relationships are stable through time. It will not be useful for very rare events. (Relationships must be stable through time for empirical probabilities to be accurate.)\n",
    "- A `priori probabilities` are deduced using logic rather than observation. Both empirical and a priori probabilities are considered objective because they are typically the same for all people.\n",
    "\n",
    "Often probabilities are stated as `odds`. If the probability of an event $E$ is $P(E)$, then the odds are:\n",
    "\n",
    "$$Odds\\ for\\ E=\\frac{P(E)}{1−P(E)}$$\n",
    "\n",
    "- `Unconditional probabilities` (a.k.a. marginal probabilities) are based on the universe of all possibilities. They can be thought of as `stand-alone probabilities`.\n",
    "- `Conditional probabilities` restrict the set of possibilities. This is how investors incorporate new information.\n",
    "- A `joint probability` reflects the probability of two or more outcomes occurring. $P(AB)$ represents the probability of A and B occurring.\n",
    "\n",
    "#### Conditional Probability: P(A | B) \n",
    "$$P(A | B) = \\frac{P(AB)}{P(B)}, P(B) ≠ 0$$ \n",
    "<div style=\"text-align: right\"> Equation 1 </div>\n",
    "The above can be manipulated algebraically to derive the **multiplication rule for probability**.\n",
    "$$P(AB)=P(A∣B)P(B)=P(BA)=P(B∣A)P(A)$$\n",
    "\n",
    "#### Joint Probability: P(AB)\n",
    "*Multiplication Rule for Probability*<br>\n",
    "The joint probability of A and B can be expressed as $P(AB) = P(A | B)P(B)$  \n",
    "<div style=\"text-align: right\"> Equation 2 </div>\n",
    "\n",
    "**Addition Rule for Probabilities**\n",
    "\n",
    "$$P(A or B) = P(A) + P(B) – P(AB)$$\n",
    "\n",
    "or \n",
    "\n",
    "$$P(A U B) = P(A) + P(B) - P(A ∩ B)$$\n",
    "<div style=\"text-align: right\"> Equation 3 </div>\n",
    "\n",
    "#### Independence\n",
    "*Independence and dependence*<br>\n",
    "An area of intense interest to investment managers and their clients is whether records of past performance are useful in identifying repeat winners and losers. Mathematically, events A and B are independent if and only if P(A∣B)=P(A) or, equivalently, P(B∣A)=P(B). For independent events, P(AB)=P(A)P(B).\n",
    "\n",
    "*Multiplication Rule for Independent Events*<br>\n",
    "$P(AB) = P(A)P(B)$  \n",
    "<div style=\"text-align: right\"> Equation 4</div>\n",
    "\n",
    "*The Total Probability Rule*<br>\n",
    "$P(A)=P(AS)+P(AS^C)=P(A|S)P(S)+P(A∣S^C)P(S^C)$\n",
    "<div style=\"text-align: right\"> Equation 5</div>\n",
    "\n",
    "$P(A)=P(AS_1)+P(AS_2)+…+P(AS_n)=P(A|S_1)P(S_1)+P(A|S_2)P(S_2)+…+P(A|S_n)P(S_n)$      \n",
    "<div style=\"text-align: right\"> Equation 6</div>\n",
    "where $S_1, S_2, …, S_n$ are mutually exclusive and exhaustive scenarios or events.\n",
    "\n",
    "The complement of $S$ is $S^C$. It represents event $S$ not occurring. It follows $P(S)+P(S^C)=1$. This notation can be used to define the total probability rule."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Expected Value (Mean) and Variance\n",
    "- The `expected value` is the **probability-weighted average of the possible outcomes of the random variable**. \n",
    "- The `variance` is the probability-weighted average of the squared deviations from the random variable's expected value. A higher variance indicates more dispersion in the random variable. \n",
    "- The standard deviation is the positive square root of the variance.\n",
    "\n",
    "$$E(X) = \\displaystyle\\sum_{i=1}^n{P(X_{i})X_{i}}$$ \n",
    "<div style=\"text-align: right\"> Equation 7</div>\n",
    "\n",
    "where $X_{i}$ is one of $n$ possible outcomes of $X$\n",
    "#### Properties of Expected Value\n",
    "Let $w_{i}$ be any constant and $R_{i}$ be a random variable.\n",
    "\n",
    "- The expected value of a constant times a random variable equals the constant times the expected value of the random variable.\n",
    "\n",
    "$E(w_{i}R_{i}) = w_{i}E(R_{i})$\n",
    "\n",
    "- The expected value of a weighted sum of random variables equals the weighted sum of the expected values, using the same weights.\n",
    "\n",
    "$$E(w_{1}R_{1}+w_{2}R_{2}+…+w_{n}R_{n})=w_{1}E(R_{1})+w_{2}E(R_{2})+…+w_{n}E(R_{n})$$\n",
    "#### Variance and Standard Deviation \n",
    "the risks associated with $E(X)$\n",
    "\n",
    "$$Var(X)=\\sigma^2(X)=E\\{[X−E(X)]^2\\}$$\n",
    "<div style=\"text-align: right\"> Equation 8 </div>\n",
    "\n",
    "(Variance of $X$ is a quantity in the squared units of $X$. For example, if the random variable is return in percent, variance of return is in units of percent squared.) \n",
    "\n",
    "Summarize the variance formula as:\n",
    "\n",
    "$$\\sigma^2(X)=P(X_{1})[X_{1}−E(X)]^2+P(X_{2})[X_{2}−E(X)]^2+…+P(X_{n})[X_{n}−E(X)]^2=\\sum_{i=1}^n{P(X_{i})[X_{i}−E(X)]^2}$$  \n",
    "<div style=\"text-align: right\"> Equation 9\n",
    "    \n",
    "where $X_{i}$ is one of $n$ possible outcomes of the random variable $X$.\n",
    "\n",
    "$σ^2(X)=E{[X−E(X)]^2}=\\sum_{i=1}^n{\\{P(X_i)[X_i−E(X)]^2\\}}=E(X^2)−[E(X)]^2$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Conditional Measures of Expected Value and Variance\n",
    "#### Conditional Expected Values\n",
    "`Conditional expected values` take into account new information or events. This can be used to define the total probability rule for the expected value.\n",
    "\n",
    "$$E(X)=E(X∣S_1)P(S_1)+E(X∣S_2)P(S_2)+…+E(X∣S_n)P(S_n)$$\n",
    "Conditional variances can also be calculated.\n",
    "\n",
    "\n",
    "$E(X|S) = P(X_{1}|S)X_{1} + P(X_{2}|S)X_{2} + … + P(X_{n}|S)X_{n}$\n",
    "<div style=\"text-align: right\"> Equation 10</div>\n",
    "\n",
    "#### The Total Probability Rule for Expected Values\n",
    "A principle for stating (unconditional) expected values in terms of conditional expected values.\n",
    "\n",
    "$E(X) = E(X|S)P(S) + E(X|S^C)P(S^C)$  \n",
    "<div style=\"text-align: right\"> Equation 11</div>\n",
    "\n",
    "$$E(X) = E(X|S_{1})P(S_{1}) + E(X|S_{2})P(S_{2}) + … + E(X|S_{n})P(S_{n})$$\n",
    "<div style=\"text-align: right\"> Equation 12</div>\n",
    "\n",
    "where $S_1, S_2, …, S_n$ are mutually exclusive and exhaustive scenarios or events.\n",
    "\n",
    "#### Conditional Variance \n",
    "1. Variance, like expected value, has a conditional counterpart to the unconditional concept.\n",
    "2. We can use conditional variance to assess risk given a particular scenario.\n",
    "\n",
    "*The relationship between conditional variance and unconditional one to be continued in later sessions.*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 290,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Probability</th>\n",
       "      <th>EPS ($)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.15</td>\n",
       "      <td>2.60</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0.45</td>\n",
       "      <td>2.45</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.24</td>\n",
       "      <td>2.20</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0.16</td>\n",
       "      <td>2.00</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Probability  EPS ($)\n",
       "0         0.15     2.60\n",
       "1         0.45     2.45\n",
       "2         0.24     2.20\n",
       "3         0.16     2.00"
      ]
     },
     "execution_count": 290,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# BankCorp’s Earnings per Share\n",
    "eps = {'Probability': [0.15, 0.45, 0.24, 0.16], \n",
    "     'EPS ($)': [2.60, 2.45, 2.20, 2.00]}\n",
    "eps_ev = pd.DataFrame(data=eps)\n",
    "eps_ev"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 291,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Probability</th>\n",
       "      <th>EPS ($)</th>\n",
       "      <th>probability_weighted_EPS</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.15</td>\n",
       "      <td>2.60</td>\n",
       "      <td>0.3900</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0.45</td>\n",
       "      <td>2.45</td>\n",
       "      <td>1.1025</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.24</td>\n",
       "      <td>2.20</td>\n",
       "      <td>0.5280</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0.16</td>\n",
       "      <td>2.00</td>\n",
       "      <td>0.3200</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Probability  EPS ($)  probability_weighted_EPS\n",
       "0         0.15     2.60                    0.3900\n",
       "1         0.45     2.45                    1.1025\n",
       "2         0.24     2.20                    0.5280\n",
       "3         0.16     2.00                    0.3200"
      ]
     },
     "execution_count": 291,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the E(X) of BankCorp’s Earnings per Share: multiplying each prediction with its probability\n",
    "eps_ev['probability_weighted_EPS'] = eps_ev['Probability']*eps_ev['EPS ($)']\n",
    "eps_ev"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 292,
   "metadata": {},
   "outputs": [],
   "source": [
    "expected_eps = eps_ev['probability_weighted_EPS'].tolist()\n",
    "ev=sum(expected_eps)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 293,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.03878475000000002\n"
     ]
    }
   ],
   "source": [
    "# Calculating the variance \n",
    "list1=[]\n",
    "for i in range(len(eps_ev['EPS ($)'])):\n",
    "    diff = eps_ev['EPS ($)'][i] - ev\n",
    "    diff_sqr = diff**2\n",
    "    eps = diff_sqr*eps_ev['Probability'][i]\n",
    "    list1.append(eps)\n",
    "\n",
    "var_x=sum(list1)\n",
    "print(var_x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Standard deviation is easier to interpret than variance, as it is in the same units as the random variable. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 294,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.19693844215896505"
      ]
     },
     "execution_count": 294,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculating the standard deviation \n",
    "var_x**0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculate the variance of variable $X$ given each scenario:<Br>\n",
    "$\\sigma^2{(X|S)} $<Br>\n",
    "$\\sigma^2{(X|S^C)}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TBC"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use case: Calculating the default risk premium for a one-period debt instrument<br>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Portfolio Expected Return and Variance of Return\n",
    "`Modern portfolio theory` is based on the principle that investment opportunities should be evaluated in the context of how they impact the tradeoff between a portfolio's expected return and the level of portfolio risk, as measured by variance.\n",
    "### Portfolio Expected Return\n",
    "A portfolio's expected return is simply the weighted average of expected returns for individual securities:\n",
    "\n",
    "$$E(R_p)=\\sum_{i=1}^n{w_iE(R_i)}$$\n",
    "\n",
    "Let $w_{i}$ be any constant and $R_{i}$ be a random variable, recall the properties of expected returns. \n",
    "$E(R_{p})=E(w_{1}R_{1}+w_{2}R_{2}+…+w_{n}R_{n})=w_{1}E(R_{1})+w_{2}E(R_{2})+…+w_{n}E(R_{n})$\n",
    "<div style=\"text-align: right\"> Equation 13</div>\n",
    "\n",
    "In an earlier reading, we calculated variance (and standard deviation) for returns on individual assets. When considering whether to add an asset to a portfolio, we must consider `covariance`, which measures the tendency for two variables to move in sync.\n",
    "### Covariance \n",
    "The `covariance` between two random variables is the probability-weighted average of the cross-products of each random variable’s deviation from its own expected value. \n",
    "#### Forward-Looking or Population Covariance\n",
    "Given two random variables $R_{i}$ and $R_{j}$, the covariance between $R_{i}$ and $R_{j}$ is\n",
    "\n",
    "$$Cov(R_{i},R_{j})=E[(R_{i}−ER_{i})(R_{j}−ER_{j})]$$\n",
    "<div style=\"text-align: right\"> Equation 14</div>\n",
    "a.k.a $\\displaystyle\\sigma_{(R_{i}, R_{j})}$ or $\\displaystyle\\sigma{_{ij}}$\n",
    "\n",
    "Based on this formula, we can draw the following conclusions:\n",
    "1. Covariance is positive if, when one asset is generating above-average returns, the other asset is as well. Both assets will also tend to generate returns below their respective averages in the same periods.\n",
    "2. Covariance is negative if one asset is generating above-average returns while the other's returns are below its average (or vice versa).\n",
    "3. The covariance of an asset's returns with itself (own covariance) is equal to its variance.\n",
    "\n",
    "#### Historical or Sample Data Covariance\n",
    "The sample covariance is the average value of the product of the deviations of observations on two random variables from their sample means.\n",
    "<br>\n",
    "$Cov(R_{i},R_{j})=\\displaystyle\\frac{\\sum_{i=1}^{n} {(R_{i,t}-\\bar{R_{i}})(R_{j,t}-\\bar{R_{j})}}}{n−1}$<br><br>\n",
    "<div style=\"text-align: right\"> Equation 15</div>\n",
    "\n",
    "**If the random variables are returns, the units of both forward-looking covariance and historical variance would be `returns squared`($\\%^2$). Because covariance statistics are also not scaled, it isn't easy to make comparisons.**\n",
    "### Variance and Standard Deviation of Portfolio Returns\n",
    "`Portfolio variance` is calculated using the following formula:\n",
    "\n",
    "$\\sigma^2(R_{p})$\n",
    "\n",
    "$=w_{1}^2\\sigma^2(R_{1})+w_{1}w_{2}Cov(R_{1},R_{2})+w_{1}w_{3}Cov(R_{1},R_{3})+w_{1}w_{2}Cov(R_{1},R_{2})$\n",
    "\n",
    "$+w_{2}^2\\sigma^2(R_{2})+w_{2}w_{3}Cov(R_{2},R_{3})+w_{1}w_{3}Cov(R_{1},R_{3})+w_{2}w_{3}Cov(R_{2},R_{3})$\n",
    "\n",
    "$+w_{3}^2\\sigma^2(R_{3})$\n",
    "\n",
    "\n",
    "$=\\displaystyle\\sum_{i=1}^{3}\\sum_{j=1}^{3}{w_{i}w_{j}Cov(R_{i},R_{j})}$\n",
    "<div style=\"text-align: right\"> Equation 16</div>\n",
    "\n",
    "*The double summation signs say: “Set $i = 1$ and let $j$ run from 1 to 3; then set $i = 2$ and let $j$ run from 1 to 3; next set $i = 3$ and let $j$ run from 1 to 3; finally, add the nine terms.”*\n",
    "\n",
    "$$\\sigma^2(R_{p})=\\displaystyle\\sum_{i=1}^{n}\\sum_{j=1}^{n}{w_{i}w_{j}Cov(R_{i},R_{j})}$$\n",
    "<div style=\"text-align: right\"> Equation 17</div>\n",
    "\n",
    "While this formula may appear relatively straightforward, it becomes increasingly complex as the more assets are added to a portfolio. If there are $n$ securities in a portfolio, there are $\\frac{n(n−1)}{2}$ covariances to estimate and $n$ variances to estimate. The variance for a two-asset portfolio is calculated as follows:\n",
    "\n",
    "$σ^2(R_p)$\n",
    "\n",
    "$=w^2_Aσ^2_A+w^2_Bσ^2_B+2w_Aw_BCov[R_A,R_B]$\n",
    "\n",
    "$=w^2_Aσ^2_A+w^2_Bσ^2_B+2w_Aw_Bρ_{AB}σ_Aσ_B$\n",
    "\n",
    "\n",
    "For a portfolio with three assets, this formula expands to:\n",
    "\n",
    "$σ^2(R_p)$\n",
    "\n",
    "$=w^2_Aσ^2_A+w^2_Bσ^2_B+w^2_Cσ^2_C$\n",
    "\n",
    "$+2w_Aw_BCov[R_A,R_B]$\n",
    "\n",
    "$+2w_Aw_CCov[R_A,R_C]$\n",
    "\n",
    "$+2w_Bw_CCov[R_B,R_C]$\n",
    "\n",
    "### Correlation\n",
    "The `correlation` between two sets of returns, $R_i$ and $R_j$, is calculated using the covariance and standard deviations of individual asset returns:\n",
    "$$\\rho(R_i,R_j) = \\rho_{ij} = Corr(R_{i},R_{j}) = \\displaystyle\\frac{Cov(R_{i},R_{j})}{[\\sigma(R_{i})\\sigma(R_{j})]}$$\n",
    "\n",
    "$Cov(R_i,R_j) = \\rho(R_i,R_j)σ(R_i)σ(R_j)$\n",
    "<div style=\"text-align: right\"> Equation 18</div>\n",
    "\n",
    "The above formula are forward-looking, we can similarly compute an historical or sample correlation by dividing historical or sample covariance between two variables by the product of sample standard deviations of the two variables.\n",
    "\n",
    "Unlike covariance, `correlation statistics` are scaled within a **range of -1 to 1**. These relationships can be shown on a scatter plot with returns for one asset on the vertical axis and returns for the other asset on the horizontal axis.\n",
    "- Positive correlations indicate a positive linear relationship. A `perfectly positive` relationship will have a correlation of 1 and will be depicted on a scatter plot as a straight upward sloping line.\n",
    "- Negative correlations indicate a negative linear relationship. A correlation of -1 indicates a `perfectly inverse` relationship.\n",
    "- A correlation of 0 indicates `no LINEAR relationship`. There will be no discernable linear pattern on a scatter plot. *(However, it is important to note that zero correlation does not mean that there is no relationship, only that there is no linear relationship.)*\n",
    "\n",
    "Limitation:\n",
    "- Correlation statistics can be sensitive to `outliers`, so judgment should be used when deciding whether to include observations that are far from the mean.\n",
    "- Correlation does not imply causation. Correlation measures may be spurious if they capture a relationship that occurs by chance or if the relationship between two variables is actually attributable to a relationship that both have with a third variable.\n",
    "\n",
    "### Covariance Given a Joint Probability Function\n",
    "The calculation of covariance in a forward-looking sense requires the specification of a joint probability function, which gives the probability of joint occurrences of values of the two random variables.\n",
    "\n",
    "$Cov(R_{A},R_{B})=\\displaystyle\\sum_i\\sum_j{P(R_{A,i},R_{B,j})(R_{A,i}−ER_{A})(R_{B,j}−ER_{B})}$\n",
    "<div style=\"text-align: right\"> Equation 19</div>\n",
    "\n",
    "Sum all possible deviation cross-products weighted by the appropriate joint probability.\n",
    "\n",
    "Covariance and correlation measures can be estimated from historical return data or a joint probability function. We can calculate the covariance and correlation between assets based on their probability-weighted returns under different market conditions.\n",
    "\n",
    "However, this calculation is much simpler if we can assume that two random variables (e.g., returns) are independent of each other. The assumption of independence is stronger than the assumption that two variables are uncorrelated because independence covers any relationship. By contrast, a zero correlation only tells us that there is no linear relationship, but it doesn't tell us anything about the possibility of a non-linear relationship.\n",
    "\n",
    "Random variables are independent if, and only if, $P(XY)=P(X)P(Y)$. **If two random variables are independent, they must also be uncorrelated. By contrast, knowing that variables are uncorrelated does not allow us to conclude that they are independent.**\n",
    "\n",
    "If two random variables are uncorrelated, then the `multiplication rule` for the expected value of the product of uncorrelated random variables tells us that $E(XY)=E(X)E(Y)$. If we cannot assume that the variables are uncorrelated, we must calculate a conditional expected value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 295,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Asset Class</th>\n",
       "      <th>Weights</th>\n",
       "      <th>Expected Return (%)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>S&amp;P 500</td>\n",
       "      <td>0.50</td>\n",
       "      <td>0.13</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>US long-term corporate bonds</td>\n",
       "      <td>0.25</td>\n",
       "      <td>0.06</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>MSCI EAFE</td>\n",
       "      <td>0.25</td>\n",
       "      <td>0.15</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                    Asset Class Weights  Expected Return (%)\n",
       "0                       S&P 500    0.50                 0.13\n",
       "1  US long-term corporate bonds    0.25                 0.06\n",
       "2                     MSCI EAFE    0.25                 0.15"
      ]
     },
     "execution_count": 295,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "portfolio = {'Asset Class': ['S&P 500', 'US long-term corporate bonds', 'MSCI EAFE'], \n",
    "     'Weights': ['0.50', '0.25', '0.25'],\n",
    "     'Expected Return (%)':[0.13, 0.06, 0.15]}\n",
    "portfolio = pd.DataFrame(data=portfolio)\n",
    "portfolio"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use case: calculating portfolio expected return "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 296,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Change data type to numeric\n",
    "portfolio['Expected Return (%)'] = pd.to_numeric(portfolio['Expected Return (%)'])\n",
    "portfolio['Weights'] = pd.to_numeric(portfolio['Weights'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 297,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Asset Class</th>\n",
       "      <th>Weights</th>\n",
       "      <th>Expected Return (%)</th>\n",
       "      <th>Weighted E(R)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>S&amp;P 500</td>\n",
       "      <td>0.50</td>\n",
       "      <td>0.13</td>\n",
       "      <td>0.0650</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>US long-term corporate bonds</td>\n",
       "      <td>0.25</td>\n",
       "      <td>0.06</td>\n",
       "      <td>0.0150</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>MSCI EAFE</td>\n",
       "      <td>0.25</td>\n",
       "      <td>0.15</td>\n",
       "      <td>0.0375</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                    Asset Class  Weights  Expected Return (%)  Weighted E(R)\n",
       "0                       S&P 500     0.50                 0.13         0.0650\n",
       "1  US long-term corporate bonds     0.25                 0.06         0.0150\n",
       "2                     MSCI EAFE     0.25                 0.15         0.0375"
      ]
     },
     "execution_count": 297,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculating each asset class's weighted expected return\n",
    "portfolio['Weighted E(R)'] = portfolio['Weights'] * portfolio['Expected Return (%)']\n",
    "portfolio"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 298,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.1175"
      ]
     },
     "execution_count": 298,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Summing up the weighted expected returns to obtain the portfolio expected return\n",
    "portfolio_er = sum(portfolio['Weighted E(R)'].tolist())\n",
    "portfolio_er "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use case: calculating portfolio variance (impletement Equation 8), introducing the concept of covariance. (Considering variance in a forward-looking sense). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 299,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Asset Class</th>\n",
       "      <th>Weights</th>\n",
       "      <th>Expected Return (%)</th>\n",
       "      <th>Weighted E(R)</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>S&amp;P 500</td>\n",
       "      <td>0.50</td>\n",
       "      <td>0.13</td>\n",
       "      <td>0.0650</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>US long-term corporate bonds</td>\n",
       "      <td>0.25</td>\n",
       "      <td>0.06</td>\n",
       "      <td>0.0150</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>MSCI EAFE</td>\n",
       "      <td>0.25</td>\n",
       "      <td>0.15</td>\n",
       "      <td>0.0375</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                    Asset Class  Weights  Expected Return (%)  Weighted E(R)\n",
       "0                       S&P 500     0.50                 0.13         0.0650\n",
       "1  US long-term corporate bonds     0.25                 0.06         0.0150\n",
       "2                     MSCI EAFE     0.25                 0.15         0.0375"
      ]
     },
     "execution_count": 299,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "portfolio"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 300,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.5 , 0.25, 0.25])"
      ]
     },
     "execution_count": 300,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Arrange the weights in to a matrix\n",
    "weights = portfolio['Weights'].to_numpy()\n",
    "weights"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 301,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Unit: percent squared</th>\n",
       "      <th>S&amp;P 500</th>\n",
       "      <th>US long-term corporate bonds</th>\n",
       "      <th>MSCI EAFE</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>S&amp;P 500</td>\n",
       "      <td>400</td>\n",
       "      <td>45</td>\n",
       "      <td>189</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>US long-term corporate bonds</td>\n",
       "      <td>45</td>\n",
       "      <td>81</td>\n",
       "      <td>38</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>MSCI EAFE</td>\n",
       "      <td>189</td>\n",
       "      <td>38</td>\n",
       "      <td>441</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          Unit: percent squared  S&P 500  US long-term corporate bonds  \\\n",
       "0                       S&P 500      400                            45   \n",
       "1  US long-term corporate bonds       45                            81   \n",
       "2                     MSCI EAFE      189                            38   \n",
       "\n",
       "   MSCI EAFE  \n",
       "0        189  \n",
       "1         38  \n",
       "2        441  "
      ]
     },
     "execution_count": 301,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Covariance Matrix\n",
    "matrix = {'Unit: percent squared': ['S&P 500', 'US long-term corporate bonds', 'MSCI EAFE'], \n",
    "     'S&P 500': [400, 45, 189], 'US long-term corporate bonds':[45, 81, 38], 'MSCI EAFE':[189, 38, 441]}\n",
    "cov_matrix = pd.DataFrame(data=matrix)\n",
    "cov_matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 302,
   "metadata": {},
   "outputs": [],
   "source": [
    "#for row in cov_matrix[['S&P 500', 'US long-term corporate bonds', 'MSCI EAFE']].to_numpy():\n",
    "#    print(row[0]*portfolio['Weights'][0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 303,
   "metadata": {},
   "outputs": [],
   "source": [
    "#row = cov_matrix.iloc[[0]][['S&P 500', 'US long-term corporate bonds', 'MSCI EAFE']].to_numpy()\n",
    "#row"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 304,
   "metadata": {},
   "outputs": [],
   "source": [
    "to_sum = []\n",
    "for i in range(len(portfolio)):\n",
    "    row = cov_matrix.iloc[[i]][['S&P 500', 'US long-term corporate bonds', 'MSCI EAFE']].to_numpy()\n",
    "    ind_var = weights[i]*row\n",
    "    # Arrange into a matrix\n",
    "    to_sum.append(ind_var.squeeze())\n",
    "matrix = np.row_stack(to_sum)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 305,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "195.875"
      ]
     },
     "execution_count": 305,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Matrix operation\n",
    "# Calculate the variance of the portfolio return\n",
    "var = (matrix*weights).sum()\n",
    "var"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 306,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "13.995535002278405"
      ]
     },
     "execution_count": 306,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the standard deviation\n",
    "var**0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looking at the first three terms in the calculation above. Their sum, 100 + 5.0625 + 27.5625 = 132.625, is the contribution of the individual variances to portfolio variance. If the returns on the three assets were independent, covariances would be 0 and the standard deviation of portfolio return would be 132.6251/2 = 11.52 percent as compared to 14 percent before. The portfolio would have less risk. Suppose the covariance terms were negative. Then a negative number would be added to 132.625, so portfolio variance and risk would be even smaller. At the same time, we have not changed expected return. For the same expected portfolio return, the portfolio has less risk. This risk reduction is a diversification benefit, meaning a risk-reduction benefit from holding a portfolio of assets. <br>\n",
    "**The diversification benefit increases with decreasing covariance.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 307,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Unit: percent squared</th>\n",
       "      <th>S&amp;P 500</th>\n",
       "      <th>US long-term corporate bonds</th>\n",
       "      <th>MSCI EAFE</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>S&amp;P 500</td>\n",
       "      <td>1.00</td>\n",
       "      <td>0.25</td>\n",
       "      <td>0.45</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>US long-term corporate bonds</td>\n",
       "      <td>0.25</td>\n",
       "      <td>1.00</td>\n",
       "      <td>0.25</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>MSCI EAFE</td>\n",
       "      <td>0.45</td>\n",
       "      <td>0.20</td>\n",
       "      <td>1.00</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          Unit: percent squared  S&P 500  US long-term corporate bonds  \\\n",
       "0                       S&P 500     1.00                          0.25   \n",
       "1  US long-term corporate bonds     0.25                          1.00   \n",
       "2                     MSCI EAFE     0.45                          0.20   \n",
       "\n",
       "   MSCI EAFE  \n",
       "0       0.45  \n",
       "1       0.25  \n",
       "2       1.00  "
      ]
     },
     "execution_count": 307,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Correlation matrix of returns\n",
    "matrix2 = {'Unit: percent squared': ['S&P 500', 'US long-term corporate bonds', 'MSCI EAFE'], \n",
    "     'S&P 500': [1.0, 0.25, 0.45], 'US long-term corporate bonds':[0.25, 1.0, 0.20], 'MSCI EAFE':[0.45, 0.25, 1.0]}\n",
    "cov_matrix2 = pd.DataFrame(data=matrix2)\n",
    "cov_matrix2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 308,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Using for loop (another approach)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example: a portfolio of two mutual funds, with 75% in Fund A. $E(R_{A}) = 20\\%$, $E(R_{B}) = 12\\%$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 309,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Fund</th>\n",
       "      <th>A</th>\n",
       "      <th>B</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>A</td>\n",
       "      <td>625</td>\n",
       "      <td>120</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>B</td>\n",
       "      <td>120</td>\n",
       "      <td>196</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  Fund    A    B\n",
       "0    A  625  120\n",
       "1    B  120  196"
      ]
     },
     "execution_count": 309,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Covariance Matrix\n",
    "funds = {'Fund':['A', 'B'], 'A':[625, 120], 'B':[120, 196]}\n",
    "cov_matrix3 = pd.DataFrame(data=funds)\n",
    "cov_matrix3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 310,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.18000000000000002"
      ]
     },
     "execution_count": 310,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the expected return \n",
    "0.75 * 0.20 + 0.25 * 0.12"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 311,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.34285714285714286"
      ]
     },
     "execution_count": 311,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the correlation matrix \n",
    "std_A = 625**0.5\n",
    "std_B = 196**0.5\n",
    "corr = 120 / (std_A*std_B)\n",
    "corr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 312,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compute portfolio standard deviation of return\n",
    "std_port = 408.81**0.5\n",
    "# Check answers in the reading"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use case: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 313,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Example\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Topics in Probability\n",
    "### Bayes' Formula\n",
    "Bayes’ formula is a method for **updating probabilities based on new information**, which uses the occurrence of the event to infer the probability of the scenario generating it. It is based on the `total probability rule` $P(A) = P(A∩B) + P(A∩B^c)$.\n",
    "- `Prior probabilities` represent the probabilities before the arrival of any new information. \n",
    "- `Posterior probability` reflects the new information.\n",
    "\n",
    "Given a set of prior probabilities for an event of interest, if you receive new information, the rule for updating your probability of the event is:\n",
    "$$P(Event | Info)= \\displaystyle\\frac{P(Info | Event)}{P(Info)}\\times P(Event)$$\n",
    "\n",
    "Use case: evaluation of mutual fund performance\n",
    "### Principles of Counting\n",
    "#### Multiplication Rule of Counting\n",
    "The multiplication rule of counting states $k$ tasks can be done $(n_1)(n_2)(n_3)…(n_k)$ ways. In this notation, the first task has $n_1$ ways of being done.\n",
    "\n",
    "$n! = n(n – 1)(n – 2)(n – 3)…1$\n",
    "#### Multinomial Formula (General Formula for Labeling Problems)\n",
    "Formulas can also be developed for `labeling`. The number of ways that $n$ objects can be labeled with $k$ different labels, with $n_{1}$ of the first type, $n_{2}$ of the second type, and so on, with $n_{1} + n_{2} + … + n_{k} = n$, is given by\n",
    "\n",
    "$$\\displaystyle\\frac{n_{!}}{n_{1}!n_{2}!…n_{k}!}$$\n",
    "\n",
    "The multinomial formula with two different labels ($k = 2$) is especially important. \n",
    "#### Combination Formula (Binomial Formula)\n",
    "If order does not matter (i.e., `combination`), then the number of ways $r$ objects can be chosen from $n$ objects is:\n",
    "\n",
    "$$\\displaystyle_{n}C_{r}=\\binom{n}r=\\frac{n!}{(n−r)!r!}$$\n",
    "\n",
    "Here $_{n}C_{r}$ and $\\displaystyle\\binom{n}r$ are shorthand notations for $\\displaystyle\\frac{n!}{(n − r)!r!}$ (read: $n$ choose $r$, or $n$ combination $r$).\n",
    "#### Permutation Formula\n",
    "If order does matter (i.e., `permutation`), then the number of ways $r$ objects can be chosen from $n$ objects is:\n",
    "\n",
    "$$\\displaystyle_{n}P_{r}=\\frac{n!}{(n−r)!}$$\n",
    "\n",
    "A useful fact can be illustrated as follows: \n",
    "\n",
    "$_{5}C_{3}$ = $\\displaystyle\\frac{5!}{(2!3!)}$ equals $_{5}C_{2}$ = $\\displaystyle\\frac{5!}{(3!2!)}$, as 3 + 2 = 5; \n",
    "\n",
    "$_{5}C_{4}$ = $\\displaystyle\\frac{5!}{(1!4!)}$ equals $_{5}C_{1}$ = $\\displaystyle\\frac{5!}{(4!1!)}$, as 4 + 1 = 5. \n",
    "\n",
    "This `symmetrical relationship` can save work when we need to calculate many possible combinations."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "User case:<br>A mutual fund guide ranked 18 bond mutual funds by total returns for the last year. The guide also assigned each fund one of five risk labels: high risk (four funds), above-average risk (four funds), average risk (three funds), below-average risk (four funds), and low risk (three funds); as 4 + 4 + 3 + 4 + 3 = 18, all the funds are accounted for. How many different ways can we take 18 mutual funds and label 4 of them high risk, 4 above-average risk, 3 average risk, 4 below-average risk, and 3 low risk, so that each fund is labeled?<br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 314,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "6402373705728000"
      ]
     },
     "execution_count": 314,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Calculate the number of possible sequences \n",
    "math.factorial(18)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 315,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "12864852000.0\n"
     ]
    }
   ],
   "source": [
    "# List of number of assignments in each category\n",
    "num_in_category = [4,4,3,4,3]\n",
    "denom = 1\n",
    "# Counting number of ways of the labeling \n",
    "for i in num_in_category:\n",
    "    n = sum(num_in_category)\n",
    "    numerator = math.factorial(n)\n",
    "    # Calculate the number of possible arrangement within each categories (because the sequence within each type does not affect the final result)\n",
    "    n_cat = math.factorial(i)\n",
    "    denom = denom * n_cat\n",
    "    result = numerator / denom\n",
    "print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### Example: binomial option pricing model<br>\n",
    "This model describes the movement of the underlying asset as a series of moves, price up (U) or price down (D). For example, two sequences of five moves containing three up moves, such as UUUDD and UDUUD, result in the same final stock price. At least for an option with a payoff dependent on final stock price, the number but not the order of up moves in a sequence matters. How many sequences of five moves belong to the group with three up moves? "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 316,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "10.0"
      ]
     },
     "execution_count": 316,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "math.factorial(5) / (math.factorial(3) * math.factorial(2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**RECAP**\n",
    "\n",
    "- Differentiate terms commonly used when valuing cash flows that are distributed over many time periods.\n",
    "- What is Deferred Ordinary Annuity, how to calculate it? (2-step method)\n",
    "- What's a harmonic mean? When would it be suitable?\n",
    "- Explain the choice of the said mean values in different context. \n",
    "- What's the relationship between geometric and arithmetic means?\n",
    "- Conditional expected value formula?\n",
    "- If the random variables are returns, the units of both forward-looking covariance and historical variance would be?\n",
    "- Portfolio covariance formula.\n",
    "- Correlation formula.\n",
    "- Portfolio expected return formula?\n",
    "- Portfolio variance formula: two assets, three assets?\n",
    "- Labeling problem?\n",
    "- Multiplication rule of counting?\n",
    "- Combination formula?\n",
    "- Permutation formula?"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}