\documentclass[a4paper,12pt]{article}
\usepackage{geometry}
\geometry{left=25mm,top=20mm,bottom=25mm,right=20mm}
%\usepackage[T1]{fontenc}
%\usepackage[latin1]{inputenc}
\usepackage{graphics}
\usepackage{epsfig}
%\renewcommand{\familydefault}{\sfdefault}
%\usepackage{helvet}
%\usepackage{sfmath}
%\usepackage[small]{caption2}
%\usepackage{eurosym}
%\usepackage{multicol}




\begin{document}
\setlength{\parindent}{0pt}
\setlength{\parskip}{6pt}

\title{Local Model Network Toolbox \\ for Nonlinear System Identification}
\author{Oliver Nelles, \\ Benjamin Hartmann, Tobias Ebert, Torsten Fischer, \\ Julian Belz, Geritt Kampmann}
\date{May 2012}
\maketitle

\section*{Abstract}
A new \textsc{Matlab} toolbox for nonlinear system identification with local model networks is introduced. The toolbox will be available as freeware (no commercial usage) and can be downloaded from \mbox{\tt www.uni-siegen.de/fb11/mrt} soon. Its goal is to provide even unskilled users with an easy possibility to generate nonlinear models from data. It trains a number of different models of different complexity with different architectures and a polynomial for comparison. At the end it gives an overview on the achieved performances and makes a recommendation for the best model. 

The toolbox can be used by the trivial function call {\tt LMNTrain(data)} where the training data {\tt data = [u\_1~u\_2~...~u\_p~y]} contains the inputs in the first columns and the output in the last column. The model with the best penalty loss function value (AIC$_{\rm c}$) is recommended \cite{akaike1973information, burnham2002model}.

\section{Introduction}
This toolbox limits itself to nonlinear models with $p$ inputs and $1$ output. For several outputs the user has to generate separate models for each. In the static case, the generated models are of the type
\begin{equation}
\hat{y} = f(u_1, u_2, \ldots, u_p)
\end{equation}
Particularly for high-dimensional input spaces, i.e., high values of $p$, this can be a very demanding task. 

Dynamic models are not considered so far. They will be included in future toolbox releases.
%Dynamic models can be constructed by using tapped delay lines, i.e.:
%\begin{eqnarray}
%\hat{y}(k+1) = f(u_1(k), \ldots, u_1(k-n_{u1}), u_2(k), \ldots, u_2(k-n_{u2}), \ldots,\\
%  u_p(k), \ldots, u_p(k-n_{up}), y(k), \ldots, y(k-n_{ny}))
%\end{eqnarray}
%However, such models performs only one-step-ahead predictions from $y(k)$ to $\hat{y}(k+1)$. For a simulation feedback is necessary. 

\section{Toolbox Methods}
The toolbox carries out several incremental learning algorithms for local model networks. It is build to run automatically and recommends the "best" model achieved. Here "best" stands for the least expected squared error on new data. Therefore, it tries to find a good bias/variance tradeoff and thus to avoid overfitting. Of course, each method can be called on each own and many parameters and options are available to fine tune the results. The primary goal of this toolbox, however, is to address the non-expert and deliver a very robust performance. It is highly encouraged and recommended to look into more details of the proposed methods. This can improve the performance even further.

The investigated methods are:
\begin{enumerate}
\item LOLIMOT with local linear models
%\item LOLIMOT with local sparse quadratic models
\item LOLIMOT with local quadratic models
\item HILOMOT with local linear models
%\item HILOMOT with local sparse quadratic models
\item HILOMOT with local quadratic models
\end{enumerate}
These methods are briefly discussed in the following. They will be extended and improved in the future. More methods will be added. Some methods can easily be commented out in order to save computation time. 

%\subsection{Polynomials} 
%As a reference, a global polynomial model is constructed. A polynomial of degree $n$ for $p$ inputs has 
%\begin{equation}
%\frac{(n+p)!}{n!\,p!}
%\label{poly_M}
%\end{equation}
%regressors and therefore as many parameters to be estimated. It is obvious, that this number becomes huge for high-dimensional input spaces and high degrees. Therefore a linear structure selection procedure is applied to select only the significant regressors and estimate the associated parameters. This is done with the function {\tt stepwise} out of the statistics toolbox from MATLAB. 


\subsection{LOLIMOT with Local Linear Models} 
LOLIMOT (LOcal LInear MOdel Tree) is extensively covered in \cite{Nelles_Book}. It is an incremental tree-structured algorithm that starts with a simple (typically linear) global model and refines its performance in each iteration by splitting the input space in an axes-orthogonal manner. In each iteration the worst local model is split into two equal halves. The associated local linear models are estimated by a local least squares approach. The algorithm is very fast and robust. Its performance can deteriorate with increasing input dimensions due to the sub-optimality of the axes-orthogonal splits.

%\subsection{LOLIMOT with Local Sparse Quadratic Models \label{sub_LOLIMOTsparse}} 
%This it the LOLIMOT algorithm described above with enhanced local models that additionally contain the terms $u_1^2, u_2^2, \ldots, u_p^2$. This is especially useful, if the model shall be utilized for optimization with respect to one or several of its inputs. 

\subsection{LOLIMOT with Local Quadratic Models \label{sub_LOLIMOTquad}} 
In this case, the LOLIMOT algorithm trains local models of full polynomial type. All quadratic local model terms including all cross-terms of type $u_1u_2$ and similar are considered. It also is recommended for optimization purposes but becomes inefficient for high-dimensional input spaces due to the huge number of cross-terms. Then the sparse quadratic models might be a superior choice.

\subsection{HILOMOT with Local Linear Models \label{sub_HILOMOT}} 
HILOMOT and its original ideas are discussed in e.g. \cite{Ernst_CDC_98, Hartmann_ECC_2009}. It realizes two major improvements over the LOLIMOT algorithm. a)~The strong limitations of axes-orthongonal splits fall and b)~the flat (parallel) model structure is transferred into a hierarchical model structure. Advantage~a) means that axis-oblique splits are carried out which overcomes the key weakness of LOLIMOT. The price to be paid is that the split has to be optimized which necessarily is a nonlinear optimization problem. This makes HILOMOT one or two orders of magnitude slower than LOLIMOT. However, it results in much better and more compact models.  Advantage~b) means that re-activation effects \cite{ShoMur97} due to the normalization can be completely avoided and only involves drawbacks on truly parallel hardware. 

%\subsection{HILOMOT with Local Sparse Quadratic Models} 
%Axes-oblique version of the method from Section \ref{sub_LOLIMOTsparse}; see Section \ref{sub_HILOMOT}. 

\subsection{HILOMOT with Local Quadratic Models} 
Axes-oblique version of the method from Section \ref{sub_LOLIMOTquad}; see Section \ref{sub_HILOMOT}. 

%\subsection{POLYMOT} 
%POLYMOT is a significant extension of the LOLIMOT algorithm that sticks to the axes-orthogonal splits of the input space. It thus shares some of the related benefits (purely linear optimization) and drawback (inefficient in high dimensions). It was first introduced in \cite{Baenfer_ICCA_09, Bänfer_Nelles_CDC_11}. In addition to the axes-othogonal splits that LOLIMOT investigates in each iteration, an increase of polynomial complexity is considered. All splits are compared with the additional alternative of increasing the polynomial complexity. The polynomial complexity is adjusted by using the {\tt stepwise} function and demand the appropriate number of regressors. 
%
%As this happens in each iteration, a trade-off between many simple local models or few complex models or some mixture of both results. Extreme outcomes could be a local linear model network as it could results from the LOLIMOT algorithm or a complex global polynomial. More realistically something in between performs best. The major tuning parameter in the POLYMOT algorithm is the penalty factor that favors or discourages splits. This toolbox uses a penalty factor of 0.9. In cases where the model generated by POLMOT is selected as best, it is strongly advised to "play" with this penalty factor.  
%
%
%\subsection{[HILOMOT/POLYMOT Hybrid]} 
%The POLYMOT idea can also be transferred to an axes-oblique partitioning strategy. The combination with HILOMOT yields the most powerful but most computationally demanding of the above described algorithms. In order to keep computation times acceptable. The structure selection with {\tt stepwise} is only carried out before and after the nonlinear split optimization, not in each iteration step. 


\section{Toolbox Inputs and Outputs}
The toolbox can be called with two inputs {\tt traindata} and \emph{optionally} {\tt valdata} and/or {\tt testdata}, respectively, and delivers two outputs {\tt LMNBest } and {\tt AllLMN }:
\[
\tt{ [LMNBest,~AllLMN]= LMNTrain(traindata, valdata, testdata)} 
\]
This is the help of the toolbox:

Each method builds its models incrementally. The choice for a good model complexity is of paramount importance. 

Commonly, no separate validation data set is available. Then the function is just called as {\tt LMNTrain(data)}. The complexity is determined for each method individually with the help of the AIC and GCV criteria, respectively \cite{akaike1973information, burnham2002model}. It is not aimed for the optimal bias/variance trade-off but other benefits for parsimonious models are also considered such as: computational demand, interpretability, and ease of handling. Furthermore, a too complex model with significant overfitting is considered more dangerous than a too simple model, because it causes the illusion of a good model fit to the unexperienced user. 

If validation data is available, the toolbox shall be called by {\tt LMNTrain(traindata, valdata)}. Then the model is assessed and the recommendation is made according to this separate validation data set which is more reliable. 
Additionally, the user can provide a separate test data set by calling {\tt LMNTrain(traindata, [~], testdata)}, which is then used only for testing the model. This ensures comparability to other models.

\section{Model Use}
The model considered best is characterized by the structure {\tt LMNBest}. The model output can be evaluated for a single data point or a whole data set by 
\begin{center} 
{\tt output = LMNBest.calculateModelOutput(input)}
\end{center}
with {\tt input} as $1 \times p$ vector ($p =$ number of inputs) or $N \times p$ matrix ($N=$ number of data points), respectively. 

For a one- or two-dimensional visualization of the model the user can type {\tt LMNBest.plotModel} into the {\sc Matlab} command window. For further information type {\tt help plotModel}. Furthermore, the partitioning can be plotted by typing {\tt LMNBest.plotPartition}.

%\newpage
\section{Toolbox Examples}

The toolbox functionality is demonstrated on four static example data sets. The user just has to run the function {\tt LMNTrainDemo}. After choosing the example the training procedure starts and after completion the user gets informed about the modeling via the LMNTool output screen (see Figs.~\ref{example} and \ref{exampleVal}).

\begin{figure}[hb]
\centering
\includegraphics[width=0.65\textwidth]{example.eps}
\caption{Output screen of the LMNTool for modeling \emph{without} test and validation data. In contrast to modeling just using the training data now the model complexity is selected on the validation data (top right). The selected model is marked with the black cross and highlighted in the title of the figure for the penalty loss function. As additional plots the model vs. test data correlation (bottom middle) and the model vs. validation data correlation (bottom rigth) are drawn.}
\label{example}
\end{figure}


\begin{figure}[h]
\centering
\includegraphics[width=\textwidth]{exampleVal.eps}
\caption{Output screen of the LMNTool for modeling \emph{with} test or validation data. The model complexity is chosen with respect to the penalty loss function value (top left). Furthermore, training errors of each modeling approach (top right), error histogram (bottom left) and data vs. model correlation (bottom right) are illustrated.}
\label{exampleVal}
\end{figure}


%\begin{figure}[htbp]
%\begin{center}
%\includegraphics[width=\textwidth]{Bilder/NumOfIn-5-Noise-0_Best_Worst.eps}
%\caption{Alle $4$ Beispielprozesse zeigten mit $5$ Eingängen und unverrauschten Trainingsdaten die größte Differenz innerhalb dieser Untersuchung.}
%\end{center}
%\label{Num-5-Noi-0.0-B-W}
%\end{figure}

\newpage
\bibliography{toolbox}
\bibliographystyle{plain}

\end{document}


