| 32 | | |
| 33 | | //! Structure of RV (used internally), i.e. expanded RVs |
| 34 | | class str { |
| 35 | | public: |
| 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 | | }; |
| | 52 | }; |
| | 53 | |
| | 54 | /*! |
| | 55 | * \brief Class representing variables, most often random variables |
| | 56 | |
| | 57 | 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. |
| | 58 | |
| | 59 | The class is implemented using global variables to assure uniqueness of description: |
| | 60 | |
| | 61 | In is a vector |
| | 62 | \dot |
| | 63 | digraph datalink { |
| | 64 | rankdir=LR; |
| | 65 | subgraph cluster0 { |
| | 66 | node [shape=record]; |
| | 67 | label = "RV_MAP \n std::map<string,int>"; |
| | 68 | map [label="{{\"a\"| \"b\" | \"c\"} | {<3> 3 |<1> 1|<2> 2}}"]; |
| | 69 | color = "white" |
| | 70 | } |
| | 71 | subgraph cluster1{ |
| | 72 | node [shape=record]; |
| | 73 | label = "RV_NAMES"; |
| | 74 | names [label="{<1> \"b\" | <2> \"c\" | <3>\"a\" }"]; |
| | 75 | color = "white" |
| | 76 | } |
| | 77 | subgraph cluster2{ |
| | 78 | node [shape=record]; |
| | 79 | label = "RV_SIZES"; |
| | 80 | labelloc = b; |
| | 81 | sizes [label="{<1>1 |<2> 4 |<3> 1}"]; |
| | 82 | color = "white" |
| | 83 | } |
| | 84 | map:1 -> names:1; |
| | 85 | map:1 -> sizes:1; |
| | 86 | map:3 -> names:3; |
| | 87 | map:3 -> sizes:3; |
| | 88 | } |
| | 89 | \enddot |
| | 90 | */ |
| | 91 | |
| | 92 | class RV :public bdmroot { |
| | 93 | protected: |
| | 94 | //! size of the data vector |
| | 95 | int dsize; |
| | 96 | //! number of individual rvs |
| | 97 | int len; |
| | 98 | //! Vector of unique IDs |
| | 99 | ivec ids; |
| | 100 | //! Vector of shifts from current time |
| | 101 | ivec times; |
| | 102 | |
| | 103 | private: |
| | 104 | //! auxiliary function used in constructor |
| | 105 | void init ( Array<std::string> in_names, ivec in_sizes, ivec in_times ); |
| | 106 | int init ( const string &name, int size ); |
| | 107 | public: |
| | 108 | //! \name Constructors |
| | 109 | //!@{ |
| | 110 | |
| | 111 | //! Full constructor |
| | 112 | RV ( Array<std::string> in_names, ivec in_sizes, ivec in_times ) {init ( in_names,in_sizes,in_times );}; |
| | 113 | //! Constructor with times=0 |
| | 114 | RV ( Array<std::string> in_names, ivec in_sizes ) {init ( in_names,in_sizes,zeros_i ( in_names.length() ) );}; |
| | 115 | //! Constructor with sizes=1, times=0 |
| | 116 | RV ( Array<std::string> in_names ) {init ( in_names,ones_i ( in_names.length() ),zeros_i ( in_names.length() ) );} |
| | 117 | //! Constructor of empty RV |
| | 118 | RV () :dsize ( 0 ),len ( 0 ),ids ( 0 ),times ( 0 ) {}; |
| | 119 | //! Constructor of a single RV with given id |
| | 120 | RV ( string name, int sz, int tm=0 ); |
| | 121 | //!@} |
| | 122 | |
| | 123 | //! \name Access functions |
| | 124 | //!@{ |
| | 125 | |
| | 126 | //! Printing output e.g. for debugging. |
| | 127 | friend std::ostream &operator<< ( std::ostream &os, const RV &rv ); |
| | 128 | int _dsize() const {return dsize;} ; |
| | 129 | //! Recount size of the corresponding data vector |
| | 130 | int countsize() const; |
| | 131 | int length() const {return len;} ; |
| | 132 | int id ( int at ) const{return ids ( at );}; |
| | 133 | int size ( int at ) const {return RV_SIZES ( at );}; |
| | 134 | int time ( int at ) const{return times ( at );}; |
| | 135 | std::string name ( int at ) const {return RV_NAMES ( at );}; |
| | 136 | void set_time ( int at, int time0 ) {times ( at ) =time0;}; |
| | 137 | //!@} |
| | 138 | |
| | 139 | //TODO why not inline and later?? |
| | 140 | |
| | 141 | //! \name Algebra on Random Variables |
| | 142 | //!@{ |
| | 143 | |
| | 144 | //! Find indices of self in another rv, \return ivec of the same size as self. |
| | 145 | ivec findself ( const RV &rv2 ) const; |
| | 146 | //! Compare if \c rv2 is identical to this \c RV |
| | 147 | bool equal ( const RV &rv2 ) const; |
| | 148 | //! Add (concat) another variable to the current one, \return true if all rv2 were added, false if rv2 is in conflict |
| | 149 | bool add ( const RV &rv2 ); |
| | 150 | //! Subtract another variable from the current one |
| | 151 | RV subt ( const RV &rv2 ) const; |
| | 152 | //! Select only variables at indeces ind |
| | 153 | RV subselect ( const ivec &ind ) const; |
| | 154 | //! Select only variables at indeces ind |
| | 155 | RV operator() ( const ivec &ind ) const {return subselect ( ind );}; |
| | 156 | //! Shift \c time shifted by delta. |
| | 157 | void t ( int delta ); |
| | 158 | //!@} |
| | 159 | |
| | 160 | //!\name Relation to vectors |
| | 161 | //!@{ |
| | 162 | |
| | 163 | //! generate \c str from rv, by expanding sizes |
| | 164 | str tostr() const; |
| | 165 | //! when this rv is a part of bigger rv, this function returns indeces of self in the data vector of the bigger crv. |
| | 166 | //! Then, data can be copied via: data_of_this = cdata(ind); |
| | 167 | ivec dataind ( const RV &crv ) const; |
| | 168 | //! generate mutual indeces when copying data betwenn self and crv. |
| | 169 | //! Data are copied via: data_of_this(selfi) = data_of_rv2(rv2i) |
| | 170 | void dataind ( const RV &rv2, ivec &selfi, ivec &rv2i ) const; |
| | 171 | //! Minimum time-offset |
| | 172 | int mint () const {return min ( times );}; |
| | 173 | //!@} |
| | 174 | |
| | 175 | }; |
| | 176 | |
| | 177 | |
| | 178 | //! Concat two random variables |
| | 179 | RV concat ( const RV &rv1, const RV &rv2 ); |
| | 180 | |
| | 181 | //!Default empty RV that can be used as default argument |
| | 182 | extern RV RV0; |
| | 183 | |
| | 184 | //! Class representing function \f$f(x)\f$ of variable \f$x\f$ represented by \c rv |
| | 185 | |
| | 186 | class fnc :public bdmroot { |
| | 187 | protected: |
| | 188 | //! Length of the output vector |
| | 189 | int dimy; |
| | 190 | public: |
| | 191 | //!default constructor |
| | 192 | fnc ( ) {}; |
| | 193 | //! function evaluates numerical value of \f$f(x)\f$ at \f$x=\f$ \c cond |
| | 194 | virtual vec eval ( const vec &cond ) { |
| | 195 | return vec ( 0 ); |
| 46 | | /*! |
| 47 | | * \brief Class representing variables, most often random variables |
| 48 | | |
| 49 | | * More?... |
| 50 | | */ |
| 51 | | |
| 52 | | class RV :public bdmroot { |
| 53 | | protected: |
| 54 | | //! size = sum of sizes |
| 55 | | int tsize; |
| 56 | | //! len = number of individual rvs |
| 57 | | int len; |
| 58 | | //! Vector of unique IDs |
| 59 | | ivec ids; |
| 60 | | //! Vector of sizes |
| 61 | | ivec sizes; |
| 62 | | //! Vector of shifts from current time |
| 63 | | ivec times; |
| 64 | | //! Array of names |
| 65 | | Array<std::string> names; |
| 66 | | |
| 67 | | private: |
| 68 | | //! auxiliary function used in constructor |
| 69 | | void init ( ivec in_ids, Array<std::string> in_names, ivec in_sizes, ivec in_times ); |
| 70 | | public: |
| 71 | | //! Full constructor |
| 72 | | RV ( Array<std::string> in_names, ivec in_sizes, ivec in_times ); |
| 73 | | //! Constructor with times=0 |
| 74 | | RV ( Array<std::string> in_names, ivec in_sizes ); |
| 75 | | //! Constructor with sizes=1, times=0 |
| 76 | | RV ( Array<std::string> in_names ); |
| 77 | | //! Constructor of empty RV |
| 78 | | RV (); |
| 79 | | //! Constructor of a single RV with given id |
| 80 | | RV (string name, int id, int sz=1, int tm=0); |
| 81 | | |
| 82 | | //! Printing output e.g. for debugging. |
| 83 | | friend std::ostream &operator<< ( std::ostream &os, const RV &rv ); |
| 84 | | |
| 85 | | //! Return number of scalars in the RV. |
| 86 | | int count() const {return tsize;} ; |
| 87 | | //! Return length (number of entries) of the RV. |
| 88 | | int length() const {return len;} ; |
| 89 | | |
| 90 | | //TODO why not inline and later?? |
| 91 | | |
| 92 | | //! Find indexes of self in another rv, \return ivec of the same size as self. |
| 93 | | ivec findself ( const RV &rv2 ) const; |
| 94 | | //! Compare if \c rv2 is identical to this \c RV |
| 95 | | bool equal ( const RV &rv2 ) const; |
| 96 | | //! Add (concat) another variable to the current one, \return true if all rv2 were added, false if rv2 is in conflict |
| 97 | | bool add ( const RV &rv2 ); |
| 98 | | //! Subtract another variable from the current one |
| 99 | | RV subt ( const RV &rv2 ) const; |
| 100 | | //! Select only variables at indeces ind |
| 101 | | RV subselect ( const ivec &ind ) const; |
| 102 | | //! Select only variables at indeces ind |
| 103 | | RV operator() ( const ivec &ind ) const; |
| 104 | | //! Shift \c time shifted by delta. |
| 105 | | void t ( int delta ); |
| 106 | | //! generate \c str from rv, by expanding sizes |
| 107 | | str tostr() const; |
| 108 | | //! when this rv is a part of bigger rv, this function returns indeces of self in the data vector of the bigger crv. |
| 109 | | //! Then, data can be copied via: data_of_this = cdata(ind); |
| 110 | | ivec dataind ( const RV &crv ) const; |
| 111 | | //! generate mutual indeces when copying data betwenn self and crv. |
| 112 | | //! Data are copied via: data_of_this(selfi) = data_of_rv2(rv2i) |
| 113 | | void dataind ( const RV &rv2, ivec &selfi, ivec &rv2i ) const; |
| 114 | | //! Minimum time-offset |
| 115 | | int mint () const {return min ( times );}; |
| 116 | | |
| 117 | | //!access function |
| 118 | | Array<std::string>& _names() {return names;}; |
| 119 | | |
| 120 | | //!access function |
| 121 | | int id ( int at ) {return ids ( at );}; |
| 122 | | //!access function |
| 123 | | int size ( int at ) {return sizes ( at );}; |
| 124 | | //!access function |
| 125 | | int time ( int at ) {return times ( at );}; |
| 126 | | //!access function |
| 127 | | std::string name ( int at ) {return names ( at );}; |
| 128 | | |
| 129 | | //!access function |
| 130 | | void set_id ( int at, int id0 ) {ids ( at ) =id0;}; |
| 131 | | //!access function |
| 132 | | void set_size ( int at, int size0 ) {sizes ( at ) =size0; tsize=sum ( sizes );}; |
| 133 | | //!access function |
| 134 | | void set_time ( int at, int time0 ) {times ( at ) =time0;}; |
| 135 | | |
| 136 | | //!Assign unused ids to this rv |
| 137 | | void newids(); |
| 138 | | }; |
| 139 | | |
| 140 | | //! Concat two random variables |
| 141 | | RV concat ( const RV &rv1, const RV &rv2 ); |
| 142 | | |
| 143 | | //!Default empty RV that can be used as default argument |
| 144 | | extern RV RV0; |
| 145 | | |
| 146 | | //! Class representing function \f$f(x)\f$ of variable \f$x\f$ represented by \c rv |
| 147 | | |
| 148 | | class fnc :public bdmroot { |
| 149 | | protected: |
| 150 | | //! Length of the output vector |
| 151 | | int dimy; |
| 152 | | public: |
| 153 | | //!default constructor |
| 154 | | fnc ( int dy ) :dimy ( dy ) {}; |
| 155 | | //! function evaluates numerical value of \f$f(x)\f$ at \f$x=\f$ \c cond |
| 156 | | virtual vec eval ( const vec &cond ) { |
| 157 | | return vec ( 0 ); |
| 158 | | }; |
| 159 | | |
| 160 | | //! function substitutes given value into an appropriate position |
| 161 | | virtual void condition ( const vec &val ) {}; |
| 162 | | |
| 163 | | //! access function |
| 164 | | int _dimy() const{return dimy;} |
| 165 | | |
| 166 | | //! Destructor for future use; |
| 167 | | virtual ~fnc() {}; |
| 168 | | }; |
| 169 | | |
| 170 | | class mpdf; |
| | 198 | //! function substitutes given value into an appropriate position |
| | 199 | virtual void condition ( const vec &val ) {}; |
| | 200 | |
| | 201 | //! access function |
| | 202 | int _dimy() const{return dimy;} |
| | 203 | }; |
| | 204 | |
| | 205 | class mpdf; |
| 174 | | class epdf :public bdmroot { |
| 175 | | protected: |
| 176 | | //! Identified of the random variable |
| 177 | | RV rv; |
| 178 | | public: |
| 179 | | //!default constructor |
| 180 | | epdf() :rv ( ) {}; |
| 181 | | |
| 182 | | //!default constructor |
| 183 | | epdf ( const RV &rv0 ) :rv ( rv0 ) {}; |
| 184 | | |
| 185 | | // //! Returns the required moment of the epdf |
| 186 | | // virtual vec moment ( const int order = 1 ); |
| 187 | | |
| 188 | | //! Returns a sample, \f$x\f$ from density \f$epdf(rv)\f$ |
| 189 | | virtual vec sample () const =0; |
| 190 | | //! Returns N samples from density \f$epdf(rv)\f$ |
| 191 | | virtual mat sample_m ( int N ) const; |
| 192 | | |
| 193 | | //! Compute log-probability of argument \c val |
| 194 | | virtual double evallog ( const vec &val ) const =0; |
| 195 | | |
| 196 | | //! Compute log-probability of multiple values argument \c val |
| 197 | | virtual vec evallog_m ( const mat &Val ) const { |
| 198 | | vec x ( Val.cols() ); |
| 199 | | for ( int i=0;i<Val.cols();i++ ) {x ( i ) =evallog ( Val.get_col ( i ) ) ;} |
| 200 | | return x; |
| 201 | | } |
| 202 | | //! Return conditional density on the given RV, the remaining rvs will be in conditioning |
| 203 | | virtual mpdf* condition ( const RV &rv ) const {it_warning ( "Not implemented" ); return NULL;} |
| 204 | | //! Return marginal density on the given RV, the remainig rvs are intergrated out |
| 205 | | virtual epdf* marginal ( const RV &rv ) const {it_warning ( "Not implemented" ); return NULL;} |
| 206 | | |
| 207 | | //! return expected value |
| 208 | | virtual vec mean() const =0; |
| 209 | | |
| 210 | | //! return expected variance (not covariance!) |
| 211 | | virtual vec variance() const = 0; |
| 212 | | |
| 213 | | //! Destructor for future use; |
| 214 | | virtual ~epdf() {}; |
| 215 | | //! access function, possibly dangerous! |
| 216 | | const RV& _rv() const {return rv;} |
| 217 | | //! modifier function - useful when copying epdfs |
| 218 | | void _renewrv ( const RV &in_rv ) {rv=in_rv;} |
| 219 | | //! |
| 220 | | }; |
| | 209 | class epdf :public bdmroot { |
| | 210 | protected: |
| | 211 | //! dimension of the random variable |
| | 212 | int dim; |
| | 213 | //! Description of the random variable |
| | 214 | RV rv; |
| | 215 | |
| | 216 | public: |
| | 217 | /*! \name Constructors |
| | 218 | Construction of each epdf should support two types of constructors: |
| | 219 | \li empty constructor, |
| | 220 | \li copy constructor, |
| | 221 | |
| | 222 | The following constructors should be supported for convenience: |
| | 223 | \li constructor followed by calling \c set_parameters() |
| | 224 | \li constructor accepting random variables calling \c set_rv() |
| | 225 | |
| | 226 | 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. |
| | 227 | @{*/ |
| | 228 | epdf() :dim(0),rv ( ) {}; |
| | 229 | epdf(const epdf &e) :dim(e.dim),rv (e.rv) {}; |
| | 230 | epdf(const RV &rv0) {set_rv(rv0);}; |
| | 231 | void set_parameters(int dim0){dim=dim0;} |
| | 232 | //!@} |
| | 233 | |
| | 234 | //! \name Matematical Operations |
| | 235 | //!@{ |
| | 236 | |
| | 237 | //! Returns a sample, \f$ x \f$ from density \f$ f_x()\f$ |
| | 238 | virtual vec sample () const {it_error("not implemneted");return vec(0);}; |
| | 239 | //! Returns N samples, \f$ [x_1 , x_2 , \ldots \ \f$ from density \f$ f_x(rv)\f$ |
| | 240 | virtual mat sample_m ( int N ) const; |
| | 241 | //! Compute log-probability of argument \c val |
| | 242 | virtual double evallog ( const vec &val ) const {it_error("not implemneted");return 0.0;}; |
| | 243 | //! Compute log-probability of multiple values argument \c val |
| | 244 | virtual vec evallog_m ( const mat &Val ) const { |
| | 245 | vec x ( Val.cols() ); |
| | 246 | for ( int i=0;i<Val.cols();i++ ) {x ( i ) =evallog ( Val.get_col ( i ) ) ;} |
| | 247 | return x; |
| | 248 | } |
| | 249 | //! Return conditional density on the given RV, the remaining rvs will be in conditioning |
| | 250 | virtual mpdf* condition ( const RV &rv ) const {it_warning ( "Not implemented" ); return NULL;} |
| | 251 | //! Return marginal density on the given RV, the remainig rvs are intergrated out |
| | 252 | virtual epdf* marginal ( const RV &rv ) const {it_warning ( "Not implemented" ); return NULL;} |
| | 253 | //! return expected value |
| | 254 | virtual vec mean() const {it_error("not implemneted");return vec(0);}; |
| | 255 | //! return expected variance (not covariance!) |
| | 256 | virtual vec variance() const {it_error("not implemneted");return vec(0);}; |
| | 257 | //!@} |
| | 258 | |
| | 259 | //! \name Connection to other classes |
| | 260 | //! Description of the random quantity via attribute \c rv is optional. |
| | 261 | //! For operations such as sampling \c rv does not need to be set. However, for \c marginalization |
| | 262 | //! and \c conditioning \c rv has to be set. NB: |
| | 263 | //! @{ |
| | 264 | |
| | 265 | //!Name its rv |
| | 266 | void set_rv ( const RV &rv0 ) {rv = rv0; }//it_assert_debug(isnamed(),""); }; |
| | 267 | //! True if rv is assigned |
| | 268 | bool isnamed() const {bool b=( dim==rv._dsize() );return b;} |
| | 269 | //! Return name (fails when isnamed is false) |
| | 270 | const RV& _rv() const {it_assert_debug ( isnamed(),"" ); return rv;} |
| | 271 | //!@} |
| | 272 | |
| | 273 | //! \name Access to attributes |
| | 274 | //! @{ |
| | 275 | |
| | 276 | //! Size of the random variable |
| | 277 | int dimension() const {return dim;} |
| | 278 | //!@} |
| | 279 | |
| | 280 | }; |
| 226 | | class mpdf : public bdmroot { |
| 227 | | protected: |
| 228 | | //! modeled random variable |
| 229 | | RV rv; |
| 230 | | //! random variable in condition |
| 231 | | RV rvc; |
| 232 | | //! pointer to internal epdf |
| 233 | | epdf* ep; |
| 234 | | public: |
| 235 | | |
| 236 | | //! 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 |
| 237 | | virtual vec samplecond ( const vec &cond) { |
| 238 | | this->condition ( cond ); |
| 239 | | vec temp= ep->sample(); |
| 240 | | return temp; |
| 241 | | }; |
| 242 | | //! 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 \param ll is a return value of log-likelihood of the sample. |
| 243 | | virtual mat samplecond_m ( const vec &cond, vec &ll, int N ) { |
| 244 | | this->condition ( cond ); |
| 245 | | mat temp ( rv.count(),N ); vec smp ( rv.count() ); |
| 246 | | for ( int i=0;i<N;i++ ) {smp=ep->sample() ;temp.set_col ( i, smp );ll ( i ) =ep->evallog ( smp );} |
| 247 | | return temp; |
| 248 | | }; |
| 249 | | //! Update \c ep so that it represents this mpdf conditioned on \c rvc = cond |
| 250 | | virtual void condition ( const vec &cond ) {it_error ( "Not implemented" );}; |
| 251 | | |
| 252 | | //! Shortcut for conditioning and evaluation of the internal epdf. In some cases, this operation can be implemented efficiently. |
| 253 | | virtual double evallogcond ( const vec &dt, const vec &cond ) { |
| 254 | | double tmp; this->condition ( cond );tmp = ep->evallog ( dt ); it_assert_debug ( std::isfinite ( tmp ),"Infinite value" ); return tmp; |
| 255 | | }; |
| 256 | | |
| 257 | | //! Matrix version of evallogcond |
| 258 | | virtual vec evallogcond_m ( const mat &Dt, const vec &cond ) {this->condition ( cond );return ep->evallog_m ( Dt );}; |
| 259 | | |
| 260 | | //! Destructor for future use; |
| 261 | | virtual ~mpdf() {}; |
| 262 | | |
| 263 | | //! Default constructor |
| 264 | | mpdf ( const RV &rv0, const RV &rvc0 ) :rv ( rv0 ),rvc ( rvc0 ) {}; |
| 265 | | //! access function |
| 266 | | RV _rvc() const {return rvc;} |
| 267 | | //! access function |
| 268 | | RV _rv() const {return rv;} |
| 269 | | //!access function |
| 270 | | epdf& _epdf() {return *ep;} |
| 271 | | //!access function |
| 272 | | epdf* _e() {return ep;} |
| | 286 | class mpdf : public bdmroot { |
| | 287 | protected: |
| | 288 | //!dimension of the condition |
| | 289 | int dimc; |
| | 290 | //! random variable in condition |
| | 291 | RV rvc; |
| | 292 | //! pointer to internal epdf |
| | 293 | epdf* ep; |
| | 294 | public: |
| | 295 | //! \name Constructors |
| | 296 | //! @{ |
| | 297 | |
| | 298 | mpdf ( ) :dimc(0),rvc ( ) {}; |
| | 299 | //! copy constructor does not set pointer \c ep - has to be done in offsprings! |
| | 300 | mpdf (const mpdf &m ) :dimc(m.dimc),rvc (m.rvc ) {}; |
| | 301 | //!@} |
| | 302 | |
| | 303 | //! \name Matematical operations |
| | 304 | //!@{ |
| | 305 | |
| | 306 | //! 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 |
| | 307 | virtual vec samplecond ( const vec &cond ) { |
| | 308 | this->condition ( cond ); |
| | 309 | vec temp= ep->sample(); |
| | 310 | return temp; |
| 274 | | |
| 275 | | /*! \brief DataLink is a connection between two data vectors Up and Down |
| 276 | | |
| 277 | | Up can be longer than Down. Down must be fully present in Up (TODO optional) |
| 278 | | See chart: |
| 279 | | \dot |
| 280 | | digraph datalink { |
| 281 | | node [shape=record]; |
| 282 | | subgraph cluster0 { |
| 283 | | label = "Up"; |
| 284 | | up [label="<1>|<2>|<3>|<4>|<5>"]; |
| 285 | | color = "white" |
| 286 | | } |
| 287 | | subgraph cluster1{ |
| 288 | | label = "Down"; |
| 289 | | labelloc = b; |
| 290 | | down [label="<1>|<2>|<3>"]; |
| 291 | | color = "white" |
| 292 | | } |
| 293 | | up:1 -> down:1; |
| 294 | | up:3 -> down:2; |
| 295 | | up:5 -> down:3; |
| 296 | | } |
| 297 | | \enddot |
| 298 | | |
| 299 | | */ |
| 300 | | class datalink_e2e { |
| 301 | | protected: |
| 302 | | //! Remember how long val should be |
| 303 | | int valsize; |
| 304 | | //! Remember how long val of "Up" should be |
| 305 | | int valupsize; |
| 306 | | //! val-to-val link, indeces of the upper val |
| 307 | | ivec v2v_up; |
| 308 | | public: |
| 309 | | //! Constructor |
| 310 | | datalink_e2e ( const RV &rv, const RV &rv_up ) : |
| 311 | | valsize ( rv.count() ), valupsize ( rv_up.count() ), v2v_up ( rv.dataind ( rv_up ) ) { |
| 312 | | it_assert_debug ( v2v_up.length() ==valsize,"rv is not fully in rv_up" ); |
| 313 | | } |
| 314 | | //! Get val for myself from val of "Up" |
| 315 | | vec get_val ( const vec &val_up ) {it_assert_debug ( valupsize==val_up.length(),"Wrong val_up" ); return get_vec ( val_up,v2v_up );} |
| 316 | | //! Fill val of "Up" by my pieces |
| 317 | | void fill_val ( vec &val_up, const vec &val ) { |
| 318 | | it_assert_debug ( valsize==val.length(),"Wrong val" ); |
| 319 | | it_assert_debug ( valupsize==val_up.length(),"Wrong val_up" ); |
| 320 | | set_subvector ( val_up, v2v_up, val ); |
| 321 | | } |
| | 312 | //! 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 \param ll is a return value of log-likelihood of the sample. |
| | 313 | virtual mat samplecond_m ( const vec &cond, int N ) { |
| | 314 | this->condition ( cond ); |
| | 315 | mat temp ( ep->dimension(),N ); vec smp ( ep->dimension() ); |
| | 316 | for ( int i=0;i<N;i++ ) {smp=ep->sample() ;temp.set_col ( i, smp );} |
| | 317 | return temp; |
| | 319 | //! Update \c ep so that it represents this mpdf conditioned on \c rvc = cond |
| | 320 | virtual void condition ( const vec &cond ) {it_error ( "Not implemented" );}; |
| | 321 | |
| | 322 | //! Shortcut for conditioning and evaluation of the internal epdf. In some cases, this operation can be implemented efficiently. |
| | 323 | virtual double evallogcond ( const vec &dt, const vec &cond ) { |
| | 324 | double tmp; this->condition ( cond );tmp = ep->evallog ( dt ); it_assert_debug ( std::isfinite ( tmp ),"Infinite value" ); return tmp; |
| | 325 | }; |
| | 326 | |
| | 327 | //! Matrix version of evallogcond |
| | 328 | virtual vec evallogcond_m ( const mat &Dt, const vec &cond ) {this->condition ( cond );return ep->evallog_m ( Dt );}; |
| | 329 | |
| | 330 | //! \name Access to attributes |
| | 331 | //! @{ |
| | 332 | |
| | 333 | RV _rv() {return ep->_rv();} |
| | 334 | RV _rvc() {it_assert_debug ( isnamed(),"" ); return rvc;} |
| | 335 | int dimension() {return ep->dimension();} |
| | 336 | int dimensionc() {return dimc;} |
| | 337 | epdf& _epdf() {return *ep;} |
| | 338 | epdf* _e() {return ep;} |
| | 339 | //!@} |
| | 340 | |
| | 341 | //! \name Connection to other objects |
| | 342 | //!@{ |
| | 343 | void set_rvc ( const RV &rvc0 ) {rvc=rvc0;} |
| | 344 | void set_rv ( const RV &rv0 ) {ep->set_rv(rv0);} |
| | 345 | bool isnamed() {return (ep->isnamed())&&(dimc==rvc._dsize());} |
| | 346 | //!@} |
| | 347 | }; |
| | 348 | |
| | 349 | /*! \brief DataLink is a connection between two data vectors Up and Down |
| | 350 | |
| | 351 | Up can be longer than Down. Down must be fully present in Up (TODO optional) |
| | 352 | See chart: |
| | 353 | \dot |
| | 354 | digraph datalink { |
| | 355 | node [shape=record]; |
| | 356 | subgraph cluster0 { |
| | 357 | label = "Up"; |
| | 358 | up [label="<1>|<2>|<3>|<4>|<5>"]; |
| | 359 | color = "white" |
| | 360 | } |
| | 361 | subgraph cluster1{ |
| | 362 | label = "Down"; |
| | 363 | labelloc = b; |
| | 364 | down [label="<1>|<2>|<3>"]; |
| | 365 | color = "white" |
| | 366 | } |
| | 367 | up:1 -> down:1; |
| | 368 | up:3 -> down:2; |
| | 369 | up:5 -> down:3; |
| | 370 | } |
| | 371 | \enddot |
| | 372 | |
| | 373 | */ |
| | 374 | class datalink { |
| | 375 | protected: |
| | 376 | //! Remember how long val should be |
| | 377 | int downsize; |
| | 378 | //! Remember how long val of "Up" should be |
| | 379 | int upsize; |
| | 380 | //! val-to-val link, indeces of the upper val |
| | 381 | ivec v2v_up; |
| | 382 | public: |
| | 383 | //! Constructor |
| | 384 | datalink ( const RV &rv, const RV &rv_up ) : |
| | 385 | downsize ( rv._dsize() ), upsize ( rv_up._dsize() ), v2v_up ( rv.dataind ( rv_up ) ) { |
| | 386 | it_assert_debug ( v2v_up.length() ==downsize,"rv is not fully in rv_up" ); |
| | 387 | } |
| | 388 | //! Get val for myself from val of "Up" |
| | 389 | vec pushdown ( const vec &val_up ) { |
| | 390 | it_assert_debug ( upsize==val_up.length(),"Wrong val_up" ); |
| | 391 | return get_vec ( val_up,v2v_up ); |
| | 392 | } |
| | 393 | //! Fill val of "Up" by my pieces |
| | 394 | void pushup ( vec &val_up, const vec &val ) { |
| | 395 | it_assert_debug ( downsize==val.length(),"Wrong val" ); |
| | 396 | it_assert_debug ( upsize==val_up.length(),"Wrong val_up" ); |
| | 397 | set_subvector ( val_up, v2v_up, val ); |
| | 398 | } |
| | 399 | }; |
| 325 | | class datalink_m2e: public datalink_e2e { |
| 326 | | protected: |
| 327 | | //! Remember how long cond should be |
| 328 | | int condsize; |
| 329 | | //!upper_val-to-local_cond link, indeces of the upper val |
| 330 | | ivec v2c_up; |
| 331 | | //!upper_val-to-local_cond link, ideces of the local cond |
| 332 | | ivec v2c_lo; |
| 333 | | |
| 334 | | public: |
| 335 | | //! Constructor |
| 336 | | datalink_m2e ( const RV &rv, const RV &rvc, const RV &rv_up ) : |
| 337 | | datalink_e2e ( rv,rv_up ), condsize ( rvc.count() ) { |
| 338 | | //establish v2c connection |
| 339 | | rvc.dataind ( rv_up, v2c_lo, v2c_up ); |
| 340 | | } |
| 341 | | //!Construct condition |
| 342 | | vec get_cond ( const vec &val_up ) { |
| 343 | | vec tmp ( condsize ); |
| 344 | | set_subvector ( tmp,v2c_lo,val_up ( v2c_up ) ); |
| 345 | | return tmp; |
| 346 | | } |
| 347 | | void fill_val_cond ( vec &val_up, const vec &val, const vec &cond ) { |
| 348 | | it_assert_debug ( valsize==val.length(),"Wrong val" ); |
| 349 | | it_assert_debug ( valupsize==val_up.length(),"Wrong val_up" ); |
| 350 | | set_subvector ( val_up, v2v_up, val ); |
| 351 | | set_subvector ( val_up, v2c_up, cond ); |
| 352 | | } |
| 353 | | }; |
| | 402 | class datalink_m2e: public datalink { |
| | 403 | protected: |
| | 404 | //! Remember how long cond should be |
| | 405 | int condsize; |
| | 406 | //!upper_val-to-local_cond link, indeces of the upper val |
| | 407 | ivec v2c_up; |
| | 408 | //!upper_val-to-local_cond link, ideces of the local cond |
| | 409 | ivec v2c_lo; |
| | 410 | |
| | 411 | public: |
| | 412 | //! Constructor |
| | 413 | datalink_m2e ( const RV &rv, const RV &rvc, const RV &rv_up ) : |
| | 414 | datalink ( rv,rv_up ), condsize ( rvc._dsize() ) { |
| | 415 | //establish v2c connection |
| | 416 | rvc.dataind ( rv_up, v2c_lo, v2c_up ); |
| | 417 | } |
| | 418 | //!Construct condition |
| | 419 | vec get_cond ( const vec &val_up ) { |
| | 420 | vec tmp ( condsize ); |
| | 421 | set_subvector ( tmp,v2c_lo,val_up ( v2c_up ) ); |
| | 422 | return tmp; |
| | 423 | } |
| | 424 | void pushup_cond ( vec &val_up, const vec &val, const vec &cond ) { |
| | 425 | it_assert_debug ( downsize==val.length(),"Wrong val" ); |
| | 426 | it_assert_debug ( upsize==val_up.length(),"Wrong val_up" ); |
| | 427 | set_subvector ( val_up, v2v_up, val ); |
| | 428 | set_subvector ( val_up, v2c_up, cond ); |
| | 429 | } |
| | 430 | }; |
| 356 | | class datalink_m2m: public datalink_m2e { |
| 357 | | protected: |
| 358 | | //!cond-to-cond link, indeces of the upper cond |
| 359 | | ivec c2c_up; |
| 360 | | //!cond-to-cond link, indeces of the local cond |
| 361 | | ivec c2c_lo; |
| 362 | | public: |
| 363 | | //! Constructor |
| 364 | | datalink_m2m ( const RV &rv, const RV &rvc, const RV &rv_up, const RV &rvc_up ) : |
| 365 | | datalink_m2e ( rv, rvc, rv_up ) { |
| 366 | | //establish c2c connection |
| 367 | | rvc.dataind ( rvc_up, c2c_lo, c2c_up ); |
| 368 | | it_assert_debug ( c2c_lo.length() +v2c_lo.length() ==condsize, "cond is not fully given" ); |
| 369 | | } |
| 370 | | //! Get cond for myself from val and cond of "Up" |
| 371 | | vec get_cond ( const vec &val_up, const vec &cond_up ) { |
| 372 | | vec tmp ( condsize ); |
| 373 | | set_subvector ( tmp,v2c_lo,val_up ( v2c_up ) ); |
| 374 | | set_subvector ( tmp,c2c_lo,cond_up ( c2c_up ) ); |
| 375 | | return tmp; |
| 376 | | } |
| 377 | | //! Fill |
| 378 | | |
| 379 | | }; |
| 380 | | |
| 381 | | /*! |
| 382 | | @brief Class for storing results (and semi-results) of an experiment |
| 383 | | |
| 384 | | 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. |
| 385 | | */ |
| 386 | | class logger : public bdmroot { |
| 387 | | protected: |
| 388 | | //! RVs of all logged variables. |
| 389 | | Array<RV> entries; |
| 390 | | //! Names of logged quantities, e.g. names of algorithm variants |
| 391 | | Array<string> names; |
| 392 | | public: |
| 393 | | //!Default constructor |
| 394 | | logger ( ) : entries ( 0 ),names ( 0 ) {} |
| 395 | | |
| 396 | | //! returns an identifier which will be later needed for calling the log() function |
| 397 | | virtual int add ( const RV &rv, string name="" ) { |
| 398 | | int id=entries.length(); |
| 399 | | names=concat ( names, name ); // diff |
| 400 | | entries.set_length ( id+1,true ); |
| 401 | | entries ( id ) = rv; |
| 402 | | return id; // identifier of the last entry |
| 403 | | } |
| 404 | | |
| 405 | | //! log this vector |
| 406 | | virtual void logit ( int id, const vec &v ) =0; |
| 407 | | |
| 408 | | //! Shifts storage position for another time step. |
| 409 | | virtual void step() =0; |
| 410 | | |
| 411 | | //! Finalize storing information |
| 412 | | virtual void finalize() {}; |
| 413 | | |
| 414 | | //! Initialize the storage |
| 415 | | virtual void init() {}; |
| 416 | | |
| 417 | | //! for future use |
| 418 | | virtual ~logger() {}; |
| 419 | | }; |
| 420 | | |
| 421 | | /*! \brief Unconditional mpdf, allows using epdf in the role of mpdf. |
| 422 | | |
| 423 | | */ |
| 424 | | class mepdf : public mpdf { |
| 425 | | public: |
| 426 | | //!Default constructor |
| 427 | | mepdf ( const epdf* em ) :mpdf ( em->_rv(),RV() ) {ep=const_cast<epdf*> ( em );}; |
| 428 | | void condition ( const vec &cond ) {} |
| 429 | | }; |
| | 433 | class datalink_m2m: public datalink_m2e { |
| | 434 | protected: |
| | 435 | //!cond-to-cond link, indeces of the upper cond |
| | 436 | ivec c2c_up; |
| | 437 | //!cond-to-cond link, indeces of the local cond |
| | 438 | ivec c2c_lo; |
| | 439 | public: |
| | 440 | //! Constructor |
| | 441 | datalink_m2m ( const RV &rv, const RV &rvc, const RV &rv_up, const RV &rvc_up ) : |
| | 442 | datalink_m2e ( rv, rvc, rv_up ) { |
| | 443 | //establish c2c connection |
| | 444 | rvc.dataind ( rvc_up, c2c_lo, c2c_up ); |
| | 445 | it_assert_debug ( c2c_lo.length() +v2c_lo.length() ==condsize, "cond is not fully given" ); |
| | 446 | } |
| | 447 | //! Get cond for myself from val and cond of "Up" |
| | 448 | vec get_cond ( const vec &val_up, const vec &cond_up ) { |
| | 449 | vec tmp ( condsize ); |
| | 450 | set_subvector ( tmp,v2c_lo,val_up ( v2c_up ) ); |
| | 451 | set_subvector ( tmp,c2c_lo,cond_up ( c2c_up ) ); |
| | 452 | return tmp; |
| | 453 | } |
| | 454 | //! Fill |
| | 455 | |
| | 456 | }; |
| | 457 | |
| | 458 | /*! |
| | 459 | @brief Class for storing results (and semi-results) of an experiment |
| | 460 | |
| | 461 | 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. |
| | 462 | */ |
| | 463 | class logger : public bdmroot { |
| | 464 | protected: |
| | 465 | //! RVs of all logged variables. |
| | 466 | Array<RV> entries; |
| | 467 | //! Names of logged quantities, e.g. names of algorithm variants |
| | 468 | Array<string> names; |
| | 469 | public: |
| | 470 | //!Default constructor |
| | 471 | logger ( ) : entries ( 0 ),names ( 0 ) {} |
| | 472 | |
| | 473 | //! returns an identifier which will be later needed for calling the log() function |
| | 474 | virtual int add ( const RV &rv, string name="" ) { |
| | 475 | int id=entries.length(); |
| | 476 | names=concat ( names, name ); // diff |
| | 477 | entries.set_length ( id+1,true ); |
| | 478 | entries ( id ) = rv; |
| | 479 | return id; // identifier of the last entry |
| | 480 | } |
| | 481 | |
| | 482 | //! log this vector |
| | 483 | virtual void logit ( int id, const vec &v ) =0; |
| | 484 | |
| | 485 | //! Shifts storage position for another time step. |
| | 486 | virtual void step() =0; |
| | 487 | |
| | 488 | //! Finalize storing information |
| | 489 | virtual void finalize() {}; |
| | 490 | |
| | 491 | //! Initialize the storage |
| | 492 | virtual void init() {}; |
| | 493 | |
| | 494 | }; |
| | 495 | |
| | 496 | /*! \brief Unconditional mpdf, allows using epdf in the role of mpdf. |
| | 497 | |
| | 498 | */ |
| | 499 | class mepdf : public mpdf { |
| | 500 | public: |
| | 501 | //!Default constructor |
| | 502 | mepdf ( const epdf* em ) :mpdf ( ) {ep=const_cast<epdf*> ( em );}; |
| | 503 | void condition ( const vec &cond ) {} |
| | 504 | }; |
| 433 | | class compositepdf { |
| 434 | | protected: |
| 435 | | //!Number of mpdfs in the composite |
| 436 | | int n; |
| 437 | | //! Elements of composition |
| 438 | | Array<mpdf*> mpdfs; |
| 439 | | public: |
| 440 | | compositepdf ( Array<mpdf*> A0 ) : n ( A0.length() ), mpdfs ( A0 ) {}; |
| 441 | | //! find common rv, flag \param checkoverlap modifies whether overlaps are acceptable |
| 442 | | RV getrv ( bool checkoverlap=false ); |
| 443 | | //! common rvc of all mpdfs is written to rvc |
| 444 | | void setrvc ( const RV &rv, RV &rvc ); |
| 445 | | }; |
| 446 | | |
| 447 | | /*! \brief Abstract class for discrete-time sources of data. |
| 448 | | |
| 449 | | 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. |
| 450 | | 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). |
| 451 | | |
| | 508 | class compositepdf { |
| | 509 | protected: |
| | 510 | //!Number of mpdfs in the composite |
| | 511 | int n; |
| | 512 | //! Elements of composition |
| | 513 | Array<mpdf*> mpdfs; |
| | 514 | public: |
| | 515 | compositepdf ( Array<mpdf*> A0 ) : n ( A0.length() ), mpdfs ( A0 ) {}; |
| | 516 | //! find common rv, flag \param checkoverlap modifies whether overlaps are acceptable |
| | 517 | RV getrv ( bool checkoverlap=false ); |
| | 518 | //! common rvc of all mpdfs is written to rvc |
| | 519 | void setrvc ( const RV &rv, RV &rvc ); |
| | 520 | }; |
| | 521 | |
| | 522 | /*! \brief Abstract class for discrete-time sources of data. |
| | 523 | |
| | 524 | 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. |
| | 525 | 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). |
| | 526 | |
| | 527 | */ |
| | 528 | |
| | 529 | class DS : public bdmroot { |
| | 530 | protected: |
| | 531 | int dtsize; |
| | 532 | int utsize; |
| | 533 | //!Description of data returned by \c getdata(). |
| | 534 | RV Drv; |
| | 535 | //!Description of data witten by by \c write(). |
| | 536 | RV Urv; // |
| | 537 | //! Remember its own index in Logger L |
| | 538 | int L_dt, L_ut; |
| | 539 | public: |
| | 540 | //! default constructors |
| | 541 | DS() :Drv ( ),Urv ( ) {}; |
| | 542 | //! Returns full vector of observed data=[output, input] |
| | 543 | virtual void getdata ( vec &dt ) {it_error ( "abstract class" );}; |
| | 544 | //! Returns data records at indeces. |
| | 545 | virtual void getdata ( vec &dt, const ivec &indeces ) {it_error ( "abstract class" );}; |
| | 546 | //! Accepts action variable and schedule it for application. |
| | 547 | virtual void write ( vec &ut ) {it_error ( "abstract class" );}; |
| | 548 | //! Accepts action variables at specific indeces |
| | 549 | virtual void write ( vec &ut, const ivec &indeces ) {it_error ( "abstract class" );}; |
| | 550 | |
| | 551 | //! Moves from \f$ t \f$ to \f$ t+1 \f$, i.e. perfroms the actions and reads response of the system. |
| | 552 | virtual void step() =0; |
| | 553 | |
| | 554 | //! Register DS for logging into logger L |
| | 555 | virtual void log_add ( logger &L ) { |
| | 556 | it_assert_debug ( dtsize==Drv._dsize(),"" ); |
| | 557 | it_assert_debug ( utsize==Urv._dsize(),"" ); |
| | 558 | |
| | 559 | L_dt=L.add ( Drv,"" ); |
| | 560 | L_ut=L.add ( Urv,"" ); |
| | 561 | } |
| | 562 | //! Register DS for logging into logger L |
| | 563 | virtual void logit ( logger &L ) { |
| | 564 | vec tmp ( Drv._dsize() +Urv._dsize() ); |
| | 565 | getdata ( tmp ); |
| | 566 | // d is first in getdata |
| | 567 | L.logit ( L_dt,tmp.left ( Drv._dsize() ) ); |
| | 568 | // u follows after d in getdata |
| | 569 | L.logit ( L_ut,tmp.mid ( Drv._dsize(), Urv._dsize() ) ); |
| | 570 | } |
| | 571 | //!access function |
| | 572 | virtual RV _drv() const {return concat ( Drv,Urv );} |
| | 573 | //!access function |
| | 574 | const RV& _urv() const {return Urv;} |
| | 575 | }; |
| | 576 | |
| | 577 | /*! \brief Bayesian Model of a system, i.e. all uncertainty is modeled by probabilities. |
| | 578 | |
| | 579 | */ |
| | 580 | |
| | 581 | class BM :public bdmroot { |
| | 582 | protected: |
| | 583 | //! Random variable of the data (optional) |
| | 584 | RV drv; |
| | 585 | //!Logarithm of marginalized data likelihood. |
| | 586 | double ll; |
| | 587 | //! 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. |
| | 588 | bool evalll; |
| | 589 | public: |
| | 590 | //! \name Constructors |
| | 591 | //! @{ |
| | 592 | |
| | 593 | BM () :ll (0),evalll ( false) {}; |
| | 594 | BM ( const BM &B ) : drv(B.drv), ll ( B.ll ), evalll ( B.evalll ) {} |
| | 595 | //! Copy function required in vectors, Arrays of BM etc. Have to be DELETED manually! |
| | 596 | //! Prototype: \code BM* _copy_(){return new BM(*this);} \endcode |
| | 597 | virtual BM* _copy_ () {return NULL;}; |
| | 598 | //!@} |
| | 599 | |
| | 600 | //! \name Mathematical operations |
| | 601 | //!@{ |
| | 602 | |
| | 603 | /*! \brief Incremental Bayes rule |
| | 604 | @param dt vector of input data |
| 453 | | |
| 454 | | class DS : public bdmroot { |
| 455 | | protected: |
| 456 | | //!Observed variables, returned by \c getdata(). |
| 457 | | RV Drv; |
| 458 | | //!Action variables, accepted by \c write(). |
| 459 | | RV Urv; // |
| 460 | | //! Remember its own index in Logger L |
| 461 | | int L_dt, L_ut; |
| 462 | | public: |
| 463 | | DS() :Drv ( RV0 ),Urv ( RV0 ) {}; |
| 464 | | DS ( const RV &Drv0, const RV &Urv0 ) :Drv ( Drv0 ),Urv ( Urv0 ) {}; |
| 465 | | //! Returns full vector of observed data=[output, input] |
| 466 | | virtual void getdata ( vec &dt ) {it_error ( "abstract class" );}; |
| 467 | | //! Returns data records at indeces. |
| 468 | | virtual void getdata ( vec &dt, const ivec &indeces ) {it_error ( "abstract class" );}; |
| 469 | | //! Accepts action variable and schedule it for application. |
| 470 | | virtual void write ( vec &ut ) {it_error ( "abstract class" );}; |
| 471 | | //! Accepts action variables at specific indeces |
| 472 | | virtual void write ( vec &ut, const ivec &indeces ) {it_error ( "abstract class" );}; |
| 473 | | |
| 474 | | //! Moves from \f$t\f$ to \f$t+1\f$, i.e. perfroms the actions and reads response of the system. |
| 475 | | virtual void step() =0; |
| 476 | | |
| 477 | | //! Register DS for logging into logger L |
| 478 | | virtual void log_add ( logger &L ) { |
| 479 | | L_dt=L.add ( Drv,"" ); |
| 480 | | L_ut=L.add ( Urv,"" ); |
| 481 | | } |
| 482 | | //! Register DS for logging into logger L |
| 483 | | virtual void logit ( logger &L ) { |
| 484 | | vec tmp(Drv.count()+Urv.count()); |
| 485 | | getdata(tmp); |
| 486 | | // d is first in getdata |
| 487 | | L.logit ( L_dt,tmp.left ( Drv.count() ) ); |
| 488 | | // u follows after d in getdata |
| 489 | | L.logit ( L_ut,tmp.mid ( Drv.count(), Urv.count() ) ); |
| 490 | | } |
| 491 | | //!access function |
| 492 | | virtual RV _drv() const {return concat(Drv,Urv);} |
| 493 | | //!access function |
| 494 | | const RV& _urv() const {return Urv;} |
| 495 | | }; |
| 496 | | |
| 497 | | /*! \brief Bayesian Model of a system, i.e. all uncertainty is modeled by probabilities. |
| 498 | | |
| 499 | | */ |
| 500 | | |
| 501 | | class BM :public bdmroot { |
| 502 | | protected: |
| 503 | | //!Random variable of the posterior |
| 504 | | RV rv; |
| 505 | | //! Random variable of the data (optional) |
| 506 | | RV drv; |
| 507 | | //!Logarithm of marginalized data likelihood. |
| 508 | | double ll; |
| 509 | | //! 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. |
| 510 | | bool evalll; |
| 511 | | public: |
| 512 | | |
| 513 | | //!Default constructor |
| 514 | | BM ( const RV &rv0, double ll0=0,bool evalll0=true ) :rv ( rv0 ), ll ( ll0 ),evalll ( evalll0 ) {//Fixme: test rv |
| 515 | | }; |
| 516 | | //!Copy constructor |
| 517 | | BM ( const BM &B ) : rv ( B.rv ), ll ( B.ll ), evalll ( B.evalll ) {} |
| 518 | | |
| 519 | | /*! \brief Incremental Bayes rule |
| 520 | | @param dt vector of input data |
| 521 | | */ |
| 522 | | virtual void bayes ( const vec &dt ) = 0; |
| 523 | | //! Batch Bayes rule (columns of Dt are observations) |
| 524 | | virtual void bayesB ( const mat &Dt ); |
| 525 | | //! Returns a reference to the epdf representing posterior density on parameters. |
| 526 | | virtual const epdf& _epdf() const =0; |
| 527 | | |
| 528 | | //! Returns a pointer to the epdf representing posterior density on parameters. Use with care! |
| 529 | | virtual const epdf* _e() const =0; |
| 530 | | |
| 531 | | //! Evaluates predictive log-likelihood of the given data record |
| 532 | | //! I.e. marginal likelihood of the data with the posterior integrated out. |
| 533 | | virtual double logpred ( const vec &dt ) const{it_error ( "Not implemented" );return 0.0;} |
| 534 | | //! Matrix version of logpred |
| 535 | | vec logpred_m ( const mat &dt ) const{vec tmp ( dt.cols() );for ( int i=0;i<dt.cols();i++ ) {tmp ( i ) =logpred ( dt.get_col ( i ) );}return tmp;} |
| 536 | | |
| 537 | | //!Constructs a predictive density (marginal density on data) |
| 538 | | virtual epdf* predictor ( const RV &rv ) const {it_error ( "Not implemented" );return NULL;}; |
| 539 | | |
| 540 | | //! Destructor for future use; |
| 541 | | virtual ~BM() {}; |
| 542 | | //!access function |
| 543 | | const RV& _rv() const {return rv;} |
| 544 | | //!access function |
| 545 | | const RV& _drv() const {return drv;} |
| 546 | | //!set drv |
| 547 | | void set_drv(const RV &rv){drv=rv;} |
| 548 | | //!access function |
| 549 | | double _ll() const {return ll;} |
| 550 | | //!access function |
| 551 | | void set_evalll ( bool evl0 ) {evalll=evl0;} |
| 552 | | |
| 553 | | //! Copy function required in vectors, Arrays of BM etc. Have to be DELETED manually! |
| 554 | | //! Prototype: BM* _copy_(){BM Tmp*=new Tmp(this*); return Tmp; } |
| 555 | | virtual BM* _copy_ ( bool changerv=false ) {it_error ( "function _copy_ not implemented for this BM" ); return NULL;}; |
| 556 | | }; |
| 557 | | |
| 558 | | /*! |
| 559 | | \brief Conditional Bayesian Filter |
| 560 | | |
| 561 | | Evaluates conditional filtering density \f$f(rv|rvc,data)\f$ for a given \c rvc which is specified in each step by calling function \c condition. |
| 562 | | |
| 563 | | This is an interface class used to assure that certain BM has operation \c condition . |
| 564 | | |
| 565 | | */ |
| 566 | | |
| 567 | | class BMcond :public bdmroot { |
| 568 | | protected: |
| 569 | | //! Identificator of the conditioning variable |
| 570 | | RV rvc; |
| 571 | | public: |
| 572 | | //! Substitute \c val for \c rvc. |
| 573 | | virtual void condition ( const vec &val ) =0; |
| 574 | | //! Default constructor |
| 575 | | BMcond ( RV &rv0 ) :rvc ( rv0 ) {}; |
| 576 | | //! Destructor for future use |
| 577 | | virtual ~BMcond() {}; |
| 578 | | //! access function |
| 579 | | const RV& _rvc() const {return rvc;} |
| 580 | | }; |
| | 606 | virtual void bayes ( const vec &dt ) = 0; |
| | 607 | //! Batch Bayes rule (columns of Dt are observations) |
| | 608 | virtual void bayesB ( const mat &Dt ); |
| | 609 | //! Evaluates predictive log-likelihood of the given data record |
| | 610 | //! I.e. marginal likelihood of the data with the posterior integrated out. |
| | 611 | virtual double logpred ( const vec &dt ) const{it_error ( "Not implemented" );return 0.0;} |
| | 612 | //! Matrix version of logpred |
| | 613 | vec logpred_m ( const mat &dt ) const{vec tmp ( dt.cols() );for ( int i=0;i<dt.cols();i++ ) {tmp ( i ) =logpred ( dt.get_col ( i ) );}return tmp;} |
| | 614 | |
| | 615 | //!Constructs a predictive density \f$ f(d_{t+1} |d_{t}, \ldots d_{0}) \f$ |
| | 616 | virtual epdf* epredictor ( ) const {it_error ( "Not implemented" );return NULL;}; |
| | 617 | //!Constructs a conditional density 1-step ahead predictor \f$ f(d_{t+1} |d_{t+h-1}, \ldots d_{t}) |
| | 618 | virtual mpdf* predictor ( ) const {it_error ( "Not implemented" );return NULL;}; |
| | 619 | //!@} |
| | 620 | |
| | 621 | //! \name Access to attributes |
| | 622 | //!@{ |
| | 623 | |
| | 624 | const RV& _drv() const {return drv;} |
| | 625 | void set_drv ( const RV &rv ) {drv=rv;} |
| | 626 | double _ll() const {return ll;} |
| | 627 | void set_evalll ( bool evl0 ) {evalll=evl0;} |
| | 628 | virtual const epdf& _epdf() const =0; |
| | 629 | virtual const epdf* _e() const =0; |
| | 630 | //!@} |
| | 631 | |
| | 632 | }; |
| | 633 | |
| | 634 | /*! |
| | 635 | \brief Conditional Bayesian Filter |
| | 636 | |
| | 637 | Evaluates conditional filtering density \f$f(rv|rvc,data)\f$ for a given \c rvc which is specified in each step by calling function \c condition. |
| | 638 | |
| | 639 | This is an interface class used to assure that certain BM has operation \c condition . |
| | 640 | |
| | 641 | */ |
| | 642 | |
| | 643 | class BMcond :public bdmroot { |
| | 644 | protected: |
| | 645 | //! Identificator of the conditioning variable |
| | 646 | RV rvc; |
| | 647 | public: |
| | 648 | //! Substitute \c val for \c rvc. |
| | 649 | virtual void condition ( const vec &val ) =0; |
| | 650 | //! Default constructor |
| | 651 | BMcond ( ) :rvc ( ) {}; |
| | 652 | //! Destructor for future use |
| | 653 | virtual ~BMcond() {}; |
| | 654 | //! access function |
| | 655 | const RV& _rvc() const {return rvc;} |
| | 656 | }; |
