# How to pass an array to function in pgfmathdeclarefunction

I am trying to plot below function as 3d surface plot, where $\theta_1,\theta_2$ are the x and y axes respectively. The Z axis is supposed to be $L(\theta_1,\theta_2)$.

$$L(\theta_1,\theta_2) = \prod_{i=1}^{m} \dfrac{1}{\sqrt{2\pi\theta_2}}{\text{exp}}{\Big[ -\dfrac{ (x_i-\theta_1)^2 }{2\theta_2} \Big]} \\ = \Big( \dfrac{1}{\sqrt{2\pi\theta_2}} \Big)^{m}{\text{exp}}{\Big[ \dfrac{-\sum_{i=1}^{m}(x_i-\theta_1)^2}{2\theta_2} \Big]} \tag{1}$$


Just in case of above latex is not rendered this is my equation:

I managed to get to an extent, but unable to find how to compute the summation by passing a series of values for \m., and eventually plot proper 3d surface plot for the function above. Kindly help.

MWE:

\documentclass{article}
\usepackage{pgfplots}
\pgfplotsset{compat=1.15}
\pgfmathdeclarefunction{joint_normal}{3}{%
\pgfmathparse{ (1/(2*pi*#3))^(#1)*exp( -(#1-#2)^2/(2*#3^2)  )}%
}
\begin{document}
\begin{tikzpicture}
\begin{axis}[
grid=both,
restrict z to domain*=0:1,
zmin=0,
colormap/hot,
%point meta min=-0.2,
%point meta max=1,
view={20}{20}  %tune here to change viewing angle
]

\def\m{5}
\addplot3[surf,domain=-30:30,domain y=0:30, samples=25] { joint_normal(\m, x, y) };
\end{axis}
\end{tikzpicture}
\end{document}


Current output:

online editor to try out: here

Note: One could assume $X_i$ values varying from any range say 0 to 10, from a normal distribution N(5, 4). Any range would do, but for not spending too much time on thinking of range, I provide this one. I could later tinker it as needed.

• You can do summations with the tikzmath library by using recursions. Other than that, I am not aware of any other way of doing the sum in an elegant way in this framework. – marmot Oct 18 '18 at 18:10
• I am new to tikz, can you please show with an example for above problem? – Parthiban Rajendran Oct 18 '18 at 18:11
• Well, to be able to do that, one would need to know what the x_i are. You sum over the x_i but as long as you do not specify what they are it is impossible to plot the function. – marmot Oct 18 '18 at 18:20
• One could assume $X_i$ values varying from any range say 0 to 10, from a normal distribution N(5, 4). Any range would do, but for not spending too much time on thinking of range, I provide this one. I could later tinker it as needed. – Parthiban Rajendran Oct 18 '18 at 18:25
• If you can make a concrete proposal for the x_i I'll be happy to give it a shot. (Sorry, was offline for a few hours and will be really online in another few hours) – marmot Oct 18 '18 at 23:50

Here is an MWE, which at the same time contains the explanation.

\documentclass{article}
\usepackage{amsmath}
\usepackage{pgfplots}
\pgfplotsset{compat=1.16}
\usetikzlibrary{math}
\begin{document}
% based on https://tex.stackexchange.com/a/307032/121799
% and https://tex.stackexchange.com/a/451326/121799
\def\xvalues{{0,1,2,4,5,7}} % notice that the 0 is the 0th entry, which is not used here
\tikzset{evaluate={
function myN(\x,\z,\k) { % \x = \theta_1 and \z=\theta_2
if \k == 1 then {
return myn(\x,\xvalues[1],\z);
} else {
return myN(\x,\z,\k-1)
+myn(\x,\xvalues[\k],\z);
};
};
},
declare function={myn(\x,\y,\z)=(-(\x-\y)*(\x-\y))/(2*\z*\z) ;
L(\x,\z,\k)=pow(2*pi*\z,-\k/2)*exp(myN(\x,\z,\k));}}

\section*{How to plot sums in Ti\emph{k}Z/pgfplots}

We define the argument of the exponential as
$$n_k(\theta_1,\theta_2)~=~-\frac{(\theta_1-x_k)^2}{2\theta_2^2}$$
and their sum as
$$N_k(\theta_1,\theta_2)~=~\sum\limits_{\ell=1}^k n_k(\theta_1,\theta_2)\;.$$
This means that $N_k$ can be defined recursively as
$$N_k(\theta_1,\theta_2)~=~N_{k-1}(\theta_1,\theta_2)+n_k(\theta_1,\theta_2)\;,$$
and this is the point where the Ti\emph{k}Z library \texttt{math} comes into
play. It allows us to do the recursive deinition. Examples are shown in
Figure~\ref{fig:N_k}.

\begin{figure}[htb]
\centering
\begin{tikzpicture}
\begin{axis}[samples=101,
use fpu=false,mark=none,
xlabel=$x$,ylabel=$y$,
xmin=0, xmax=10,
domain=0:10,legend pos=south west
]
\addlegendentry{$N_1$}
\addlegendentry{$N_2$}
\addlegendentry{$N_3$}
\end{axis}
\end{tikzpicture}
\caption{$N_1$, $N_2$ and $N_3$ for $\theta_2=1$ and $\{x_k\}=\{1,2,4\}$.}
\label{fig:N_k}
\end{figure}

\clearpage

Of course, one can then define functions of these sums,
$$L_k(\theta_1,\theta_2)~=~ \Big( \dfrac{1}{\sqrt{2\pi\theta_2}} \Big)^{m}\,\exp\Bigl[ \dfrac{-\sum_{i=1}^{k}(x_i-\theta_1)^2}{2\theta_2} \Bigr]\;.$$
Examples are shown in Figure~\ref{fig:L_k}.

\begin{figure}[htb]
\centering
\begin{tikzpicture}
\begin{axis}[samples=101,
use fpu=false,mark=none,
xlabel=$x$,ylabel=$y$,
xmin=0, xmax=10,
domain=0:10,legend pos=north east
]
\addlegendentry{$L_1$}
\addlegendentry{$L_2$}
\addlegendentry{$L_3$}
\end{axis}
\end{tikzpicture}
\caption{$L_1$, $L_2$ and $L_3$ for $\theta_2=1$ and $\{x_k\}=\{1,2,4\}$.}
\label{fig:L_k}
\end{figure}

\end{document}


The second page contains (hopefully) what you are seeking for.

I'd also like to urge you not to confuse TikZ/pgfplots with a computer algegra system. You can do these things, but should not be too surprised if the performance is below the one of, say, Mathematica.

And here is a 3D example, similar to what you do in your MWE.

\documentclass[tikz,border=3.14mm]{standalone}
\usepackage{pgfplots}
\pgfplotsset{compat=1.16}
\usetikzlibrary{math}
\begin{document}
% based on https://tex.stackexchange.com/a/307032/121799
% and https://tex.stackexchange.com/a/451326/121799
\def\xvalues{{0,1,2,4,5,7}} % notice that the 0 is the 0th entry, which is not used here
\tikzset{evaluate={
function myN(\x,\z,\k) { % \x = \theta_1 and \z=\theta_2
if \k == 1 then {
return myn(\x,\xvalues[1],\z);
} else {
return myN(\x,\z,\k-1)
+myn(\x,\xvalues[\k],\z);
};
};
},
declare function={myn(\x,\y,\z)=(-(\x-\y)*(\x-\y))/(2*\z*\z) ;
L(\x,\z,\k)=pow(2*pi*\z,-\k/2)*exp(myN(\x,\z,\k));}}

\begin{tikzpicture}
\begin{axis}[use fpu=false,
grid=both,
restrict z to domain*=0:1,
zmin=0,
colormap/hot,
%point meta min=-0.2,
%point meta max=1,
view={20}{20}  %tune here to change viewing angle
]

\addplot3[surf,domain=-1:9,domain y=1:4, samples=25] { L(x, y,4) };
\end{axis}
\end{tikzpicture}
\end{document}


• @manooooh Good catch, I have just copied that from the OP. – marmot Oct 19 '18 at 3:54
• wow, wow.. I am so grateful for you for this MWE with explanation. yes, as you said beyond a point, performance bothers me. So I am also trying via python in parallel. – Parthiban Rajendran Oct 19 '18 at 18:02