Probability density function

Absolutely continuous univariate distributions

A probability density function is most commonly associated with absolutely continuous univariate distributions. A random variable ${\displaystyle X}$ has density ${\displaystyle f_X}$, where ${\displaystyle f_X}$ is a non-negative Lebesgue-integrable function, if: $${\displaystyle Pr[a\le X\le b]={\int }_a^b{f}_X(x)dx.}$$

Hence, if ${\displaystyle F_X}$ is the cumulative distribution function of ${\displaystyle X}$, then: $${\displaystyle F_X(x)={\int }_{-\mathrm{\infty }}^x{f}_X(u)du,}$$ and (if ${\displaystyle F_X}$ is differentiable at ${\displaystyle x}$) $${\displaystyle f_X(x)=\frac{d}{dx}{F}_X(x).}$$

Intuitively, one can think of ${\displaystyle f_X(x)dx}$ as being the probability of ${\displaystyle X}$ falling within the infinitesimal interval ${\displaystyle [x,x+dx]}$.

Formal definition

(This definition may be extended to any probability distribution using the measure-theoretic definition of probability.)

A random variable ${\displaystyle X}$ with values in a measurable space ${\displaystyle (\mathcal{X},\mathcal{A})}$ (usually ${\displaystyle {\mathbb{R}}^n}$ with the Borel sets as measurable subsets) has as probability distribution the pushforward measure X_∗P on ${\displaystyle (\mathcal{X},\mathcal{A})}$: the density of ${\displaystyle X}$ with respect to a reference measure ${\displaystyle \mu }$ on ${\displaystyle (\mathcal{X},\mathcal{A})}$ is the Radon–Nikodym derivative: $${\displaystyle f=\frac{d{X}_{\ast }P}{d\mu }.}$$

That is, f is any measurable function with the property that: $${\displaystyle Pr[X\in A]={\int }_{X^{-1}A}dP={\int }_{A}fd\mu }$$ for any measurable set ${\displaystyle A\in \mathcal{A}.}$

Discussion

In the continuous univariate case above, the reference measure is the Lebesgue measure. The probability mass function of a discrete random variable is the density with respect to the counting measure over the sample space (usually the set of integers, or some subset thereof).

It is not possible to define a density with reference to an arbitrary measure (e.g. one can not choose the counting measure as a reference for a continuous random variable). Furthermore, when it does exist, the density is almost unique, meaning that any two such densities coincide almost everywhere.

Further details

Unlike a probability, a probability density function can take on values greater than one; for example, the continuous uniform distribution on the interval [0, 1/2] has probability density f(x) = 2 for 0 ≤ x ≤ 1/2 and f(x) = 0 elsewhere.

The standard normal distribution has probability density$${\displaystyle f(x)=\frac{1}{\sqrt{2\pi }}{e}^{-x^2/2}.}$$

If a random variable X is given and its distribution admits a probability density function f, then the expected value of X (if the expected value exists) can be calculated as$${\displaystyle \mathrm{E}[X]={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}xf(x)dx.}$$

Not every probability distribution has a density function: the distributions of discrete random variables do not; nor does the Cantor distribution, even though it has no discrete component, i.e., does not assign positive probability to any individual point.

A distribution has a density function if its cumulative distribution function F(x) is absolutely continuous. In this case: F is almost everywhere differentiable, and its derivative can be used as probability density:$${\displaystyle \frac{d}{dx}F(x)=f(x).}$$

If a probability distribution admits a density, then the probability of every one-point set {a} is zero; the same holds for finite and countable sets.

Two probability densities f and g represent the same probability distribution precisely if they differ only on a set of Lebesgue measure zero.

In the field of statistical physics, a non-formal reformulation of the relation above between the derivative of the cumulative distribution function and the probability density function is generally used as the definition of the probability density function. This alternate definition is the following:

If dt is an infinitely small number, the probability that X is included within the interval (t, t + dt) is equal to f(t) dt, or: $${\displaystyle Pr(t<X<t+dt)=f(t)dt.}$$

Link between discrete and continuous distributions

It is possible to represent certain discrete random variables as well as random variables involving both a continuous and a discrete part with a generalized probability density function using the Dirac delta function. (This is not possible with a probability density function in the sense defined above, it may be done with a distribution.) For example, consider a binary discrete random variable having the Rademacher distribution—that is, taking −1 or 1 for values, with probability 1⁄2 each. The density of probability associated with this variable is: $${\displaystyle f(t)=\frac{1}{2}(\delta (t+1)+\delta (t-1)).}$$

More generally, if a discrete variable can take n different values among real numbers, then the associated probability density function is: $${\displaystyle f(t)=\sum _{i=1}^n{p}_i\delta (t-x_i),}$$ where ${\displaystyle x_1,\dots ,x_n}$ are the discrete values accessible to the variable and ${\displaystyle p_1,\dots ,p_n}$ are the probabilities associated with these values.

This substantially unifies the treatment of discrete and continuous probability distributions. The above expression allows for determining statistical characteristics of such a discrete variable (such as the mean, variance, and kurtosis), starting from the formulas given for a continuous distribution of the probability.

Families of densities

It is common for probability density functions (and probability mass functions) to be parametrized—that is, to be characterized by unspecified parameters. For example, the normal distribution is parametrized in terms of the mean and the variance, denoted by ${\displaystyle \mu }$ and ${\displaystyle {\sigma }^2}$ respectively, giving the family of densities $${\displaystyle f(x;\mu ,{\sigma }^2)=\frac{1}{\sigma \sqrt{2\pi }}{e}^{-\frac{1}{2}{\left(\frac{x-\mu }{\sigma }\right)}^2}.}$$ Different values of the parameters describe different distributions of different random variables on the same sample space (the same set of all possible values of the variable); this sample space is the domain of the family of random variables that this family of distributions describes. A given set of parameters describes a single distribution within the family sharing the functional form of the density. From the perspective of a given distribution, the parameters are constants, and terms in a density function that contain only parameters, but not variables, are part of the normalization factor of a distribution (the multiplicative factor that ensures that the area under the density—the probability of something in the domain occurring— equals 1). This normalization factor is outside the kernel of the distribution.

Since the parameters are constants, reparametrizing a density in terms of different parameters to give a characterization of a different random variable in the family, means simply substituting the new parameter values into the formula in place of the old ones.

Densities associated with multiple variables

For continuous random variables X₁, ..., X_n, it is also possible to define a probability density function associated to the set as a whole, often called joint probability density function. This density function is defined as a function of the n variables, such that, for any domain D in the n-dimensional space of the values of the variables X₁, ..., X_n, the probability that a realisation of the set variables falls inside the domain D is $${\displaystyle Pr\left(X_1,\dots ,X_n\in D\right)={\int }_D{f}_{X_1,\dots ,X_n}(x_1,\dots ,x_n)d{x}_1\cdots d{x}_n.}$$

If F(x₁, ..., x_n) = Pr(X₁ ≤ x₁, ..., X_n ≤ x_n) is the cumulative distribution function of the vector (X₁, ..., X_n), then the joint probability density function can be computed as a partial derivative $${\displaystyle f(x)={\frac{{\mathrm{\partial }}^{n}F}{\mathrm{\partial }{x}_1\cdots \mathrm{\partial }{x}_n}|}_x}$$

Marginal densities

For i = 1, 2, ..., n, let f_Xi(x_i) be the probability density function associated with variable X_i alone. This is called the marginal density function, and can be deduced from the probability density associated with the random variables X₁, ..., X_n by integrating over all values of the other n − 1 variables: $${\displaystyle f_{X_i}(x_i)=\int f(x_1,\dots ,x_n)d{x}_1\cdots d{x}_{i-1}d{x}_{i+1}\cdots d{x}_n.}$$

Independence

Continuous random variables X₁, ..., X_n admitting a joint density are all independent from each other if $${\displaystyle f_{X_1,\dots ,X_n}(x_1,\dots ,x_n)=f_{X_1}(x_1)\cdots f_{X_n}(x_n).}$$

Corollary

If the joint probability density function of a vector of n random variables can be factored into a product of n functions of one variable $${\displaystyle f_{X_1,\dots ,X_n}(x_1,\dots ,x_n)=f_1(x_1)\cdots f_n(x_n),}$$ (where each f_i is not necessarily a density) then the n variables in the set are all independent from each other, and the marginal probability density function of each of them is given by $${\displaystyle f_{X_i}(x_i)=\frac{f_i(x_i)}{\int f_i(x)dx}.}$$

Example

This elementary example illustrates the above definition of multidimensional probability density functions in the simple case of a function of a set of two variables. Let us call ${\displaystyle \overrightarrow{R}}$ a 2-dimensional random vector of coordinates (X, Y): the probability to obtain ${\displaystyle \overrightarrow{R}}$ in the quarter plane of positive x and y is $${\displaystyle Pr\left(X>0,Y>0\right)={\int }_0^{\mathrm{\infty }}{\int }_0^{\mathrm{\infty }}{f}_{X,Y}(x,y)dxdy.}$$

Function of random variables and change of variables in the probability density function

If the probability density function of a random variable (or vector) X is given as f_X(x), it is possible (but often not necessary; see below) to calculate the probability density function of some variable Y = g(X). This is also called a "change of variable" and is in practice used to generate a random variable of arbitrary shape f_g(X) = f_Y using a known (for instance, uniform) random number generator.

It is tempting to think that in order to find the expected value E(g(X)), one must first find the probability density f_g(X) of the new random variable Y = g(X). However, rather than computing $${\displaystyle \mathrm{E}(g(X))={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}y{f}_{g(X)}(y)dy,}$$ one may find instead $${\displaystyle \mathrm{E}(g(X))={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}g(x)f_X(x)dx.}$$

The values of the two integrals are the same in all cases in which both X and g(X) actually have probability density functions. It is not necessary that g be a one-to-one function. In some cases the latter integral is computed much more easily than the former. See Law of the unconscious statistician.

Scalar to scalar

Let ${\displaystyle g:\mathbb{R}\to \mathbb{R}}$ be a monotonic function, then the resulting density function is $${\displaystyle f_Y(y)=f_X(g^{-1}(y))\left|\frac{d}{dy}(g^{-1}(y))\right|.}$$

Here g⁻¹ denotes the inverse function.

This follows from the fact that the probability contained in a differential area must be invariant under change of variables. That is, $${\displaystyle \left|f_Y(y)dy\right|=\left|f_X(x)dx\right|,}$$ or $${\displaystyle f_Y(y)=\left|\frac{dx}{dy}\right|f_X(x)=\left|\frac{d}{dy}(x)\right|f_X(x)=\left|\frac{d}{dy}(g^{-1}(y))\right|f_X(g^{-1}(y))=\left|{\left(g^{-1}\right)}^{\prime }(y)\right|\cdot f_X(g^{-1}(y)).}$$

For functions that are not monotonic, the probability density function for y is $${\displaystyle \sum _{k=1}^{n(y)}\left|\frac{d}{dy}{g}_k^{-1}(y)\right|\cdot f_X(g_k^{-1}(y)),}$$ where n(y) is the number of solutions in x for the equation ${\displaystyle g(x)=y}$, and ${\displaystyle g_k^{-1}(y)}$ are these solutions.

Vector to vector

Suppose x is an n-dimensional random variable with joint density f. If y = G(x), where G is a bijective, differentiable function, then y has density p_Y: $${\displaystyle p_Y(\mathbf{y})=f(G^{-1}(\mathbf{y}))\left|det\left[{\frac{d{G}^{-1}(\mathbf{z})}{d\mathbf{z}}|}_{\mathbf{z}=\mathbf{y}}\right]\right|}$$ with the differential regarded as the Jacobian of the inverse of G(⋅), evaluated at y.

For example, in the 2-dimensional case x = (x₁, x₂), suppose the transform G is given as y₁ = G₁(x₁, x₂), y₂ = G₂(x₁, x₂) with inverses x₁ = G₁⁻¹(y₁, y₂), x₂ = G₂⁻¹(y₁, y₂). The joint distribution for y = (y₁, y₂) has density $${\displaystyle p_{Y_1,Y_2}(y_1,y_2)=f_{X_1,X_2}(G_1^{-1}(y_1,y_2),G_2^{-1}(y_1,y_2))\left|\frac{\mathrm{\partial }{G}_1^{-1}}{\mathrm{\partial }{y}_1}\frac{\mathrm{\partial }{G}_2^{-1}}{\mathrm{\partial }{y}_2}-\frac{\mathrm{\partial }{G}_1^{-1}}{\mathrm{\partial }{y}_2}\frac{\mathrm{\partial }{G}_2^{-1}}{\mathrm{\partial }{y}_1}\right|.}$$

Vector to scalar

Let ${\displaystyle V:{\mathbb{R}}^n\to \mathbb{R}}$ be a differentiable function and ${\displaystyle X}$ be a random vector taking values in ${\displaystyle {\mathbb{R}}^n}$, ${\displaystyle f_X}$ be the probability density function of ${\displaystyle X}$ and ${\displaystyle \delta (\cdot )}$ be the Dirac delta function. It is possible to use the formulas above to determine ${\displaystyle f_Y}$, the probability density function of ${\displaystyle Y=V(X)}$, which will be given by $${\displaystyle f_Y(y)={\int }_{{\mathbb{R}}^n}{f}_X(\mathbf{x})\delta (y-V(\mathbf{x}))d\mathbf{x}.}$$

This result leads to the law of the unconscious statistician: $${\displaystyle \begin{array}{cc}{\mathrm{E}}_Y[Y]& ={\int }_{\mathbb{R}}y{f}_Y(y)dy\\ & ={\int }_{\mathbb{R}}y{\int }_{{\mathbb{R}}^n}{f}_X(\mathbf{x})\delta (y-V(\mathbf{x}))d\mathbf{x}dy\\ & ={\int }_{{\mathbb{R}}^n}{\int }_{\mathbb{R}}y{f}_X(\mathbf{x})\delta (y-V(\mathbf{x}))dyd\mathbf{x}\\ & ={\int }_{{\mathbb{R}}^n}V(\mathbf{x})f_X(\mathbf{x})d\mathbf{x}={\mathrm{E}}_X[V(X)].\end{array}}$$

Proof:

Let ${\displaystyle Z}$ be a collapsed random variable with probability density function ${\displaystyle p_Z(z)=\delta (z)}$ (i.e., a constant equal to zero). Let the random vector ${\displaystyle \tilde{X}}$ and the transform ${\displaystyle H}$ be defined as $${\displaystyle H(Z,X)=\left[\begin{array}{c}Z+V(X)\\ X\end{array}\right]=\left[\begin{array}{c}Y\\ \tilde{X}\end{array}\right].}$$

It is clear that ${\displaystyle H}$ is a bijective mapping, and the Jacobian of ${\displaystyle H^{-1}}$ is given by: $${\displaystyle \frac{d{H}^{-1}(y,\tilde{\mathbf{x}})}{dyd\tilde{\mathbf{x}}}=\left[\begin{array}{cc}1& -\frac{dV(\tilde{\mathbf{x}})}{d\tilde{\mathbf{x}}}\\ {\mathbf{0}}_{n\times 1}& {\mathbf{I}}_{n\times n}\end{array}\right],}$$ which is an upper triangular matrix with ones on the main diagonal, therefore its determinant is 1. Applying the change of variable theorem from the previous section we obtain that $${\displaystyle f_{Y,X}(y,x)=f_X(\mathbf{x})\delta (y-V(\mathbf{x})),}$$ which if marginalized over ${\displaystyle x}$ leads to the desired probability density function.

Sums of independent random variables

The probability density function of the sum of two independent random variables U and V, each of which has a probability density function, is the convolution of their separate density functions: $${\displaystyle f_{U+V}(x)={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{f}_U(y)f_V(x-y)dy=\left(f_U\ast f_V\right)(x)}$$

It is possible to generalize the previous relation to a sum of N independent random variables, with densities U₁, ..., U_N: $${\displaystyle f_{U_1+\cdots +U}(x)=\left(f_{U_1}\ast \cdots \ast f_{U_N}\right)(x)}$$

This can be derived from a two-way change of variables involving Y = U + V and Z = V, similarly to the example below for the quotient of independent random variables.

Products and quotients of independent random variables

Given two independent random variables U and V, each of which has a probability density function, the density of the product Y = UV and quotient Y = U/V can be computed by a change of variables.

Example: Quotient distribution

To compute the quotient Y = U/V of two independent random variables U and V, define the following transformation: $${\displaystyle \begin{array}{cc}Y& =U/V\\ Z& =V\end{array}}$$

Then, the joint density p(y,z) can be computed by a change of variables from U,V to Y,Z, and Y can be derived by marginalizing out Z from the joint density.

The inverse transformation is $${\displaystyle \begin{array}{cc}U& =YZ\\ V& =Z\end{array}}$$

The absolute value of the Jacobian matrix determinant ${\displaystyle J(U,V\mid Y,Z)}$ of this transformation is: $${\displaystyle \left|det\left[\begin{array}{cc}\frac{\mathrm{\partial }u}{\mathrm{\partial }y}& \frac{\mathrm{\partial }u}{\mathrm{\partial }z}\\ \frac{\mathrm{\partial }v}{\mathrm{\partial }y}& \frac{\mathrm{\partial }v}{\mathrm{\partial }z}\end{array}\right]\right|=\left|det\left[\begin{array}{cc}z& y\\ 0& 1\end{array}\right]\right|=|z|.}$$

Thus: $${\displaystyle p(y,z)=p(u,v)J(u,v\mid y,z)=p(u)p(v)J(u,v\mid y,z)=p_U(yz)p_V(z)|z|.}$$

And the distribution of Y can be computed by marginalizing out Z: $${\displaystyle p(y)={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{p}_U(yz)p_V(z)|z|dz}$$

This method crucially requires that the transformation from U,V to Y,Z be bijective. The above transformation meets this because Z can be mapped directly back to V, and for a given V the quotient U/V is monotonic. This is similarly the case for the sum U + V, difference U − V and product UV.

Exactly the same method can be used to compute the distribution of other functions of multiple independent random variables.

Example: Quotient of two standard normals

Given two standard normal variables U and V, the quotient can be computed as follows. First, the variables have the following density functions: $${\displaystyle \begin{array}{cc}p(u)& =\frac{1}{\sqrt{2\pi }}{e}^{-u^2/2}\\ p(v)& =\frac{1}{\sqrt{2\pi }}{e}^{-v^2/2}\end{array}}$$

We transform as described above: $${\displaystyle \begin{array}{cc}Y& =U/V\\ Z& =V\end{array}}$$

This leads to: $${\displaystyle \begin{array}{cc}p(y)& ={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{p}_U(yz)p_V(z)|z|dz\\ & ={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{1}{\sqrt{2\pi }}{e}^{-\frac{1}{2}{y}^2{z}^2}\frac{1}{\sqrt{2\pi }}{e}^{-\frac{1}{2}{z}^2}|z|dz\\ & ={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{1}{2\pi }{e}^{-\frac{1}{2}\left(y^2+1\right)z^2}|z|dz\\ & =2{\int }_0^{\mathrm{\infty }}\frac{1}{2\pi }{e}^{-\frac{1}{2}\left(y^2+1\right)z^2}zdz\\ & ={\int }_0^{\mathrm{\infty }}\frac{1}{\pi }{e}^{-\left(y^2+1\right)u}du& & u=\frac{1}{2}{z}^2\\ & ={-\frac{1}{\pi \left(y^2+1\right)}{e}^{-\left(y^2+1\right)u}|}_{u=0}^{\mathrm{\infty }}\\ & =\frac{1}{\pi \left(y^2+1\right)}\end{array}}$$

This is the density of a standard Cauchy distribution.