41

For a small class I want to draw a bivariate normal distribution and show where the means from the two variables meet in the space. Has anyone done this in TikZ?

2 Answers2

86

Christian was faster, but because I invested the time, here's my take on showing how the bivariate distribution results from the two univariate ones.

There are a couple of things in the code that might be useful for you:

You can define mathematical functions using declare function={<name>(<argument macros>)=<function>;}, which will help to keep the code clean and avoid repetitions.

You can define a new colormap using \pgfplotsset{ colormap={<name>}{<color model>(<distance>)=(<value1>); <color model>(<distance 2>)=(<value2>)} }. This is a very powerful feature, so you should definitely read up on it in the pgfplots manual.

The legend created by colorbar is a whole new plot, so you can configure it with all the usual axis options.

There are different ways for defining 3D functions: \addplot3 {<function>}; will evaluate <function> at every point on a grid and assume the result to be a z-value. \addplot3 ({<x>},{<y>},{<z>}); defines a parametric function in 3D space, which allows you to (among other things) draw three-dimensional lines.

\documentclass{standalone}

\usepackage{pgfplots}

\begin{document}

\pgfplotsset{
colormap={whitered}{color(0cm)=(white); color(1cm)=(orange!75!red)}
}

\begin{tikzpicture}[
    declare function={mu1=1;},
    declare function={mu2=2;},
    declare function={sigma1=0.5;},
    declare function={sigma2=1;},
    declare function={normal(\m,\s)=1/(2*\s*sqrt(pi))*exp(-(x-\m)^2/(2*\s^2));},
    declare function={bivar(\ma,\sa,\mb,\sb)=
        1/(2*pi*\sa*\sb) * exp(-((x-\ma)^2/\sa^2 + (y-\mb)^2/\sb^2))/2;}]
\begin{axis}[
    colormap name=whitered,
    width=15cm,
    view={45}{65},
    enlargelimits=false,
    grid=major,
    domain=-1:4,
    y domain=-1:4,
    samples=26,
    xlabel=$x_1$,
    ylabel=$x_2$,
    zlabel={$P$},
    colorbar,
    colorbar style={
        at={(1,0)},
        anchor=south west,
        height=0.25*\pgfkeysvalueof{/pgfplots/parent axis height},
        title={$P(x_1,x_2)$}
    }
]
\addplot3 [surf] {bivar(mu1,sigma1,mu2,sigma2)};
\addplot3 [domain=-1:4,samples=31, samples y=0, thick, smooth] (x,4,{normal(mu1,sigma1)});
\addplot3 [domain=-1:4,samples=31, samples y=0, thick, smooth] (-1,x,{normal(mu2,sigma2)});

\draw [black!50] (axis cs:-1,0,0) -- (axis cs:4,0,0);
\draw [black!50] (axis cs:0,-1,0) -- (axis cs:0,4,0);

\node at (axis cs:-1,1,0.18) [pin=165:$P(x_1)$] {};
\node at (axis cs:1.5,4,0.32) [pin=-15:$P(x_2)$] {};
\end{axis}
\end{tikzpicture}
\end{document}
Christian
  • 19,238
Jake
  • 232,450
  • 1
    Looks pretty cool :) – Christian Feuersänger Oct 16 '11 at 11:41
  • 5
    It's a pretty cool package! – Jake Oct 16 '11 at 11:42
  • 5
    Note that due to the symmetry of the gaussian, you may want to choose an odd number of samples for the 1d curves to avoid the flat top. – Christian Feuersänger Oct 16 '11 at 11:52
  • Ah yes, good point. I've changed it now. – Jake Oct 16 '11 at 21:07
  • 2
    It never stops to impress me. – Harold Cavendish Oct 16 '11 at 22:24
  • 1
    Wauv - very beautiful!! – Thomas Jensen Oct 17 '11 at 08:44
  • @Jake Awesome answer and awesome package; Thank you both – Alain Matthes Mar 23 '12 at 17:30
  • 2
    Just wanted to let everybody who might already have used this code know that there was a tiny but significant typo in the code. 1/(2*pi*\sa*\sa) should have been 1/(2*pi*\sa*\sb). It's fixed now but you might want to check any documents in which you used this beautiful plot. – Christian May 16 '13 at 11:48
  • 1
    I will use this in my Prob and Stats class when we cover joint distributions. – R. Schumacher Jul 10 '13 at 03:19
  • Seems to me there are still some typos here. The normal pdf should be normal(\m,\s)=1/(\s*sqrt(2*pi))*exp(-(x-\m)^2/(2*\s^2)); and the bivariate normal bivar(\ma,\sa,\mb,\sb)= 1/(2*pi*\sa*\sb) * exp( -((x-\ma)^2/\sa^2 + (y-\mb)^2/\sb^2)/2 ); (note that this is true for independent random variables only). Also, the univariate normal should be multiplied by a factor to be exactly equal to the projection of the bivariate to the x- or y- axis. That factor is 1/(sigma2*sqrt(2*pi)) for P(x_1) and 1/(sigma1*sqrt(2*pi)) for P(x_2). But still: awesome picture! ^^. – Bibi Jul 07 '16 at 10:29
  • Was this resolved? – Alex Lenail Oct 24 '17 at 00:16
  • Any chance you could extend this to the general case of nonzero correlation? https://en.wikipedia.org/wiki/Multivariate_normal_distribution#Bivariate_case I tried, but I am running into incomprehensible error messages :/ – jmb Nov 21 '18 at 18:16
46

What you need is to sample a bivariate function and to choose a proper plot type. Typical choices are mesh, surface, or contour plots. The most common is a surface plot, I think.

In any case, TikZ does not support bivariate functions - but you can use pgfplots which is build on top of it and provides simple integration into tikz. I should note that I am author of pgfplots. Here is how it could help here:

enter image description here

\documentclass[a4paper]{article}
\usepackage{pgfplots}

\begin{document}

\def\centerx{2}
\def\centery{-1}
\begin{tikzpicture}
    \begin{axis}
    \addplot3[surf,domain=-2:6,domain y=-5:3] 
        {exp(-( (x-\centerx)^2 + (y-\centery)^2)/3 )};
    \node[circle,inner sep=1pt,fill=blue,pin=90:$\mu$] 
        at (axis cs:\centerx,\centery,1) {};
    \end{axis}
\end{tikzpicture}
\end{document}

You will see that the example is NO real normal distribution because it is not normalized to unit integral (it is normalized to unit max norm).

But the example shows the key aspects:

  1. choosing surf
  2. specifying the plotting expression
  3. adjust domain and domain y (consider also samples and samples y)
  4. a \node which describes the center in "some way".

Pgfplots also supports mesh plots (use mesh instead of surf) or contour plots (using contour gnuplot instead of surf; in this case, you need gnuplot and the shell-escape feature like pdflatex -shell-escape):

enter image description here

\def\centerx{2}
\def\centery{-1}
\begin{tikzpicture}
    \begin{axis}[view={0}{90},axis equal]
    \addplot3[contour gnuplot,domain=-2:6,domain y=-5:3] 
        {exp(-( (x-\centerx)^2 + (y-\centery)^2)/3 )};
    \node[circle,inner sep=1pt,fill=blue,pin=90:$\mu$] 
        at (axis cs:\centerx,\centery,1) {};
    \end{axis}
\end{tikzpicture}

This second example also shows that view can be used to change the view angles (horizontal and vertical angle) and that axis equal might be of interest to get the axis ratios comparable (I must admit to my shame that the combination view={0}{90},axis equal appears to produce unexpected compilation errors - I will fix that).

Of interest may also be shader=interp or various choices for colormap which are listed in the manual.