From da67ffa00cae00afd5bdce0f90892c56de35e4d3 Mon Sep 17 00:00:00 2001 From: Chris Fonnesbeck Date: Wed, 18 Jan 2017 10:20:56 -0600 Subject: [PATCH 1/3] Changed case of init arguments in sample to upper case, but kept backward-compatibility with lower case --- pymc3/examples/arbitrary_stochastic.py | 2 +- pymc3/examples/factor_potential.py | 2 +- pymc3/sampling.py | 31 +++++++++++++++----------- pymc3/tests/test_examples.py | 6 ++--- pymc3/tests/test_glm.py | 2 +- pymc3/tests/test_models_linear.py | 8 +++---- pymc3/tests/test_plots.py | 4 ++-- 7 files changed, 30 insertions(+), 25 deletions(-) diff --git a/pymc3/examples/arbitrary_stochastic.py b/pymc3/examples/arbitrary_stochastic.py index 9739c56dac..a342a453e5 100644 --- a/pymc3/examples/arbitrary_stochastic.py +++ b/pymc3/examples/arbitrary_stochastic.py @@ -20,7 +20,7 @@ def run(n_samples=3000): start = model.test_point h = pm.find_hessian(start, model=model) step = pm.Metropolis(model.vars, h, blocked=True, model=model) - trace = pm.sample(n_samples, step, start, model=model) + trace = pm.sample(n_samples, step=step, start=start, model=model) return trace if __name__ == "__main__": diff --git a/pymc3/examples/factor_potential.py b/pymc3/examples/factor_potential.py index 17b1d6c743..a189d729e2 100644 --- a/pymc3/examples/factor_potential.py +++ b/pymc3/examples/factor_potential.py @@ -13,7 +13,7 @@ def run(n=3000): if n == "short": n = 50 with model: - trace = sample(n, step, start) + trace = sample(n, step=step, start=start) if __name__ == '__main__': run() diff --git a/pymc3/sampling.py b/pymc3/sampling.py index 69931ca21c..95f5660e89 100644 --- a/pymc3/sampling.py +++ b/pymc3/sampling.py @@ -98,20 +98,20 @@ def sample(draws, step=None, init='advi', n_init=200000, start=None, A step function or collection of functions. If no step methods are specified, or are partially specified, they will be assigned automatically (defaults to None). - init : str {'advi', 'advi_map', 'map', 'nuts', None} + init : str {'ADVI', 'ADVI_MAP', 'MAP', 'NUTS', None} Initialization method to use. - * advi : Run ADVI to estimate starting points and diagonal covariance + * ADVI : Run ADVI to estimate starting points and diagonal covariance matrix. If njobs > 1 it will sample starting points from the estimated posterior, otherwise it will use the estimated posterior mean. - * advi_map: Initialize ADVI with MAP and use MAP as starting point. - * map : Use the MAP as starting point. - * nuts : Run NUTS to estimate starting points and covariance matrix. If + * ADVI_MAP: Initialize ADVI with MAP and use MAP as starting point. + * MAP : Use the MAP as starting point. + * NUTS : Run NUTS to estimate starting points and covariance matrix. If njobs > 1 it will sample starting points from the estimated posterior, otherwise it will use the estimated posterior mean. * None : Do not initialize. n_init : int Number of iterations of initializer - If 'advi', number of iterations, if 'nuts', number of draws. + If 'ADVI', number of iterations, if 'nuts', number of draws. start : dict Starting point in parameter space (or partial point) Defaults to trace.point(-1)) if there is a trace provided and @@ -147,6 +147,8 @@ def sample(draws, step=None, init='advi', n_init=200000, start=None, """ model = modelcontext(model) + if init is not None: + init = init.lower() if step is None and init is not None and pm.model.all_continuous(model.vars): # By default, use NUTS sampler @@ -398,7 +400,7 @@ def sample_ppc(trace, samples=None, model=None, vars=None, size=None, random_see return {k: np.asarray(v) for k, v in ppc.items()} -def init_nuts(init='advi', njobs=1, n_init=500000, model=None, +def init_nuts(init='ADVI', njobs=1, n_init=500000, model=None, random_seed=-1, **kwargs): """Initialize and sample from posterior of a continuous model. @@ -409,17 +411,17 @@ def init_nuts(init='advi', njobs=1, n_init=500000, model=None, Parameters ---------- - init : str {'advi', 'advi_map', 'map', 'nuts'} + init : str {'ADVI', 'ADVI_MAP', 'MAP', 'NUTS'} Initialization method to use. - * advi : Run ADVI to estimate posterior mean and diagonal covariance matrix. - * advi_map: Initialize ADVI with MAP and use MAP as starting point. - * map : Use the MAP as starting point. - * nuts : Run NUTS and estimate posterior mean and covariance matrix. + * ADVI : Run ADVI to estimate posterior mean and diagonal covariance matrix. + * ADVI_MAP: Initialize ADVI with MAP and use MAP as starting point. + * MAP : Use the MAP as starting point. + * NUTS : Run NUTS and estimate posterior mean and covariance matrix. njobs : int Number of parallel jobs to start. n_init : int Number of iterations of initializer - If 'advi', number of iterations, if 'metropolis', number of draws. + If 'ADVI', number of iterations, if 'metropolis', number of draws. model : Model (optional if in `with` context) **kwargs : keyword arguments Extra keyword arguments are forwarded to pymc3.NUTS. @@ -439,6 +441,9 @@ def init_nuts(init='advi', njobs=1, n_init=500000, model=None, pm._log.info('Initializing NUTS using {}...'.format(init)) random_seed = int(np.atleast_1d(random_seed)[0]) + + if init is not None: + init = init.lower() if init == 'advi': v_params = pm.variational.advi(n=n_init, random_seed=random_seed) diff --git a/pymc3/tests/test_examples.py b/pymc3/tests/test_examples.py index 2ce8f1a793..c2da6c1c2b 100644 --- a/pymc3/tests/test_examples.py +++ b/pymc3/tests/test_examples.py @@ -54,7 +54,7 @@ def test_run(self): h = H(start) step = pm.HamiltonianMC(model.vars, h) - pm.sample(50, step, start) + pm.sample(50, step=step, start=start) class TestARM12_6(SeededTest): @@ -88,7 +88,7 @@ def too_slow(self): start = pm.find_MAP(start=start, vars=[model['groupmean'], model['sd_interval_'], model['floor_m']]) step = pm.NUTS(model.vars, scaling=start) - pm.sample(50, step, start) + pm.sample(50, step=step, start=start) class TestARM12_6Uranium(SeededTest): @@ -129,7 +129,7 @@ def too_slow(self): h = np.diag(H(start)) step = pm.HamiltonianMC(model.vars, h) - pm.sample(50, step, start) + pm.sample(50, step=step, start=start) def build_disaster_model(masked=False): diff --git a/pymc3/tests/test_glm.py b/pymc3/tests/test_glm.py index 259ffdb48b..b4cc107ed7 100644 --- a/pymc3/tests/test_glm.py +++ b/pymc3/tests/test_glm.py @@ -34,7 +34,7 @@ def test_linear_component(self): Normal('y_obs', mu=y_est, sd=sigma, observed=self.y_linear) start = find_MAP(vars=[sigma]) step = Slice(model.vars) - trace = sample(500, step, start, progressbar=False, random_seed=self.random_seed) + trace = sample(500, step=step, start=start, progressbar=False, random_seed=self.random_seed) self.assertAlmostEqual(np.mean(trace['Intercept']), self.intercept, 1) self.assertAlmostEqual(np.mean(trace['x']), self.slope, 1) diff --git a/pymc3/tests/test_models_linear.py b/pymc3/tests/test_models_linear.py index fac66fc0de..b0a00953e7 100644 --- a/pymc3/tests/test_models_linear.py +++ b/pymc3/tests/test_models_linear.py @@ -44,7 +44,7 @@ def test_linear_component(self): Normal('y_obs', mu=lm.y_est, sd=sigma, observed=self.y_linear) # yields y_obs start = find_MAP(vars=[sigma]) step = Slice(model.vars) - trace = sample(500, step, start, progressbar=False, random_seed=self.random_seed) + trace = sample(500, step=step, start=start, progressbar=False, random_seed=self.random_seed) self.assertAlmostEqual(np.mean(trace['lm_Intercept']), self.intercept, 1) self.assertAlmostEqual(np.mean(trace['lm_x0']), self.slope, 1) @@ -58,7 +58,7 @@ def test_linear_component_from_formula(self): Normal('y_obs', mu=lm.y_est, sd=sigma, observed=self.y_linear) start = find_MAP(vars=[sigma]) step = Slice(model.vars) - trace = sample(500, step, start, progressbar=False, random_seed=self.random_seed) + trace = sample(500, step=step, start=start, progressbar=False, random_seed=self.random_seed) self.assertAlmostEqual(np.mean(trace['Intercept']), self.intercept, 1) self.assertAlmostEqual(np.mean(trace['x']), self.slope, 1) @@ -79,7 +79,7 @@ def test_glm(self): ) start = find_MAP() step = Slice(model.vars) - trace = sample(500, step, start, progressbar=False, random_seed=self.random_seed) + trace = sample(500, step=step, start=start, progressbar=False, random_seed=self.random_seed) self.assertAlmostEqual(np.mean(trace['glm_Intercept']), self.intercept, 1) self.assertAlmostEqual(np.mean(trace['glm_x0']), self.slope, 1) self.assertAlmostEqual(np.mean(trace['glm_sd']), self.sd, 1) @@ -91,7 +91,7 @@ def test_glm_from_formula(self): Glm.from_formula('y ~ x', self.data_linear, name=NAME) start = find_MAP() step = Slice(model.vars) - trace = sample(500, step, start, progressbar=False, random_seed=self.random_seed) + trace = sample(500, step=step, start=start, progressbar=False, random_seed=self.random_seed) self.assertAlmostEqual(np.mean(trace['%s_Intercept' % NAME]), self.intercept, 1) self.assertAlmostEqual(np.mean(trace['%s_x' % NAME]), self.slope, 1) diff --git a/pymc3/tests/test_plots.py b/pymc3/tests/test_plots.py index 05dc8d5a46..dc2ea22e44 100644 --- a/pymc3/tests/test_plots.py +++ b/pymc3/tests/test_plots.py @@ -20,7 +20,7 @@ def test_plots(): start = model.test_point h = find_hessian(start) step = Metropolis(model.vars, h) - trace = sample(3000, step, start) + trace = sample(3000, step=step, start=start) traceplot(trace) forestplot(trace) @@ -34,7 +34,7 @@ def test_plots_multidimensional(): with model as model: h = np.diag(find_hessian(start)) step = Metropolis(model.vars, h) - trace = sample(3000, step, start) + trace = sample(3000, step=step, start=start) traceplot(trace) plot_posterior(trace) From 016e80784d4a4e51291f0f6e1b160224e3e169c1 Mon Sep 17 00:00:00 2001 From: Chris Fonnesbeck Date: Wed, 18 Jan 2017 12:29:00 -0600 Subject: [PATCH 2/3] Changed positional arguments to keyword arguments in tests --- pymc3/tests/test_diagnostics.py | 14 +++++++------- pymc3/tests/test_plots.py | 2 +- pymc3/tests/test_sampling.py | 2 +- 3 files changed, 9 insertions(+), 9 deletions(-) diff --git a/pymc3/tests/test_diagnostics.py b/pymc3/tests/test_diagnostics.py index d5e6987385..770ca2ce7d 100644 --- a/pymc3/tests/test_diagnostics.py +++ b/pymc3/tests/test_diagnostics.py @@ -20,9 +20,9 @@ def get_ptrace(self, n_samples): # Run sampler step1 = Slice([model.early_mean_log_, model.late_mean_log_]) step2 = Metropolis([model.switchpoint]) - start = {'early_mean': 2., 'late_mean': 3., 'switchpoint': 90} - ptrace = sample(n_samples, [step1, step2], start, njobs=2, progressbar=False, - random_seed=[1, 4]) + start = {'early_mean': 7., 'late_mean': 1., 'switchpoint': 90} + ptrace = sample(n_samples, step=[step1, step2], start=start, njobs=2, + progressbar=False, random_seed=[1, 4]) return ptrace def test_good(self): @@ -53,7 +53,7 @@ def test_right_shape_python_float(self, shape=None, test_shape=None): # start sampling at the MAP start = find_MAP() step = NUTS(scaling=start) - ptrace = sample(n_samples, step, start, + ptrace = sample(n_samples, step=step, start=start, njobs=n_jobs, random_seed=42) rhat = gelman_rubin(ptrace)['x'] @@ -92,7 +92,7 @@ def get_switchpoint(self, n_samples): # Run sampler step1 = Slice([model.early_mean_log_, model.late_mean_log_]) step2 = Metropolis([model.switchpoint]) - trace = sample(n_samples, [step1, step2], progressbar=False, random_seed=1) + trace = sample(n_samples, step=[step1, step2], progressbar=False, random_seed=1) return trace['switchpoint'] def test_geweke_negative(self): @@ -153,7 +153,7 @@ def test_effective_n(self): # start sampling at the MAP start = find_MAP() step = NUTS(scaling=start) - ptrace = sample(n_samples, step, start, + ptrace = sample(n_samples, step=step, start=start, njobs=n_jobs, random_seed=42) n_effective = effective_n(ptrace)['x'] @@ -174,7 +174,7 @@ def test_effective_n_right_shape_python_float(self, # start sampling at the MAP start = find_MAP() step = NUTS(scaling=start) - ptrace = sample(n_samples, step, start, + ptrace = sample(n_samples, step=step, start=start, njobs=n_jobs, random_seed=42) n_effective = effective_n(ptrace)['x'] diff --git a/pymc3/tests/test_plots.py b/pymc3/tests/test_plots.py index dc2ea22e44..a0d23464da 100644 --- a/pymc3/tests/test_plots.py +++ b/pymc3/tests/test_plots.py @@ -47,7 +47,7 @@ def test_multichain_plots(): step1 = Slice([model.early_mean_log_, model.late_mean_log_]) step2 = Metropolis([model.switchpoint]) start = {'early_mean': 2., 'late_mean': 3., 'switchpoint': 50} - ptrace = sample(1000, [step1, step2], start, njobs=2) + ptrace = sample(1000, step=[step1, step2], start=start, njobs=2) forestplot(ptrace, varnames=['early_mean', 'late_mean']) autocorrplot(ptrace, varnames=['switchpoint']) diff --git a/pymc3/tests/test_sampling.py b/pymc3/tests/test_sampling.py index 5c81be566c..7d85d2284d 100644 --- a/pymc3/tests/test_sampling.py +++ b/pymc3/tests/test_sampling.py @@ -59,7 +59,7 @@ def test_sample(self): with self.model: for njobs in test_njobs: for steps in [1, 10, 300]: - pm.sample(steps, self.step, {}, None, njobs=njobs, random_seed=self.random_seed) + pm.sample(steps, step=self.step, njobs=njobs, random_seed=self.random_seed) def test_sample_init(self): with self.model: From faab2e1f6b2679909d15bec76c2a6cd803657ac9 Mon Sep 17 00:00:00 2001 From: Chris Fonnesbeck Date: Wed, 18 Jan 2017 16:29:35 -0600 Subject: [PATCH 3/3] Fix failures in test_examples and test_diagnostics --- pymc3/tests/test_diagnostics.py | 4 ++-- pymc3/tests/test_examples.py | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/pymc3/tests/test_diagnostics.py b/pymc3/tests/test_diagnostics.py index 770ca2ce7d..d85738012d 100644 --- a/pymc3/tests/test_diagnostics.py +++ b/pymc3/tests/test_diagnostics.py @@ -20,9 +20,9 @@ def get_ptrace(self, n_samples): # Run sampler step1 = Slice([model.early_mean_log_, model.late_mean_log_]) step2 = Metropolis([model.switchpoint]) - start = {'early_mean': 7., 'late_mean': 1., 'switchpoint': 90} + start = {'early_mean': 7., 'late_mean': 1., 'switchpoint': 100} ptrace = sample(n_samples, step=[step1, step2], start=start, njobs=2, - progressbar=False, random_seed=[1, 4]) + progressbar=False, random_seed=[20090425, 19700903]) return ptrace def test_good(self): diff --git a/pymc3/tests/test_examples.py b/pymc3/tests/test_examples.py index c2da6c1c2b..d071826cce 100644 --- a/pymc3/tests/test_examples.py +++ b/pymc3/tests/test_examples.py @@ -263,7 +263,7 @@ def test_run(self): start = {'psi': 0.5, 'z': (self.y > 0).astype(int), 'theta': 5} step_one = pm.Metropolis([model.theta_interval_, model.psi_logodds_]) step_two = pm.BinaryMetropolis([model.z]) - pm.sample(50, [step_one, step_two], start) + pm.sample(50, step=[step_one, step_two], start=start) class TestRSV(SeededTest):