Hey guys, in this video I'm going to go over the statistical technique called the Box-Muller
transform and this is a transformation that's used to sample from the normal distribution. In order
to generate random numbers from a distribution, it's common to use another type of transformation
called the inverse transform sampling method but to use that method you have to know the CDF and
then the inverse CDF of the distribution you want to sample from. So for the normal distribution,
there's no well-defined CDF. So there's no closed form expression for the normal distribution,
the CDF of the normal distribution. In that case you can't use the inverse sampling method
to sample from the normal distribution. So the Box-Muller transform is one way of sampling from
the normal distribution and I'm going to go over the algorithm in this video and then at the end
I'll go over some R code as well. The Box-Muller transform uses polar coordinates and also
uniform random variables in order to transform into samples from the normal distribution. It
uses polar coordinates and uniformly distributed variables as part of the transform to get to
samples from the normal distribution. Let's go ahead and get started with the video on how the
transform works. The end result is going to be two independent random variables that are distributed
normally and so for the sake of this video let's say we want to sample from the standard normal
distribution that has mean 0 and variance 1. We have two normally distributed random variables
that are independent that I'm labeling X and Y. To sample from X and Y, we're going to use polar
coordinates. Let's say every random sample from X is labeled lowercase x and let's say every random
sample we generate from the random variable capital Y we're going to label as lowercase y.
So x and y are two independent normal samples from the normal distribution. We can consider
these random samples as a coordinate are in the Cartesian plane. These samples can be considered x
comma y and we can plot that on a Cartesian plane as shown here. Then from this consideration of
the samples as points, we can then convert the points or the samples into polar coordinates.
To do that, we need to know r, which is the distance from the origin to the point x, y. So
that's r labeled here. And we need to know theta, which is the angle from the x-axis to the point or
to the line from the point. The formula for r is calculated just from a simple distance from the
point to the origin. So r is equal to the square root of x squared plus y squared. And then we can
square both sides to get that r squared is equal to x squared plus y squared. Now if we go back to
thinking of x and y as random variables that are normally distributed, then r squared, which is x
squared plus y squared, you can think of this as being the square of a normal distribution for x
squared. That's going to be chi squared with degrees of freedom one and y squared is also
chi squared with degrees of freedom one. And then the entire expression is the sum of two chi-square
distributions. That would be distributed chi squared with degrees of freedom two.
So basically you square a normal distribution, so that gives you a chi squared, and then you
take the sum of two chi-square distributed distributions with degrees of freedom one
and that would give you a chi squared with degrees of freedom 2. Written out another way, r squared
is equal to x squared plus y squared and that is chi squared distributed with degrees of freedom 2.
Now this chi-squared distribution with degrees of freedom two can also be transformed equivalently
to other distributions. This has been proven in the research literature and I will put a link
in the description box below with an article that has all the proofs for this part but basically a
chi-squared degrees of freedom two distribution is equal to the exponential distribution with
lambda equal to one-half. Again I'm not proving this part but this has been shown to be true in
the literature and you can rewrite or express an exponential distribution in terms of the negative
log of a uniform distribution. For exponential one-half, that is equivalent to -2 natural log
of a uniformly distributed random variable U. I'm going to label that as U1 where U1 is uniformly
distributed from 0,1. And so r squared is distributed -2 log U. And so this--again
this transformation from chi-squared to uniform has been proven to be true in the literature.
I'll put a link in the description box below but basically we have our first uniform random
variable that we can sample from and then we can transform the uniform into a normal variable.
Basically we're working backwards. We want to be able to get samples from these two normal
random variables and the Box-Muller does that by transforming them into uniform random variables,
which we can sample from directly. We did that for r and we have to do the same thing for theta. For
theta, that's going to be able to be distributed anywhere in the unit circle. That's just going to
be 2 pi U2, where U2 is uniformly distributed 0,1 as well. Now we know r2 and theta and we're able
to transform back to the Cartesian coordinates. We know what the final polar coordinates are. For r,
that's going to be equal to the square root of -2 natural log U1 and for theta, that's going to be
equal to 2 pi U2 and then we can convert to the Cartesian coordinates to get the final equations
that we're actually going to use to sample from the normal distribution. For x, which is the first
random variable, that's going to be equal to r cosine theta. And again that just comes from the
definitions from polar coordinates. And so that's going to be equal to the square root of -2 log U1
times cosine 2 pi U2 and y is going to be equal to r sine theta, which is equal to the square root of
-2 log U1 times the sine of 2 pi U2. And so these equations can be used to sample from the normal
distribution and you will get two independently normally distributed normal 0,1 samples
by using and plugging in two uniformly distributed random variables. Basically we have random samples
that we can collect from U1 and U2 and we can use the Box-Muller transformation to get X, Y,
which are normally distributed 0,1. And so that is the basic idea with the algorithm and then in
the next part of the video I will go over an R demo. All right so in this part of the video, I
just have a short snippet of code to show visually that the samples that we are generating from the
Box-Muller transform are normally distributed. In this case I'm going to get two variables, which
I'm labeling samples1 and samples2, which are each going to be independent and normally distributed
0,1 and they're going to have a thousand samples each. For n, that's going to be my sample size. So
for each variable I want a thousand samples from the normal distribution. And then the next lines
of code I'm just going to fill the two vectors or two variables with NAs. And then this is the line
of code that runs the Box-Muller transform. So for the first two lines of code in the for loop, I'm
just sampling from the uniform distribution for my random variables U1 and U2. Then in this line
of code, I'm transforming it for the r distance equation. I'm taking the square root of -2 log
U1 and then for theta I'm doing 2 pi U2 and then I'm plugging in r times cosine theta for X, or
samples1, and then for samples2 or y, I'm plugging in r sine theta. I'm going to run that for loop
and then my samples1 and samples2 variables should be normally distributed. I'm going to compare that
to the normally distributed random samples that are generated in R from the code rnorm.
To do that I'm just running these next couple lines of code and basically what I want to do is
just create a histogram that I want to show. This final line of code will run a histogram in ggplot2
and when that finishes running we see the visualization here. And also by the
way I will put the link to this code in the description box below but you can see that the
blue values is the histogram for the samples that were--this is like the X random variable--and
so those were the samples that were generated from the Box-Muller transform and then the red samples
were a thousand samples generated from rnorm, which is the way you generate random samples from
the normal distribution in R. You can see that I have both those histograms overlapping and they
roughly both look to be normally distributed. You can see that there's some random variation between
the two histograms but they both look roughly normally distributed and if you plotted samples2,
you would get the same sort of result. And that's it for this video and thank you guys for watching!
*chiptune music*