FHWS/Fusion2018
Archived
4
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

Merge branch 'master' of https://git.frank-ebner.de/FHWS/Fusion2018

Browse Source
This commit is contained in:
toni
2018-02-20 15:52:47 +01:00
parent f04f38de30 0b842a2f62
commit 4c72c4f15f
6 changed files with 108 additions and 68 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
tex/bare_conf.tex
View File

@@ -114,6 +114,8 @@
\newcommand{\Lideal} {\ensuremath{ L_{\text{ideal}} }}
\newcommand{\floor} [1]{\ensuremath{ \lfloor #1 \rfloor }}
\newcommand{\etal} [1]{#1~et~al.}
\newcommand{\G} [2]{\ensuremath{ \mathcal{N} \left(#1,#2\right) }}
\newcommand{\VecTwo}[2]{\ensuremath{\left[\begin{smallmatrix} #1 \\ #2 \end{smallmatrix}\right] }}
\newcommand{\qq} [1]{``#1''}

40
tex/chapters/experiments.tex
View File

@@ -1,21 +1,34 @@
\section{Experiments}
We now empirically evaluate the accuracy of our method and compare its runtime performance with other state of the art approaches.
To conclude our findings we present a real world example from a indoor localisation system.
All tests are performed on a Intel Core \mbox{i5-7600K} CPU with a frequency of $4.5 \text{GHz}$, which supports the AVX2 instruction set, hence 256-bit wide SIMD registers are available.
We compare our C++ implementation of the box filter based KDE to the KernSmooth R package and the \qq{FastKDE} implementation \cite{oBrien2016fast}.
The KernSmooth packages provides a FFT-based BKDE implementation based on optimized C functions at its core.
\subsection{Error}
In order to quantity the accuracy of our method the mean integrated squared error (MISE) is used.
The ground truth is given as a synthetic data set drawn from a mixture normal density.
Clearly, the choice of the ground truth distribution affects the resulting error.
However, as our method approximates the KDE it is only of interest to evaluate the closeness to the KDE and not to the ground truth itself.
\subsection{Mean Integrated Squared Error}
We now empirically evaluate the accuracy of our method, using the mean integrated squared error (MISE).
The ground truth is given as $N=1000$ synthetic samples drawn from a bivariate mixture normal density $f$
\begin{equation}
\begin{split}
\bm{X} \sim &\G{\VecTwo{0}{0}}{0.5\bm{I}} + \G{\VecTwo{3}{0}}{\bm{I}} \\
&+ \G{\VecTwo{0}{3}}{\bm{I}} + \G{\VecTwo{-3}{0} }{\bm{I}} + \G{\VecTwo{0}{-3}}{\bm{I}}
\end{split}
\end{equation}
where the majority of the probability mass lies in the range $[-6; 6]^2$.
Clearly, the structure of the ground truth affects the error in the estimate, but as our method approximates the KDE only the closeness to the KDE is of interest.
Therefore, the particular choice of the ground truth is only of minor importance here.
At first we evaluate the accuracy of our method as a function of the bandwidth $h$ in comparison to the exact KDE and the BKDE.
Both the BKDE and the extended box filter estimate resemble the error curve of the KDE quite well and stable.
They are rather close to each other, with a tendency to diverge for larger $h$.
In contrast, the error curve of the box filter estimate has noticeable jumps at $h=(0.4; 0.252; 0.675; 0.825)$.
These jumps are caused by the rounding of the integer-valued box width given by \eqref{eq:boxidealwidth}.
As the extend box filter is able to approximate an exact $\sigma$, it lacks these discontinues.
The exact KDE, evaluated at $50^2$ points, is compared to the BKDE, box filter, and extended box filter approximation, which are evaluated at a smaller grid with $30^2$ points.
The MISE between $f$ and the estimates as a function of $h$ are evaluated, and the resulting plot is given in figure~\ref{fig:evalBandwidth}.
\begin{figure}
\label{fig:evalBandwidth}
\end{figure}
Other test cases of theoretical relevance are error as a function of the grid size $G$ and the sample size $N$.
However, both cases do not give a deeper insight of the error behaviour of our method, as it closely mimics the error curve of the KDE and only confirm the theoretical expectations.
% kde, box filter, exbox in abhänigkeit von h (bild)
% sample size und grid size text
@@ -23,5 +36,8 @@ At first we evaluate the accuracy of our method as a function of the bandwidth $
\subsection{Performance}
All tests are performed on a Intel Core \mbox{i5-7600K} CPU with a frequency of $4.5 \text{GHz}$, which supports the AVX2 instruction set, hence 256-bit wide SIMD registers are available.
We compare our C++ implementation of the box filter based KDE to the KernSmooth R package and the \qq{FastKDE} implementation \cite{oBrien2016fast}.
The KernSmooth packages provides a FFT-based BKDE implementation based on optimized C functions at its core.
\input{chapters/realworld}

1
tex/chapters/kde.tex
View File

@@ -85,6 +85,7 @@ However, for many applications it is recommend to use the simple binning rule
\end{align}
or the common linear binning rule which divides the sample into two fractional weights shared by the nearest grid points
\begin{align}
\label{eq:linearBinning}
b_j(x,\delta) &=
\begin{cases}
w_j(1-|\delta^{-1}x-j|) & \text{if } |\delta^{-1}x-j|\le1 \\

27
tex/chapters/multivariate.tex
View File

@@ -1,5 +1,4 @@
\section{Extension to multi-dimensional data}
\todo{WIP}
So far only univariate sample sets were considered.
This is due to the fact, that the equations of the KDE \eqref{eq:kde}, BKDE \eqref{eq:binKde}, Gaussian filter \eqref{eq:gausFilt}, and the box filter \eqref{eq:boxFilt} are quite easily extended to multi-dimensional input.
Each method can be seen as several one-dimensional problems combined to a multi-dimensional result.
@@ -11,22 +10,7 @@ Multivariate kernel functions can be constructed in various ways, however, a pop
Such a kernel is constructed by combining several univariate kernels into a product, where each kernel is applied in each dimension with a possibly different bandwidth.
Given a multivariate random variable $X=(x_1,\dots ,x_d)$ in $d$ dimensions.
The sample $\bm{X}$ is a $n\times d$ matrix defined as \cite[162]{scott2015}
\begin{equation}
\bm{X}=
\begin{pmatrix}
X_1 \\
\vdots \\
X_n \\
\end{pmatrix}
=
\begin{pmatrix}
x_{11} & \dots & x_{1d} \\
\vdots & \ddots & \vdots \\
x_{n1} & \dots & x_{nd}
\end{pmatrix} \text{.}
\end{equation}
The sample $\bm{X}$ is a $n\times d$ matrix defined as \cite[162]{scott2015}.
The multivariate KDE $\hat{f}$ which defines the estimate pointwise at $\bm{x}=(x_1, \dots, x_d)^T$ is given as \cite[162]{scott2015}
\begin{equation}
\label{eq:mvKDE}
@@ -41,16 +25,17 @@ In addition, only smoothing in the direction of the axes are possible.
If smoothing in other directions is necessary, the computation needs to be done on a prerotated sample set and the estimate needs to be rotated back to fit the original coordinate system \cite{wand1994fast}.
For the multivariate BKDE, in addition to the kernel function the grid and the binning rules need to be extended to multivariate data.
\todo{Reicht hier text oder müssen Formeln her?}
Their extensions are rather straightforward, as the grid is easily defined on many dimensions.
Likewise, the ideas of common and linear binning rule scale with the dimensionality \cite{wand1994fast}.
In general multi-dimensional filters are multi-dimensional convolution operations.
However, by utilizing the separability property of convolution a straightforward and a more efficient implementation can be found.
Convolution is separable if the filter kernel is separable, i.e. it can be split into successive convolutions of several kernels.
In example, the Gaussian filter is separable, because of $e^{x^2+y^2} = e^{x^2}\cdot e^{y^2}$.
Likewise digital filters based on such kernels are called separable filters.
They are easily applied to multi-dimensional signals, because the input signal can be filtered in each dimension separately by an one-dimensional filter.
They are easily applied to multi-dimensional signals, because the input signal can be filtered in each dimension individually by an one-dimensional filter.
The Gaussian filter is separable, because of $e^{x^2+y^2} = e^{x^2}\cdot e^{y^2}$.
% KDE:

37
tex/chapters/mvg.tex
View File

@@ -9,7 +9,7 @@ Consequently, the filter kernel of a Gaussian filter is a Gaussian with finite s
Assuming a finite-support Gaussian filter kernel of size $M$ and a input signal $x$, discrete convolution produces the smoothed output signal
\begin{equation}
\label{eq:gausFilt}
y[n] = \frac{1}{\sigma\sqrt{2\pi}} \sum_{k=0}^{M-1} x[k]\expp{-\frac{(n-k)^2}{2\sigma^2}} \text{,}
(G_{\sigma}*x)(i) = \frac{1}{\sigma\sqrt{2\pi}} \sum_{k=0}^{M-1} x(k) \expp{-\frac{(i-k)^2}{2\sigma^2}} \text{,}
\end{equation}
where $\sigma$ is a smoothing parameter called standard deviation.
@@ -31,14 +31,14 @@ The box filter is a simplistic filter defined as convolution of the input signal
A discrete box filter of \qq{radius} $l\in N_0$ is given as
\begin{equation}
\label{eq:boxFilt}
y[n] = \sum_{k=0}^{L-1} x[k] B_L(n-k)
(B_L*x)(i) = \sum_{k=0}^{L-1} x(k) B_L(i-k)
\end{equation}
where $L=2l+1$ is the width of the filter and
\begin{equation}
\label{eq:boxFx}
B_L(n)=
B_L(i)=
\begin{cases}
L^{-1} & \text{if } -l \le n \le l \\
L^{-1} & \text{if } -l \le i \le l \\
0 & \text{else }
\end{cases}
\end{equation}
@@ -55,11 +55,26 @@ Opposed to the Gaussian filter where several evaluations of the exponential func
The most efficient way to compute a single output value of the box filter is:
\begin{equation}
\label{eq:recMovAvg}
y[i] = y[i-1] + x[i+r] - x[i-r-1]
y(i) = y(i-1) + x(i+l) - x(i-l-1)
\end{equation}
This recursive calculation scheme further reduces the time complexity of the box filter to $\landau{N}$, by reusing already calculated values.
Furthermore, only one addition and subtraction is required to calculate a single output value.
The overall algorithm to efficiently compute \eqref{eq:boxFilt} is listed in Algorithm~\ref{alg:naiveboxalgo}.
\begin{algorithm}[ht]
\caption{Recursive 1D box filter}
\label{alg:naiveboxalgo}
\begin{algorithmic}[1]
\Statex \textbf{Input:} $f$ defined on $[1,N]$ and box radius $l$
\Statex \textbf{Output:} $u \gets B_L*f$
\State $a := \sum_{i=-l}^{l} f(i) $
\State $u(1) := a/(2l+1)$
\For{$i=2 \textbf{ to } N$}
\State $a:= a + f(i+l) - f(i-l-1) $ \Comment{\eqref{eq:recMovAvg}}
\State $u(i) := a/(2l+1)$
\EndFor
\end{algorithmic}
\end{algorithm}
Given a fast approximation scheme, it is necessary to construct a box filter analogous to a given Gaussian filter.
As seen in \eqref{eq:gausFilt}, the solely parameter of the Gaussian kernel is the standard deviation $\sigma$.
@@ -100,18 +115,6 @@ However, the overall shape will not change.
Just like the conventional box filter, the extended version has a uniform value in the range $[-l; l]$, but unlike the conventional the extended box filter has different values at its edges.
This extension introduces only marginal computational overhead over conventional box filtering.
\commentByToni{Warum benutzen wir den extended box filter nicht? oder tun wir das? Liest sich so, als wäre er der heilige Gral.
Ansonsten: Kapitel find ich gut. Vielleicht kann man hier und da noch paar sätze fusionieren um etwas Platz zu sparen.}
\commentByToni{Aber irgendwie fehlt mir hier noch so ein Absatz oder Kapitel "Zusammenstecken" der Mathematik. Die Usage kommt wieder so Hart. Und nochmal eine kleine Diskussion zum Zusammenstecken. Vielleicht kann man das mit dem 2D verbinden? Weil aktuell haben wir die beiden Verfahren besprochen, aber der eigene Anteil ist irgendwie nicht ersichtlich. Was war jetzt deine tolle Leistung hier? Das muss durch eine gute Diskussion klar werden.}

69
tex/chapters/usage.tex
View File

@@ -1,22 +1,60 @@
\section{Usage}
The objective of our method is to allow a reliable recover of the most probable state from a time-sequential Monte Carlo sensor fusion system.
Assuming a sample based representation, our method allows to estimate the density of the unknown distribution of the state space in a narrow time frame.
Such systems are often used to obtain an estimation of the most probable state in near real time.
As the density estimation poses only a single step in the whole process, its computation needs to be as fast as possible.
%The objective of our method is to allow a reliable recover of the most probable state from a time-sequential Monte Carlo sensor fusion system.
%Assuming a sample based representation, our method allows to estimate the density of the unknown distribution of the state space in a narrow time frame.
%Such systems are often used to obtain an estimation of the most probable state in near real time.
%As the density estimation poses only a single step in the whole process, its computation needs to be as fast as possible.
% not taking to much time from the frame
%Consider a set of two-dimensional samples, presumably generated from e.g. particle filter system.
Assuming that the generated samples are often stored in a sequential list, the first step is to create a grid representation.
In order to efficiently compute the grid and to allocate the required memory the extrema of the samples need to be known in advance.
\begin{algorithm}[ht]
\caption{Bivariate \textsc{boxKDE}}
\label{alg:boxKDE}
\begin{algorithmic}[1]
\Statex \textbf{Input:} Samples $\bm{X}_1, \dots, \bm{X}_N$ and weights $w_1, \dots, w_N$
\Statex \textbf{Output:} Approximative density estimate $\hat{f}$ on $G_1 \times G_2$
\Statex
\For{$i=1 \textbf{ to } N$} \Comment{Data binning}
\State Find the $4$ nearest grid points to $\bm{X}_i$
\State Compute bin count $C_{i,j}$ as recommended by \cite{wand1994fast}
\EndFor
\Statex
\State $\tilde{\bm{h}} := \bm{\delta}^{-1} \bm{h}$ \Comment{Scaled bandwidth}
\State $\bm{L} := \floor{\sqrt{12\tilde{\bm{h}}^2n^{-1}+\bm{1}}}$ \Comment{\eqref{eq:boxidealwidth}}
% \State $l := \floor{(L-1)/2}$
\Statex
%\For{$1 \textbf{ to } n$}
\Loop{ $n$ \textbf{times}} \Comment{$n$ box filter iterations}
\For{$ i=1 \textbf{ to } G_1$}
\State Compute $\hat{f}_{i,1:G_2} \gets B_{L_2} * C_{i,1:G_2}$ \Comment{Alg. \ref{alg:naiveboxalgo}}
\EndFor
\For{$ j=1 \textbf{ to } G_2$}
\State Compute $\hat{f}_{1:G_1,j} \gets B_{L_1} * C_{1:G_1,j}$ \Comment{Alg. \ref{alg:naiveboxalgo}}
\EndFor
\EndLoop
\end{algorithmic}
\end{algorithm}
Consider a set of two-dimensional samples with associated weights, e.g. presumably generated from a particle filter system.
The overall process for bivariate data is described in Algorithm~\ref{alg:boxKDE}.
Assuming that the given $N$ samples are stored in a sequential list, the first step is to create a grid representation.
In order to efficiently construct the grid and to allocate the required memory the extrema of the samples need to be known in advance.
These limits might be given by the application, for example, the position of a pedestrian within a building is limited by the physical dimensions of the building.
Such knowledge should be integrated into the system to avoid a linear search over the sample set, naturally reducing the computation time.
The second parameter to be defined by the application is the size of the grid, which can be set directly or defined in terms of bin sizes.
As the number of grid points directly affects both computation time and accuracy, a suitable grid should be as coarse as possible but at the same time narrow enough to produce an estimate sufficiently fast with an acceptable approximation error.
Given the extreme values of the samples and grid sizes $G_1$ and $G_2$ defined by the user, a $G_1\times G_2$ grid can be constructed, using a binning rule from \eqref{eq:simpleBinning} or \eqref{eq:linearBinning}.
As the number of grid points directly affects both computation time and accuracy, a suitable grid should be as coarse as possible, but at the same time narrow enough to produce an estimate sufficiently fast with an acceptable approximation error.
Given the extreme values of the samples and the number of grid points $G$, the computation of the grid has a linear complexity of \landau{N} where $N$ is the number of samples.
If the extreme values are unknown, an additional $\landau{N}$ search is required.
The grid is stored as an linear array in memory, thus its space complexity is $\landau{G}$.
If the extreme values are known in advanced, the computation of the grid is $\landau{N}$, otherwise an additional $\landau{N}$ search is required.
The grid is stored as an linear array in memory, thus its space complexity is $\landau{G_1\cdot G_2}$.
Next, the binned data is filtered with a Gaussian using the box filter approximation.
The box filter width is derived from the standard deviation of the approximated Gaussian, which is in turn equal to the bandwidth of the KDE.
@@ -28,7 +66,7 @@ For this reason, $h$ needs to be divided by the bin size to account the discrepa
Given the scaled bandwidth the required box filter width can be computed. % as in \eqref{label}
Due to its best runtime performance the recursive box filter implementation is used.
If multivariate data is processed, the algorithm is easily extended due to its separability.
Each filter pass is computed in $\landau{G}$ operations, however an additional memory buffer is required.
Each filter pass is computed in $\landau{G}$ operations, however, an additional memory buffer is required.
While the integer-sized box filter requires fewest operations, it causes a larger approximation error due to rounding errors.
Depending on the required accuracy the extended box filter algorithm can further improve the estimation results, with only a small additional overhead.
@@ -36,8 +74,3 @@ Due to its simple indexing scheme, the recursive box filter can easily be comput
Finally, the most likely state can be obtained from the filtered data, i.e. from the estimated discrete density, by searching filtered data for its maximum value.
\commentByToni{An sich ganz cooles Kapitel, aber wir müssen den Bezug nach oben stärker aufbauen. Also die Formeln zitieren. irgendwie halt nach oben referenzieren, damit niemand abgehängt wird.
Würde es Sinn machen das obere irgendwie Algorithmisch darzustellen? Also mit Pseudocode? Weil irgendwie/wo müssen wir ja "DAS IST UNSER APPROACH" stehen haben}.