Develop metatree

bayesml/__init__.py

@@ -6,6 +6,7 @@
 from . import multivariate_normal
 from . import normal
 from . import poisson
+from . import metatree
 
 __all__ = ['bernoulli',
            'categorical',
@@ -14,5 +15,6 @@
            'linearregression',
            'multivariate_normal',
            'normal',
-           'poisson'
+           'poisson',
+           'metatree'
            ]

bayesml/_check.py

@@ -65,6 +65,7 @@ def nonneg_int_vecs(val,val_name,exception_class):
             return val
     raise(exception_class(val_name + " must be a numpy.ndarray whose ndim >= 1 and dtype is int. Its values must be non-negative (including 0)."))
 
+
 def nonneg_float_vec(val,val_name,exception_class):
     if type(val) is np.ndarray:
         if np.issubdtype(val.dtype,np.floating) and val.ndim == 1 and np.all(val>=0):
@@ -242,8 +243,14 @@ def onehot_vecs(val,val_name,exception_class):
             return val
     raise(exception_class(val_name + " must be a numpy.ndarray whose dtype is int and whose last axis constitutes one-hot vectors."))
 
+def int_vecs(val,val_name,exception_class):
+    if type(val) is np.ndarray:
+        if np.issubdtype(val.dtype,np.integer) and val.ndim >= 1:
+            return val
+    raise(exception_class(val_name + " must be a numpy.ndarray whose dtype is int and ndim >= 1."))
+
 def shape_consistency(val: int, val_name: str, correct: int, correct_name: str, exception_class):
     if val != correct:
         message = (f"{val_name} must coincide with {correct_name}: "
                    + f"{val_name} = {val}, {correct_name} = {correct}")
-        raise(exception_class(message))
+        raise(exception_class(message))

bayesml/bernoulli/_bernoulli.py

@@ -272,18 +272,20 @@ def update_posterior(self,x):
         self.hn_alpha += np.sum(x==1)
         self.hn_beta += np.sum(x==0)
 
-    def estimate_params(self,loss="squared"):
+    def estimate_params(self,loss="squared",dict_out=False):
         """Estimate the parameter of the stochastic data generative model under the given criterion.
 
         Parameters
         ----------
         loss : str, optional
             Loss function underlying the Bayes risk function, by default \"squared\".
             This function supports \"squared\", \"0-1\", \"abs\", and \"KL\".
+        dict_out : bool, optional
+            If ``True``, output will be a dict, by default ``False``.
 
         Returns
         -------
-        Estimator : {float, None, rv_frozen}
+        estimator : {float, None, rv_frozen} or dict of {str : float, None}
             The estimated values under the given loss function. If it is not exist, `None` will be returned.
             If the loss function is \"KL\", the posterior distribution itself will be returned
             as rv_frozen object of scipy.stats.
@@ -294,19 +296,37 @@ def estimate_params(self,loss="squared"):
         scipy.stats.rv_discrete
         """
         if loss == "squared":
-            return self.hn_alpha / (self.hn_alpha + self.hn_beta)
+            if dict_out:
+                return {'theta':self.hn_alpha / (self.hn_alpha + self.hn_beta)}
+            else:
+                return self.hn_alpha / (self.hn_alpha + self.hn_beta)
         elif loss == "0-1":
             if self.hn_alpha > 1.0 and self.hn_beta > 1.0:
-                return (self.hn_alpha - 1.0) / (self.hn_alpha + self.hn_beta - 2.0)
+                if dict_out:
+                    return {'theta':(self.hn_alpha - 1.0) / (self.hn_alpha + self.hn_beta - 2.0)}
+                else:
+                    return (self.hn_alpha - 1.0) / (self.hn_alpha + self.hn_beta - 2.0)
             elif self.hn_alpha > 1.0:
-                return 1.0
+                if dict_out:
+                    return {'theta':1.0}
+                else:
+                    return 1.0
             elif self.hn_beta > 1.0:
-                return 0.0
+                if dict_out:
+                    return {'theta':0.0}
+                else:
+                    return 0.0
             else:
                 warnings.warn("MAP estimate doesn't exist for the current hn_alpha and hn_beta.",ResultWarning)
-                return None
+                if dict_out:
+                    return {'theta':None}
+                else:
+                    return None
         elif loss == "abs":
-            return ss_beta.median(self.hn_alpha,self.hn_beta)
+            if dict_out:
+                return {'theta':ss_beta.median(self.hn_alpha,self.hn_beta)}
+            else:
+                return ss_beta.median(self.hn_alpha,self.hn_beta)
         elif loss == "KL":
             return ss_beta(self.hn_alpha,self.hn_beta)
         else:

bayesml/categorical/_categorical.py

@@ -356,18 +356,20 @@ def update_posterior(self, x):
         for k in range(self.degree):
             self.hn_alpha_vec[k] += x[:,k].sum()
 
-    def estimate_params(self, loss="squared"):
+    def estimate_params(self, loss="squared",dict_out=False):
         """Estimate the parameter of the stochastic data generative model under the given criterion.
 
         Parameters
         ----------
         loss : str, optional
             Loss function underlying the Bayes risk function, by default \"squared\".
             This function supports \"squared\", \"0-1\", and \"KL\".
+        dict_out : bool, optional
+            If ``True``, output will be a dict, by default ``False``.
 
         Returns
         -------
-        Estimates : {numpy ndarray, float, None, or rv_frozen}
+        estimates : {numpy ndarray, float, None, or rv_frozen}
             The estimated values under the given loss function. If it is not exist, `None` will be returned.
             If the loss function is \"KL\", the posterior distribution itself will be returned
             as rv_frozen object of scipy.stats.
@@ -378,10 +380,16 @@ def estimate_params(self, loss="squared"):
         scipy.stats.rv_discrete
         """
         if loss == "squared":
-            return self.hn_alpha_vec / np.sum(self.hn_alpha_vec)
+            if dict_out:
+                return {'theta_vec':self.hn_alpha_vec / np.sum(self.hn_alpha_vec)}
+            else:
+                return self.hn_alpha_vec / np.sum(self.hn_alpha_vec)
         elif loss == "0-1":
             if np.all(self.hn_alpha_vec > 1):
-                return (self.hn_alpha_vec - 1) / (np.sum(self.hn_alpha_vec) - self.degree)
+                if dict_out:
+                    return {'theta_vec':(self.hn_alpha_vec - 1) / (np.sum(self.hn_alpha_vec) - self.degree)}
+                else:
+                    return (self.hn_alpha_vec - 1) / (np.sum(self.hn_alpha_vec) - self.degree)
             else:
                 warnings.warn("MAP estimate of lambda_mat doesn't exist for the current hn_alpha_vec.",ResultWarning)
                 return None

bayesml/exponential/_exponential.py

@@ -66,7 +66,7 @@ def gen_params(self):
 
         The generated vaule is set at ``self.lambda_``.
         """
-        self.lambda_ = self.rng.gamma(self.h_alpha,1.0/self.h_beta, 1)
+        self.lambda_ = self.rng.gamma(self.h_alpha,1.0/self.h_beta)
 
     def set_params(self,lambda_):
         """Set the parameter of the sthocastic data generative model.
@@ -277,18 +277,20 @@ def update_posterior(self,x):
         self.hn_alpha += x.size
         self.hn_beta += np.sum(x)
 
-    def estimate_params(self,loss="squared"):
+    def estimate_params(self,loss="squared",dict_out=False):
         """Estimate the parameter of the stochastic data generative model under the given criterion.
 
         Parameters
         ----------
         loss : str, optional
             Loss function underlying the Bayes risk function, by default \"squared\".
             This function supports \"squared\", \"0-1\", \"abs\", and \"KL\".
+        dict_out : bool, optional
+            If ``True``, output will be a dict, by default ``False``.
 
         Returns
         -------
-        Estimator : {float, None, rv_frozen}
+        estimator : {float, None, rv_frozen}
             The estimated values under the given loss function. If it is not exist, `None` will be returned.
             If the loss function is \"KL\", the posterior distribution itself will be returned
             as rv_frozen object of scipy.stats.
@@ -299,14 +301,26 @@ def estimate_params(self,loss="squared"):
         scipy.stats.rv_discrete
         """
         if loss == "squared":
-            return self.hn_alpha / self.hn_beta
+            if dict_out:
+                return {'lambda_':self.hn_alpha / self.hn_beta}
+            else:
+                return self.hn_alpha / self.hn_beta
         elif loss == "0-1":
             if self.hn_alpha > 1.0 :
-                return (self.hn_alpha - 1.0) / self.hn_beta
+                if dict_out:
+                    return {'lambda_':(self.hn_alpha - 1.0) / self.hn_beta}
+                else:
+                    return (self.hn_alpha - 1.0) / self.hn_beta
             else:
-                return 0.0
+                if dict_out:
+                    return {'lambda_':0.0}
+                else:
+                    return 0.0
         elif loss == "abs":
-            return ss_gamma.median(a=self.hn_alpha,scale=1/self.hn_beta)
+            if dict_out:
+                return {'lambda_':ss_gamma.median(a=self.hn_alpha,scale=1/self.hn_beta)}
+            else:
+                return ss_gamma.median(a=self.hn_alpha,scale=1/self.hn_beta)
         elif loss == "KL":
             return ss_gamma(a=self.hn_alpha,scale=1/self.hn_beta)
         else:

bayesml/linearregression/_linearregression.py

@@ -25,6 +25,10 @@ class GenModel(base.Generative):
         a value consistent with ``theta_vec``, ``h_mu_vec``, 
         and ``h_lambda_mat`` is used. If all of them are not given,
         degree is assumed to be 1.
+    theta_vec : numpy ndarray, optional
+        a vector of real numbers, by default [0.0, 0.0, ... , 0.0]
+    tau : float, optional
+        a positive real number, by default 1.0
     h_mu_vec : numpy ndarray, optional
         a vector of real numbers, by default [0.0, 0.0, ... , 0.0]
     h_lambda_mat : numpy ndarray, optional
@@ -558,7 +562,7 @@ def update_posterior(self, x, y):
         self.hn_beta += (-self.hn_mu_vec[np.newaxis,:] @ self.hn_lambda_mat @ self.hn_mu_vec[:,np.newaxis]
                          + y @ y + hn1_mu[np.newaxis,:] @ hn1_Lambda @ hn1_mu[:,np.newaxis])[0,0] /2.0
 
-    def estimate_params(self,loss="squared"):
+    def estimate_params(self,loss="squared",dict_out=False):
         """Estimate the parameter of the stochastic data generative model under the given criterion.
 
         Note that the criterion is applied to estimating ``theta_vec`` and ``tau`` independently.
@@ -569,10 +573,12 @@ def estimate_params(self,loss="squared"):
         loss : str, optional
             Loss function underlying the Bayes risk function, by default \"squared\".
             This function supports \"squared\", \"0-1\", \"abs\", and \"KL\".
+        dict_out : bool, optional
+            If ``True``, output will be a dict, by default ``False``.
 
         Returns
         -------
-        Estimates : tuple of {numpy ndarray, float, None, or rv_frozen}
+        estimates : tuple of {numpy ndarray, float, None, or rv_frozen}
             * ``theta_vec`` : the estimate for w
             * ``tau_hat`` : the estimate for tau
             The estimated values under the given loss function. If it is not exist, `None` will be returned.
@@ -584,15 +590,27 @@ def estimate_params(self,loss="squared"):
         scipy.stats.rv_continuous
         scipy.stats.rv_discrete
         """
-        if loss == "squared": 
-            return  self.hn_mu_vec, self.hn_alpha/self.hn_beta
+        if loss == "squared":
+            if dict_out:
+                return {'theta_vec':self.hn_mu_vec,'tau':self.hn_alpha/self.hn_beta}
+            else:
+                return  self.hn_mu_vec, self.hn_alpha/self.hn_beta
         elif loss == "0-1":
             if self.hn_alpha >= 1.0:
-                return self.hn_mu_vec, (self.hn_alpha - 1.0) / self.hn_beta
+                if dict_out:
+                    return {'theta_vec':self.hn_mu_vec,'tau':(self.hn_alpha - 1.0) / self.hn_beta}
+                else:
+                    return self.hn_mu_vec, (self.hn_alpha - 1.0) / self.hn_beta
             else:
-                return self.hn_mu_vec, 0
+                if dict_out:
+                    return {'theta_vec':self.hn_mu_vec,'tau':0.0}
+                else:
+                    return self.hn_mu_vec, 0.0
         elif loss == "abs":
-            return self.hn_mu_vec, ss_gamma.median(a=self.hn_alpha,scale=1.0/self.hn_beta)
+            if dict_out:
+                return {'theta_vec':self.hn_mu_vec,'tau':ss_gamma.median(a=self.hn_alpha,scale=1.0/self.hn_beta)}
+            else:
+                return self.hn_mu_vec, ss_gamma.median(a=self.hn_alpha,scale=1.0/self.hn_beta)
         elif loss == "KL":
             return (ss_multivariate_t(loc=self.hn_mu_vec,
                                         shape=np.linalg.inv(self.hn_alpha / self.hn_beta * self.hn_lambda_mat),