update distributions, changing standardisation for subspaces

equadratures/distributions/beta.py

@@ -22,22 +22,40 @@ class Beta(Distribution):
         Upper bound of the support of the beta distribution.
     """
     def __init__(self, lower=None, upper=None, shape_A=None, shape_B=None):
-        self.shape_A = shape_A
-        self.shape_B = shape_B
-        self.lower = lower
-        self.upper = upper
-        if (self.shape_A is not None) and (self.shape_B is not None):
-            if self.shape_A >= 1. and self.shape_B >= 1.0:
-                self.mean = (self.shape_A) / (self.shape_A + self.shape_B)
-                self.variance = (self.shape_A * self.shape_B) / ( (self.shape_A + self.shape_B)**2 * (self.shape_A + self.shape_B + 1.0) )
-                self.skewness = 2.0 * (self.shape_B - self.shape_A) * np.sqrt(self.shape_A + self.shape_B + 1.0) / ( (self.shape_A + self.shape_B + 2.0) * np.sqrt(self.shape_A * self.shape_B) )
-                self.kurtosis = 6.0 * ((self.shape_A - self.shape_B)**2 * (self.shape_A + self.shape_B + 1.0) - self.shape_A * self.shape_B * (self.shape_A + self.shape_B + 2.0)  ) /( (self.shape_A * self.shape_B) * (self.shape_A + self.shape_B + 2.0) * (self.shape_A + self.shape_B + 3.0)) + 3.0
-                self.bounds = np.array([0, 1])
-                self.shape_parameter_A = self.shape_B - 1.0
-                self.shape_parameter_B = self.shape_A - 1.0
-                self.parent = beta(self.shape_A, self.shape_B)
-        if (self.lower is not None) and (self.upper is not None):
-            self.x_range_for_pdf = np.linspace(self.lower, self.upper, RECURRENCE_PDF_SAMPLES)
+        if shape_A is None:
+            self.shape_A = 2.0
+        else:
+            self.shape_A = shape_A
+        if shape_B is None:
+            self.shape_B = 2.0
+        else:
+            self.shape_B = shape_B
+        if self.shape_A <= 0 or self.shape_B <= 0:
+            raise ValueError('Invalid Beta distribution parameters. Alpha and beta should be positive.')
+
+        self.shape_parameter_A = self.shape_B
+        self.shape_parameter_B = self.shape_A
+
+        if lower is None:
+            self.lower = 0.0
+        else:
+            self.lower = lower
+        if upper is None:
+            self.upper = 1.0
+        else:
+            self.upper = upper
+
+        if self.lower > self.upper:
+            raise ValueError('Invalid Beta distribution parameters. Lower should be smaller than upper.')
+
+        self.bounds = np.array([self.lower, self.upper])
+        loc = self.lower
+        scale = self.upper - self.lower
+        self.parent = beta(self.shape_A, self.shape_B, loc=loc, scale=scale)
+        self.mean, self.variance, self.skewness, self.kurtosis = beta.stats(self.shape_A, self.shape_B, loc=loc,
+                                                                            scale=scale, moments='mvsk')
+        self.x_range_for_pdf = np.linspace(self.lower, self.upper, RECURRENCE_PDF_SAMPLES)
+
     def get_description(self):
         """
         A description of the beta distribution.
@@ -91,7 +109,9 @@ def get_recurrence_coefficients(self, order):
         :return:
             Recurrence coefficients associated with the beta distribution.
         """
-        ab =  jacobi_recurrence_coefficients(self.shape_parameter_A, self.shape_parameter_B, self.lower, self.upper, order)
+        ab = jacobi_recurrence_coefficients(self.shape_parameter_A - 1.0
+                                            , self.shape_parameter_B - 1.0
+                                            , self.lower, self.upper, order)
         return ab
 
     def get_icdf(self, xx):

equadratures/distributions/cauchy.py

@@ -16,16 +16,22 @@ class Cauchy(Distribution):
     """
     def __init__(self, location=None, scale=None):
         self.location = location
-        self.scale = scale
+        if scale is None:
+            self.scale = 1.0
+        else:
+            self.scale = scale
+
         self.bounds = np.array([-np.inf, np.inf])
         self.mean = np.nan
+        self.variance = np.nan
         self.skewness = np.nan
         self.kurtosis = np.nan
-        if self.scale is not None:
-            self.x_range_for_pdf = np.linspace(-15*self.scale, 15*self.scale, RECURRENCE_PDF_SAMPLES)
-            self.parent = cauchy(loc=self.location, scale=self.scale)
-            self.mean = np.mean(self.get_samples(m=1000))
-            self.variance = np.var(self.get_samples(m=1000))
+
+        self.x_range_for_pdf = np.linspace(-15*self.scale, 15*self.scale, RECURRENCE_PDF_SAMPLES)
+        self.parent = cauchy(loc=self.location, scale=self.scale)
+        # self.mean = np.mean(self.get_samples(m=1000))
+        # self.variance = np.var(self.get_samples(m=1000))
+
     def get_description(self):
         """
         A description of the Cauchy distribution.

equadratures/distributions/chebyshev.py

@@ -14,16 +14,28 @@ class Chebyshev(Distribution):
 		Upper bound of the support of the Chebyshev (arcsine) distribution.
     """
     def __init__(self, lower, upper):
-        self.lower = lower
-        self.upper = upper
-        self.bounds = np.array([0.0, 1.0])
-        if ( self.lower is not None ) and (self.upper is not None) :
-         	self.x_range_for_pdf = np.linspace(self.lower, self.upper, RECURRENCE_PDF_SAMPLES)
-        self.mean = 0.5
-        self.variance = 1.0/8.0
-        self.skewness = 0.0
+        if lower is None:
+            self.lower = 0.0
+        else:
+            self.lower = lower
+        if upper is None:
+            self.upper = 1.0
+        else:
+            self.upper = upper
+
+        if self.lower > self.upper:
+            raise ValueError('Invalid Beta distribution parameters. Lower should be smaller than upper.')
+
+        self.bounds = np.array([self.lower, self.upper])
+        self.x_range_for_pdf = np.linspace(self.lower, self.upper, RECURRENCE_PDF_SAMPLES)
+        loc = self.lower
+        scale = self.upper - self.lower
+
+        self.parent = arcsine(loc=loc, scale=scale)
+        self.mean, self.variance, self.skewness, self.kurtosis = self.parent.stats(moments='mvsk')
         self.shape_parameter_A = -0.5
         self.shape_parameter_B = -0.5
+
     def get_description(self):
         """
         A description of the Chebyshev (arcsine) distribution.
@@ -35,6 +47,7 @@ def get_description(self):
         """
         text = "is a Chebyshev or arcsine distribution that is characterised by its lower bound, which is"+str(self.lower)+" and its upper bound, which is"+str(self.upper)+"."
         return text
+
     def get_pdf(self, points=None):
         """
         A Chebyshev probability density function.
@@ -49,7 +62,7 @@ def get_pdf(self, points=None):
             Probability density values along the support of the Chebyshev (arcsine) distribution.
         """
         if points is not None:
-            return arcsine.pdf(points)
+            return self.parent.pdf(points)
         else:
             raise ValueError( 'Please digit an input for getPDF method')
     def get_cdf(self, points=None):
@@ -66,7 +79,7 @@ def get_cdf(self, points=None):
             Cumulative density values along the support of the Chebyshev (arcsine) distribution.
         """
         if points is not None:
-             return arcsine.cdf(points)
+            return self.parent.cdf(points)
         else:
             raise ValueError( 'Please digit an input for getCDF method')
     def get_recurrence_coefficients(self, order):
@@ -80,7 +93,7 @@ def get_recurrence_coefficients(self, order):
         :return:
             Recurrence coefficients associated with the Chebyshev distribution.
         """
-        ab =  jacobi_recurrence_coefficients(self.shape_parameter_A, self.shape_parameter_B, self.lower, self.upper, order)
+        ab = jacobi_recurrence_coefficients(self.shape_parameter_A, self.shape_parameter_B, self.lower, self.upper, order)
         return ab
     def get_icdf(self, xx):
         """
@@ -93,7 +106,7 @@ def get_icdf(self, xx):
         :return:
             Inverse cumulative density function values of the Arcisine distribution.
         """
-        return arcsine.ppf(xx)
+        return self.parent.ppf(xx)
     def get_samples(self, m=None):
         """
         Generates samples from the Arcsine distribution.
@@ -108,4 +121,4 @@ def get_samples(self, m=None):
             number = m
         else:
             number = 500000
-        return arcsine.rvs(size=number)
+        return self.parent.rvs(size=number)

equadratures/distributions/chi.py

@@ -12,21 +12,27 @@ class Chi(Distribution):
 		Degrees of freedom for the chi-squared distribution.
     """
     def __init__(self, dofs):
-        self.dofs = dofs
-        if self.dofs is not None:
-            if self.dofs == 1:
-                self.bounds = np.array([1e-15, np.inf])
-            else:
-                self.bounds = np.array([0.0, np.inf])
-            if self.dofs >= 1:
-                mean, var, skew, kurt = chi.stats(dofs, moments='mvsk')
-                self.parent = chi(dofs)
-                self.mean = mean
-                self.variance = var
-                self.skewness = skew
-                self.kurtosis = kurt
-                self.x_range_for_pdf = np.linspace(0.0, 10.0*self.mean,RECURRENCE_PDF_SAMPLES)
-                self.parent = chi(self.dofs)
+        if dofs is None:
+            self.dofs = 1
+        else:
+            self.dofs = dofs
+
+        if self.dofs < 0:
+            raise ValueError('Invalid parameter in chi distribution: dofs must be positive.')
+
+        if self.dofs == 1:
+            self.bounds = np.array([1e-15, np.inf])
+        else:
+            self.bounds = np.array([0.0, np.inf])
+
+        mean, var, skew, kurt = chi.stats(dofs, moments='mvsk')
+        self.mean = mean
+        self.variance = var
+        self.skewness = skew
+        self.kurtosis = kurt
+        self.x_range_for_pdf = np.linspace(0.0, 10.0*self.mean,RECURRENCE_PDF_SAMPLES)
+        self.parent = chi(self.dofs)
+
     def get_description(self):
         """
         A description of the Chi-squared distribution.

equadratures/distributions/chisquared.py

@@ -12,19 +12,23 @@ class Chisquared(Distribution):
 		Degrees of freedom for the chi-squared distribution.
     """
     def __init__(self, dofs):
-        self.dofs = dofs
-        if self.dofs is not None:
-            if self.dofs == 1:
-                self.bounds = np.array([1e-15, np.inf])
-            else:
-                self.bounds = np.array([0.0, np.inf])
-            if self.dofs >= 1:
-                self.mean = float(self.dofs)
-                self.variance = 2 * self.mean
-                self.skewness = np.sqrt(8.0 / self.mean)
-                self.kurtosis = 12.0/self.mean + 3.0
-                self.x_range_for_pdf = np.linspace(0.0, 10.0*self.mean,RECURRENCE_PDF_SAMPLES)
-                self.parent = chi2(self.dofs)
+        if dofs is None:
+            self.dofs = 1
+        else:
+            self.dofs = int(dofs)
+
+        if not isinstance(self.dofs, int) or self.dofs < 1:
+            raise ValueError('Invalid parameter in chi2 distribution: dofs must be positive integer.')
+
+        if self.dofs == 1:
+            self.bounds = np.array([1e-15, np.inf])
+        else:
+            self.bounds = np.array([0.0, np.inf])
+
+        self.parent = chi2(self.dofs)
+        self.mean, self.variance, self.skewness, self.kurtosis = self.parent.stats(moments='mvsk')
+        self.x_range_for_pdf = np.linspace(0.0, 10.0*self.mean,RECURRENCE_PDF_SAMPLES)
+
     def get_description(self):
         """
         A description of the Chi-squared distribution.

equadratures/distributions/exponential.py

@@ -13,17 +13,19 @@ class Exponential(Distribution):
 		Rate parameter of the Exponential distribution.
     """
     def __init__(self, rate=None):
-        self.rate = rate
+        if rate is None:
+            self.rate = 1.0
+        else:
+            self.rate = rate
+
         if (self.rate is not None) and (self.rate > 0.0):
-            #self.mean = 1. / self.rate
-            #self.variance = 1./(self.rate)**2
-            self.skewness = 2.0
-            self.kurtosis = 6.0
+            self.parent = expon(scale=1.0 / rate)
+            self.mean, self.variance, self.skewness, self.kurtosis = self.parent.stats(moments='mvsk')
             self.bounds = np.array([0.0, np.inf])
             self.x_range_for_pdf = np.linspace(0.0, 20*self.rate, RECURRENCE_PDF_SAMPLES)
-            self.parent = expon(scale=1.0/rate)
-            self.mean = self.parent.mean()
-            self.variance = self.parent.var()
+        else:
+            raise ValueError('Invalid parameters in exponential distribution. Rate should be positive.')
+
     def get_description(self):
         """
         A description of the Exponential distribution.

equadratures/distributions/gamma.py

@@ -1,7 +1,6 @@
 """The Gamma distribution."""
 from equadratures.distributions.template import Distribution
 import numpy as np
-from scipy.special import erf, erfinv, gamma, beta, betainc, gammainc
 from scipy.stats import gamma
 RECURRENCE_PDF_SAMPLES = 8000
 class Gamma(Distribution):
@@ -14,16 +13,22 @@ class Gamma(Distribution):
 		Scale parameter of the gamma distribution.
     """
     def __init__(self, shape=None, scale=None):
-        self.shape = shape
-        self.scale = scale
+        if shape is None:
+            self.shape = 1.0
+        else:
+            self.shape = shape
+        if scale is None:
+            self.scale = 1.0
+        else:
+            self.scale = scale
+
         self.bounds = np.array([0.0, np.inf])
-        if (self.shape is not None) and (self.scale is not None) and (self.shape > 0.0) :
-            self.mean = self.shape * self.scale
-            self.variance = self.shape * self.scale**2
-            self.skewness = 2.0 / np.sqrt(self.shape)
-            self.kurtosis = 6.0 / self.shape # double-check!
-            self.x_range_for_pdf = np.linspace(0, self.shape*self.scale*10, RECURRENCE_PDF_SAMPLES)
-            self.parent = gamma(a=self.shape, scale=self.scale)
+        if self.shape < 0 or self.scale < 0:
+            raise ValueError('Invalid parameters in Gamma distribution. Shape and Scale should be positive.')
+        self.parent = gamma(a=self.shape, scale=self.scale)
+        self.mean, self.variance, self.skewness, self.kurtosis = self.parent.stats(moments='mvsk')
+        self.x_range_for_pdf = np.linspace(0, self.scale*10, RECURRENCE_PDF_SAMPLES)
+
     def get_description(self):
         """
         A description of the gamma distribution.
@@ -35,6 +40,7 @@ def get_description(self):
         """
         text = "is a gamma distribution with a shape parameter of "+str(self.shape)+", and a scale parameter of "+str(self.scale)+"."
         return text
+
     def get_pdf(self, points=None):
         """
         A gamma probability density function.

equadratures/distributions/gaussian.py

@@ -18,12 +18,21 @@ class Gaussian(Distribution):
 		Variance of the Gaussian distribution.
     """
     def __init__(self, mean, variance):
-        self.mean = mean
-        self.variance = variance
-        if self.variance is not None:
-            self.sigma = np.sqrt(self.variance)
-            self.x_range_for_pdf = np.linspace(-15.0 * self.sigma, 15.0*self.sigma, RECURRENCE_PDF_SAMPLES) + self.mean
-            self.parent = norm(loc=self.mean, scale=self.sigma)
+        if mean is None:
+            self.mean = 0.0
+        else:
+            self.mean = mean
+        if variance is None:
+            self.variance = 1.0
+        else:
+            self.variance = variance
+
+        if self.variance <= 0:
+            raise ValueError('Invalid Gaussian distribution parameters. Variance should be positive.')
+
+        self.sigma = np.sqrt(self.variance)
+        self.x_range_for_pdf = np.linspace(-15.0 * self.sigma, 15.0*self.sigma, RECURRENCE_PDF_SAMPLES) + self.mean
+        self.parent = norm(loc=self.mean, scale=self.sigma)
         self.skewness = 0.0
         self.kurtosis = 0.0
         self.bounds = np.array([-np.inf, np.inf])