# Plot the probability density function of the gamma distribution

I want to draw gamma distribution with different shape and scale parameters. The following MWE displayed error. The error information is :

Here is MWE:

\documentclass[border=5mm]{standalone}
\usepackage{pgfplots}
\begin{document}

\pgfmathdeclarefunction{gammaPDF}{2}{
\pgfmathparse{1/(#2^#1*gamma(#1))*x^(#1-1)*exp(-x/#2)}
}

\begin{tikzpicture}
\begin{axis}[
width=4cm,
height=3cm,
axis x line=bottom,
axis y line=left
]
blue!50!white,
mark=none,
domain=0:20,
samples=100,
smooth] {\gammaPDF(9,0.5)};
\end{axis}
\end{tikzpicture}

-
remove the backslash from gammaPDF and define gamma –  percusse Jun 21 '13 at 19:38
@percusse: From the link:tex.stackexchange.com/questions/80345/…, where the gamma function is used but no predefined. So why I can't use it directly here. –  user22986 Jun 21 '13 at 20:07
@user22986: In the answer you linked to, the function is evaluated using the gnuplot backend. The gamma function isn't defined in the PGF math parser. –  Jake Jun 21 '13 at 20:21

Here's the gamma PDF using the gamma function approximation that Wolfram Alpha spits out:

\documentclass[border=5mm]{standalone}
\usepackage{pgfplots}
\begin{document}

\begin{tikzpicture}[
declare function={gamma(\z)=
(2.506628274631*sqrt(1/\z) + 0.20888568*(1/\z)^(1.5) + 0.00870357*(1/\z)^(2.5) - (174.2106599*(1/\z)^(3.5))/25920 - (715.6423511*(1/\z)^(4.5))/1244160)*exp((-ln(1/\z)-1)*\z);},
declare function={gammapdf(\x,\k,\theta) = \x^(\k-1)*exp(-\x/\theta) / (\theta^\k*gamma(\k));}
]

\begin{axis}[
axis lines=left,
enlargelimits=upper,
samples=50,
legend entries={$k=1\quad \theta=2$,$k=2\quad \theta=2$, $k=9\quad \theta=0.5$}
]