From d3064e8be17aebcee87cad0365058b59f93423ac Mon Sep 17 00:00:00 2001 From: Gerald Teschl Date: Wed, 5 Jul 2023 13:39:18 +0200 Subject: [PATCH 1/7] Added exponential and gamma distribution --- .../probability/probability_distribution.pyx | 74 +++++++++++++++++-- 1 file changed, 68 insertions(+), 6 deletions(-) diff --git a/src/sage/probability/probability_distribution.pyx b/src/sage/probability/probability_distribution.pyx index 4e27a861f6e..1ddc1b00ce2 100644 --- a/src/sage/probability/probability_distribution.pyx +++ b/src/sage/probability/probability_distribution.pyx @@ -46,8 +46,7 @@ import sage.rings.real_double from sage.modules.free_module_element import vector # TODO: Add more distributions available in gsl -# available but not currently wrapped are exponential, laplace, cauchy, landau, gamma, -# gamma, beta logistic. +# available but not currently wrapped are laplace, cauchy, landau, logistic. cdef enum: uniform @@ -61,6 +60,8 @@ cdef enum: exppow weibull beta + exponential + gamma cdef class ProbabilityDistribution: """ @@ -495,6 +496,33 @@ cdef class RealDistribution(ProbabilityDistribution): sage: T.cum_distribution_function(1) 1.0 + The exponential distribution has one parameter ``mu``:: + + sage: mu = 2 + sage: T = RealDistribution('exponential', mu) + sage: s = T.get_random_element() + sage: 0 <= s + True + sage: s.parent() + Real Double Field + sage: T.distribution_function(0) + 0.0 + + The gamma distribution has two parameters ``a`` and ``b``:: + + sage: a = 2 + sage: b = 2 + sage: T = RealDistribution('gamma', [a, b]) + sage: s = T.get_random_element() + sage: 0 <= s + True + sage: s.parent() + Real Double Field + sage: T.distribution_function(0) + 0.0 + sage: T.cum_distribution_function(1) + 1.0 + The weibull distribution has two parameters ``a`` and ``b``:: sage: a = 1 @@ -676,7 +704,11 @@ cdef class RealDistribution(ProbabilityDistribution): result = gsl_ran_weibull(self.r, self.parameters[0], self.parameters[1]) elif self.distribution_type == beta: result = gsl_ran_beta(self.r, self.parameters[0], self.parameters[1]) - else: + elif self.distribution_type == exponential: + result = gsl_ran_exponential(self.r, self.parameters[0]) + elif self.distribution_type == gamma: + result = gsl_ran_gamma(self.r, self.parameters[0], self.parameters[1]) + else: raise TypeError("Not a supported probability distribution") return sage.rings.real_double.RDF(result) @@ -723,7 +755,6 @@ cdef class RealDistribution(ProbabilityDistribution): self.parameters[1] = float(parameters[1]) self.distribution_type = pareto elif name == 'rayleigh': - self.distribution_type = rayleigh try: float(parameters) except Exception: @@ -804,7 +835,26 @@ cdef class RealDistribution(ProbabilityDistribution): self.parameters[0] = float(parameters[0]) self.parameters[1] = float(parameters[1]) self.distribution_type = beta - else: + elif name == 'exponential': + try: + float(parameters) + except Exception: + raise TypeError("exponential distribution requires parameter mu coercible to float") + self.parameters = sig_malloc(sizeof(double)) + self.parameters[0] = float(parameters) + self.distribution_type = exponential + elif name == 'gamma': + if len(parameters) != 2: + raise TypeError("gamma distribution requires two real parameters") + try: + map(float, parameters) + except Exception: + raise TypeError("gamma distribution requires real parameters") + self.parameters = sig_malloc(sizeof(double)*2) + self.parameters[0] = float(parameters[0]) + self.parameters[1] = float(parameters[1]) + self.distribution_type = gamma + else: raise TypeError("Not a supported probability distribution") self.name = name @@ -869,6 +919,10 @@ cdef class RealDistribution(ProbabilityDistribution): return sage.rings.real_double.RDF(gsl_ran_weibull_pdf(x, self.parameters[0], self.parameters[1])) elif self.distribution_type == beta: return sage.rings.real_double.RDF(gsl_ran_beta_pdf(x, self.parameters[0], self.parameters[1])) + elif self.distribution_type == exponential: + return sage.rings.real_double.RDF(gsl_ran_exponential_pdf(x, self.parameters[0])) + elif self.distribution_type == gamma: + return sage.rings.real_double.RDF(gsl_ran_gamma_pdf(x, self.parameters[0], self.parameters[1])) else: raise TypeError("Not a supported probability distribution") @@ -905,7 +959,11 @@ cdef class RealDistribution(ProbabilityDistribution): return sage.rings.real_double.RDF(gsl_cdf_weibull_P(x, self.parameters[0], self.parameters[1])) elif self.distribution_type == beta: return sage.rings.real_double.RDF(gsl_cdf_beta_P(x, self.parameters[0], self.parameters[1])) - else: + elif self.distribution_type == exponential: + return sage.rings.real_double.RDF(gsl_cdf_exponential_P(x, self.parameters[0]])) + elif self.distribution_type == gamma: + return sage.rings.real_double.RDF(gsl_cdf_gamma_P(x, self.parameters[0], self.parameters[1])) + else: raise TypeError("Not a supported probability distribution") def cum_distribution_function_inv(self, x): @@ -942,6 +1000,10 @@ cdef class RealDistribution(ProbabilityDistribution): return sage.rings.real_double.RDF(gsl_cdf_weibull_Pinv(x, self.parameters[0], self.parameters[1])) elif self.distribution_type == beta: return sage.rings.real_double.RDF(gsl_cdf_beta_Pinv(x, self.parameters[0], self.parameters[1])) + elif self.distribution_type == exponential: + return sage.rings.real_double.RDF(gsl_cdf_exponential_Pinv(x, self.parameters[0])) + elif self.distribution_type == gamma: + return sage.rings.real_double.RDF(gsl_cdf_gamma_Pinv(x, self.parameters[0], self.parameters[1])) else: raise TypeError("Not a supported probability distribution") From e5edc6473e44c183170e6be941e123faa34d81f6 Mon Sep 17 00:00:00 2001 From: Gerald Teschl Date: Wed, 5 Jul 2023 13:42:25 +0200 Subject: [PATCH 2/7] Update probability_distribution.pyx --- src/sage/probability/probability_distribution.pyx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/sage/probability/probability_distribution.pyx b/src/sage/probability/probability_distribution.pyx index 1ddc1b00ce2..8fa55b44dfb 100644 --- a/src/sage/probability/probability_distribution.pyx +++ b/src/sage/probability/probability_distribution.pyx @@ -706,9 +706,9 @@ cdef class RealDistribution(ProbabilityDistribution): result = gsl_ran_beta(self.r, self.parameters[0], self.parameters[1]) elif self.distribution_type == exponential: result = gsl_ran_exponential(self.r, self.parameters[0]) - elif self.distribution_type == gamma: + elif self.distribution_type == gamma: result = gsl_ran_gamma(self.r, self.parameters[0], self.parameters[1]) - else: + else: raise TypeError("Not a supported probability distribution") return sage.rings.real_double.RDF(result) @@ -854,7 +854,7 @@ cdef class RealDistribution(ProbabilityDistribution): self.parameters[0] = float(parameters[0]) self.parameters[1] = float(parameters[1]) self.distribution_type = gamma - else: + else: raise TypeError("Not a supported probability distribution") self.name = name @@ -963,7 +963,7 @@ cdef class RealDistribution(ProbabilityDistribution): return sage.rings.real_double.RDF(gsl_cdf_exponential_P(x, self.parameters[0]])) elif self.distribution_type == gamma: return sage.rings.real_double.RDF(gsl_cdf_gamma_P(x, self.parameters[0], self.parameters[1])) - else: + else: raise TypeError("Not a supported probability distribution") def cum_distribution_function_inv(self, x): From 90603da019638dc32f3a463792d55983143e9159 Mon Sep 17 00:00:00 2001 From: Gerald Teschl Date: Wed, 5 Jul 2023 13:44:28 +0200 Subject: [PATCH 3/7] Update probability_distribution.pyx --- src/sage/probability/probability_distribution.pyx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/sage/probability/probability_distribution.pyx b/src/sage/probability/probability_distribution.pyx index 8fa55b44dfb..fe08e9a5f16 100644 --- a/src/sage/probability/probability_distribution.pyx +++ b/src/sage/probability/probability_distribution.pyx @@ -843,7 +843,7 @@ cdef class RealDistribution(ProbabilityDistribution): self.parameters = sig_malloc(sizeof(double)) self.parameters[0] = float(parameters) self.distribution_type = exponential - elif name == 'gamma': + elif name == 'gamma': if len(parameters) != 2: raise TypeError("gamma distribution requires two real parameters") try: @@ -854,7 +854,7 @@ cdef class RealDistribution(ProbabilityDistribution): self.parameters[0] = float(parameters[0]) self.parameters[1] = float(parameters[1]) self.distribution_type = gamma - else: + else: raise TypeError("Not a supported probability distribution") self.name = name From 4c8c0d0dcc3bacc5d515d1f6f4b5017d281d0ae7 Mon Sep 17 00:00:00 2001 From: Gerald Teschl Date: Wed, 5 Jul 2023 13:46:47 +0200 Subject: [PATCH 4/7] Update probability_distribution.pyx --- src/sage/probability/probability_distribution.pyx | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/sage/probability/probability_distribution.pyx b/src/sage/probability/probability_distribution.pyx index fe08e9a5f16..ed45990cedb 100644 --- a/src/sage/probability/probability_distribution.pyx +++ b/src/sage/probability/probability_distribution.pyx @@ -520,8 +520,6 @@ cdef class RealDistribution(ProbabilityDistribution): Real Double Field sage: T.distribution_function(0) 0.0 - sage: T.cum_distribution_function(1) - 1.0 The weibull distribution has two parameters ``a`` and ``b``:: From 089da87d9f332046effeeed7add919da4c20aa04 Mon Sep 17 00:00:00 2001 From: Gerald Teschl Date: Wed, 12 Jul 2023 11:35:47 +0200 Subject: [PATCH 5/7] Update probability_distribution.pyx --- src/sage/probability/probability_distribution.pyx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/sage/probability/probability_distribution.pyx b/src/sage/probability/probability_distribution.pyx index ed45990cedb..81966188dc3 100644 --- a/src/sage/probability/probability_distribution.pyx +++ b/src/sage/probability/probability_distribution.pyx @@ -958,7 +958,7 @@ cdef class RealDistribution(ProbabilityDistribution): elif self.distribution_type == beta: return sage.rings.real_double.RDF(gsl_cdf_beta_P(x, self.parameters[0], self.parameters[1])) elif self.distribution_type == exponential: - return sage.rings.real_double.RDF(gsl_cdf_exponential_P(x, self.parameters[0]])) + return sage.rings.real_double.RDF(gsl_cdf_exponential_P(x, self.parameters[0])) elif self.distribution_type == gamma: return sage.rings.real_double.RDF(gsl_cdf_gamma_P(x, self.parameters[0], self.parameters[1])) else: From 5239b7029a0e9092fceb83273809e7f80dc8a87b Mon Sep 17 00:00:00 2001 From: Gerald Teschl Date: Thu, 24 Aug 2023 09:12:50 +0200 Subject: [PATCH 6/7] Update probability_distribution.pyx --- src/sage/probability/probability_distribution.pyx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/sage/probability/probability_distribution.pyx b/src/sage/probability/probability_distribution.pyx index f0792385e12..005f50b992f 100644 --- a/src/sage/probability/probability_distribution.pyx +++ b/src/sage/probability/probability_distribution.pyx @@ -507,7 +507,7 @@ cdef class RealDistribution(ProbabilityDistribution): Real Double Field sage: T.distribution_function(0) 0.0 - + The gamma distribution has two parameters ``a`` and ``b``:: sage: a = 2 From 3be7da857f5032d7ce3c451748530a72e594f25d Mon Sep 17 00:00:00 2001 From: Gerald Teschl Date: Thu, 24 Aug 2023 12:45:29 +0200 Subject: [PATCH 7/7] Update probability_distribution.pyx --- src/sage/probability/probability_distribution.pyx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/sage/probability/probability_distribution.pyx b/src/sage/probability/probability_distribution.pyx index 005f50b992f..bfc91d94ec0 100644 --- a/src/sage/probability/probability_distribution.pyx +++ b/src/sage/probability/probability_distribution.pyx @@ -506,7 +506,7 @@ cdef class RealDistribution(ProbabilityDistribution): sage: s.parent() Real Double Field sage: T.distribution_function(0) - 0.0 + 0.5 The gamma distribution has two parameters ``a`` and ``b``::