[2] | 1 | /*! |
---|
[5] | 2 | \file |
---|
| 3 | \brief Bayesian Models (bm) that use Bayes rule to learn from observations |
---|
| 4 | \author Vaclav Smidl. |
---|
[2] | 5 | |
---|
[5] | 6 | ----------------------------------- |
---|
| 7 | BDM++ - C++ library for Bayesian Decision Making under Uncertainty |
---|
| 8 | |
---|
| 9 | Using IT++ for numerical operations |
---|
| 10 | ----------------------------------- |
---|
| 11 | */ |
---|
| 12 | |
---|
[384] | 13 | #ifndef BDMBASE_H |
---|
| 14 | #define BDMBASE_H |
---|
[2] | 15 | |
---|
[351] | 16 | #include <map> |
---|
[263] | 17 | |
---|
[190] | 18 | #include "../itpp_ext.h" |
---|
[357] | 19 | #include "../bdmroot.h" |
---|
[461] | 20 | #include "../shared_ptr.h" |
---|
[384] | 21 | #include "user_info.h" |
---|
[2] | 22 | |
---|
[340] | 23 | using namespace libconfig; |
---|
[270] | 24 | using namespace itpp; |
---|
| 25 | using namespace std; |
---|
[2] | 26 | |
---|
[477] | 27 | namespace bdm { |
---|
[340] | 28 | |
---|
[270] | 29 | typedef std::map<string, int> RVmap; |
---|
| 30 | extern ivec RV_SIZES; |
---|
| 31 | extern Array<string> RV_NAMES; |
---|
| 32 | |
---|
[422] | 33 | //! Structure of RV, i.e. RVs expanded into a flat list of IDs, used for debugging. |
---|
[477] | 34 | class str { |
---|
[270] | 35 | public: |
---|
[477] | 36 | //! vector id ids (non-unique!) |
---|
| 37 | ivec ids; |
---|
| 38 | //! vector of times |
---|
| 39 | ivec times; |
---|
| 40 | //!Default constructor |
---|
| 41 | str ( ivec ids0, ivec times0 ) : ids ( ids0 ), times ( times0 ) { |
---|
| 42 | it_assert_debug ( times0.length() == ids0.length(), "Incompatible input" ); |
---|
| 43 | }; |
---|
[270] | 44 | }; |
---|
[145] | 45 | |
---|
[270] | 46 | /*! |
---|
| 47 | * \brief Class representing variables, most often random variables |
---|
[5] | 48 | |
---|
[270] | 49 | The purpose of this class is to decribe a vector of data. Such description is used for connecting various vectors between each other, see class datalink. |
---|
[32] | 50 | |
---|
[270] | 51 | The class is implemented using global variables to assure uniqueness of description: |
---|
[2] | 52 | |
---|
[270] | 53 | In is a vector |
---|
| 54 | \dot |
---|
| 55 | digraph datalink { |
---|
| 56 | rankdir=LR; |
---|
| 57 | subgraph cluster0 { |
---|
| 58 | node [shape=record]; |
---|
| 59 | label = "RV_MAP \n std::map<string,int>"; |
---|
| 60 | map [label="{{\"a\"| \"b\" | \"c\"} | {<3> 3 |<1> 1|<2> 2}}"]; |
---|
| 61 | color = "white" |
---|
| 62 | } |
---|
| 63 | subgraph cluster1{ |
---|
| 64 | node [shape=record]; |
---|
| 65 | label = "RV_NAMES"; |
---|
| 66 | names [label="{<1> \"b\" | <2> \"c\" | <3>\"a\" }"]; |
---|
| 67 | color = "white" |
---|
| 68 | } |
---|
| 69 | subgraph cluster2{ |
---|
| 70 | node [shape=record]; |
---|
| 71 | label = "RV_SIZES"; |
---|
| 72 | labelloc = b; |
---|
| 73 | sizes [label="{<1>1 |<2> 4 |<3> 1}"]; |
---|
| 74 | color = "white" |
---|
| 75 | } |
---|
| 76 | map:1 -> names:1; |
---|
| 77 | map:1 -> sizes:1; |
---|
| 78 | map:3 -> names:3; |
---|
| 79 | map:3 -> sizes:3; |
---|
| 80 | } |
---|
| 81 | \enddot |
---|
| 82 | */ |
---|
[32] | 83 | |
---|
[477] | 84 | class RV : public root { |
---|
[270] | 85 | protected: |
---|
[477] | 86 | //! size of the data vector |
---|
| 87 | int dsize; |
---|
| 88 | //! number of individual rvs |
---|
| 89 | int len; |
---|
| 90 | //! Vector of unique IDs |
---|
| 91 | ivec ids; |
---|
| 92 | //! Vector of shifts from current time |
---|
| 93 | ivec times; |
---|
[5] | 94 | |
---|
[270] | 95 | private: |
---|
[477] | 96 | //! auxiliary function used in constructor |
---|
| 97 | void init ( const Array<std::string> &in_names, const ivec &in_sizes, const ivec &in_times ); |
---|
| 98 | int init ( const string &name, int size ); |
---|
[270] | 99 | public: |
---|
[477] | 100 | //! \name Constructors |
---|
| 101 | //!@{ |
---|
[271] | 102 | |
---|
[477] | 103 | //! Full constructor |
---|
| 104 | RV ( const Array<std::string> &in_names, const ivec &in_sizes, const ivec &in_times ) { |
---|
| 105 | init ( in_names, in_sizes, in_times ); |
---|
| 106 | } |
---|
[422] | 107 | |
---|
[477] | 108 | //! Constructor with times=0 |
---|
| 109 | RV ( const Array<std::string> &in_names, const ivec &in_sizes ) { |
---|
| 110 | init ( in_names, in_sizes, zeros_i ( in_names.length() ) ); |
---|
| 111 | } |
---|
[422] | 112 | |
---|
[477] | 113 | //! Constructor with sizes=1, times=0 |
---|
| 114 | RV ( const Array<std::string> &in_names ) { |
---|
| 115 | init ( in_names, ones_i ( in_names.length() ), zeros_i ( in_names.length() ) ); |
---|
| 116 | } |
---|
[422] | 117 | |
---|
[477] | 118 | //! Constructor of empty RV |
---|
| 119 | RV() : dsize ( 0 ), len ( 0 ), ids ( 0 ), times ( 0 ) {} |
---|
| 120 | //! Constructor of a single RV with given id |
---|
| 121 | RV ( string name, int sz, int tm = 0 ); |
---|
| 122 | //!@} |
---|
[271] | 123 | |
---|
[477] | 124 | //! \name Access functions |
---|
| 125 | //!@{ |
---|
[271] | 126 | |
---|
[477] | 127 | //! State output, e.g. for debugging. |
---|
| 128 | friend std::ostream &operator<< ( std::ostream &os, const RV &rv ); |
---|
[422] | 129 | |
---|
[477] | 130 | int _dsize() const { |
---|
| 131 | return dsize; |
---|
| 132 | } |
---|
[422] | 133 | |
---|
[477] | 134 | //! Recount size of the corresponding data vector |
---|
| 135 | int countsize() const; |
---|
| 136 | ivec cumsizes() const; |
---|
| 137 | int length() const { |
---|
| 138 | return len; |
---|
| 139 | } |
---|
| 140 | int id ( int at ) const { |
---|
| 141 | return ids ( at ); |
---|
| 142 | } |
---|
| 143 | int size ( int at ) const { |
---|
| 144 | return RV_SIZES ( ids ( at ) ); |
---|
| 145 | } |
---|
| 146 | int time ( int at ) const { |
---|
| 147 | return times ( at ); |
---|
| 148 | } |
---|
| 149 | std::string name ( int at ) const { |
---|
| 150 | return RV_NAMES ( ids ( at ) ); |
---|
| 151 | } |
---|
| 152 | void set_time ( int at, int time0 ) { |
---|
| 153 | times ( at ) = time0; |
---|
| 154 | } |
---|
| 155 | //!@} |
---|
[271] | 156 | |
---|
[477] | 157 | //TODO why not inline and later?? |
---|
[32] | 158 | |
---|
[477] | 159 | //! \name Algebra on Random Variables |
---|
| 160 | //!@{ |
---|
[271] | 161 | |
---|
[477] | 162 | //! Find indices of self in another rv, \return ivec of the same size as self. |
---|
| 163 | ivec findself ( const RV &rv2 ) const; |
---|
| 164 | //! Compare if \c rv2 is identical to this \c RV |
---|
| 165 | bool equal ( const RV &rv2 ) const; |
---|
| 166 | //! Add (concat) another variable to the current one, \return true if all rv2 were added, false if rv2 is in conflict |
---|
| 167 | bool add ( const RV &rv2 ); |
---|
| 168 | //! Subtract another variable from the current one |
---|
| 169 | RV subt ( const RV &rv2 ) const; |
---|
| 170 | //! Select only variables at indices ind |
---|
| 171 | RV subselect ( const ivec &ind ) const; |
---|
[422] | 172 | |
---|
[477] | 173 | //! Select only variables at indices ind |
---|
| 174 | RV operator() ( const ivec &ind ) const { |
---|
| 175 | return subselect ( ind ); |
---|
| 176 | } |
---|
[422] | 177 | |
---|
[477] | 178 | //! Select from data vector starting at di1 to di2 |
---|
| 179 | RV operator() ( int di1, int di2 ) const; |
---|
[422] | 180 | |
---|
[477] | 181 | //! Shift \c time by delta. |
---|
| 182 | void t ( int delta ); |
---|
| 183 | //!@} |
---|
[271] | 184 | |
---|
[477] | 185 | //!\name Relation to vectors |
---|
| 186 | //!@{ |
---|
[271] | 187 | |
---|
[477] | 188 | //! generate \c str from rv, by expanding sizes TODO to_string.. |
---|
| 189 | str tostr() const; |
---|
| 190 | //! when this rv is a part of bigger rv, this function returns indices of self in the data vector of the bigger crv. |
---|
| 191 | //! Then, data can be copied via: data_of_this = cdata(ind); |
---|
| 192 | ivec dataind ( const RV &crv ) const; |
---|
| 193 | //! generate mutual indices when copying data between self and crv. |
---|
| 194 | //! Data are copied via: data_of_this(selfi) = data_of_rv2(rv2i) |
---|
| 195 | void dataind ( const RV &rv2, ivec &selfi, ivec &rv2i ) const; |
---|
| 196 | //! Minimum time-offset |
---|
| 197 | int mint() const { |
---|
| 198 | return min ( times ); |
---|
| 199 | } |
---|
| 200 | //!@} |
---|
[271] | 201 | |
---|
[477] | 202 | // TODO aktualizovat dle soucasneho UI |
---|
| 203 | /*! \brief UI for class RV (description of data vectors) |
---|
[357] | 204 | |
---|
[477] | 205 | \code |
---|
| 206 | rv = { |
---|
| 207 | type = "rv"; //identifier of the description |
---|
| 208 | // UNIQUE IDENTIFIER same names = same variable |
---|
| 209 | names = ["a", "b", "c", ...]; // which will be used e.g. in loggers |
---|
[357] | 210 | |
---|
[477] | 211 | //optional arguments |
---|
| 212 | sizes = [1, 2, 3, ...]; // (optional) default = ones() |
---|
| 213 | times = [-1, -2, 0, ...]; // time shifts with respect to current time (optional) default = zeros() |
---|
| 214 | } |
---|
| 215 | \endcode |
---|
| 216 | */ |
---|
| 217 | void from_setting ( const Setting &set ); |
---|
[357] | 218 | |
---|
[477] | 219 | // TODO dodelat void to_setting( Setting &set ) const; |
---|
[436] | 220 | |
---|
[477] | 221 | //! Invalidate all named RVs. Use before initializing any RV instances, with care... |
---|
| 222 | static void clear_all(); |
---|
[270] | 223 | }; |
---|
[477] | 224 | UIREGISTER ( RV ); |
---|
[32] | 225 | |
---|
[145] | 226 | //! Concat two random variables |
---|
[477] | 227 | RV concat ( const RV &rv1, const RV &rv2 ); |
---|
[2] | 228 | |
---|
[211] | 229 | //!Default empty RV that can be used as default argument |
---|
[270] | 230 | extern RV RV0; |
---|
[145] | 231 | |
---|
[85] | 232 | //! Class representing function \f$f(x)\f$ of variable \f$x\f$ represented by \c rv |
---|
[2] | 233 | |
---|
[477] | 234 | class fnc : public root { |
---|
[270] | 235 | protected: |
---|
[477] | 236 | //! Length of the output vector |
---|
| 237 | int dimy; |
---|
[270] | 238 | public: |
---|
[477] | 239 | //!default constructor |
---|
| 240 | fnc() {}; |
---|
| 241 | //! function evaluates numerical value of \f$f(x)\f$ at \f$x=\f$ \c cond |
---|
| 242 | virtual vec eval ( const vec &cond ) { |
---|
| 243 | return vec ( 0 ); |
---|
| 244 | }; |
---|
[27] | 245 | |
---|
[477] | 246 | //! function substitutes given value into an appropriate position |
---|
| 247 | virtual void condition ( const vec &val ) {}; |
---|
[28] | 248 | |
---|
[477] | 249 | //! access function |
---|
| 250 | int dimension() const { |
---|
| 251 | return dimy; |
---|
| 252 | } |
---|
[270] | 253 | }; |
---|
[2] | 254 | |
---|
[270] | 255 | class mpdf; |
---|
[7] | 256 | |
---|
[4] | 257 | //! Probability density function with numerical statistics, e.g. posterior density. |
---|
[32] | 258 | |
---|
[477] | 259 | class epdf : public root { |
---|
[270] | 260 | protected: |
---|
[477] | 261 | //! dimension of the random variable |
---|
| 262 | int dim; |
---|
| 263 | //! Description of the random variable |
---|
| 264 | RV rv; |
---|
[32] | 265 | |
---|
[270] | 266 | public: |
---|
[477] | 267 | /*! \name Constructors |
---|
| 268 | Construction of each epdf should support two types of constructors: |
---|
| 269 | \li empty constructor, |
---|
| 270 | \li copy constructor, |
---|
[271] | 271 | |
---|
[477] | 272 | The following constructors should be supported for convenience: |
---|
| 273 | \li constructor followed by calling \c set_parameters() |
---|
| 274 | \li constructor accepting random variables calling \c set_rv() |
---|
[271] | 275 | |
---|
[477] | 276 | All internal data structures are constructed as empty. Their values (including sizes) will be set by method \c set_parameters(). This way references can be initialized in constructors. |
---|
| 277 | @{*/ |
---|
| 278 | epdf() : dim ( 0 ), rv() {}; |
---|
| 279 | epdf ( const epdf &e ) : dim ( e.dim ), rv ( e.rv ) {}; |
---|
| 280 | epdf ( const RV &rv0 ) : dim ( rv0._dsize() ) { |
---|
| 281 | set_rv ( rv0 ); |
---|
| 282 | }; |
---|
| 283 | void set_parameters ( int dim0 ) { |
---|
| 284 | dim = dim0; |
---|
| 285 | } |
---|
| 286 | //!@} |
---|
[271] | 287 | |
---|
[477] | 288 | //! \name Matematical Operations |
---|
| 289 | //!@{ |
---|
[271] | 290 | |
---|
[477] | 291 | //! Returns a sample, \f$ x \f$ from density \f$ f_x()\f$ |
---|
| 292 | virtual vec sample() const { |
---|
| 293 | it_error ( "not implemneted" ); |
---|
| 294 | return vec ( 0 ); |
---|
| 295 | }; |
---|
| 296 | //! Returns N samples, \f$ [x_1 , x_2 , \ldots \ \f$ from density \f$ f_x(rv)\f$ |
---|
| 297 | virtual mat sample_m ( int N ) const; |
---|
| 298 | //! Compute log-probability of argument \c val |
---|
| 299 | //! In case the argument is out of suport return -Infinity |
---|
| 300 | virtual double evallog ( const vec &val ) const { |
---|
| 301 | it_error ( "not implemneted" ); |
---|
| 302 | return 0.0; |
---|
| 303 | }; |
---|
| 304 | //! Compute log-probability of multiple values argument \c val |
---|
| 305 | virtual vec evallog_m ( const mat &Val ) const { |
---|
| 306 | vec x ( Val.cols() ); |
---|
| 307 | for ( int i = 0; i < Val.cols(); i++ ) { |
---|
| 308 | x ( i ) = evallog ( Val.get_col ( i ) ) ; |
---|
| 309 | } |
---|
| 310 | return x; |
---|
| 311 | } |
---|
| 312 | //! Compute log-probability of multiple values argument \c val |
---|
| 313 | virtual vec evallog_m ( const Array<vec> &Avec ) const { |
---|
| 314 | vec x ( Avec.size() ); |
---|
| 315 | for ( int i = 0; i < Avec.size(); i++ ) { |
---|
| 316 | x ( i ) = evallog ( Avec ( i ) ) ; |
---|
| 317 | } |
---|
| 318 | return x; |
---|
| 319 | } |
---|
| 320 | //! Return conditional density on the given RV, the remaining rvs will be in conditioning |
---|
| 321 | virtual mpdf* condition ( const RV &rv ) const { |
---|
| 322 | it_warning ( "Not implemented" ); |
---|
| 323 | return NULL; |
---|
| 324 | } |
---|
[271] | 325 | |
---|
[477] | 326 | //! Return marginal density on the given RV, the remainig rvs are intergrated out |
---|
| 327 | virtual epdf* marginal ( const RV &rv ) const { |
---|
| 328 | it_warning ( "Not implemented" ); |
---|
| 329 | return NULL; |
---|
| 330 | } |
---|
[271] | 331 | |
---|
[477] | 332 | //! return expected value |
---|
| 333 | virtual vec mean() const { |
---|
| 334 | it_error ( "not implemneted" ); |
---|
| 335 | return vec ( 0 ); |
---|
| 336 | }; |
---|
[271] | 337 | |
---|
[477] | 338 | //! return expected variance (not covariance!) |
---|
| 339 | virtual vec variance() const { |
---|
| 340 | it_error ( "not implemneted" ); |
---|
| 341 | return vec ( 0 ); |
---|
| 342 | }; |
---|
| 343 | //! Lower and upper bounds of \c percentage % quantile, returns mean-2*sigma as default |
---|
| 344 | virtual void qbounds ( vec &lb, vec &ub, double percentage = 0.95 ) const { |
---|
| 345 | vec mea = mean(); |
---|
| 346 | vec std = sqrt ( variance() ); |
---|
| 347 | lb = mea - 2 * std; |
---|
| 348 | ub = mea + 2 * std; |
---|
| 349 | }; |
---|
| 350 | //!@} |
---|
[271] | 351 | |
---|
[477] | 352 | //! \name Connection to other classes |
---|
| 353 | //! Description of the random quantity via attribute \c rv is optional. |
---|
| 354 | //! For operations such as sampling \c rv does not need to be set. However, for \c marginalization |
---|
| 355 | //! and \c conditioning \c rv has to be set. NB: |
---|
| 356 | //! @{ |
---|
[271] | 357 | |
---|
[477] | 358 | //!Name its rv |
---|
| 359 | void set_rv ( const RV &rv0 ) { |
---|
| 360 | rv = rv0; |
---|
| 361 | } //it_assert_debug(isnamed(),""); }; |
---|
| 362 | //! True if rv is assigned |
---|
| 363 | bool isnamed() const { |
---|
| 364 | bool b = ( dim == rv._dsize() ); |
---|
| 365 | return b; |
---|
| 366 | } |
---|
| 367 | //! Return name (fails when isnamed is false) |
---|
| 368 | const RV& _rv() const { |
---|
| 369 | it_assert_debug ( isnamed(), "" ); |
---|
| 370 | return rv; |
---|
| 371 | } |
---|
| 372 | //!@} |
---|
[377] | 373 | |
---|
[477] | 374 | //! \name Access to attributes |
---|
| 375 | //! @{ |
---|
[377] | 376 | |
---|
[477] | 377 | //! Size of the random variable |
---|
| 378 | int dimension() const { |
---|
| 379 | return dim; |
---|
| 380 | } |
---|
| 381 | //! Load from structure with elements: |
---|
| 382 | //! \code |
---|
| 383 | //! { rv = {class="RV", names=(...),}; // RV describing meaning of random variable |
---|
| 384 | //! // elements of offsprings |
---|
| 385 | //! } |
---|
| 386 | //! \endcode |
---|
| 387 | //!@} |
---|
| 388 | void from_setting ( const Setting &set ) { |
---|
| 389 | RV* r = UI::build<RV> ( set, "rv" ); |
---|
| 390 | if ( r ) { |
---|
| 391 | set_rv ( *r ); |
---|
| 392 | delete r; |
---|
| 393 | } |
---|
| 394 | } |
---|
[377] | 395 | |
---|
[270] | 396 | }; |
---|
[32] | 397 | |
---|
[190] | 398 | |
---|
[5] | 399 | //! Conditional probability density, e.g. modeling some dependencies. |
---|
[32] | 400 | //TODO Samplecond can be generalized |
---|
[477] | 401 | class mpdf : public root { |
---|
[270] | 402 | protected: |
---|
[477] | 403 | //!dimension of the condition |
---|
| 404 | int dimc; |
---|
| 405 | //! random variable in condition |
---|
| 406 | RV rvc; |
---|
[461] | 407 | private: |
---|
[487] | 408 | //! internal epdf, used only as cache to avoid virtual calls of \c _rv() and \c _dimension() |
---|
| 409 | epdf* ep; |
---|
| 410 | |
---|
| 411 | protected: |
---|
| 412 | void set_ep(epdf &iepdf) { |
---|
| 413 | ep = &iepdf; |
---|
| 414 | } |
---|
[461] | 415 | |
---|
[270] | 416 | public: |
---|
[477] | 417 | //! \name Constructors |
---|
| 418 | //! @{ |
---|
[271] | 419 | |
---|
[487] | 420 | mpdf() : dimc ( 0 ), rvc(), ep(NULL) { } |
---|
[2] | 421 | |
---|
[487] | 422 | mpdf ( const mpdf &m ) : dimc ( m.dimc ), rvc ( m.rvc ), ep ( m.ep ) { } |
---|
[477] | 423 | //!@} |
---|
[271] | 424 | |
---|
[477] | 425 | //! \name Matematical operations |
---|
| 426 | //!@{ |
---|
[102] | 427 | |
---|
[477] | 428 | //! Returns a sample from the density conditioned on \c cond, \f$x \sim epdf(rv|cond)\f$. \param cond is numeric value of \c rv |
---|
[487] | 429 | virtual vec samplecond ( const vec &cond ){it_error("Not implemented");return vec(0);}; |
---|
[461] | 430 | |
---|
[477] | 431 | //! Returns \param N samples from the density conditioned on \c cond, \f$x \sim epdf(rv|cond)\f$. \param cond is numeric value of \c rv |
---|
[487] | 432 | virtual mat samplecond_m ( const vec &cond, int N ){it_error("Not implemented");return mat();} |
---|
[461] | 433 | |
---|
| 434 | |
---|
[477] | 435 | //! Shortcut for conditioning and evaluation of the internal epdf. In some cases, this operation can be implemented efficiently. |
---|
[487] | 436 | virtual double evallogcond ( const vec &dt, const vec &cond ){it_error("Not implemented");return 0.0;} |
---|
[32] | 437 | |
---|
[477] | 438 | //! Matrix version of evallogcond |
---|
[487] | 439 | virtual vec evallogcond_m ( const mat &Dt, const vec &cond ){it_error("Not implemented");return vec(0);} |
---|
[201] | 440 | |
---|
[477] | 441 | //! Array<vec> version of evallogcond |
---|
[487] | 442 | virtual vec evallogcond_m ( const Array<vec> &Dt, const vec &cond ){it_error("Not implemented");return vec(0);} |
---|
[271] | 443 | |
---|
[477] | 444 | //! \name Access to attributes |
---|
| 445 | //! @{ |
---|
[271] | 446 | |
---|
[477] | 447 | RV _rv() { |
---|
[487] | 448 | return ep->_rv(); |
---|
[477] | 449 | } |
---|
| 450 | RV _rvc() { |
---|
| 451 | it_assert_debug ( isnamed(), "" ); |
---|
| 452 | return rvc; |
---|
| 453 | } |
---|
| 454 | int dimension() { |
---|
[487] | 455 | return ep->dimension(); |
---|
[477] | 456 | } |
---|
| 457 | int dimensionc() { |
---|
| 458 | return dimc; |
---|
| 459 | } |
---|
[461] | 460 | |
---|
[477] | 461 | //! Load from structure with elements: |
---|
| 462 | //! \code |
---|
[487] | 463 | //! { class = "mpdf_offspring", |
---|
| 464 | //! rv = {class="RV", names=(...),}; // RV describing meaning of random variable |
---|
[477] | 465 | //! rvc= {class="RV", names=(...),}; // RV describing meaning of random variable in condition |
---|
| 466 | //! // elements of offsprings |
---|
| 467 | //! } |
---|
| 468 | //! \endcode |
---|
| 469 | //!@} |
---|
| 470 | void from_setting ( const Setting &set ); |
---|
| 471 | //!@} |
---|
[461] | 472 | |
---|
[477] | 473 | //! \name Connection to other objects |
---|
| 474 | //!@{ |
---|
| 475 | void set_rvc ( const RV &rvc0 ) { |
---|
| 476 | rvc = rvc0; |
---|
| 477 | } |
---|
| 478 | void set_rv ( const RV &rv0 ) { |
---|
[487] | 479 | ep->set_rv ( rv0 ); |
---|
[477] | 480 | } |
---|
| 481 | bool isnamed() { |
---|
[487] | 482 | return ( ep->isnamed() ) && ( dimc == rvc._dsize() ); |
---|
[477] | 483 | } |
---|
| 484 | //!@} |
---|
[270] | 485 | }; |
---|
[32] | 486 | |
---|
[487] | 487 | template <class EPDF> |
---|
| 488 | class mpdf_internal: public mpdf{ |
---|
| 489 | protected : |
---|
| 490 | EPDF iepdf; |
---|
| 491 | public: |
---|
| 492 | //! constructor |
---|
| 493 | mpdf_internal(): mpdf(),iepdf(){set_ep(iepdf);} |
---|
| 494 | //! Update \c iepdf so that it represents this mpdf conditioned on \c rvc = cond |
---|
| 495 | //! This function provides convenient reimplementation in offsprings |
---|
| 496 | virtual void condition ( const vec &cond ) { |
---|
| 497 | it_error ( "Not implemented" ); |
---|
| 498 | }; |
---|
| 499 | //!access function to iepdf |
---|
| 500 | EPDF& e(){return iepdf;} |
---|
| 501 | |
---|
| 502 | //! Reimplements samplecond using \c condition() |
---|
| 503 | vec samplecond ( const vec &cond ); |
---|
| 504 | //! Reimplements evallogcond using \c condition() |
---|
| 505 | double evallogcond ( const vec &val, const vec &cond ); |
---|
| 506 | //! Efficient version of evallogcond for matrices |
---|
| 507 | virtual vec evallogcond_m ( const mat &Dt, const vec &cond ); |
---|
| 508 | //! Efficient version of evallogcond for Array<vec> |
---|
| 509 | virtual vec evallogcond_m ( const Array<vec> &Dt, const vec &cond ); |
---|
| 510 | //! Efficient version of samplecond |
---|
| 511 | virtual mat samplecond_m ( const vec &cond, int N ); |
---|
| 512 | }; |
---|
| 513 | |
---|
[270] | 514 | /*! \brief DataLink is a connection between two data vectors Up and Down |
---|
[2] | 515 | |
---|
[270] | 516 | Up can be longer than Down. Down must be fully present in Up (TODO optional) |
---|
| 517 | See chart: |
---|
| 518 | \dot |
---|
| 519 | digraph datalink { |
---|
[377] | 520 | node [shape=record]; |
---|
| 521 | subgraph cluster0 { |
---|
| 522 | label = "Up"; |
---|
| 523 | up [label="<1>|<2>|<3>|<4>|<5>"]; |
---|
| 524 | color = "white" |
---|
[270] | 525 | } |
---|
[377] | 526 | subgraph cluster1{ |
---|
| 527 | label = "Down"; |
---|
| 528 | labelloc = b; |
---|
| 529 | down [label="<1>|<2>|<3>"]; |
---|
| 530 | color = "white" |
---|
[270] | 531 | } |
---|
| 532 | up:1 -> down:1; |
---|
| 533 | up:3 -> down:2; |
---|
| 534 | up:5 -> down:3; |
---|
| 535 | } |
---|
| 536 | \enddot |
---|
[263] | 537 | |
---|
[270] | 538 | */ |
---|
[477] | 539 | class datalink { |
---|
[270] | 540 | protected: |
---|
[477] | 541 | //! Remember how long val should be |
---|
| 542 | int downsize; |
---|
[424] | 543 | |
---|
[477] | 544 | //! Remember how long val of "Up" should be |
---|
| 545 | int upsize; |
---|
[424] | 546 | |
---|
[477] | 547 | //! val-to-val link, indices of the upper val |
---|
| 548 | ivec v2v_up; |
---|
[424] | 549 | |
---|
[270] | 550 | public: |
---|
[477] | 551 | //! Constructor |
---|
| 552 | datalink() : downsize ( 0 ), upsize ( 0 ) { } |
---|
| 553 | datalink ( const RV &rv, const RV &rv_up ) { |
---|
| 554 | set_connection ( rv, rv_up ); |
---|
| 555 | } |
---|
[424] | 556 | |
---|
[477] | 557 | //! set connection, rv must be fully present in rv_up |
---|
| 558 | void set_connection ( const RV &rv, const RV &rv_up ) { |
---|
| 559 | downsize = rv._dsize(); |
---|
| 560 | upsize = rv_up._dsize(); |
---|
| 561 | v2v_up = rv.dataind ( rv_up ); |
---|
[271] | 562 | |
---|
[477] | 563 | it_assert_debug ( v2v_up.length() == downsize, "rv is not fully in rv_up" ); |
---|
| 564 | } |
---|
[424] | 565 | |
---|
[477] | 566 | //! set connection using indices |
---|
| 567 | void set_connection ( int ds, int us, const ivec &upind ) { |
---|
| 568 | downsize = ds; |
---|
| 569 | upsize = us; |
---|
| 570 | v2v_up = upind; |
---|
[286] | 571 | |
---|
[477] | 572 | it_assert_debug ( v2v_up.length() == downsize, "rv is not fully in rv_up" ); |
---|
| 573 | } |
---|
[424] | 574 | |
---|
[477] | 575 | //! Get val for myself from val of "Up" |
---|
| 576 | vec pushdown ( const vec &val_up ) { |
---|
| 577 | it_assert_debug ( upsize == val_up.length(), "Wrong val_up" ); |
---|
| 578 | return get_vec ( val_up, v2v_up ); |
---|
| 579 | } |
---|
[424] | 580 | |
---|
[477] | 581 | //! Fill val of "Up" by my pieces |
---|
| 582 | void pushup ( vec &val_up, const vec &val ) { |
---|
| 583 | it_assert_debug ( downsize == val.length(), "Wrong val" ); |
---|
| 584 | it_assert_debug ( upsize == val_up.length(), "Wrong val_up" ); |
---|
| 585 | set_subvector ( val_up, v2v_up, val ); |
---|
| 586 | } |
---|
[270] | 587 | }; |
---|
[115] | 588 | |
---|
[424] | 589 | //! Data link with a condition. |
---|
[477] | 590 | class datalink_m2e: public datalink { |
---|
[270] | 591 | protected: |
---|
[477] | 592 | //! Remember how long cond should be |
---|
| 593 | int condsize; |
---|
[424] | 594 | |
---|
[477] | 595 | //!upper_val-to-local_cond link, indices of the upper val |
---|
| 596 | ivec v2c_up; |
---|
[424] | 597 | |
---|
[477] | 598 | //!upper_val-to-local_cond link, indices of the local cond |
---|
| 599 | ivec v2c_lo; |
---|
[192] | 600 | |
---|
[270] | 601 | public: |
---|
[477] | 602 | //! Constructor |
---|
| 603 | datalink_m2e() : condsize ( 0 ) { } |
---|
[424] | 604 | |
---|
[477] | 605 | void set_connection ( const RV &rv, const RV &rvc, const RV &rv_up ) { |
---|
| 606 | datalink::set_connection ( rv, rv_up ); |
---|
| 607 | condsize = rvc._dsize(); |
---|
| 608 | //establish v2c connection |
---|
| 609 | rvc.dataind ( rv_up, v2c_lo, v2c_up ); |
---|
| 610 | } |
---|
[424] | 611 | |
---|
[477] | 612 | //!Construct condition |
---|
| 613 | vec get_cond ( const vec &val_up ) { |
---|
| 614 | vec tmp ( condsize ); |
---|
| 615 | set_subvector ( tmp, v2c_lo, val_up ( v2c_up ) ); |
---|
| 616 | return tmp; |
---|
| 617 | } |
---|
[424] | 618 | |
---|
[477] | 619 | void pushup_cond ( vec &val_up, const vec &val, const vec &cond ) { |
---|
| 620 | it_assert_debug ( downsize == val.length(), "Wrong val" ); |
---|
| 621 | it_assert_debug ( upsize == val_up.length(), "Wrong val_up" ); |
---|
| 622 | set_subvector ( val_up, v2v_up, val ); |
---|
| 623 | set_subvector ( val_up, v2c_up, cond ); |
---|
| 624 | } |
---|
[270] | 625 | }; |
---|
[424] | 626 | |
---|
[192] | 627 | //!DataLink is a connection between mpdf and its superordinate (Up) |
---|
| 628 | //! This class links |
---|
[477] | 629 | class datalink_m2m: public datalink_m2e { |
---|
[270] | 630 | protected: |
---|
[477] | 631 | //!cond-to-cond link, indices of the upper cond |
---|
| 632 | ivec c2c_up; |
---|
| 633 | //!cond-to-cond link, indices of the local cond |
---|
| 634 | ivec c2c_lo; |
---|
[424] | 635 | |
---|
[270] | 636 | public: |
---|
[477] | 637 | //! Constructor |
---|
| 638 | datalink_m2m() {}; |
---|
| 639 | void set_connection ( const RV &rv, const RV &rvc, const RV &rv_up, const RV &rvc_up ) { |
---|
| 640 | datalink_m2e::set_connection ( rv, rvc, rv_up ); |
---|
| 641 | //establish c2c connection |
---|
| 642 | rvc.dataind ( rvc_up, c2c_lo, c2c_up ); |
---|
| 643 | it_assert_debug ( c2c_lo.length() + v2c_lo.length() == condsize, "cond is not fully given" ); |
---|
| 644 | } |
---|
[424] | 645 | |
---|
[477] | 646 | //! Get cond for myself from val and cond of "Up" |
---|
| 647 | vec get_cond ( const vec &val_up, const vec &cond_up ) { |
---|
| 648 | vec tmp ( condsize ); |
---|
| 649 | set_subvector ( tmp, v2c_lo, val_up ( v2c_up ) ); |
---|
| 650 | set_subvector ( tmp, c2c_lo, cond_up ( c2c_up ) ); |
---|
| 651 | return tmp; |
---|
| 652 | } |
---|
| 653 | //! Fill |
---|
[190] | 654 | |
---|
[270] | 655 | }; |
---|
[190] | 656 | |
---|
[270] | 657 | /*! |
---|
| 658 | @brief Class for storing results (and semi-results) of an experiment |
---|
[267] | 659 | |
---|
[270] | 660 | This class abstracts logging of results from implementation. This class replaces direct logging of results (e.g. to files or to global variables) by calling methods of a logger. Specializations of this abstract class for specific storage method are designed. |
---|
| 661 | */ |
---|
[477] | 662 | class logger : public root { |
---|
[270] | 663 | protected: |
---|
[477] | 664 | //! RVs of all logged variables. |
---|
| 665 | Array<RV> entries; |
---|
| 666 | //! Names of logged quantities, e.g. names of algorithm variants |
---|
| 667 | Array<string> names; |
---|
[270] | 668 | public: |
---|
[477] | 669 | //!Default constructor |
---|
| 670 | logger() : entries ( 0 ), names ( 0 ) {} |
---|
[267] | 671 | |
---|
[477] | 672 | //! returns an identifier which will be later needed for calling the \c logit() function |
---|
| 673 | //! For empty RV it returns -1, this entry will be ignored by \c logit(). |
---|
| 674 | virtual int add ( const RV &rv, string prefix = "" ) { |
---|
| 675 | int id; |
---|
| 676 | if ( rv._dsize() > 0 ) { |
---|
| 677 | id = entries.length(); |
---|
| 678 | names = concat ( names, prefix ); // diff |
---|
| 679 | entries.set_length ( id + 1, true ); |
---|
| 680 | entries ( id ) = rv; |
---|
| 681 | } else { |
---|
| 682 | id = -1; |
---|
| 683 | } |
---|
| 684 | return id; // identifier of the last entry |
---|
| 685 | } |
---|
[267] | 686 | |
---|
[477] | 687 | //! log this vector |
---|
| 688 | virtual void logit ( int id, const vec &v ) = 0; |
---|
| 689 | //! log this double |
---|
| 690 | virtual void logit ( int id, const double &d ) = 0; |
---|
[267] | 691 | |
---|
[477] | 692 | //! Shifts storage position for another time step. |
---|
| 693 | virtual void step() = 0; |
---|
[267] | 694 | |
---|
[477] | 695 | //! Finalize storing information |
---|
| 696 | virtual void finalize() {}; |
---|
[267] | 697 | |
---|
[477] | 698 | //! Initialize the storage |
---|
| 699 | virtual void init() {}; |
---|
[267] | 700 | |
---|
[270] | 701 | }; |
---|
[267] | 702 | |
---|
[270] | 703 | /*! \brief Unconditional mpdf, allows using epdf in the role of mpdf. |
---|
[190] | 704 | |
---|
[270] | 705 | */ |
---|
[384] | 706 | class mepdf : public mpdf { |
---|
[487] | 707 | |
---|
| 708 | shared_ptr<epdf> ipdf; |
---|
[270] | 709 | public: |
---|
[477] | 710 | //!Default constructor |
---|
| 711 | mepdf() { } |
---|
[461] | 712 | |
---|
[477] | 713 | mepdf ( shared_ptr<epdf> em ) { |
---|
[487] | 714 | ipdf = em; |
---|
| 715 | set_ep(*ipdf.get()); |
---|
[477] | 716 | dimc = 0; |
---|
| 717 | } |
---|
[461] | 718 | |
---|
[477] | 719 | //! empty |
---|
| 720 | void condition ( const vec &cond ); |
---|
[461] | 721 | |
---|
[477] | 722 | //! Load from structure with elements: |
---|
| 723 | //! \code |
---|
| 724 | //! { class = "mepdf", |
---|
[487] | 725 | //! epdf = {class="epdf_offspring",...} |
---|
[477] | 726 | //! } |
---|
| 727 | //! \endcode |
---|
| 728 | //!@} |
---|
| 729 | void from_setting ( const Setting &set ); |
---|
[270] | 730 | }; |
---|
[477] | 731 | UIREGISTER ( mepdf ); |
---|
[115] | 732 | |
---|
[477] | 733 | //!\brief Chain rule of pdfs - abstract part common for mprod and merger. |
---|
[192] | 734 | //!this abstract class is common to epdf and mpdf |
---|
[384] | 735 | //!\todo Think of better design - global functions rv=get_rv(Array<mpdf*> mpdfs); ?? |
---|
| 736 | class compositepdf { |
---|
[270] | 737 | protected: |
---|
[477] | 738 | //! Elements of composition |
---|
| 739 | Array<mpdf*> mpdfs; |
---|
| 740 | bool owning_mpdfs; |
---|
[270] | 741 | public: |
---|
[477] | 742 | compositepdf() : mpdfs ( 0 ) {}; |
---|
| 743 | compositepdf ( Array<mpdf*> A0, bool own = false ) { |
---|
| 744 | set_elements ( A0, own ); |
---|
| 745 | }; |
---|
| 746 | void set_elements ( Array<mpdf*> A0, bool own = false ) { |
---|
| 747 | mpdfs = A0; |
---|
| 748 | owning_mpdfs = own; |
---|
| 749 | }; |
---|
| 750 | //! find common rv, flag \param checkoverlap modifies whether overlaps are acceptable |
---|
| 751 | RV getrv ( bool checkoverlap = false ); |
---|
| 752 | //! common rvc of all mpdfs is written to rvc |
---|
| 753 | void setrvc ( const RV &rv, RV &rvc ); |
---|
| 754 | ~compositepdf() { |
---|
| 755 | if ( owning_mpdfs ) for ( int i = 0; i < mpdfs.length(); i++ ) { |
---|
| 756 | delete mpdfs ( i ); |
---|
| 757 | } |
---|
| 758 | }; |
---|
[270] | 759 | }; |
---|
[175] | 760 | |
---|
[270] | 761 | /*! \brief Abstract class for discrete-time sources of data. |
---|
[12] | 762 | |
---|
[270] | 763 | The class abstracts operations of: (i) data aquisition, (ii) data-preprocessing, (iii) scaling of data, and (iv) data resampling from the task of estimation and control. |
---|
| 764 | Moreover, for controlled systems, it is able to receive the desired control action and perform it in the next step. (Or as soon as possible). |
---|
[12] | 765 | |
---|
[270] | 766 | */ |
---|
[32] | 767 | |
---|
[477] | 768 | class DS : public root { |
---|
[270] | 769 | protected: |
---|
[477] | 770 | int dtsize; |
---|
| 771 | int utsize; |
---|
| 772 | //!Description of data returned by \c getdata(). |
---|
| 773 | RV Drv; |
---|
| 774 | //!Description of data witten by by \c write(). |
---|
| 775 | RV Urv; // |
---|
| 776 | //! Remember its own index in Logger L |
---|
| 777 | int L_dt, L_ut; |
---|
[270] | 778 | public: |
---|
[477] | 779 | //! default constructors |
---|
| 780 | DS() : Drv(), Urv() {}; |
---|
| 781 | //! Returns full vector of observed data=[output, input] |
---|
| 782 | virtual void getdata ( vec &dt ) { |
---|
| 783 | it_error ( "abstract class" ); |
---|
| 784 | }; |
---|
| 785 | //! Returns data records at indeces. |
---|
| 786 | virtual void getdata ( vec &dt, const ivec &indeces ) { |
---|
| 787 | it_error ( "abstract class" ); |
---|
| 788 | }; |
---|
| 789 | //! Accepts action variable and schedule it for application. |
---|
| 790 | virtual void write ( vec &ut ) { |
---|
| 791 | it_error ( "abstract class" ); |
---|
| 792 | }; |
---|
| 793 | //! Accepts action variables at specific indeces |
---|
| 794 | virtual void write ( vec &ut, const ivec &indeces ) { |
---|
| 795 | it_error ( "abstract class" ); |
---|
| 796 | }; |
---|
[32] | 797 | |
---|
[477] | 798 | //! Moves from \f$ t \f$ to \f$ t+1 \f$, i.e. perfroms the actions and reads response of the system. |
---|
| 799 | virtual void step() = 0; |
---|
[32] | 800 | |
---|
[477] | 801 | //! Register DS for logging into logger L |
---|
| 802 | virtual void log_add ( logger &L ) { |
---|
| 803 | it_assert_debug ( dtsize == Drv._dsize(), "" ); |
---|
| 804 | it_assert_debug ( utsize == Urv._dsize(), "" ); |
---|
[32] | 805 | |
---|
[477] | 806 | L_dt = L.add ( Drv, "" ); |
---|
| 807 | L_ut = L.add ( Urv, "" ); |
---|
| 808 | } |
---|
| 809 | //! Register DS for logging into logger L |
---|
| 810 | virtual void logit ( logger &L ) { |
---|
| 811 | vec tmp ( Drv._dsize() + Urv._dsize() ); |
---|
| 812 | getdata ( tmp ); |
---|
| 813 | // d is first in getdata |
---|
| 814 | L.logit ( L_dt, tmp.left ( Drv._dsize() ) ); |
---|
| 815 | // u follows after d in getdata |
---|
| 816 | L.logit ( L_ut, tmp.mid ( Drv._dsize(), Urv._dsize() ) ); |
---|
| 817 | } |
---|
| 818 | //!access function |
---|
| 819 | virtual RV _drv() const { |
---|
| 820 | return concat ( Drv, Urv ); |
---|
| 821 | } |
---|
| 822 | //!access function |
---|
| 823 | const RV& _urv() const { |
---|
| 824 | return Urv; |
---|
| 825 | } |
---|
| 826 | //! set random rvariables |
---|
| 827 | virtual void set_drv ( const RV &drv, const RV &urv ) { |
---|
| 828 | Drv = drv; |
---|
| 829 | Urv = urv; |
---|
| 830 | } |
---|
[270] | 831 | }; |
---|
[18] | 832 | |
---|
[270] | 833 | /*! \brief Bayesian Model of a system, i.e. all uncertainty is modeled by probabilities. |
---|
[32] | 834 | |
---|
[283] | 835 | This object represents exact or approximate evaluation of the Bayes rule: |
---|
| 836 | \f[ |
---|
| 837 | f(\theta_t | d_1,\ldots,d_t) = \frac{f(y_t|\theta_t,\cdot) f(\theta_t|d_1,\ldots,d_{t-1})}{f(y_t|d_1,\ldots,d_{t-1})} |
---|
| 838 | \f] |
---|
| 839 | |
---|
| 840 | Access to the resulting posterior density is via function \c posterior(). |
---|
| 841 | |
---|
| 842 | As a "side-effect" it also evaluates log-likelihood of the data, which can be accessed via function _ll(). |
---|
| 843 | It can also evaluate predictors of future values of \f$y_t\f$, see functions epredictor() and predictor(). |
---|
| 844 | |
---|
| 845 | Alternatively, it can evaluate posterior density conditioned by a known constant, \f$ c_t \f$: |
---|
| 846 | \f[ |
---|
| 847 | f(\theta_t | c_t, d_1,\ldots,d_t) \propto f(y_t,\theta_t|c_t,\cdot, d_1,\ldots,d_{t-1}) |
---|
| 848 | \f] |
---|
| 849 | |
---|
| 850 | The value of \f$ c_t \f$ is set by function condition(). |
---|
| 851 | |
---|
[270] | 852 | */ |
---|
[32] | 853 | |
---|
[477] | 854 | class BM : public root { |
---|
[270] | 855 | protected: |
---|
[477] | 856 | //! Random variable of the data (optional) |
---|
| 857 | RV drv; |
---|
| 858 | //!Logarithm of marginalized data likelihood. |
---|
| 859 | double ll; |
---|
| 860 | //! If true, the filter will compute likelihood of the data record and store it in \c ll . Set to false if you want to save computational time. |
---|
| 861 | bool evalll; |
---|
[270] | 862 | public: |
---|
[477] | 863 | //! \name Constructors |
---|
| 864 | //! @{ |
---|
[271] | 865 | |
---|
[477] | 866 | BM() : ll ( 0 ), evalll ( true ), LIDs ( 4 ), LFlags ( 4 ) { |
---|
| 867 | LIDs = -1;/*empty IDs*/ |
---|
| 868 | LFlags = 0; |
---|
| 869 | LFlags ( 0 ) = 1;/*log only mean*/ |
---|
| 870 | }; |
---|
| 871 | BM ( const BM &B ) : drv ( B.drv ), ll ( B.ll ), evalll ( B.evalll ) {} |
---|
| 872 | //! Copy function required in vectors, Arrays of BM etc. Have to be DELETED manually! |
---|
| 873 | //! Prototype: \code BM* _copy_() const {return new BM(*this);} \endcode |
---|
| 874 | virtual BM* _copy_() const { |
---|
| 875 | return NULL; |
---|
| 876 | }; |
---|
| 877 | //!@} |
---|
[18] | 878 | |
---|
[477] | 879 | //! \name Mathematical operations |
---|
| 880 | //!@{ |
---|
[271] | 881 | |
---|
[477] | 882 | /*! \brief Incremental Bayes rule |
---|
| 883 | @param dt vector of input data |
---|
| 884 | */ |
---|
| 885 | virtual void bayes ( const vec &dt ) = 0; |
---|
| 886 | //! Batch Bayes rule (columns of Dt are observations) |
---|
| 887 | virtual void bayesB ( const mat &Dt ); |
---|
| 888 | //! Evaluates predictive log-likelihood of the given data record |
---|
| 889 | //! I.e. marginal likelihood of the data with the posterior integrated out. |
---|
| 890 | virtual double logpred ( const vec &dt ) const { |
---|
| 891 | it_error ( "Not implemented" ); |
---|
| 892 | return 0.0; |
---|
| 893 | } |
---|
| 894 | //! Matrix version of logpred |
---|
| 895 | vec logpred_m ( const mat &dt ) const { |
---|
| 896 | vec tmp ( dt.cols() ); |
---|
| 897 | for ( int i = 0; i < dt.cols(); i++ ) { |
---|
| 898 | tmp ( i ) = logpred ( dt.get_col ( i ) ); |
---|
| 899 | } |
---|
| 900 | return tmp; |
---|
| 901 | } |
---|
[32] | 902 | |
---|
[477] | 903 | //!Constructs a predictive density \f$ f(d_{t+1} |d_{t}, \ldots d_{0}) \f$ |
---|
| 904 | virtual epdf* epredictor() const { |
---|
| 905 | it_error ( "Not implemented" ); |
---|
| 906 | return NULL; |
---|
| 907 | }; |
---|
| 908 | //!Constructs a conditional density 1-step ahead predictor \f$ f(d_{t+1} |d_{t+h-1}, \ldots d_{t}) |
---|
| 909 | virtual mpdf* predictor() const { |
---|
| 910 | it_error ( "Not implemented" ); |
---|
| 911 | return NULL; |
---|
| 912 | }; |
---|
| 913 | //!@} |
---|
[271] | 914 | |
---|
[477] | 915 | //! \name Extension to conditional BM |
---|
| 916 | //! This extension is useful e.g. in Marginalized Particle Filter (\ref bdm::MPF). |
---|
| 917 | //! Alternatively, it can be used for automated connection to DS when the condition is observed |
---|
| 918 | //!@{ |
---|
[283] | 919 | |
---|
[477] | 920 | //! Name of extension variable |
---|
| 921 | RV rvc; |
---|
| 922 | //! access function |
---|
| 923 | const RV& _rvc() const { |
---|
| 924 | return rvc; |
---|
| 925 | } |
---|
[283] | 926 | |
---|
[477] | 927 | //! Substitute \c val for \c rvc. |
---|
| 928 | virtual void condition ( const vec &val ) { |
---|
| 929 | it_error ( "Not implemented!" ); |
---|
| 930 | }; |
---|
[283] | 931 | |
---|
[477] | 932 | //!@} |
---|
[283] | 933 | |
---|
| 934 | |
---|
[477] | 935 | //! \name Access to attributes |
---|
| 936 | //!@{ |
---|
[271] | 937 | |
---|
[477] | 938 | const RV& _drv() const { |
---|
| 939 | return drv; |
---|
| 940 | } |
---|
| 941 | void set_drv ( const RV &rv ) { |
---|
| 942 | drv = rv; |
---|
| 943 | } |
---|
| 944 | void set_rv ( const RV &rv ) { |
---|
| 945 | const_cast<epdf&> ( posterior() ).set_rv ( rv ); |
---|
| 946 | } |
---|
| 947 | double _ll() const { |
---|
| 948 | return ll; |
---|
| 949 | } |
---|
| 950 | void set_evalll ( bool evl0 ) { |
---|
| 951 | evalll = evl0; |
---|
| 952 | } |
---|
| 953 | virtual const epdf& posterior() const = 0; |
---|
| 954 | virtual const epdf* _e() const = 0; |
---|
| 955 | //!@} |
---|
[28] | 956 | |
---|
[477] | 957 | //! \name Logging of results |
---|
| 958 | //!@{ |
---|
[200] | 959 | |
---|
[477] | 960 | //! Set boolean options from a string, recognized are: "logbounds,logll" |
---|
| 961 | virtual void set_options ( const string &opt ) { |
---|
| 962 | LFlags ( 0 ) = 1; |
---|
| 963 | if ( opt.find ( "logbounds" ) != string::npos ) { |
---|
| 964 | LFlags ( 1 ) = 1; |
---|
| 965 | LFlags ( 2 ) = 1; |
---|
| 966 | } |
---|
| 967 | if ( opt.find ( "logll" ) != string::npos ) { |
---|
| 968 | LFlags ( 3 ) = 1; |
---|
| 969 | } |
---|
| 970 | } |
---|
| 971 | //! IDs of storages in loggers 4:[1=mean,2=lb,3=ub,4=ll] |
---|
| 972 | ivec LIDs; |
---|
[190] | 973 | |
---|
[477] | 974 | //! Flags for logging - same size as LIDs, each entry correspond to the same in LIDs |
---|
| 975 | ivec LFlags; |
---|
| 976 | //! Add all logged variables to a logger |
---|
| 977 | virtual void log_add ( logger &L, const string &name = "" ) { |
---|
| 978 | // internal |
---|
| 979 | RV r; |
---|
| 980 | if ( posterior().isnamed() ) { |
---|
| 981 | r = posterior()._rv(); |
---|
| 982 | } else { |
---|
| 983 | r = RV ( "est", posterior().dimension() ); |
---|
| 984 | }; |
---|
[190] | 985 | |
---|
[477] | 986 | // Add mean value |
---|
| 987 | if ( LFlags ( 0 ) ) LIDs ( 0 ) = L.add ( r, name + "mean_" ); |
---|
| 988 | if ( LFlags ( 1 ) ) LIDs ( 1 ) = L.add ( r, name + "lb_" ); |
---|
| 989 | if ( LFlags ( 2 ) ) LIDs ( 2 ) = L.add ( r, name + "ub_" ); |
---|
| 990 | if ( LFlags ( 3 ) ) LIDs ( 3 ) = L.add ( RV ( "ll", 1 ), name ); //TODO: "local" RV |
---|
| 991 | } |
---|
| 992 | virtual void logit ( logger &L ) { |
---|
| 993 | L.logit ( LIDs ( 0 ), posterior().mean() ); |
---|
| 994 | if ( LFlags ( 1 ) || LFlags ( 2 ) ) { //if one of them is off, its LID==-1 and will not be stored |
---|
| 995 | vec ub, lb; |
---|
| 996 | posterior().qbounds ( lb, ub ); |
---|
| 997 | L.logit ( LIDs ( 1 ), lb ); |
---|
| 998 | L.logit ( LIDs ( 2 ), ub ); |
---|
| 999 | } |
---|
| 1000 | if ( LFlags ( 3 ) ) L.logit ( LIDs ( 3 ), ll ); |
---|
| 1001 | } |
---|
| 1002 | //!@} |
---|
[270] | 1003 | }; |
---|
[32] | 1004 | |
---|
[339] | 1005 | |
---|
[254] | 1006 | }; //namespace |
---|
[384] | 1007 | #endif // BDMBASE_H |
---|