Probability

Data as vectors

The usual way you would recieve data is a matrix with rows being observations and columns being features. Every observation can be thought of as a vector, and these vectors are the empirical values for the “random variable” that is a general vector which is a list of such features. For example, for this matrix $\bold{X}$ :

The random variable $\bold{x} \in \mathbb{R}^2$ , which is a list of features (scalar random variables), and has the observed values of $[1\;2]^T,[3\;4]^T,[4 \;5]^T$ and $[5\;6]^T$ .

Note : Any random vector $\bold{x}$ is a column vector, unless stated otherwise.

The probability of of $\bold{x}$ taking a certain value $\bold{x_1}$ (in its full domain, not only in the observed values) is denoted as $p_\bold{x}(\bold{x_1})$ , or for our convenience, simply as $p(\bold{x})$ .

If the case that the domain of $\bold{x}$ contains infinite elements and the probability for one value is infinitesimal, $p(\bold{x})$ refers to the probability distribution function, but I’ll still refer to it as the “probability” as if we are always taking about a discrete dirtribution.

For an event described by multiple random variables, the probability is written as $p(\bold{x,y,z,\dots})$ . In such a setting, if there are less random variables inside the bracket than required to describe the event, then we are talking about the probability of any of the events with the specified random variable taking the values in the brackets, regardless of what values the other random variable take. For example, if an event is described by three random variables $\bold{x,y,z}$ , then $p(\bold{x,y})=\sum_\bold{z} p(\bold{x,y,z})$ .

Since we can always concatenate multiple random variables describing an event in a single big vector, so we only need to talk about 2 random variables at max.

We use ‘ $|$ ‘ to mean “given” it is used to denote conditional probability, defined (yes, I’m using this equation to define the probability, rather than the usual set theory based definition) as

This is the probability of and event having a certain $\bold{x}$ , given that it has $\bold{y}$ as specified.

From here it’s easy to derive Baye’s rule.

Bayes’ rule

One good example where it’s used is the Gaussian Naive Bayes’ Classifier.

A simpler example would be the $\beta$-distribution, which is probability distribution of the probability $\theta$ of a person winning a match, if we know that he has won $\alpha$ matches and lost $\beta$ assuming that this probability of winning stays the same in every match. Here we are trying to find $p(\theta\;|\;\alpha,\beta)$ . we already know that $p(\alpha,\beta \;|\;\theta) =\; ^{\alpha+\beta}C_\alpha\, \theta^\alpha \,(1-\theta)^\beta$ and we can find $p(\alpha,\beta)$ by integrating $p(\alpha,\beta\;|\;\theta)$ over $\theta$ , which gives us $^{\alpha+\beta}C_\alpha(^{\alpha+\beta}C_\beta(\alpha+\beta+1))^{-1} = (\alpha + \beta - 1)^{-1}$ and we assume that $p(\theta)=1 \;\;\forall\; 0\leq\theta\leq1$ and $0$ otherwise. This gives us that :

This is called the beta distribution.

Using Bayes’ rule if you have some prior pdf for $\bold{x}$ , given more data in the form of vectors $\bold{y}_i$ used as row vectors in the matrix $\bold{Y}$ , you can update $\bold{x}$ . We want to find the posterior pdf for $\bold{x}$ , that is $p(\bold{x\;|\;Y})$ . This is calculated as $\frac{p(\bold{Y\;|\;x})\,p(\bold{x})}{p(\bold{Y})}$ . Since the observations are all independent to each other, thus we can break this as $\prod_i\frac{p(\bold{y_i|x})\,p(\bold{x})}{p(\bold{y_i})}$ . This means we need to know $p(\bold{x\;|\;y})=\frac{p(\bold{y\;|\;x})\,p(\bold{x})}{p(\bold{y})}$ first. If $\bold{y}$ was in the same domain as $\bold{x}$ , then we could just calculate $p(\bold{y\;|\;x})$ to be the delta function in $\bold{x}$ centred around $\bold{y}$ , and then we would get the same for $p(\bold{x\;|\;y})$ . So it makes no sense to take them to be in the same set of features. Although, we can certainly take them to be in a set of dependent features.

Statistically Independent variables :

Two random variables $\bold{x,y}$ are statistically independent iff $p(\bold{x,y}) = p(\bold{x})\,p(\bold{y})$ .

Expected values

For a function $f(\bold{x})$ of the random variable $\bold{x}$ , the expected value is given by :

It’s easy to see that $\text{E}(\bold{x})$ is the mean of $\bold{x}$ .

Note : The definition is only valid for discrete distributions

Covariance and variance

Suppose you you want to know how much are 2 feature $x,y$ similar (linearly dependent), given some data in the form of a matrix $\bold{X}$ with these two features present in the observations, then you can just take the dot product of the column vectors of these two features after subtracting the mean from each entry (This is important to do as if these quantities don’t vary a lot but have high values all the time, then the covariance of such quantities must intuitively be small). Then, divide by the length of the column vectors since you don’t want this value to depend on the size of the data-set. This quantity is called the covariance. Basically,

Here $[i]$ means “in the $i^\text{th}$ observation” where there are $n$ observations.

And the covariance of a feature with itself is its variance.

The last equality is obtained by expanding the square inside the summation and using the fact that $\text{E}(x)$ is a constant.

The covariance can be used as a definition of an inner product between two random variables. Then the induced norm would be the square root of the variance. And the angle $\theta$ between two random variables $x,y$ would be given by :

This value is called the “Correlation” .

When the random variable $\bold{x}$ is a vector, we can express the covariances of all pairs of features $x_i,x_j$ in a covariance matrix, which is basically the “Variance” of $\bold{x}$ .

Again, $\bold{X}$ is the matrix representing the data-set, with rows being $\bold{x}[i]^T$

Now, for taking the covariance of two vector valued random variables, $\bold{x,y}$ , we do ;

Suppose two such random variables are linearly related as $\bold{y = Ax + b}$, then the covariance is :

Here $\text{Var}(\bold{x}) = \bold{\Sigma}$ is the covariance matrix of the features that $\bold{x}$ is composed of.

And thus ,

The first equality is easy to prove because $\text{Cov}(\bold{x,y})$ is linear in $\bold{x}$ for a fixed $\bold{y}$ .

Sum of variables

For any two random variables, $\bold{x,y}$ , the variable $\bold{z=x+y}$ is a random variable with $\text{E}(\bold{z}) = \text{E}(\bold{x}) + \text{E}(\bold{y})$ and

$\text{Var}(\bold{z}) = \text{Var}(\bold{x}) + \text{Var}(\bold{y}) + \text{Cov}(\bold{x,y}) + \text{Cov}(\bold{y,x})$ .

Statistically independent variables are linearly independent

For statistically independent variables $\bold{x_i,x_j}$

$\text{Cov}(\bold{x_i,x_j}) = \text{E}((\bold{x_i}-\text{E}(\bold{x_i}))((\bold{x_j}-\text{E}(\bold{x_j}))^T) = \sum_\bold{x_i,x_j}(\bold{x_i}-\text{E}(\bold{x_i}))(\bold{x_j}-\text{E}(\bold{x_j}))^Tp(\bold{x_i,x_j}) = \sum_\bold{x_i,x_j}(\bold{x_i}-\text{E}(\bold{x_i}))(\bold{x_j}-\text{E}(\bold{x_j}))^Tp(\bold{x_i})p(\bold{x_j}) = \sum_\bold{x_i}(\bold{x_i}-\text{E}(\bold{x_i}))p(\bold{x_i})\sum_\bold{x_j}(\bold{x_j}-\text{E}(\bold{x_j}))^T p(\bold{x_j}) = \text{E}(\bold{x_i}-\text{E}(\bold{x_i})) \text{E}(\bold{x_j}-\text{E}(\bold{x_j}))^T =\bold{ 0 0^T } = \bold{0}$

I.I.D. random variables

Random variables $\bold{x_1,x_2,x_2,\dots,x_n}$ are independent and identically distributed if :

• They are statistically independent , i.e. $p(\bold{x_1,x_2,\dots,x_k}) = p(\bold{x_1})\,p(\bold{x_2})\dots\,p(\bold{x_k}) \forall k \leq n$

• All of their distributions are the same function. This can be summarised as : $\bold{x_i=x_j} \;\implies\;p(\bold{x_i})=p(\bold{x_j})$ .

It’s easy to see that in a case like this, the variables are also linearly independent .

This is useful when you want you take $n$ observations with each observation coming from the same distribution and predict beforehand the properties of these observations. For example, the expected value of the computed mean of observations $\bold{x_i}$ is the same as the expected value for a single observation, say $\boldsymbol{\mu}$, and the variance is thus $\text{E}((\frac{1}{n}\bold{\sum_ix_i} - \boldsymbol{\mu})(\frac{1}{n}\sum_i\bold{x_i} - \boldsymbol{\mu})^T) = \frac{1}{n^2}\text{E}((\sum_i\bold{(x_i-\boldsymbol{\mu_i})})(\sum_i\bold{(x_i-\boldsymbol{\mu_i})})^T) = \frac{1}{n^2}\sum_i\sum_j \text{E}(\bold{(x_i-\boldsymbol{\mu}_i)(x_j-\boldsymbol{\mu}_j)^T}) = \frac{1}{n^2} \sum_{i,j} \text{Cov}(\bold{x_i,x_j}) = \frac{1}{n^2} \sum_{i} \text{Cov}(\bold{x_i,x_i}) = \frac{1}{n^2} \sum_{i} \text{Var}(\bold{x_i}) = \frac{1}{n^2}(n\bold{\Sigma}) = \frac{1}{n}\boldsymbol\Sigma$

The 4th equality is because of the linear independence of the variables.

Gaussian distribution

This distribution is important to the next theorem. Although it lacks physical meaning, it has immense mathematical meaning.

For a random variable $x \in \mathbb{R}$ , if it follows the Gaussian distribution

then it has mean $\mu$ and variance $\sigma^2$ , where $\sigma = \sqrt{\text{Var}(x)}$ is the standard deviation .

To make life simpler, every time we want to say that a variable $x$ follows a Gaussian distribution with mean and variance $\mu,\sigma^2$, instead of saying it in words, we write $x \backsim N(\mu,\sigma^2)$.

For a vector valued random variable $\bold{x}\in \mathbb{R}^n$ , made of statistically independent scalar Gaussian random variables $x_1,x_2,\dots,x_n$ with $x_i \backsim N(\mu_i,\sigma^2_i)$ , we have :

then if $\text{E}(\bold{x}) = \boldsymbol{\mu}$ and $\text{Var}(\bold{x}) = \bold{\Sigma}$ , we have $\prod_i \sigma_i^2 = |\bold{\Sigma}|$ and $\sum_i\dfrac{(x_i-{\mu}_i)^2}{\sigma^2_i} = \bold{(x-\boldsymbol{\mu})^T\;\Sigma^{-1}\;(x-\boldsymbol{\mu})}$ , giving us :

It turns out that this equation is valid even when $x_i,x_j$ are NOT { independent OR Gaussian } , but $\bold{x}$ overall is Gaussian .

Marginals and conditionals of Gaussians

Suppose $\bold{x,y}$ are two random variables, which together (concatenated) form a random variable $\bold{z} = [\bold{x^T\;y^T]^T}$ which we know to be Gaussian, and suppose we want to find out the distribution of $\bold{x}$ given the value of $\bold{y}$ (Perhaps you want to predict something about a point with a specific $\bold{y}$ ). We do that by taking the conditional distribution $p(\bold{x\;|\;y})$ . This method isn’t specific to Gaussian distributions, but here we are interested in the result obtained when we consider Gaussian distribution specifically. Basically, to find $p(\bold{x\;|\;y})$, you need to find $p(\bold{x})$ , which we’ll prove to be a Gaussian distribution $N(\bold{\boldsymbol\mu_x,\Sigma_{x,x}})$ .

(Say that $\bold{\Sigma_{a,b} = \text{Cov}(a,b)}$ and $\boldsymbol{\mu}_\bold{a} = \text{E}(\bold{a})$ for any two random variables $\bold{a,b}$ .)

Rather than brute forcing our way through complex linear algebra, we’ll prove this inductively, by proving this for a scalar valued $y$, say the last coordinate of $\bold{z}$ to get a reduced set of coordinates being in a guassian distribution, and then using the fact that we can do this over an over and over untill we eventually delete an arbitary number of arbitary coordinates, making this true for any general $\bold{x,y}$ .

First let’s write the covariance and mean of $\bold{z}$ in terms of $\bold{x}$ and $y$

You are free to verify that the inverse of any positive definite symmetric matrix $S = \begin{bmatrix}A & B \\ B^T & D\end{bmatrix}$ where $A,B,D$ are sub-matrices is also positive definite symmetric, say $\begin{bmatrix}A' & B' \\ B'^T & D'\end{bmatrix}$ .

In our case $A = \bold{\Sigma_{x,x}},\; B = \bold{\Sigma_{x,y}},\; B^T = \bold{\Sigma_{y,x}},\; D = {\sigma_{y}^2}.$ In favour of my sanity, I’ll use $A,B,D$ for the linear algebra rather than bold greek symbols with subscripts all the time. Also, let’s call $p=\bold{x-\boldsymbol\mu_x}$ and $q={y-\mu_y}$ .

So now we have :

It’s easy to see that the stuff inside of the exp function (except the -1/2) evaluates to

Now since $\int_{-\infty}^\infty\exp(-\frac{1}{2}(q+\frac{p^TB'}{D'})^2)dq = \sqrt{2\pi}$ and other stuff is constant wrt $q$, thus we only have to worry about $p^T(A'-\frac{B'^TB'}{D'})p$ inside the exp function (except the -1/2) after the integration is done. Namely, $p(\bold{x}) \;\alpha \; \exp(-\frac{1}{2}(\bold{x}-\boldsymbol\mu_\bold{x})^T\bold{M}(\bold{x}-\boldsymbol\mu_\bold{x}))$ where $\bold{M} = A'-\frac{B'^TB'}{D'}$ is some constant square matrix. Now after normalising the expression to get a pdf, we’ll just end up getting a normal distribution $\bold{x} \backsim N(\boldsymbol\mu_\bold{x},\bold{M})$ . What this means is that $\bold{M=\Sigma_{x,x}^{-1}}$ , which isn’t surprising because the way it was defined, it’s basically $A^{-1}$ . And thus, using the inductive method I proposed, we can say that for any random variable formed by deleting some coordinates (more than one) of $\bold{z}$ , say $\bold{x}$ , has a normal distribution .

Now, moving towards the conditional:

it’s easy to see that the thing inside the exp function (except the -1/2) is a symmetric quadratic form in $\bold{x}$, and can thus be expressed as $\bold{(x-a)^TM(x-a)}$ for some symmetric matrix $\bold{M}$ . Thus $p(\bold{x|y})$ is also a normal distribution. To find the parameters of this distribution, you have to get your hands dirty and do a lot of linear algebra involving block matrices. I am going to skip all of that and just give the results to you :

Change of variables

For probability density functions $p_\bold{y}(\bold{y})$ and $p_\bold{x}(\bold{x})$ for variables $\bold{x,y}$ related as $\bold{y} = f(\bold{x})$ , we have $p_\bold{y}(\bold{y})\,dy_1\,dy_2\dots \,dy_n = p_\bold{x}(\bold{x})\,dx_1\,dx_2\,dx_3\dots \,dx_n$

We know that the determinant of the Jacobian of $\bold{y}$ wrt $\bold{x}$ gives us the factor by which the volume $dx_1\,dx_2\,dx_3\dots \,dx_n$ scales to become the volume $dy_1\,dy_2\dots \,dy_n$ ..Thus, we can write $p_\bold{y}(\bold{y})|\bold{J}_f(\bold{x})| = p_\bold{x}(\bold{x})$ .

Affine Transform of a Gaussian variable

Consider $\bold{x} \backsim N\bold{(\boldsymbol\mu_x,\boldsymbol\Sigma_x)}$ . Then, for a variable $\bold{y=Ax+b}$ , we have

• $\bold{J_y(x)=A}$

• $\text{E}\bold{(y)=\boldsymbol\mu_y=A\boldsymbol\mu_x+b \implies \boldsymbol\mu_x = A^{-1}(\boldsymbol\mu_y-b)}$

• $\text{Var}\bold{(y)=\Sigma_y = A\Sigma_xA^T \implies \Sigma_x = A^{-1}\Sigma_y (A^{-1})^T}$

and thus

Thus $\bold{y}$ also follows the Gaussian distribution,

namely $\bold{y}\backsim N(\bold{A\boldsymbol\mu_x+b,A\Sigma_xA^T})$ .

Standard Gaussian

Consider $\bold{y}\backsim N(\bold{\boldsymbol\mu,\Sigma})$ . This is the result of an affine transformation on $\bold{x}\backsim N(\bold{0,I})$ given by $\bold{y = Ax+\boldsymbol\mu}$ for some $\bold{A}$ such that $\bold{AA^T=\Sigma}$ . We know how to find a possible $\bold{A}$ using the Cholesky Decomposition .

Thus all Gaussian variables are just Affine Transformations of the Standard Gaussian $\bold{x}\backsim N(\bold{0,I})$ .