Git Product home page Git Product logo

libdirectional's Introduction

libDirectional

libDirectional is a library for directional statistics as well as recursive estimation on directional manifolds. The supported manifolds include

  • unit circle
  • unit hypersphere
  • complex unit hypersphere
  • torus
  • SE(2)

For each of these manifolds, the library contains certain probability distributions as well as recursive filtering algorithms.

Installation

Requirements:

  • Matlab 2022a or later (earlier versions of libDirectional are compatible with older versions of Matlab)
  • a suitable compiler (Visual Studio 2015, Visual Studio 2017, MinGW64, gcc 4.7 or later, XCode)

To use libDirectional, add the entire lib-folder including subdirectories to Matlab's search path. Then, change to the lib-folder and run the compileAll.m script. This script should compile all mex-files used by libDirectional.

If you experience any issues, run mex -setup and mex -setup C++ to ensure that you have selected the correct compiler. In case you have trouble compiling with gcc, make sure that you are using a version that is officially supported by MATLAB.

The following toolboxes are recommended for libDirectional:

  • image_toolbox
  • optimization_toolbox
  • statistics_toolbox
  • symbolic_toolbox

The statistics toolbox and the optimization toolbox are fairly widely used, but the other toolboxes are only required for certain very specific features.

Example: Plotting Probability Density Functions

For example, we can generate a two-dimensional plot of the pdf of a wrapped normal distribution with parameters mu = 2 and sigma = 1.3 simply by typing the following two commands.

>> wn = WNDistribution(2, 1.3);
>> wn.plot2d();

We can then set the labels and axis using the following code:

>> setupAxisCircular('x');
>> xlabel('x'); ylabel('f(x)');

Similarly, we can create plots of other distributions. A three-dimensional plot of the pdf of a von Mises distribution with parameters mu = 6 and kappa=0.5 could be generated using the following code.

>> vm = VMDistribution(6, 0.5);
>> vm.plot3d('color','red');
>> hold on; vm.plotCircle('color','black'); hold off;
>> xlabel('cos(x)'); ylabel('sin(x)'); zlabel('f(x)');

Example: Numerical and Analytical Calculation

Let us again consider the wrapped normal distribution defined in the previous example. Suppose we want to calculate the first trigonometric moment, i.e., E(exp(ix)), of this distribution. For this purpose, we simply call the corresponding function:

>> wn.trigonometricMoment(1)

This produces the output

ans =
  -0.1788 + 0.3906i

In the case of the wrapped normal distribution, trigonometricMoment is a function inside the class WNDistribution that implements an analytic calculation of the trigonometric moment. If no analytic solution was implemented, the function trigonometricMoment in the base class AbstractCircularDistribution would automatically fall back to an algorithm based on numerical integration. Even though an analytical solution is available for the wrapped normal distribution, we can still call the numerical algorithm as follows.

>> wn.trigonometricMomentNumerical(1)

We obtain the result

ans =
  -0.1788 + 0.3906i

This can, for example, be used to compare the numerical and analytical results in order to ensure correctness of the analytical implementation. In this case, both results match up to the displayed number of digits, but in certain cases, analytical and numerical solutions may differ more significantly.

Example: Nonlinear Circular Filtering

Let us consider a system with circular state x_k in [0, 2pi) and dynamics

x_{k+1} = a(x_k) + w_k 
a_k(x_k) = x + 0.1 cos(x_k) mod 2 pi 

where w_k is distributed according to WN(x; 0, 0.4) If we assume that the current state is distributed according to WN(x; 2, 0.5), we can perform the prediction step with the WN-assumed filter using the following commands.

>> filter = WNFilter();
>> filter.setState(WNDistribution(2,0.5));
>> a = @(x) mod(x + 0.1*cos(x),2*pi);
>> filter.predictNonlinear(a, WNDistribution(0, 0.4));
>> filter.getEstimate()

This produces the output

ans = 
  WNDistribution with properties:
       mu: 1.9623
    sigma: 0.6092

As you can see, the predicted density is returned as a wrapped normal distribution. Now we consider the measurement model

z_k = h_k(x_k) + v_k

with

h_k: [0,2 pi) -> R,  h_k(x) = sin(x)

where v_k is additive noise distributed according to N(x; 0, 0.7). As you can see, we have a circular state, but a real-valued measurement, in this case. However, a circular measurement (or a measurement on a completely different manifold) would be possible as well. If we obtain a measurement, say z = 0.3, we can perform the measurement update as follows.

>> h = @(x) sin(x);
>> measurementNoise = GaussianDistribution(0, 0.7);
>> likelihood = LikelihoodFactory.additiveNoiseLikelihood(h, measurementNoise);
>> filter.updateNonlinearProgressive(likelihood, 0.3)
>> filter.getEstimate()

This produces the output

ans = 
  WNDistribution with properties:
       mu: 2.0030
    sigma: 0.6414

Once again, we obtain the result as a wrapped normal distribution.

Unit Tests

The unit tests for libDirectional are located in the tests subfolder. You can automatically run all unit tests by executing runLibDirectionalUnitTests. Tests that are very computationally expensive are skipped by default. You can enable computationally expensive tests by running runLibDirectionalUnitTests(true), which is more thorough but takes much more time.

Externals

libDrectional relies on the following external libraries, which are also included in the externals folder.

Furthermore, we use the script circVMcdf by Shai Revzen (GPLv3). We also use some code from libBingham (BSD), but this library is not in the externals folder as only small parts are used.

License

libDirectional is licensed under the GPLv3 license.

Citation

If you use libDirectional in your research, please cite the library using as follows.

@Article{libdirectional,
author = {Gerhard Kurz and Igor Gilitschenski and Florian Pfaff and
          Lukas Drude and Uwe D. Hanebeck and Reinhold Haeb-Umbach
          and Roland Y. Siegwart},
title = {Directional Statistics and Filtering Using {libDirectional}},
year = {2019},
journal = {Journal of Statistical Software},
volume = {89},
number = {4},
pages = {1--31},
doi = {10.18637/jss.v089.i04},
}

Contact

Lead author: Gerhard Kurz

Mail: kurz.gerhard (at) gmail.com

Web: https://www.gerhardkurz.de

Contributors:

  • Igor Gilitschenski
  • Florian Pfaff
  • Lukas Drude

libdirectional's People

Contributors

boeddeker avatar daniel-frisch-kit avatar florianpfaff avatar gerhardkurz avatar igilitschenski avatar libdirectional avatar lukasdrude avatar

Stargazers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

Watchers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

libdirectional's Issues

Symbolic Toolbox Dependency

The list of requirements only specifies MATLAB r2014a or later and a supported C++ compiler. However, when attempting to compile, I got an error in line 18 of lib/util/cBinghamNormSymbolic.m because it uses a command from the Symbolic Toolbox, which I do not have.

The requirements should either be updated to give a list of required toolboxes, or else the outputs from the symbolic calculations should be included in the project without having to regenerate them on each user's machine.

Generate helper functions for d = 2
To use 'sym', the following product must be licensed, installed, and enabled:
  Symbolic Math Toolbox

Error in cBinghamNormSymbolic (line 18)
    X = sym('x', [D, 1]);

Error in generateSymbolic (line 4)
    [c, X] = cBinghamNormSymbolic(D);

Error in compileUtil (line 52)
    generateSymbolic(d);

Error in compileAll (line 16)
compileUtil

Compile error with MinGW64 compiler

I am trying to compile libDirectional on MATLAB R2018a running on Windows 7 using the MinGW64 C++ compiler (which is at the top of the supported compiler list). When I try, I encounter the following error:

>> addpath(genpath('lib'))
>> compileAll
Building with 'MinGW64 Compiler (C++)'.
MEX completed successfully.
Building with 'MinGW64 Compiler (C++)'.
MEX completed successfully.
Building with 'MinGW64 Compiler (C++)'.
Error using mex
E:\Users\icraig\Documents\Code\libDirectional\lib\util\numericalBinghamMLE.cpp: In function 'int binghamMLE(int, double*, double*, double*)':
E:\Users\icraig\Documents\Code\libDirectional\lib\util\numericalBinghamMLE.cpp:135:39: error: 'M_PI' was not declared in this scope
             ncDeriv(j) = -tmpNc[2]/(2*M_PI);
                                       ^


Error in compileUtil (line 40)
mex(options{:}, fullfile(utilDir,'numericalBinghamMLE.cpp'),fullfile(utilDir,'binghamNormalizationConstant.cpp'))

Error in compileAll (line 16)
compileUtil

If it helps, here is my compiler configuration:

K>> cc = mex.getCompilerConfigurations;
K>> cc(1)

ans = 

  CompilerConfiguration with properties:

             Name: 'MinGW64 Compiler (C++)'
     Manufacturer: 'GNU'
         Language: 'C++'
          Version: '5.3.0'
         Location: 'C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset'
        ShortName: 'mingw64-g++'
         Priority: 'E'
          Details: [1×1 mex.CompilerConfigurationDetails]
       LinkerName: 'C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin\g++'
    LinkerVersion: ''
           MexOpt: 'C:\Users\icraig\AppData\Roaming\MathWorks\MATLAB\R2018a\mex_C++_win64.xml'

K>> cc(1).Details

ans = 

  CompilerConfigurationDetails with properties:

         CompilerExecutable: 'C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin\g++'
              CompilerFlags: '-fexceptions -fno-omit-frame-pointer -std=c++11'
          OptimizationFlags: '-O -DNDEBUG'
                 DebugFlags: '-g'
           LinkerExecutable: 'C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin\g++'
                LinkerFlags: '-m64 -Wl,--no-undefined -shared -static'
    LinkerOptimizationFlags: '-s'
           LinkerDebugFlags: '-g'
                     SetEnv: 'set COMPILER=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin\gcc↵                set CXXCOMPILER=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin\g++↵                set COMPFLAGS=-c -fexceptions -fno-omit-frame-pointer -m64 -DMATLAB_MEX_FILE  -DMATLAB_MEX_FILE ↵                set CXXCOMPFLAGS=-c -fexceptions -fno-omit-frame-pointer -std=c++11 -m64 -DMATLAB_MEX_FILE  -DMATLAB_MEX_FILE ↵                set OPTIMFLAGS=-O -DNDEBUG↵                set DEBUGFLAGS=-g↵                set LINKER=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin\gcc↵                set CXXLINKER=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin\g++↵                set LINKFLAGS=-m64 -Wl,--no-undefined -shared -static -L"E:\Program Files\MATLAB\R2018a\extern\lib\win64\mingw64" -llibmx -llibmex -llibmat -lm -llibmwlapack -llibmwblas -llibMatlabDataArray -llibMatlabEngine -Wl,"E:\Program Files\MATLAB\R2018a/extern/lib/win64/mingw64/mexFunction.def"↵                set LINKDEBUGFLAGS=-g↵                set NAME_OUTPUT= -o "%OUTDIR%%MEX_NAME%%MEX_EXT%"↵set PATH=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\bin;E:\Program Files\MATLAB\R2018a\extern\include\win64;E:\Program Files\MATLAB\R2018a\extern\include;E:\Program Files\MATLAB\R2018a\simulink\include;E:\Program Files\MATLAB\R2018a\lib\win64;%MATLAB_BIN%;%PATH%↵set INCLUDE=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\include;C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\x86_64-w64-mingw32\include\c++;C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\x86_64-w64-mingw32\include;C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\x86_64-w64-mingw32\include\c++\x86_64-w64-mingw32;C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\x86_64-w64-mingw32\include\c++\backward;%INCLUDE%↵set LIB=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset\lib;;%LIB%↵set LIBPATH=E:\Program Files\MATLAB\R2018a\extern\lib\win64;%LIBPATH%↵'
           CommandLineShell: 'set MINGW_ROOT_PATH=C:\ProgramData\MATLAB\SupportPackages\R2018a\3P.instrset\mingw_w64.instrset'
        CommandLineShellArg: ''

Square operation missing

Hi ! From the "Random sampling from the watson distribution" paper, it seems that U_1 term in step 3 is squared, while your implementation just uses U_1 as is.

Is this a typo in your implementation or in the paper ?

Cheers,
David

ToroidalMixture parameter estimation

Hi,

I have a dataset which seems to have a ToroidalMixture distribution or other mixture-like circular distributions.

Would it be possible to estimate the parameters of a ToroidalMixture (or other circular-mixture distributions) for this data?
I can see in your package examples for estimating parameters of WNDistribution from data but for mixtures it seems that there is no estimation routine. Is that correct?

Thanks in advance,
Hamidreza

Recommend Projects

  • React photo React

    A declarative, efficient, and flexible JavaScript library for building user interfaces.

  • Vue.js photo Vue.js

    🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.

  • Typescript photo Typescript

    TypeScript is a superset of JavaScript that compiles to clean JavaScript output.

  • TensorFlow photo TensorFlow

    An Open Source Machine Learning Framework for Everyone

  • Django photo Django

    The Web framework for perfectionists with deadlines.

  • D3 photo D3

    Bring data to life with SVG, Canvas and HTML. 📊📈🎉

Recommend Topics

  • javascript

    JavaScript (JS) is a lightweight interpreted programming language with first-class functions.

  • web

    Some thing interesting about web. New door for the world.

  • server

    A server is a program made to process requests and deliver data to clients.

  • Machine learning

    Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.

  • Game

    Some thing interesting about game, make everyone happy.

Recommend Org

  • Facebook photo Facebook

    We are working to build community through open source technology. NB: members must have two-factor auth.

  • Microsoft photo Microsoft

    Open source projects and samples from Microsoft.

  • Google photo Google

    Google ❤️ Open Source for everyone.

  • D3 photo D3

    Data-Driven Documents codes.