\section{Our Approach with Cubic Splines}

Our introduction to the topic of activation functions was the article on the fitness
of periodic vs. monotonic activation functions (\cite{sopena}).
When studying two classic benchmark problems, the XOR and Two-Quardants Problem,
using a sine, and in the Two-Quadrants Problem a heaviside activation function, we
figured that for an analytical approach, it would be good to have an activation 
function that can be designed as needed.
The most reasonable way of doing this is interpolating a function from a set of points. 
Furthermore, we aim to make the activation function itself evolvable, using NetGen. 
NetGen is using the SNNS for ANN training. Since the SNNS needs the first and second 
derivative, the function has to be contiuous.\\
The goal of cubic spline interpolation is to get an interpolation 
formula which is smooth in the first derivative, and continuous in the second 
derivative, both within an interval and at its boundaries.

\subsection{Mathematical Theory of Cubic Splines}

Given are $n$ distinct points with $P_i = (x_i,y_i)$, $i=1,\dots, n$. We assume that the points are ordered and unique with respect to the first coordinate. Thus these $n$ points define $n-1$ intervalls on the $x$-axes, denoted by $x \in [x_i,x_{i+1}]$. In each of these intervalls a function $f_i(x), i=1,\dots , n-1$ is defined. Since we assume that $f_i$ is a cubic spline we may write
\begin{eqnarray*}
f_i(x) & = & f_i(x_i) + a_i\times(x - x_i) + b_i\times(x - x_i)^2 + c_i\times(x - x_i)^3, \\
&&x_i \le x \le x_{i+1}, i=1,\dots ,n-1
\end{eqnarray*}
For $x_{i+1}$
\begin{eqnarray*}
f_i(x_{i+1}) & = & f_i(x_i) + a_i\times(x_{i+1} - x_i) + b_i\times(x_{i+1} - x_i)^2 + c_i\times(x_{i+1} - x_i)^3
\end{eqnarray*}
if $\Delta x_i = x_{i+1} - x_i, i=1,\dots ,n-1$ then
\begin{eqnarray*}
f_i(x_{i+1}) & = & f_i(x_i) + a_i\times \Delta x_i + b_i\times\Delta x_i^2 + c_i\times\Delta x_i^3
\end{eqnarray*}
With another substitution we achieve
\begin{eqnarray*}
D_i & = & \frac{f_i(x_{i+1}) - f_i(x_i)}{\Delta x_i} =  a_i + b_i\times\Delta x_i + c_i\times\Delta x_i^2
\end{eqnarray*}
Or first requirements states that the first derivatives have to be equal in $x_i, i=1,\dots,n-1$, which yields
\begin{eqnarray*}
f_i'(x_{i+1}) = f_{i+1}'(x_{i+1})
\end{eqnarray*}
which means
\begin{eqnarray*}
a_i + 2b_i\times\Delta x_i + 3c_i\times\Delta x_i^2 = a_{i+1}
\end{eqnarray*}
Substituting equation in the equation above gives
\begin{eqnarray*}
c_i & = & \frac{a_i + a_{i+1} - 2D_i}{\Delta x_i^2} 
\end{eqnarray*}
Our second requirements states that the second derivatives have to be equal in $x_i , i=1,\dots,n-1$, which yields
\begin{eqnarray*}
f_i''(x_{i+1}) = f_{i+1}''(x_{i+1})
\end{eqnarray*}
which means
\begin{eqnarray*}
b_i + 3c_i\times\Delta x_i = b_{i+1}
\end{eqnarray*}
Substituting equation in the equation above gives
\begin{eqnarray*}
\frac{D_i - a_i}{\Delta x_i} + 2c_i\times\Delta x_i & = & \frac{D_{i+1} - a_{i+1}}{\Delta x_{i+1}} - c_{i+1}\times \Delta x_{i+1}
\end{eqnarray*}
Combining equation and the equation above we get
\begin{eqnarray*}
a_i\times \Delta x_{i+1} + 2a_{i+1}(\Delta x_{i} + \Delta x_{i+1}) + a_{i+2}\times \Delta x_{i} & = & 3(D_{i}\times \Delta x_{i+1} + D_{i+1}\times \Delta x_{i}) \\&&i=1,\dots, n-2 
\end{eqnarray*}

\subsection{Implementation of the Cubic Spline Activation Function}
The cubic coefficients for each interval are precalculated and stored in a file.
The fastest way would be to load these coefficients when the training begins, but
this would have been a very drastic change in the SNNS' source code and would
have required deep knowledge of the entire code. \\
Since this was not possible, we decided to read the coefficients from the file 
when the value of the function is calculated. The calling unit is identified by
its name, which is a direct reference of the filename of the function parameter
file. \\
The function we added to the SNNS first calculates the sum of the weighted inputs,
and then simply steps through the file, looking for the correct interval. If the
unit's input is outside of the interval of the function, it returns the function
value at the respective boundary. If within the interval, the corresponding cubic 
coefficients are loaded and the function value, or the value of the requested 
derivative, is calculated and returned. \\
The appendix contains a detailed description of the implementation, along with the 
source code.


\subsection{Results}
We tested the cubic spline activation function with the XOR and Two-Quadrants 
Problem.

\subsubsection{XOR-Problem}
For a detailed description of this problem see section \ref{intr:XOR}.\\
The ANNs we constructed for the problem had two input units, one output and no
hidden layer. We will now compare three possible solutions (activation function
of the input neurons is always the identity function). The activation functions in
the output unit are:
\begin{enumerate}
	\item the SNNS' sine function
	\item interpolated sine
	\item customized cubic spline
\end{enumerate}
The performance in supervised learning was quite different. The SNNS' sine function
was not able to learn the correct setting of the weights. The interpolated sine 
function (interpolated through $(\pi,0), (\frac{\pi}{2},-1), (\frac{\pi}{2},1), ...$) was able to learn only with 
positive initial weights, but a standard squared error of 0 was reached after 
only 40 epochs. The best performance was reached with a customized function depicted in figure~\ref{figure:xor}.

\begin{figure}[h]
\centering
\includegraphics{xor_snns.ps}\\
\caption{\label{figure:xor}Customized XOR activation function}
\end{figure}

Using this function, a squared standard error of 0.001 was reached after only 15
epochs, the SSE for these two approaches are compared in figure~\ref{figure:xor:sine}

\begin{figure}[h]
\centering
\includegraphics{xor_sinus.ps}\\
\caption{\label{figure:xor:sine}
\small
Comparison of learing the XOR Problem with sinus and cubic spline.
Both used the same simple architecture network and
the same learning pattern, only $f1$ used an interpolated sine as activation function 
in the output neuron, and $f2$
used the customized cubic spline depicted in figure~\ref{figure:xor}.
}
\end{figure}

\newpage
\subsubsection{Two-Quadrants Problem}
\begin{figure}[!h]

\parbox{5cm}{
\setlength{\unitlength}{1cm}
\begin{picture}(5,5)
\put(1,1){\framebox(1.5,1.5){0}}
\put(1,2.5){\makebox(1.5,1.5){1}}
\put(2.5,2.5){\framebox(1.5,1.5){0}}
\put(2.5,1){\makebox(1.5,1.5){1}}
\put(1,0){\vector(0,1){5}}
\put(0.6, 4.7){y}
\put(4.7,0.6){x}
\put(0,1){\vector(1,0){5}}
\put(0.9, 4){\line(1,0){0.2}}
\put(4, 0.9){\line(0,1){0.2}}
\put(0.6, 3.85){1}
\put(3.945, 0.6){1}
\end{picture}
\caption{\label{figure:quad:layout}Target output of the Quadrants Problem.}
}
\hfill
\parbox{7cm}{
Any point within $([0,1],[0,1])$ should be assigned to either 0 or 1 by an ANN as shown 
in figure~\ref{figure:quad:layout}.\\
The activation functions in the hidden layer are interpolated heaviside functions, 
reducing the problem to the XOR-Problem, which is again solved by the customized 
cubic spline function in the output layer.


Even though the net works correctly with the interpolated heaviside function when the
weights are set to the correct values, learning is hardly ever possible (due to the 
mass of local minima). 
}
\end{figure}
We then experimented with a sine activation
function in the hidden layer. With this activation function, learning is possible, but
getting the correct weights is still a game of chance.


\begin{figure}[!h]
\centering
\xymatrix{
\bigcirc_{i_1} \ar[d]_<<<<<{1}& & \bigcirc_{i_2} \ar[d]^<<<<<{1} \ar@{..}[r] & input\ layer(2) \\
\bigcirc_{h_1} \ar[dr]_<<<<<{1}& & \bigcirc_{h_2} \ar[dl]^<<<<<{1} \ar@{..}[r] & hidden\ layer(2) \\
 & \bigcirc_{o} \ar@{..}[rr] &  & output\ layer(1)
}
\caption{\label{figure:quad:architecture}Architecture of the ANN for the Quadrants Problem.}
\end{figure}


We obtained the best results
with a modification of the interpolated heaviside function. We inclined the two 
surfaces to give the heaviside a ''monotonic touch``, which is better for the 
backpropagation algorithm.


There is still a high chance for the ANN to be trapped in a
local minimum, but the overall performance is better than with the two other functions.
The SNNS' standard way of initializing a net is not very good here, since there is only 
one interval for both the weights and biases.

\begin{figure}[h]
\centering
\includegraphics{quad_snns.ps}
\caption{\label{figure:quad1} Error surface of quadrants problem using the sine function in the hidden layer.}
\end{figure}

\newpage

\begin{figure}[h]
\centering
\includegraphics{quad2_snns.ps}
\caption{\label{figure:quad2} Error surface for quadrants problem using adapted heaviside functions in the hidden layer.}
\end{figure}

\begin{figure}[!bh]
\centering
\includegraphics{quad_projection_col.ps}
\caption{\label{figure:projection}Actual output of the Quadrants Problem using the same architecture as in 
figure~\ref{figure:quad:architecture} with adapted heaviside as activation fuction.}
\end{figure}

\newpage
