Tag Archives: math

The Maximum Likelihood Estimate of an unknown parameter

 

The Maximum Likelihood Estimate is one estimate of a parameter’s value.  It basically answers the question: “Given my data, what’s the most likely value of an unknown parameter”.  To use it, we need to compute the likelihood function by taking a sample of n observations from the same distribution, evaluating each of the n observations, (x1, …, xn), at the pdf pX(x), and taking the product of the result for all n.

Likelihood Function

Note: The expression pX(x; θ) is just the distribution of the population from which we’re sampling – usually the parameter is added with a semicolon to emphasize the fact that the distribution is characterized by the unknown parameter

The likelihood function is just the joint density of the random sample.  Since samples are independent and identically distributed (iid), the joint pdf of all n sample observations is the product of the individual densities.  This is the same principle that lets us multiply the P(A), P(B), and P(C) together to find P(A intersect B intersect C) when events A, B, and C are independent.  Suppose we take a sample of size n = 3 from the distribution, and the resulting values are x1, x2, and x3.  What’s the probability associated with the three sample values? That is, what’s the joint density of the three sample values, pX(x1, x2, x3)?

Joint

Generalizing this case for all n clearly gives us the result that the joint density of n randomly drawn sample values is the product of individual densities and the likelihood function is nothing more than the joint pdf of the sample – a multivariate probability density function of the values taken on my the n random variables in the sample.

The likelihood function is a function of the unknown parameter.  The Maximum Likelihood Estimate for the unknown parameter is the parameter value that maximizes the likelihood function:

We use calculus to find this value, by first taking the derivative of the likelihood function with respect to the unknown parameter, then setting it equal to 0 and solving for the parameter.  Don’t forget to verify conditions so as to make sure you are indeed finding a maximum.

Derivative

This will usually involve complicated, messy math.  To mitigate this, we sometimes work with the logarithm of the likelihood function and use properties of logs to simplify computations.  This won’t chance our answer – taking the logarithm of some function won’t change the point at which the maximum value is achieved.

LogDerivative

The value of the parameter that you end up with maximizes the probability of your sample values x1, …,xn.  You could say it’s the value “most consistent” with the observed sample – the Maximum Likelihood Estimate.


Deriving the Present Value and Future Value of an Annuity Immediate

Below is the derivation of the present and future value of a unit annuity immediate, or a series of $1 cash flows that occur at equal intervals of time at the end of each period.  I originally wrote this document as a review for myself in preparation for actuary exam FM/2.  The majority of questions on the exam, despite the wide array of topics covered, come down to solving for the value of some annuity.  Granted, it likely won’t be a case as simple as the one below, but many problems about loans, bonds, yield rates, and even financial derivatives biol down to an annuity problem.

Annuity_Derivation


Binomial Distribution Basics

The Binomial Distribution is a discrete probability distribution of the number of successes in n independent trials with a binary outcome and probability of success p.  There are three important conditions that must be met for the binomial distribution to be used correctly:

  • Trials are independent
  • The probability of success p is constant for all trials
  • Each trial has only two possible outcomes (binary)

One way to describe a Binomial random variable is in terms of Bernoulli Trials.  A Bernoulli trial is a random experiment in which there are only two possible outcomes – success and failure.  If we let success = 1 and failure = 0, then a Bernoulli Random Variable looks like this:

P(X = 1) = p

P(X = 0) = 1 – p

So really a Binomial experiment is just a sequence of n repeated Bernoulli trials.  Formally, a random variable X that counts the number of successes, k, in n trials, has a Binomial Distribution with parameters n and p, written X~Bin(n, p).  Let p = Probability that success occurs at any given trial.  Then the probability of k successes is:

Screen shot 2013-11-14 at 2.10.11 AM

  • The first part of the equation is “n choose k” or a combination of k objects from a population of n objects.  This expression is the number of ways we can arrange k successes in n trials, where the order of the k successes is unimportant.  This is called a Binomial Coefficient.
  • The second part of the equation is the expression p raised to the power k, where p is the probability of success in any given trial.  This expression takes advantage of the independence assumption – independent probabilities are multiplicative, so the probability of p happening k times is p*p*p*…*p (k times) = p raised to the k power.
  • The third part of the equation is the probability of realizing n-k failures where trials are independent and the probability of failure is (1-p).  Intuitively, if there are n trials and k of those are successes, the remaining trials, of which there are n-k, must result in failure.

It might help to think backwards and ignore the binomial coefficient for a moment.  If k independent events occur with probability p, and n – k independent events occur with probability (1 – p), then the probability that all mentioned events occur is obviously the product (pk)*(1-p)n-k.  While what we’ve just stated is accurate for one specific sequence of k successes and n – k failures, we must consider that there are a number of combinations of outcomes that yield k successes and n – k failures.  If we flip a coin three times and define success as landing on heads, there is more than one outcome for k = 2.  We could have HHT, HTH, or THH.   The exact number of outcomes that yield k successes and n – k failures is the binomial coefficient “n choose k”.  That is why it is multiplied by probability of success and failure.

Ex) Three coins are tossed, each with probability p that the coin lands on heads (it may not be a fair coin, so p may not be 0.5).  Since coin tosses are independent, we can write:

Screen shot 2013-11-14 at 2.11.36 AM

The table below shows the probabilities associated with each sequence:

Screen shot 2013-11-14 at 2.12.20 AM

The chart above gives the probability for each outcome – that is a sequence of heads and tails.  It would be incorrect to look at row 6 and then assume that the probability of obtaining 2 heads is the probability in column 4.  The probability there is the probability of obtaining exactly the sequence TTH, not the more general case of obtaining two heads.

But if our main interest is the number of heads that occur, the sequence HHT is the same as HTH.  Note that there are three outcomes with exactly two heads, for example, each with probability p2(1 – p).  This is where the “n choose k” term becomes necessary.  Exactly k successes may be achieved in a number of different ways.  In the context of this question, the binomial probability is:

P(k heads) = (# ways to arrange k heads and n-k tails)*(prob. of any specific sequence w/ k heads and n-k tails)

P(k heads) = (# ways to arrange k heads and n-k tails)*pk (1-p)n-k

Implicit in the operations just performed is that a random variable was used to map outcomes in a sample space to the real numbers.  Basically, we used a random variable to assign a real number (0, 1, 2, or 3 in this case) to the outcomes (sequences of H and T).  A random variable is not random nor is it a variable – it is a function that maps sample outcomes to the real numbers.  This is a very simplified definition that will not suffice in a mathematical statistics course most of the time, but that’s basically what’s going on here.  The outcomes were grouped in such a way to facilitate problem solving – in this case, it was used to give us the number of heads rather than a specific sequence.

Just as the probability function assigns probabilities to sample outcomes in the sample space, a probability density function assigns probabilities to values that a random variable takes on.  Let X be the number of heads that occur in a sequence of 3 coin tosses.  Then the probability density function is:

Screen shot 2013-11-14 at 2.13.51 AM

For example:

Screen shot 2013-11-14 at 2.14.30 AM

And here is a chart of the actual probabilities:

Screen shot 2013-11-14 at 2.15.01 AM

It should be clear now why independence of trials is such an important criterion in using the Binomial distribution.  If trials are not independent, i.e. the outcome of one trial affects the probability of subsequent outcomes, then the binomial model is not applicable.  In this situation, the Hypergeometric Distribution should be used instead.


Paradoxical expected values

We tend to think of the expected value of a random variable as a standard, practical, and simple measure of the variable’s central tendency.  However, there are situations in which the expected value is misleading or even nonexistent.  One famous case, commonly called the “St. Petersburg Paradox”, illustrates perfectly the limitations of the expected value as a measure of central tendency.

Suppose a fair coin is flipped until the first tail appears.  You win $2 if a tail appears on the first toss, $4 if a tail appears on the second toss, $8 if a tail appears on the third toss, etc.  The sample space for this experiment is infinite: S = {T, HT, HHT, HHHT, …}.  A game is called fair if the ante, or amount you must pay to play the game, is exactly equal to the expected value E(x) of the game.  Casino games are obviously never fair, because casinos stand to earn a profit (on average) from the games they offer – otherwise, they would not have a viable business model and they couldn’t afford to give you free drinks and erect fancy statues (though, while we’re on the subject, craps is the “fairest” as far as casino games go in the sense that the odds are skewed least heavily in the dealer’s favor).

Let the random variable W denote the winnings from playing in the game described above.   We want to categorize W such that we can mathematically evaluate its expected value.  Once we find this expected value, we know how much we’d have to ante for this game to be considered fair.  We know how much we’ll win on the first, second, and third toss, but let’s generalize for k tosses:

  • T on 1st toss -> $2
  • T on 2nd toss -> $4
  • T on 3rd toss -> $8
  • T on kth toss -> $2^k

So we win $2^k if the first tail appears on the kth toss.  But what is the probability of each of those outcomes?  For toss one, there is a ½ chance of getting a tail.  Since trials are independent, the probability of the first tail appearing on the second toss is (1/2)*(1/2) = 1/4.  For the third toss, the probability is 1/8.  So the probability of winning $2^k is the probability of getting the first tail on the kth toss which is:

Screen shot 2013-11-01 at 7.01.14 PM

To find the expected value of the random variable W, we need to find the sum over all values that k (remember that k = the number of tosses before the first tail) of k times the probability of k.  It’s evident from the sample space that k can take on an infinite number of values.  The expected value is:

Screen shot 2013-11-01 at 7.01.21 PM

So the expected value of the random variable W (winnings) is the sum over all k of the payoff (2^k dollars) times the probability of realizing that payoff (p(2^k)).  So the first term is the actual payoff, the second is the probability of that payoff.  Plugging in what we know about the probability of k and the infinite nature of the sample space:

Screen shot 2013-11-01 at 7.01.30 PM

The sum diverges!  There is no finite expected value.  That is, the expected value of the game is infinite, which means we’d have to pay an infinite amount of money for the game to be fair.

The point is not that you should wager an infinite amount of money in order to play the aforementioned game.  The point is that the expected value is often an inappropriate measure of central tendency that leads to an inaccurate characterization of a distribution.


Criminology and Combinatorial Probability

I just read a really interesting application of combinatorics to criminology in my stats book, An introduction to Mathematical Statistics and Its Applications, by Larsen and Marx.  Alphonse Bertillon was a French criminologist who developed a system to identify criminals.  He chose 11 anatomical characteristics that supposedly remain somewhat unchanged throughout a person’s adult life – such as ear length – and divided each characteristic into three classifications: small, medium, and large.  So each of the 11 anatomical variables are categorized according to their respective size, and the ordered sequence of 11 letters (s for ‘small’, m for ‘medium’, and l for ‘large’) comprised a “Bertillon Configuration”.

This seems rudimentary and imprecise, perhaps, but keep in mind this was before the discovery of fingerprinting and so forth.  One glaring issue with this system is the obvious subjectivity involved with assigning a classification to a given anatomical variable.  Also, I’d imagine it’s kind of hard to get all of this information about a person, especially if they’re suspected of having committed a crime but haven’t ever been in custody.

Issues aside, the Bertillon configuration is an interesting idea I think.  But the obvious question to ask is how likely is it that two people have the exact same Bertillon configuration?  The implications are two people who, from a legal classification standpoint, are exactly the same (and you can use your imagination from there).  So, how many people must live in a municipality before two of them necessarily share a Bertillon configuration?

The question is actually very simple to answer via the multiplication rule:

If Operation A can be performed in m different ways, and Operation B can be performed in n different ways, the sequence {Operation A, Operation B} can be performed in (m)(n) different ways.

Corollary: If Operation Ai, i = 1, 2, …, k can be performed in ni ways, i = 1, 2, …, k, respectively, then the ordered sequence, {Operation A1, Operation A2, …, Operation Ak) can be performed in (n1)(n2)***(nk) ways.

Think of Operation A as the anatomical variable and Operation B as assigning it a value (s,m,l).

Screen shot 2013-10-17 at 10.08.05 AM

Each Operation can be performed in three different ways, namely s, m, and l.  Since there are 11 anatomical variables, by the multiplication rule there are (3)*(3)*(3)*(3)*(3)*(3)*(3)*(3) *(3)*(3)*(3) = 311 ways to assign arrange the Bertillon configuration.  Therefore, if 311 = 177,147 people live in a town, we’d expect two of them to have the same Bertillon configuration.


Birthday Statistics

Today is my roommate’s birthday and three of my friends’ birthdays.  There are 365 days in a year, so what are the chances that so many birthdays are on the same day?  More importantly, what are the chances that I have to spend so much time and money on gifts?

The probability is higher than you’d expect – so high, in fact, it has a special name, the birthday paradox.  Informally, the birthday paradox is the surprising probabilistic phenomenon that for a relatively small group of people, the probability that at least two people share the same birthday is surprisingly high.  For example, if you have 23 people in your stats class, there’s roughly a 50% chance that at least two people have the same birthday.

One way to think about this result is that, at first, it seems like there are only 22 comparisons to be made – essentially 22 ‘chances’ that one person shares his or her birthday with someone else in the room.  This is true for the first person, who we’ll call Bob, but what happens when Bob is done comparing his birthday to everyone else’s?  Well, then the second person, who we’ll call Mary, compares her birthday to everyone in the room except for Bob, who she has already compared birthdays with.  So Mary has 21 people to compare with.  The third person has 20 potential matches, the fourth has 19, and so forth.  So what you really have is 22+21+…+1 = 253 comparisons, not 22.

Each individual comparison has only a 1/365 chance of being a match – equivalently, each trial has a 99.726% chance of not being a match.  This is probably why we’re so surprised by this paradox – we’re conditioned to think it’s rare to have the same birthday as someone, because when you ask any given person when their birthday is there is only a .274% chance you’ll have the same birthday as him or her.  But each time we compare, there is one less possibility – i.e. there is one more birth date that isn’t available, so to speak.  Below is the probability that no two people in the room share a birthday.  Notice that the first term is 365/365 because the first person has 365 birthdays to choose from without ‘taking’ someone else’s.  The second person, however, has one less choice, because the first person’s birthday has already been taken.  The third person has one less choice than the second person.  This continues until the last person, person number 23, has 22 fewer choices (365-22=343) of birthdays if he is to avoid having the same birthday as another person in the room.

Screen shot 2013-09-17 at 8.21.26 AM

The probability that two people in the room do share a birthday is obviously the complement of this, so it’s equal to 1 – 0.4927 = 50.73%.  You can write the expression above with factorials etc., but I just wanted to sort of explain why (though not actually justify in a rigorously mathematical sense) the birthday paradox is a thing.  For a more hardcore explanation, you can visit Wolfram Mathworld here.


Continuously Compounded Returns

I’m in a financial “math” class that I really shouldn’t have taken (it’s in the business school) and our lecture notes the other day included some incredibly basic properties about interest calculations.  Anyway, I told one of my friends that it was super fucking easy to derive the formula for continuous interest by taking the limit as the number of compounding periods approaches infinity.  Upon further review, I stand by the fact that it is conceptually easy to do this, but there are some computational issues that one could run into – for example, I had to use a relatively simple substitution in order to avoid using L’Hopital’s rule or something else messy like the definition of a difference quotient, etc.  But, regardless, below is the derivation, starting with the definition of the effective annual rate for an interest rate, r, compounded n periods per year for t years.  Suppose there are n compounding periods per year and r is the interest rate.

Screen shot 2013-09-17 at 2.46.18 AM

How do we show that the effective annual rate under continuously compounded interest (i.e. effective interest with arbitrarily large n) is er – 1?  We need to show the equation below, which equates the limit of the EAR as the compounding periods approach infinity to the formula we claimed represents the continuously compounded rate (Exp(rt)-1):

Screen shot 2013-09-17 at 2.45.36 AM

If I were being tested on this in a math class for whatever reason, I might not show the continuously compounded interest rate this way.  This is more of an intuitive justification than a real proof.  To see a real proof, you can check out notes from Wharton here – it’s not really any “harder” per se, but some of the steps might be less obvious.  For example, the proof involves taking the natural logarithm of both sides of the equation (as shown below) which serves a purpose, though it may not seem like it at first.

Screen shot 2013-09-17 at 2.59.49 AM

This was done in order to take advantage of a useful property of logarithms: the log of a product is the sum of the logs.  Using this, the right hand side is split up into a sum and the first term no longer depends on the limit because it doesn’t have an n in it.  You’re left with a sort of messy equation that involves logs on both sides (it’s implicitly defined), so to rectify that you just need to remember that e and the natural log are inverses, so e raised to the power ln is 1 and the ln of e is 1.

Note: just so you don’t sound like an asshole, e is Euler’s number, and it is pronounced like “oiler”.

Another note:  There are tons of definitions of e, but the one I think of first when I think of e is the one below:

Screen shot 2013-09-17 at 3.06.44 AM

In words, e is the positive real number such that the integral from 1 to e of the function 1/t is equal to 1.  Equivalently, the area of the integral of 1/t from 1 to e is 1.  ln(x) crosses the x-axis when x=1 (since ln(1)=0), so the area in between the function ln(x) and the x-axis from its x-intercept to x = e is 1.  Obviously, ln(x) takes on the value 1 when x=e since ln(e)=1.  I’m not good at making fancy graphs online, but below is a heinous picture of the point I’m trying to convey.  Sorry Euler 😦

Screen shot 2013-09-17 at 3.16.06 AM


Intro to Value at Risk

Value at risk (VaR) is a way to statistically measure the worst expected loss on an investment given certain parameters.  Specifically, the level of confidence (called the alpha level, α) and time horizon must be specified to compute it.  We must also assume that returns on the asset in question are approximately normally distributed, because the value at risk model calculates the worst expected loss over a given horizon at a given level of confidence under normal market conditions1.  Although it’s a well-established fact that, in general, stock prices are log normally distributed while stock returns are normally distributed.  It may be a good idea to verify this first, using a goodness of fit test.  The Anderson-Darling, Shapiro-Wilk, and  Kolmogorov-Smirnov tests are all good examples.  When we say normal, we mean that the data is normally distributed; we are making a judgment regarding the distribution of returns, not the likelihood of realizing a return, as the dictionary definition of normal would suggest.  As we’ll see, the VaR turns out to be a relatively unusual result in most cases, and its rarity is inversely related to the level of confidence we choose.  It follows that a formal definition of the value at risk is a measure of the worst expected loss over a given time frame and level of confidence under normal market conditions.

The actual VaR number should be the worst expected loss on an asset in terms of currency (not probability).  If we let f represent the probability density function of the return2 of the asset and c the confidence interval, we have the following formula:

Screen shot 2013-08-13 at 4.16.32 PM

In words, this formula equates the confidence level, or degree of certainty with which we can draw a conclusion, with the integral from -∞ to the negative of our VaR value (some real number measured as a return in dollars) of the probability density function (pdf).  The Fundamental Theorem of Calculus tells us that integrating over a range gives us the area between the x-axis and the function.  In the context of a probability density function, whose y-axis is probability, this says that by integrating from A to B, we’re finding the probability that the random variable takes on a value between A and B.  In the context of VaR, when we integrate from negative infinity to –VaR, we’re finding the probability of realizing a value equal between -∞ and –VaR, or, more simply, the probability of realizing a value less than –VaR.  So for a given confidence level c, we’re looking for the value –VaR such that the integral from -∞ to –VaR is equal to (1 – c).

It might be confusing that the negative of the VaR is used in this computation.  This arises because of the symmetry property of the normal distribution.  Each value, however extreme (extreme meaning far from the mean in terms of standard deviations) has a complementary value in the opposite direction.  If we’re talking about the standard normal distribution, we’re just saying that the distribution is normally distributed with a mean of 0 and a standard deviation of 1.  So a relatively extreme value, say 3, has a complementary value, -3, that is just as extreme in the opposite direction.  3 is a value 3 standard deviations greater than the mean, and -3 is 3 standard deviations below the mean.  But we’re only interested in a potential loss in VaR models, as risk is by definition a quantitative measure of realizing a loss, we’re only interested in the –VaR.  The corresponding positive VaR will be a value just as extreme but positive – in this case, a good thing since it means higher returns.

Every potential portfolio value corresponds to a z score, or a standardized measure of how many standard units (standard deviations) a value lies from the mean of the distribution, as was mentioned above.  Ordinarily, we transform a random variable X into a standardized Z score with the formula Z = (x – μ)/σ, where, in this case, x = VaR and Z = Z value corresponding to the confidence level.  However, since we’re only interested in the negative values, we substitute x = -VaR and Z = -Z.  We can rewrite the equation more simply now: VaR = Zσ – μ.  This equation will standardize our VaR value.

So if you’d like to be 95% certain you won’t lose more than the VaR, you set the left hand side of the equation equal to (1 – 95%) = (1 – .95) = .05 and solve for the number VaR, which can then be interpreted as the 95% value at risk.  What you’ve really found here is the portfolio value that represents the spot on the x axis where the total area above the axis and below the pdf on the left hand side all the way to x = negative infinity is equal to .05.

If we use the 99% confidence level, we have c = (1 – .99) = .01.  We have to find a z score that divides the normal distribution into a section with 1% of the area to the left and 99% of the area to the right.  The picture below has 1% of the area in each tail of the distribution, and thus 98% shaded in the middle.

Screen shot 2013-08-13 at 1.48.48 PM

We’re interested in the z score that divides the lower white region from the rest of the graph, because that will correspond to the value that 99% of the x values are greater than.  But since the normal curve is symmetrical, that value is just going to be the negative of z value, which represents the x value with 99% of observations below it.

Screen shot 2013-08-13 at 2.24.26 PM

The pictures reveal that the z value in question is z = -2.33.  To reiterate, standardizeing the value of any random variable x (i.e. assigning it a z value), uses the relationship Z = (x – mean)/(StanDev).  This Z value is the VaR in standardized form.  Therefore, we can substitute –VaR for Z in computations because what we’re really finding when we find the VaR is the left tail (Z score) of the distribution.

If the VaR is 200 at the 99% confidence level, you could say that largest possible loss on your asset at the 99% confidence level is 200, given normal market conditions.  Alternatively, you could be 99% certain your losses won’t exceed 200.

Since the cumulative density function (F) is the indefinite integral of the probability density function (f), we can use it to simplify this expression (sorry for the bullshit above the equation – I had to screenshot it because wordpress hates equations):

Screen shot 2013-08-13 at 4.20.37 PM

An example should make this clearer.  Suppose the daily return on a portfolio follows a normal distribution with a mean of $1000 and a standard deviation of $500.  At the 99% confidence level, we have alpha = .01 which corresponds to a z-score of 2.33, i.e. the CDF of the normal distribution is equal to .99 when the z value is 2.33; F(2.33) = 0.99.  But we want the lower end of the distribution, -2.33.  Therefore:

VaR = ασ – μ = 2.33(500) – 1000 = 165

So with 99% confidence we can assume that our worst possible expected loss in one day is $165.  The value of -$165, i.e. a $165 loss, is 2.33 standard deviations (z scores) below our mean value of $1000.  So although we still expect (on average) to earn $1000 with a standard deviation of $500, we will only lose more than $165 in a single day less than 1% of the time.  The PDF and CDF are below.

pdf

cdf

So the value -165 is exactly equivalent to the standardized z score of -2.33 as they both represent the value that 1% of observations fall below and 99% exceed.  The only difference is that with the number -2.33 we’re talking about a standardized random variable with a mean of 0 and a standard deviation of 1, whereas the VaR number of -165 represents a portfolio return with a mean of $1000 and a standard deviation of $500.  Alternatively, -$165 is the value 2.33 standard deviations to the left of the $1000 mean portfolio return.


www.openeuroscience.com/

Open source projects for neuroscience!

Systematic Investor

Systematic Investor Blog

Introduction to Data Science, Columbia University

Blog to document and reflect on Columbia Data Science Class

Heuristic Andrew

Good-enough solutions for an imperfect world

r4stats.com

"History doesn't repeat itself but it does rhyme"

My Blog

take a minute, have a seat, look around

Data Until I Die!

Data for Life :)

R Statistics and Programming

Resources and Information About R Statistics and Programming

Models are illuminating and wrong

A data scientist discussing his journey in the analytics profession

Xi'an's Og

an attempt at bloggin, nothing more...

Practical Vision Science

Vision science, open science and data analysis

Big Data Econometrics

Small posts about Big Data.

Simon Ouderkirk

Remote Work, Small Data, Digital Hospitality. Work from home, see the world.

rbresearch

Quantitative research, trading strategy ideas, and backtesting for the FX and equity markets

Statisfaction

I can't get no

The Optimal Casserole

No Line Is Ever Pointless

SOA Exam P / CAS Exam 1

Preparing for Exam P / Exam 1 thru Problem Solving

schapshow

Mathematical statistics for the layman.