Refactor distribution.draw_values (#2902)

* Fix for #2900. Changed the way in which draw_values handles the named node-inputs. Now the tree dependence is constructed to set the givens dict.

* Fixed conflicts

* Fixed more conflicts

* Fixed typo

* Changed test_dep_vars to test for successful draws even in the cases of dependent variables.

* Removed comments from test_random.py, and distribution.py. Added content to RELEASE-NOTES. Fixed bug in the interaction between draw_values, _draw_value and _compile_theano_function. In some cases, draw_values would set an item of the givens dictionary to a theano.tensor.TensorConstant. In _draw_value(param, ...), if param was a theano.tensor.TensorVariable without a random method, and not set in point, _compile_theano_function would be called, using as one of its variables, a theano.tensor.TensorConstant. This lead to TypeError: ("Constants not allowed in param list",...) exceptions being raised. The fix was to skip the inclusion into the givens dictionary of named nodes that were instances of theano.tensor.TensorConstant, because their value would already be available for theano during the function compilation.

* Fixed another bug which was similar to the theano.tensor.TensorConstant, but it occurred on theano.tensor.sharedvar.SharedVariable instances. The error that was raised was similar, SharedVariables cannot be supplied as raw input to theano.function. The fix was the same as for TensorConstants, skip them when constructing the givens dictionary.

* Guarded against a potencial bug. In draw_values, when skipping for TensorConstant and SharedVariable types, these nodes could be added to the stack again later because their names would not be in givens.keys(). To counter that, a separate set, `stored`, with the names of nodes that are either stored in givens or whos values should be available to theano.function, is used to chose which nodes to add to the stack.

* Syntax change based on twiecki's comment.

* Extended RELEASE-NOTES.md to also mention the sharedvar.ShareVariable fix.

Author: lucianopaz

RELEASE-NOTES.md

@@ -20,12 +20,15 @@
 - Add `offset` kwarg to `.glm`.
 - Changed the `compare` function to accept a dictionary of model-trace pairs instead of two separate lists of models and traces.
 - add test and support for creating multivariate mixture and mixture of mixtures
+- `distribution.draw_values`, now is also able to draw values from conditionally dependent RVs, such as autotransformed RVs (Refer to PR #2902).
 
 ### Fixes
 
 - `VonMises` does not overflow for large values of kappa. i0 and i1 have been removed and we now use log_i0 to compute the logp.
 - The bandwidth for KDE plots is computed using a modified version of Scott's rule. The new version uses entropy instead of standard deviation. This works better for multimodal distributions. Functions using KDE plots has a new argument `bw` controlling the bandwidth.
 - fix PyMC3 variable is not replaced if provided in more_replacements (#2890)
+- Fix for issue #2900. For many situations, named node-inputs do not have a `random` method, while some intermediate node may have it. This meant that if the named node-input at the leaf of the graph did not have a fixed value, `theano` would try to compile it and fail to find inputs, raising a `theano.gof.fg.MissingInputError`. This was fixed by going through the theano variable's owner inputs graph, trying to get intermediate named-nodes values if the leafs had failed.
+- In `distribution.draw_values`, some named nodes could be `theano.tensor.TensorConstant`s or `theano.tensor.sharedvar.SharedVariable`s. Nevertheless, in `distribution._draw_value`, these would be passed to `distribution._compile_theano_function` as if they were `theano.tensor.TensorVariable`s. This could lead to the following exceptions `TypeError: ('Constants not allowed in param list', ...)` or `TypeError: Cannot use a shared variable (...)`. The fix was to not add `theano.tensor.TensorConstant` or `theano.tensor.sharedvar.SharedVariable` named nodes into the `givens` dict that could be used in `distribution._compile_theano_function`.
 
 ### Deprecations
 

pymc3/distributions/distribution.py

@@ -4,7 +4,7 @@
 from theano import function
 import theano
 from ..memoize import memoize
-from ..model import Model, get_named_nodes, FreeRV, ObservedRV
+from ..model import Model, get_named_nodes_and_relations, FreeRV, ObservedRV
 from ..vartypes import string_types
 
 __all__ = ['DensityDist', 'Distribution', 'Continuous', 'Discrete',
@@ -229,17 +229,72 @@ def draw_values(params, point=None):
 
     """
     # Distribution parameters may be nodes which have named node-inputs
-    # specified in the point. Need to find the node-inputs to replace them.
-    givens = {}
+    # specified in the point. Need to find the node-inputs, their
+    # parents and children to replace them.
+    leaf_nodes = {}
+    named_nodes_parents = {}
+    named_nodes_children = {}
     for param in params:
         if hasattr(param, 'name'):
-            named_nodes = get_named_nodes(param)
-            if param.name in named_nodes:
-                named_nodes.pop(param.name)
-            for name, node in named_nodes.items():
-                if not isinstance(node, (tt.sharedvar.SharedVariable,
-                                         tt.TensorConstant)):
-                    givens[name] = (node, _draw_value(node, point=point))
+            # Get the named nodes under the `param` node
+            nn, nnp, nnc = get_named_nodes_and_relations(param)
+            leaf_nodes.update(nn)
+            # Update the discovered parental relationships
+            for k in nnp.keys():
+                if k not in named_nodes_parents.keys():
+                    named_nodes_parents[k] = nnp[k]
+                else:
+                    named_nodes_parents[k].update(nnp[k])
+            # Update the discovered child relationships
+            for k in nnc.keys():
+                if k not in named_nodes_children.keys():
+                    named_nodes_children[k] = nnc[k]
+                else:
+                    named_nodes_children[k].update(nnc[k])
+    
+    # Init givens and the stack of nodes to try to `_draw_value` from
+    givens = {}
+    stored = set([])  # Some nodes 
+    stack = list(leaf_nodes.values())  # A queue would be more appropriate
+    while stack:
+        next_ = stack.pop(0)
+        if next_ in stored:
+            # If the node already has a givens value, skip it
+            continue
+        elif isinstance(next_, (tt.TensorConstant,
+                                tt.sharedvar.SharedVariable)):
+            # If the node is a theano.tensor.TensorConstant or a
+            # theano.tensor.sharedvar.SharedVariable, its value will be
+            # available automatically in _compile_theano_function so
+            # we can skip it. Furthermore, if this node was treated as a
+            # TensorVariable that should be compiled by theano in
+            # _compile_theano_function, it would raise a `TypeError:
+            # ('Constants not allowed in param list', ...)` for 
+            # TensorConstant, and a `TypeError: Cannot use a shared
+            # variable (...) as explicit input` for SharedVariable.
+            stored.add(next_.name)
+            continue
+        else:
+            # If the node does not have a givens value, try to draw it.
+            # The named node's children givens values must also be taken
+            # into account.
+            children = named_nodes_children[next_]
+            temp_givens = [givens[k] for k in givens.keys() if k in children]
+            try:
+                # This may fail for autotransformed RVs, which don't
+                # have the random method
+                givens[next_.name] = (next_, _draw_value(next_,
+                                                         point=point,
+                                                         givens=temp_givens))
+                stored.add(next_.name)
+            except theano.gof.fg.MissingInputError:
+                # The node failed, so we must add the node's parents to
+                # the stack of nodes to try to draw from. We exclude the
+                # nodes in the `params` list.
+                stack.extend([node for node in named_nodes_parents[next_]
+                              if node is not None and
+                              node.name not in stored and
+                              node not in params])
     values = []
     for param in params:
         values.append(_draw_value(param, point=point, givens=givens.values()))

pymc3/model.py

@@ -78,29 +78,71 @@ def incorporate_methods(source, destination, methods, default=None,
         else:
             setattr(destination, method, None)
 
-
-def get_named_nodes(graph):
-    """Get the named nodes in a theano graph
-    (i.e., nodes whose name attribute is not None).
+def get_named_nodes_and_relations(graph):
+    """Get the named nodes in a theano graph (i.e., nodes whose name
+    attribute is not None) along with their relationships (i.e., the
+    node's named parents, and named children, while skipping unnamed
+    intermediate nodes)
 
     Parameters
     ----------
     graph - a theano node
 
     Returns:
-    A dictionary of name:node pairs.
+    leaf_nodes: A dictionary of name:node pairs, of the named nodes that
+        are also leafs of the graph
+    node_parents: A dictionary of node:set([parents]) pairs. Each key is
+        a theano named node, and the corresponding value is the set of
+        theano named nodes that are parents of the node. These parental
+        relations skip unnamed intermediate nodes.
+    node_children: A dictionary of node:set([children]) pairs. Each key
+        is a theano named node, and the corresponding value is the set
+        of theano named nodes that are children of the node. These child
+        relations skip unnamed intermediate nodes.
+    
     """
-    return _get_named_nodes(graph, {})
-
-
-def _get_named_nodes(graph, nodes):
-    if graph.owner is None:
-        if graph.name is not None:
-            nodes.update({graph.name: graph})
+    if graph.name is not None:
+        node_parents = {graph: set()}
+        node_children = {graph: set()}
     else:
+        node_parents = {}
+        node_children = {}
+    return _get_named_nodes_and_relations(graph, None, {}, node_parents, node_children)
+
+def _get_named_nodes_and_relations(graph, parent, leaf_nodes,
+                                        node_parents, node_children):
+    if graph.owner is None:  # Leaf node
+        if graph.name is not None:  # Named leaf node
+            leaf_nodes.update({graph.name: graph})
+            if parent is not None:  # Is None for the root node
+                try:
+                    node_parents[graph].add(parent)
+                except KeyError:
+                    node_parents[graph] = set([parent])
+                node_children[parent].add(graph)
+            # Flag that the leaf node has no children
+            node_children[graph] = set()
+    else:  # Intermediate node
+        if graph.name is not None:  # Intermediate named node
+            if parent is not None:  # Is only None for the root node
+                try:
+                    node_parents[graph].add(parent)
+                except KeyError:
+                    node_parents[graph] = set([parent])
+                node_children[parent].add(graph)
+            # The current node will be set as the parent of the next
+            # nodes only if it is a named node
+            parent = graph
+            # Init the nodes children to an empty set
+            node_children[graph] = set()
         for i in graph.owner.inputs:
-            nodes.update(_get_named_nodes(i, nodes))
-    return nodes
+            temp_nodes, temp_inter, temp_tree = \
+                _get_named_nodes_and_relations(i, parent, leaf_nodes,
+                                               node_parents, node_children)
+            leaf_nodes.update(temp_nodes)
+            node_parents.update(temp_inter)
+            node_children.update(temp_tree)
+    return leaf_nodes, node_parents, node_children
 
 
 class Context(object):

pymc3/tests/test_random.py

@@ -80,13 +80,11 @@ def test_dep_vars(self):
         point = {'a': np.array([1., 2.])}
         npt.assert_equal(draw_values([a], point=point), [point['a']])
 
-        with pytest.raises(theano.gof.MissingInputError):
-            draw_values([a])
-
-        # We need the untransformed vars
-        with pytest.raises(theano.gof.MissingInputError):
-            draw_values([a], point={'sd': np.array([2., 3.])})
-
-        val1 = draw_values([a], point={'sd_log__': np.array([2., 3.])})[0]
-        val2 = draw_values([a], point={'sd_log__': np.array([2., 3.])})[0]
-        assert np.all(val1 != val2)
+        val1 = draw_values([a])[0]
+        val2 = draw_values([a], point={'sd': np.array([2., 3.])})[0]
+        val3 = draw_values([a], point={'sd_log__': np.array([2., 3.])})[0]
+        val4 = draw_values([a], point={'sd_log__': np.array([2., 3.])})[0]
+        
+        assert all([np.all(val1 != val2), np.all(val1 != val3),
+                    np.all(val1 != val4), np.all(val2 != val3),
+                    np.all(val2 != val4), np.all(val3 != val4)])