Index: library/bdm/stat/exp_family.h
===================================================================
--- library/bdm/stat/exp_family.h (revision 1020)
+++ library/bdm/stat/exp_family.h (revision 1033)
@@ -552,4 +552,90 @@
 UIREGISTER ( eDirich );
 
+/*! Product of Beta distributions
+
+Continuous Dirichlet density of \f$n\f$-dimensional variable \f$x\f$
+\f[
+f(x|\alpha,\beta)  = \frac{\Gamma(\alpha+\beta)}{\Gamma(\alpha)\Gamma(\beta)}\, x^{\alpha-1}(1-x)^{\beta-1} 
+\f]
+is a simplification of Dirichlet to univariate case.
+*/
+class eBeta: public eEF {
+	public:
+		//!sufficient statistics
+		vec alpha;
+		//!sufficient statistics
+		vec beta;
+	public:
+		//!\name Constructors
+		//!@{
+			
+			eBeta () : eEF () {};
+			eBeta ( const eBeta &B0 ) : eEF (), alpha(B0.alpha),beta(B0.beta) {
+				validate();
+			};
+			//!@}
+			
+			//! using sampling procedure from wikipedia
+			vec sample() const {
+				vec y ( beta.length() ); // all vectors
+				for ( int i = 0; i < beta.length(); i++ ) {
+					GamRNG.setup ( alpha ( i ), 1 );
+					#pragma omp critical
+					double Ga = GamRNG();
+					
+					GamRNG.setup ( beta ( i ), 1 );
+					#pragma omp critical
+					double Gb = GamRNG();
+					
+					y ( i ) = Ga/(Ga+Gb);
+				}
+				return y;
+			}
+			
+			vec mean() const {
+				return elem_div(alpha, alpha + beta); // dot-division
+			};
+			vec variance() const {
+				vec apb=alpha+beta;
+				return elem_div (elem_mult ( alpha, beta) ,
+								 elem_mult ( elem_mult(apb,apb), apb+1 ) );
+			}
+			//! In this instance, val is ...
+			double evallog_nn ( const vec &val ) const {
+				double tmp;
+				tmp = ( alpha - 1 ) * log ( val ) + (beta-1)*log(1-val);
+				return tmp;
+			}
+			
+			double lognc () const {
+				double lgb = 0.0;
+				for ( int i = 0; i < beta.length(); i++ ) {
+					lgb += -lgamma ( alpha(i)+beta(i) ) + lgamma(alpha(i)) + lgamma(beta(i));
+				}
+				return lgb;
+			}
+			
+			/*! configuration structure
+			\code
+			class = 'eBeta';
+			beta  = [];           //parametr beta
+			alpha = [];           //parametr alpha
+			\endcode
+			*/
+			void from_setting ( const Setting &set ){
+				UI::get(alpha, set, "alpha", UI::compulsory);
+				UI::get(beta, set, "beta", UI::compulsory);
+			}
+			void validate(){
+				bdm_assert(alpha.length()==beta.length(), "eBeta:: alpha and beta length do not match");
+				dim = alpha.length();
+			};
+			void to_setting ( Setting &set ) const{
+				UI::save(alpha, set, "alpha");
+				UI::save(beta, set, "beta");
+			}
+};
+UIREGISTER ( eBeta );
+
 /*! Random Walk on Dirichlet
 Using simple assignment
@@ -594,4 +680,53 @@
 UIREGISTER ( mDirich );
 
+/*! \brief Random Walk with vector Beta distribution
+Using simple assignment
+\f{align} 
+\alpha &= rvc / k + \beta_c \\
+\beta &= (1-rvc) / k + \beta_c \\
+\f{align}
+for each element of alpha and beta, mean value = rvc, variance = (k+1)*mean*mean;
+
+The greater \f$ k \f$ is, the greater is the variance of the random walk;
+
+\f$ \beta_c \f$ is used as regularizing element to avoid corner cases, i.e. when one element of rvc is zero.
+By default is it set to 0.1;
+*/
+
+class mBeta: public pdf_internal<eBeta>{
+	//! vector of constants \f$ k \f$ of the random walk 
+	vec k;
+	//! stabilizing coefficient \f$ \beta_c \f$
+	vec betac;
+
+	public:
+		void condition ( const vec &val ) {
+			this->iepdf.alpha =  elem_div(val , k) + betac;
+			this->iepdf.beta =  elem_div (1-val , k) + betac;
+		};
+		/*! Create Beta random walk
+		\f[ f(rv|rvc) = \prod Beta(rvc,k) \f]
+		from structure
+		\code
+		class = 'mBeta';
+		k = 1;                      // multiplicative constant k
+		--- optional ---
+		rv = RV({'name'},size)      // description of RV
+		beta  = [];                 // initial value of beta
+		betac = [];                 // initial value of beta stabilizing constant
+		\endcode
+		*/
+		void from_setting ( const Setting &set );
+		void 	to_setting  (Setting  &set) const;
+		void validate(){
+			
+			pdf_internal<eBeta>::validate();
+			bdm_assert(betac.length()==dimension(),"Incomaptible betac");
+			bdm_assert(k.length()==dimension(),"Incomaptible k");
+			dimc = iepdf.dimension();
+		}
+		//! 
+};
+UIREGISTER(mBeta);
 
 //! \brief Estimator for Multinomial density
Index: library/bdm/stat/exp_family.cpp
===================================================================
--- library/bdm/stat/exp_family.cpp (revision 1015)
+++ library/bdm/stat/exp_family.cpp (revision 1033)
@@ -775,3 +775,30 @@
 		dimc = _beta.length();
 	}
-}
+
+
+void mBeta::from_setting ( const Setting &set ) {
+	pdf::from_setting ( set ); // reads rv and rvc
+	if ( _rv()._dsize() > 0 ) {
+		rvc = _rv().copy_t ( -1 );
+	}
+	if ( !UI::get ( iepdf.beta, set, "beta", UI::optional ) ) {
+		iepdf.beta = ones ( _rv()._dsize() );
+	}
+	if ( !UI::get ( iepdf.alpha, set, "alpha", UI::optional ) ) {
+		iepdf.alpha = ones ( _rv()._dsize() );
+	}
+	if ( !UI::get ( betac, set, "betac", UI::optional ) ) {
+		betac = 0.1 * ones ( _rv()._dsize() );
+	}
+	
+	UI::get ( k, set, "k", UI::compulsory );
+}
+
+void mBeta::to_setting  (Setting  &set) const {
+	pdf::to_setting(set);
+	UI::save( iepdf.beta, set, "beta");
+	UI::save( betac, set, "betac");
+	UI::save ( k, set, "k" );
+}
+
+}
