Merge pull request #2046 from pymc-devs/transformed_init_fix

Ensure auto-transformed variables use specified starting values

Author: fonnesbeck

pymc3/distributions/transforms.py

@@ -6,8 +6,8 @@
 from ..math import logit, invlogit
 import numpy as np
 
-__all__ = ['transform', 'stick_breaking', 'logodds',
-           'log', 'sum_to_1', 't_stick_breaking']
+__all__ = ['transform', 'stick_breaking', 'logodds', 'interval',
+          'lowerbound', 'upperbound', 'log', 'sum_to_1', 't_stick_breaking']
 
 
 class Transform(object):

pymc3/model.py

@@ -989,11 +989,14 @@ def __init__(self, type=None, owner=None, index=None, name=None,
         if type is None:
             type = distribution.type
         super(TransformedRV, self).__init__(type, owner, index, name)
+        
+        self.transformation = transform
 
         if distribution is not None:
             self.model = model
 
             transformed_name = "{}_{}_".format(name, transform.name)
+
             self.transformed = model.Var(
                 transformed_name, transform.apply(distribution), total_size=total_size)
 

pymc3/sampling.py

@@ -348,9 +348,9 @@ def _iter_sample(draws, step, start=None, trace=None, chain=0, tune=None,
     strace = _choose_backend(trace, chain, model=model)
 
     if len(strace) > 0:
-        _soft_update(start, strace.point(-1))
+        _update_start_vals(start, strace.point(-1), model)
     else:
-        _soft_update(start, model.test_point)
+        _update_start_vals(start, model.test_point, model)
 
     try:
         step = CompoundStep(step)
@@ -457,13 +457,19 @@ def stop_tuning(step):
 
     return step
 
-
-def _soft_update(a, b):
-    """As opposed to dict.update, don't overwrite keys if present.
+def _update_start_vals(a, b, model):
+    """Update a with b, without overwriting existing keys. Values specified for
+    transformed variables on the original scale are also transformed and inserted.
     """
+    
+    for name in a:
+        for tname in b:
+            if tname.startswith(name) and tname!=name:
+                transform_func = [d.transformation for d in model.deterministics if d.name==name][0]
+                b[tname] = transform_func.forward(a[name]).eval()
+    
     a.update({k: v for k, v in b.items() if k not in a})
 
-
 def sample_ppc(trace, samples=None, model=None, vars=None, size=None,
                random_seed=None, progressbar=True):
     """Generate posterior predictive samples from a model given a trace.

pymc3/tests/test_diagnostics.py

@@ -21,7 +21,7 @@ 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': 100}
+            start = {'early_mean': 7., 'late_mean': 5., 'switchpoint': 10}
             ptrace = sample(n_samples, step=[step1, step2], start=start, njobs=2, 
                     progressbar=False, random_seed=[20090425, 19700903])
         return ptrace
@@ -35,15 +35,15 @@ def test_good(self):
 
     def test_bad(self):
         """Confirm Gelman-Rubin statistic is far from 1 for a small number of samples."""
-        n_samples = 10
+        n_samples = 5
         rhat = gelman_rubin(self.get_ptrace(n_samples))
         assert not all(1 / self.good_ratio < r <
                              self.good_ratio for r in rhat.values())
 
     def test_right_shape_python_float(self, shape=None, test_shape=None):
         """Check Gelman-Rubin statistic shape is correct w/ python float"""
         n_jobs = 3
-        n_samples = 10
+        n_samples = 5
 
         with Model():
             if shape is not None:

pymc3/tests/test_sampling.py

@@ -1,5 +1,7 @@
 from itertools import combinations
 import numpy as np
+from numpy.testing import assert_almost_equal
+
 try:
     import unittest.mock as mock  # py3
 except ImportError:
@@ -117,6 +119,12 @@ def test_soft_update_empty(self):
         test_point = {'a': 3, 'b': 4}
         pm.sampling._soft_update(start, test_point)
         assert start == test_point
+        
+    def test_soft_update_transformed(self):
+        start = {'a': 2}
+        test_point = {'a_log_': 0}
+        pm.sampling._soft_update(start, test_point)
+        assert assert_almost_equal(start['a_log_'], np.log(start['a']))
 
 
 class TestNamedSampling(SeededTest):