def ln_prior(p, joker_params):
P, phi0, ecc, omega, s, K, *v_terms = p
lnp = 0.
# TODO: more repeated code here and hard-coded priors
if ecc < 0 or ecc > 1:
return -np.inf
lnp += beta.logpdf(ecc, 0.867, 3.03) # Kipping et al. 2013
# TODO: do we need P_min, P_max here?
if not joker_params._fixed_jitter:
# DFM's idea: wide, Gaussian prior in log(s^2)
# lnp += norm.logpdf(np.log(s), ) # TODO: put in hyper-parameters
# TODO:
pass
return lnp
python类logpdf()的实例源码
def logpdf(self, rowid, targets, constraints=None, inputs=None):
constraints = self.populate_constraints(rowid, targets, constraints)
# XXX Disable logpdf queries without constraints.
if inputs:
raise ValueError('Prohibited inputs: %s' % (inputs,))
if not constraints:
raise ValueError('Provide at least one constraint: %s'
% (constraints,))
self._validate_simulate_logpdf(rowid, targets, constraints)
# Retrieve the dataset and neighborhoods.
dataset, neighborhoods = self._find_neighborhoods(targets, constraints)
models = [self._create_local_model_joint(targets, dataset[n])
for n in neighborhoods]
# Compute logpdf in each neighborhood and simple average.
lp = [m.logpdf(targets) for m in models]
return gu.logsumexp(lp) - np.log(len(models))
def test_aic_fail_no_data_points():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1)
aics = c.comparison.aic()
assert len(aics) == 1
assert aics[0] is None
def test_aic_fail_no_num_params():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_eff_data_points=1000)
aics = c.comparison.aic()
assert len(aics) == 1
assert aics[0] is None
def test_aic_0():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1, num_eff_data_points=1000)
aics = c.comparison.aic()
assert len(aics) == 1
assert aics[0] == 0
def test_aic_posterior_dependence():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
p2 = norm.logpdf(d, scale=2)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1, num_eff_data_points=1000)
c.add_chain(d, posterior=p2, num_free_params=1, num_eff_data_points=1000)
aics = c.comparison.aic()
assert len(aics) == 2
assert aics[0] == 0
expected = 2 * np.log(2)
assert np.isclose(aics[1], expected, atol=1e-3)
def test_aic_parameter_dependence():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1, num_eff_data_points=1000)
c.add_chain(d, posterior=p, num_free_params=2, num_eff_data_points=1000)
aics = c.comparison.aic()
assert len(aics) == 2
assert aics[0] == 0
expected = 1 * 2 + (2.0 * 2 * 3 / (1000 - 2 - 1)) - (2.0 * 1 * 2 / (1000 - 1 - 1))
assert np.isclose(aics[1], expected, atol=1e-3)
def test_bic_fail_no_data_points():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1)
bics = c.comparison.bic()
assert len(bics) == 1
assert bics[0] is None
def test_bic_fail_no_num_params():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_eff_data_points=1000)
bics = c.comparison.bic()
assert len(bics) == 1
assert bics[0] is None
def test_bic_0():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1, num_eff_data_points=1000)
bics = c.comparison.bic()
assert len(bics) == 1
assert bics[0] == 0
def test_bic_posterior_dependence():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
p2 = norm.logpdf(d, scale=2)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1, num_eff_data_points=1000)
c.add_chain(d, posterior=p2, num_free_params=1, num_eff_data_points=1000)
bics = c.comparison.bic()
assert len(bics) == 2
assert bics[0] == 0
expected = 2 * np.log(2)
assert np.isclose(bics[1], expected, atol=1e-3)
def test_bic_parameter_dependence():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=1, num_eff_data_points=1000)
c.add_chain(d, posterior=p, num_free_params=2, num_eff_data_points=1000)
bics = c.comparison.bic()
assert len(bics) == 2
assert bics[0] == 0
expected = np.log(1000)
assert np.isclose(bics[1], expected, atol=1e-3)
def test_bic_data_dependence2():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p, num_free_params=2, num_eff_data_points=1000)
c.add_chain(d, posterior=p, num_free_params=3, num_eff_data_points=500)
bics = c.comparison.bic()
assert len(bics) == 2
assert bics[0] == 0
expected = 3 * np.log(500) - 2 * np.log(1000)
assert np.isclose(bics[1], expected, atol=1e-3)
def test_dic_0():
d = norm.rvs(size=1000)
p = norm.logpdf(d)
c = ChainConsumer()
c.add_chain(d, posterior=p)
dics = c.comparison.dic()
assert len(dics) == 1
assert dics[0] == 0
def test_dic_posterior_dependence():
d = norm.rvs(size=1000000)
p = norm.logpdf(d)
p2 = norm.logpdf(d, scale=2)
c = ChainConsumer()
c.add_chain(d, posterior=p)
c.add_chain(d, posterior=p2)
bics = c.comparison.dic()
assert len(bics) == 2
assert bics[1] == 0
dic1 = 2 * np.mean(-2 * p) + 2 * norm.logpdf(0)
dic2 = 2 * np.mean(-2 * p2) + 2 * norm.logpdf(0, scale=2)
assert np.isclose(bics[0], dic1 - dic2, atol=1e-3)
def test_lpdf(x, loc, scale):
aae(lpdf(x, loc, scale), norm.logpdf(x, loc, scale))
def test_lpdf_1d(x, loc, scale):
aae(lpdf_1d(x, loc, scale), norm.logpdf(x, loc, scale))
def test_lpdf_3d(x, loc, scale):
aae(lpdf_3d(x, loc, scale), norm.logpdf(x, loc, scale))
def test_lpdf_std(x):
aae(lpdf_std(x), norm.logpdf(x))
def test_preds_ll(alpha, mu, gamma, err, num, w):
current_impl = Lvm.preds_ll(alpha, mu, gamma, err, num, w)
simple_impl = np.nansum(w * norm.logpdf(num, mu+gamma*alpha, err))
simple_impl += np.sum(norm.logpdf(alpha))
assert_approx_equal(current_impl, simple_impl)
def phoDurDist(x,dur_mean,proportion_std = 0.35):
'''
build gaussian distribution which has mean = dur_mean, std = dur_mean*proportion_std
:param dur_mean: estimated from the centroid duration
:param proportion_std:
:return:
'''
prob = norm.logpdf(x,dur_mean,dur_mean*proportion_std)
return prob
def _create_local_model_joint(self, targets, dataset):
assert all(q in self.outputs for q in targets)
assert dataset.shape[1] == len(targets)
lookup = {
'numerical': self._create_local_model_numerical,
'categorical': self._create_local_model_categorical,
'nominal': self._create_local_model_categorical,
}
models = {
q: lookup[self.stattypes[self.outputs.index(q)]](q, dataset[:,i])
for i, q in enumerate(targets)}
simulate = lambda q, N=None: {c: models[c].simulate(N) for c in q}
logpdf = lambda q: sum(models[c].logpdf(x) for c,x in q.iteritems())
return LocalGpm(simulate, logpdf)
def _create_local_model_numerical(self, q, locality):
assert q not in self.levels
(mu, std) = (np.mean(locality), max(np.std(locality), .01))
simulate = lambda N=None: self.rng.normal(mu, std, size=N)
logpdf = lambda x: norm.logpdf(x, mu, std)
return LocalGpm(simulate, logpdf)
def _create_local_model_categorical(self, q, locality):
assert q in self.levels
assert all(0 <= l < self.levels[q] for l in locality)
counts = np.bincount(locality.astype(int), minlength=self.levels[q])
p = counts / np.sum(counts, dtype=float)
simulate = lambda N: self.rng.choice(self.levels[q], p=p, size=N)
logpdf = lambda x: np.log(p[x])
return LocalGpm(simulate, logpdf)
def logpdf_joint(self, x, y):
return gu.logsumexp([np.log(.25)
+ norm.logpdf(x, loc=mx, scale=self.noise)
+ norm.logpdf(y, loc=my, scale=self.noise)
for (mx,my) in zip(self.mx, self.my)])
def logpdf_joint(self, x, y):
return multivariate_normal.logpdf(
np.array([x,y]), np.array([0,0]),
np.array([[1,1-self.noise],[1-self.noise,1]]))
def logpdf_marginal(self, z):
return norm.logpdf(z, scale=1)
def logpdf_conditional(self, w, z):
mean = self.conditional_mean(z)
var = self.conditional_variance(z)
return norm.logpdf(w, loc=mean, scale=np.sqrt(var))
def logpdf(self, rowid, targets, constraints=None, inputs=None):
DistributionGpm.logpdf(self, rowid, targets, constraints, inputs)
x = targets[self.outputs[0]]
if not (self.l <= x <= self.h):
return -float('inf')
logpdf_unorm = NormalTrunc.calc_predictive_logp(
x, self.mu, self.sigma, self.l, self.h)
logcdf_norm = NormalTrunc.calc_log_normalizer(
self.mu, self.sigma, self.l, self.h)
return logpdf_unorm - logcdf_norm
def calc_predictive_logp(x, mu, sigma, l, h):
return norm.logpdf(x, loc=mu, scale=sigma)
