/*!
\page devguide_mat BDM Development - your first bdm class in Matlab
bdmtoolbox provides a prodefined points of extensions in directory applications/bdmtoolbox/mex/mex_classes.
The easiest way to extend bdmtoolbox is to use clasess starting with "mex". The simplies class of this type is mexFnc, which just contains name of an arbitrary matlab function to call.
See \ref ug_pdf_fnc example.
This tutorial assumes you are familiar with basic classes of BDM: random variables \ref userguide_pdf and estimators \ref userguide_estim.
If no premade class is suitable for your mneeds, you can create your own class in Matlab via predefined interfaces to class epdf (mexEpdf) and BM (mexBM).
These classes are treated as regular C++ classes in prepared mex files. Creation of a new estimator (offspring of mexBM) can be achieved by inheriting from class mexBM.
In practice this means creation of a new file, e.g. my_bm.m, with header:
\code
class my_bm < mexBM
\endcode
This file must be located in directories known to matlab (must be in its path). This can be any directory if you use addpath command.
See relevant documentation in Matlab.
\section dev_mat_pdf Creating your own probability density
Creation of your own density is achieved by inheriting from mexEpdf:
\code
class my_pdf < mexEpdf
\endcode
in a file called my_pdf.m.
Calling methods from mexEpdf will result in an error that tells you which functionality is missing.
It is necessary to give them a meaning.
Epdf is a class representing probability density function of a multivariate (vector) variable.
The easiest example of a density is dirac Delta function:
\code
classdef mexDirac < mexEpdf
% Dirac delta probability distribution
properties % DATA structures
point % point of the support
end
methods
function m=mean(obj) % compute mean values of the density
m = obj.point;
end
function obj=validate(obj) % check if data structures are consistent
% point should be a column
if (size(obj.point,2)>1)
if (size(obj.point,1)==1) % it is row
obj.point = obj.point';
end
else
error('Point in mexDirac is not a vector');
end
end
function dim=dimension(obj) % inform others about your dimensions
dim = size(obj.point,1);
end
function v=variance(obj) % compute variance
v=zeros(size(obj.point));
end
function l=evallog(obj,x) % return logarithm of your density at point x
if obj.point==x
l = inf;
else
l=0;
end
end
function s=sample(obj); % return random sample from your density
s = obj.point;
end
end
end
\endcode
Once you define these functions, this object can be plugged into (almost) an arbitrary algoritm, such as simulator, etc.
For seamless use with other objects, you nedd to define its random variables - via attribute rv, \ref ug_pdf_create.
For inspiration see example class mexLaplace.m and its use in simulation in tutorial/userguide/epdfds_example.m
\section dev_mat_bm Creating your own Bayesian Model (BM)
If you look into the mexBM.m file you see the basic attributes and default functions. mexBM as such will not work as a valid estimator.
In order to make it work, you need to re-define at least functions dimensions and bayes, and attribute apost_pdf:
\code
function dims=dimensions(p)
%please fill
%dims = [size_of_posterior size_of_data size_of_condition]
dims = [0,0,0] %
end
\endcode
This function announces sizes of random variables to other objects.
These sizes can be fixed for a specific problem, or derived from size of internal objects. E.g. for a Kalman filter, the function would be:
\code
function dims=dimensions(p)
%please fill
%dims = [size_of_posterior size_of_data size_of_condition]
dims = [size(A,1), size(C,1), size(B,2)] %
end
\endcode
Which indicates that dimension of \f$x_t\f$ is the rows of matrix A, size of \f$y_t\f$ is the rows of matrix C and the size of \f$u_t\f$ is given by columns of matrix B.
The argument apost_pdf should be a valid object of class mexEpdf, see previous section how to achieve it.
The second function to extend is
\code
function obj=bayes(obj,dt,cond)
% transform old estimate into new estimate
end
\endcode
Here, you need to define an update (both time- and data-update) of your Bayesian filter.
Note, that statistics of the posterior are stored in attribute apost_pdf. So, the task of this function is to update fileds in apost_pdf.
Any additional sttributes needed for this task shoul be added to properties of the object.
An example calss of BM for Laplace pdf is implemented in class mexLaplaceBM.m. This class represents maximum likelihood estimator of two parameters of the density.
One parameter is a median of data on a sliding window.
This can be interpreted as a coarse approximation of Bayesian filtering with posterior reduced to Dirac delta.
Hence, we need:
- apost_pdf to be of type mexDirac,
- store information about window length in the properties
- store the windown of data in properties
- the method bayes adds new data to the window and copies median of the window to the point in apost_pdf.
This is implemented as follows:
\code
classdef mexLaplaceBM < mexBM
% Approximate Bayesian estimator of parameters of Laplace distributed observation.
% Maximum likelihood approximation of the Bayes rule is used, posterior is in the form of dirac.
properties
max_window_length = 10; % max window length (default = 10)
data_window =[]; % sliding window of data
end
methods
function obj=validate(obj) % prepare all internal objects for use
obj.apost_pdf = mexDirac;
obj.apost_pdf.point = [0;0];
obj.log_evidence = 0; % evidence is not computed!
disp('val');
end
function dims=dimensions(obj)
%please fill: dims = [size_of_posterior size_of_data size_of_condition]
dims = [2,1,0] % we have: [2d parameters, 1d observations, 0d condition]
end
function obj=bayes(obj,dt,cond) % approximate bayes rule
if size(obj.data_window,2)>=obj.max_window_length
obj.data_window = [dt obj.data_window(1:end-1)];
else
obj.data_window = [dt obj.data_window];
end
% transform old estimate into new estimate
m_hat = mean(obj.data_window);
b_hat = sum(abs(obj.data_window-m_hat))/ size(obj.data_window,2);
obj.apost_pdf.point = [m_hat; b_hat]; % store result in psoterior pdf
end
function p=epredictor(obj,cond) % when predictive density is needed approximate it by Laplace with point estimates
% return predictive density (max likelihood)
p = mexLaplace;
p.mu = obj.apost_pdf.point(1);
p.b = obj.apost_pdf.point(2);
p=p.validate;
end
end
end
\endcode
Similar extensions to mexEpdf and mexBM can be made to other classes. See library/mex/mex_Epdf.h and library/mex/mex_BM.h how it is done.
As always, feel free to contact the maintainer for help.
*/