Skip to content

Commit 3aaef18

Browse files
author
Nick
committed
0.90.3
* Fixed error in C-box class * Manually created student-t distribution function
1 parent c26da58 commit 3aaef18

3 files changed

Lines changed: 72 additions & 43 deletions

File tree

pba/pbox.py

Lines changed: 11 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@
22
import itertools
33
from typing import *
44
from warnings import *
5-
5+
from math import isclose
66
import numpy as np
77
from matplotlib import pyplot as plt
88
from scipy import interpolate as spi
@@ -86,6 +86,12 @@ def _check_steps(a, b):
8686

8787
def _check_moments(left, right, steps, mean, var):
8888

89+
def nearly(x: Interval, y: Interval):
90+
if isclose(x.left,y.left) and isclose(x.right, y.right):
91+
return True
92+
else:
93+
return False
94+
8995
def _sideVariance(w, mu):
9096
if not isinstance(w, np.ndarray):
9197
w = np.array(w)
@@ -94,7 +100,7 @@ def _sideVariance(w, mu):
94100
return max(0, np.mean((w - mu) ** 2))
95101

96102
cmean = Interval(np.mean(left), np.mean(right))
97-
if not mean.equiv(cmean) and not mean.equiv(Interval(-np.inf, np.inf)):
103+
if not nearly(mean, cmean) and not mean.equiv(Interval(-np.inf, np.inf)):
98104
warn("Mean specified does not match calculated mean. Using calculated mean")
99105
print(f"Specified mean: {mean}, calculated mean: {cmean}")
100106
mean = cmean
@@ -142,7 +148,7 @@ def _sideVariance(w, mu):
142148

143149
cvar = Interval(vl, vr)
144150

145-
if not var.equiv(cvar) and not var.equiv(Interval(0, np.inf)):
151+
if not nearly(var, cvar) and not var.equiv(Interval(0, np.inf)):
146152
warn(
147153
"Variance specified does not match calculated variance. Using calculated variance"
148154
)
@@ -157,10 +163,10 @@ def _sideVariance(w, mu):
157163
def _arithmetic(a, b, method, op, enforce_steps=True, interpolation_method="linear"):
158164
# * If enforce_steps is True, the number of steps in the returned p-box is the maximum of the number of steps in a and b.
159165

160-
if b.__class__.__name__ == "Interval":
166+
if isinstance(b, Interval):
161167
other = Pbox(other, steps=a.steps)
162168

163-
if b.__class__.__name__ == "Pbox":
169+
if isinstance(b, Pbox) or issubclass(b.__class__,Pbox):
164170

165171
a, b = _check_steps(a, b)
166172

pba/pbox_constructors/distributions.py

Lines changed: 60 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -171,44 +171,67 @@ def foldnorm(mu, s, steps=200):
171171
)
172172

173173

174-
# def frechet_r(*args, steps = 200):
175-
# args = list(args)
176-
# for i in range(0,len(args)):
177-
# if args[i].__class__.__name__ != 'Interval':
178-
# args[i] = Interval(args[i])
179-
180-
# Left, Right, mean, var = __get_bounds('frechet_r',steps,*args)
181-
182-
# return Pbox(
183-
# Left,
184-
# Right,
185-
# steps = steps,
186-
# shape = 'frechet_r',
187-
# mean_left = mean.left,
188-
# mean_right = mean.right,
189-
# var_left = var.left,
190-
# var_right = var.right
191-
# )
192-
193-
# def frechet_l(*args, steps = 200):
194-
# args = list(args)
195-
# for i in range(0,len(args)):
196-
# if args[i].__class__.__name__ != 'Interval':
197-
# args[i] = Interval(args[i])
198-
199-
# Left, Right, mean, var = __get_bounds('frechet_l',steps,*args)
200-
201-
# return Pbox(
202-
# Left,
203-
# Right,
204-
# steps = steps,
205-
# shape = 'frechet_l',
206-
# mean_left = mean.left,
207-
# mean_right = mean.right,
208-
# var_left = var.left,
209-
# var_right = var.right
210-
# )
174+
def t(mean, std, df, steps=200):
175+
176+
def t_from_mean_std_df(mean: float, std: float, df: float):
211177

178+
scale = std * np.sqrt((df - 2.0) / df)
179+
return stats.t(df=df, loc=mean, scale=scale)
180+
181+
if mean.__class__.__name__ != "Interval":
182+
mean = Interval(mean)
183+
if std.__class__.__name__ != "Interval":
184+
std = Interval(std)
185+
186+
if df <= 2:
187+
raise ValueError(f"df must be > 2 to have a finite std (got df={df}).")
188+
189+
x = np.linspace(0.0001, 0.9999, steps)
190+
191+
new_args = [
192+
[mean.lo() / std.lo()],
193+
[mean.hi() / std.lo()],
194+
[mean.lo() / std.hi()],
195+
[mean.hi() / std.hi()],
196+
]
197+
198+
bounds = []
199+
200+
mean_hi = -np.inf
201+
mean_lo = np.inf
202+
var_lo = np.inf
203+
var_hi = 0
204+
205+
for mu,sigma in new_args:
206+
207+
t_ = t_from_mean_std_df(mu,sigma,df)
208+
bounds.append(t_.ppf(x))
209+
bmean, bvar = t_.stats(moments="mv")
210+
211+
if bmean < mean_lo:
212+
mean_lo = bmean
213+
if bmean > mean_hi:
214+
mean_hi = bmean
215+
if bvar > var_hi:
216+
var_hi = bvar
217+
if bvar < var_lo:
218+
var_lo = bvar
219+
220+
mean = Interval(np.float64(mean_lo), np.float64(mean_hi))
221+
222+
Left = [min([b[i] for b in bounds]) for i in range(steps)]
223+
Right = [max([b[i] for b in bounds]) for i in range(steps)]
224+
225+
return Pbox(
226+
Left,
227+
Right,
228+
steps=steps,
229+
shape="foldnorm",
230+
mean_left=mean.left,
231+
mean_right=mean.right,
232+
var_left=var.left,
233+
var_right=var.right,
234+
)
212235

213236
def trapz(a, b, c, d, steps=200):
214237
if a.__class__.__name__ != "Interval":

pba/pbox_constructors/parametric.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -104,7 +104,7 @@
104104
"recipinvgauss",
105105
"semicircular",
106106
"skewnorm",
107-
"t",
107+
# "t",
108108
# "trapz",
109109
"triang",
110110
"truncexpon",

0 commit comments

Comments
 (0)