Python源码示例：scipy.stats.beta.ppf()
示例1
def _clopper_pearson_confint(self, counts: int, shots: int, alpha: float
) -> Tuple[float, float]:
"""Compute the Clopper-Pearson confidence interval for `shots` i.i.d. Bernoulli trials.
Args:
counts: The number of positive counts.
shots: The number of shots.
alpha: The confidence level for the confidence interval.
Returns:
The Clopper-Pearson confidence interval.
"""
lower, upper = 0, 1
# if counts == 0, the beta quantile returns nan
if counts != 0:
lower = beta.ppf(alpha / 2, counts, shots - counts + 1)
# if counts == shots, the beta quantile returns nan
if counts != shots:
upper = beta.ppf(1 - alpha / 2, counts + 1, shots - counts)
return lower, upper
示例2
def plot(self, data, size, newdata=None):
data = np.array(data)
numsample = len(data)
colmean = np.mean(data, axis=0)
matcov = np.cov(data.T)
matinv = np.linalg.inv(matcov)
values = []
for sample in data:
dif = sample - colmean
value = matinv.dot(dif.T).dot(dif)
values.append(value)
cl = ((numsample - 1)**2) / numsample
lcl = cl * beta.ppf(0.00135, size / 2, (numsample - size - 1) / 2)
center = cl * beta.ppf(0.5, size / 2, (numsample - size - 1) / 2)
ucl = cl * beta.ppf(0.99865, size / 2, (numsample - size - 1) / 2)
return (values, center, lcl, ucl, self._title)
示例3
def get_probility(x, N, a=0.05):
'''
计算 N 此重复实验，事件 A 发生 x 次时，事件 A 的发生概率
Args:
x : 事件发生次数
N : 总实验次数
a : 设置显著性水平 1-a，默认 0.05
Returns:
P, E, std, (low, high)
事件 A 发生概率 P 的最概然取值、期望、以及其置信区间
'''
P = 1.0 * x / N
E = (x + 1.0) / (N + 2.0)
std = np.sqrt(E * (1 - E) / (N + 3))
low = beta.ppf(0.5 * a, x + 1, N - x + 1)
high = beta.ppf(1 - 0.5 * a, x + 1, N - x + 1)
return P, E, std, (low, high)
示例4
def _clopper_pearson_interval(self):
r"""
Computes the Clopper-Pearson 'exact' confidence interval.
References
----------
Wikipedia contributors. (2018, July 14). Binomial proportion confidence interval.
In Wikipedia, The Free Encyclopedia. Retrieved 00:40, August 15, 2018,
from https://en.wikipedia.org/w/index.php?title=Binomial_proportion_confidence_interval&oldid=850256725
"""
p = self.x / self.n
if self.alternative == 'less':
lower_bound = 0.0
upper_bound = beta.ppf(1 - self.alpha, self.x + 1, self.n - self.x)
elif self.alternative == 'greater':
upper_bound = 1.0
lower_bound = beta.ppf(self.alpha, self.x, self.n - self.x + 1)
else:
lower_bound = beta.ppf(self.alpha / 2, self.x, self.n - self.x + 1)
upper_bound = beta.ppf(1 - self.alpha / 2, self.x + 1, self.n - self.x)
clopper_pearson_interval = {
'probability of success': p,
'conf level': 1 - self.alpha,
'interval': (lower_bound, upper_bound)
}
return clopper_pearson_interval
示例5
def transform_normal(x,hyperparameters):
mu,sigma = hyperparameters
return norm.ppf(x,loc=mu,scale=sigma)
示例6
def transform_beta(x,hyperparameters):
a,b = hyperparameters
return beta.ppf(x,a,b)
示例7
def transform_exponential(x,hyperparameters):
a = hyperparameters
return gamma.ppf(x, a)
示例8
def transform_truncated_normal(x,hyperparameters):
mu, sigma, a, b = hyperparameters
ar, br = (a - mu) / sigma, (b - mu) / sigma
return truncnorm.ppf(x,ar,br,loc=mu,scale=sigma)
示例9
def transform_normal(x,mu,sigma):
return norm.ppf(x,loc=mu,scale=sigma)
示例10
def transform_exponential(x,a=1.):
return gamma.ppf(x, a)
示例11
def transform_truncated_normal(x,mu,sigma,a=0.,b=1.):
ar, br = (a - mu) / sigma, (b - mu) / sigma
return truncnorm.ppf(x,ar,br,loc=mu,scale=sigma)
示例12
def binom_interval(success, total, confint=0.95):
"""
Compute two-sided binomial confidence interval in Python. Based on R's binom.test.
"""
from scipy.stats import beta
quantile = (1 - confint) / 2.
lower = beta.ppf(quantile, success, total - success + 1)
upper = beta.ppf(1 - quantile, success + 1, total - success)
return (lower, upper)
示例13
def __init__(self, n, x, p=0.5, alternative='two-sided', alpha=0.05, continuity=True):
if x > n:
raise ValueError('number of successes cannot be greater than number of trials.')
if p > 1.0:
raise ValueError('expected probability of success cannot be greater than 1.')
if alternative not in ('two-sided', 'greater', 'less'):
raise ValueError("'alternative must be one of 'two-sided' (default), 'greater', or 'less'.")
self.n = n
self.x = x
self.p = p
self.q = 1.0 - self.p
self.alpha = alpha
self.alternative = alternative
self.continuity = continuity
self.p_value = self._p_value()
if self.alternative == 'greater':
self.z = norm.ppf(self.alpha)
elif self.alternative == 'less':
self.z = norm.ppf(1 - self.alpha)
else:
self.z = norm.ppf(1 - self.alpha / 2)
self.clopper_pearson_interval = self._clopper_pearson_interval()
self.wilson_score_interval = self._wilson_score_interval()
self.agresti_coull_interval = self._agresti_coull_interval()
self.arcsine_transform_interval = self._arcsine_transform_interval()
self.test_summary = {
'Number of Successes': self.x,
'Number of Trials': self.n,
'p-value': self.p_value,
'alpha': self.alpha,
'intervals': {
'Clopper-Pearson': self.clopper_pearson_interval,
'Wilson Score': self.wilson_score_interval,
'Agresti-Coull': self.agresti_coull_interval,
'Arcsine Transform': self.arcsine_transform_interval
}
}
示例14
def _conf_int(self):
r"""
Computes the confidence interval.
Returns
-------
intervals : tuple
Tuple containing the low and high confidence interval.
Notes
-----
Confidence intervals are reported as a proportion, denoted by :math:`1 - \alpha`, which
represents the ratio of intervals that would contain the population parameter if samples
were redrawn and tested with the same procedure. A confidence level is the interval reported
as a percentage, :math:`(1 - \alpha) * 100\%`. The width of the :math:`(1 - \alpha) * 100\%`
interval has several dependencies:
1. The confidence level. As :math:`1 - \alpha` increases, so does the width of the interval.
2. As the sample size :math:`n` increases, the smaller the standard error and thus a narrower interval.
3. If the standard deviation is large, then the interval will also be large.
The computation of the intervals uses Welch's t-interval which extends the two-sample pooled
t-interval for unequal population variances and sample sizes.
.. math::
\left(\bar{X_1} - \bar{X_2}\right) \pm t_{\alpha / 2, r} \sqrt{\frac{s_{x_1}^2}{n_1} + \frac{s_{x_2}^2}{n_2}}
References
----------
Rencher, A. C., & Christensen, W. F. (2012). Methods of multivariate analysis (3rd Edition).
"""
xn, xvar, xbar = self.sample_statistics[self._y1_summary_stat_name]['obs'], \
self.sample_statistics[self._y1_summary_stat_name]['variance'], \
self.sample_statistics[self._y1_summary_stat_name]['mean']
if self.y2 is not None:
yn, yvar, ybar = self.sample_statistics['Sample 2']['obs'], \
self.sample_statistics['Sample 2']['variance'], \
self.sample_statistics['Sample 2']['mean']
low_interval = (xbar - ybar) + t.ppf(self.alpha / 2., self.parameter) * np.sqrt(xvar / xn + yvar / yn)
high_interval = (xbar - ybar) - t.ppf(self.alpha / 2., self.parameter) * np.sqrt(xvar / xn + yvar / yn)
else:
low_interval = xbar + 1.96 * np.sqrt(xvar / xn)
high_interval = xbar - 1.96 * np.sqrt(xvar / xn)
if self.alternative == 'greater':
intervals = np.inf, float(high_interval)
elif self.alternative == 'less':
intervals = -np.inf, float(low_interval)
else:
intervals = float(low_interval), float(high_interval)
return intervals
示例15
def beta_confidence_intervals(ci_X, ntrials, ci=0.95):
"""
Compute confidence intervals of beta distributions.
Parameters
----------
ci_X : numpy.array
Computed confidence interval estimate from `ntrials` experiments
ntrials : int
The number of trials that were run.
ci : float, optional, default=0.95
Confidence interval to report (e.g. 0.95 for 95% confidence interval)
Returns
-------
Plow : float
The lower bound of the symmetric confidence interval.
Phigh : float
The upper bound of the symmetric confidence interval.
Examples
--------
>>> ci_X = np.random.rand(10,10)
>>> ntrials = 100
>>> [Plow, Phigh] = beta_confidence_intervals(ci_X, ntrials)
"""
# Compute low and high confidence interval for symmetric CI about mean.
ci_low = 0.5 - ci/2;
ci_high = 0.5 + ci/2;
# Compute for every element of ci_X.
from scipy.stats import beta
Plow = ci_X * 0.0;
Phigh = ci_X * 0.0;
for i in range(ci_X.shape[0]):
for j in range(ci_X.shape[1]):
Plow[i,j] = beta.ppf(ci_low, a = ci_X[i,j] * ntrials, b = (1-ci_X[i,j]) * ntrials);
Phigh[i,j] = beta.ppf(ci_high, a = ci_X[i,j] * ntrials, b = (1-ci_X[i,j]) * ntrials);
return [Plow, Phigh]
示例16
def plot(self, data, size, newdata=None):
sizes = data[:, 0]
sample = data[:, 1:]
samples = dict()
for n, value in zip(sizes, sample):
if n in samples:
samples[n] = np.vstack([samples[n], value])
else:
samples[n] = value
m = len(samples.keys())
n = len(samples[1])
p = len(samples[1].T)
variance, S = [], []
for i in range(m):
mat = np.cov(samples[i + 1].T, ddof=1)
variance.append(mat.diagonal())
S.append(cova(mat))
variance, S = np.array(variance), np.array(S)
means = np.array([samples[xs + 1].mean(axis=0) for xs in range(m)])
means_total = means.mean(axis=0)
Smat = var_cov(variance, S)
Smat_inv = np.linalg.inv(Smat)
values = []
for i in range(m):
a = means[i] - means_total
values.append(5 * a @ Smat_inv @ a.T)
p1 = (p * (m - 1) * (n - 1))
p2 = (m * n - m - p + 1)
lcl = (p1 / p2) * f.ppf(0.00135, p, p2)
center = (p1 / p2) * f.ppf(0.50, p, p2)
ucl = (p1 / p2) * f.ppf(0.99865, p, p2)
return (values, center, lcl, ucl, self._title)
