implemented advanced slicing. grad implemented in only 1 level.

theano/tensor/basic.py

@@ -1160,11 +1160,8 @@ class _tensor_py_operators:
                 break
 
         if advanced:
-            if config.experimental.advanced_indexing:
-                if len(args) == 1:
-                    return AdvancedSubtensor1()(self, *args)
-                else:
-                    return AdvancedSubtensor(args)(self, *args)
+            if len(args) == 1:
+                return advanced_subtensor1(self, *args)
             else:
                 return AdvancedSubtensor(args)(self, *args)
         else:
@@ -3973,18 +3970,13 @@ def inverse_permutation(perm):
 # Should reproduce numpy's behaviour:
 # http://docs.scipy.org/doc/numpy/reference/arrays.indexing.html#advanced-indexing
 
-AddConfigVar('experimental.advanced_indexing',
-        "enable not-well-tested advanced indexing functionality",
-        BoolParam(False))
-
-
 class AdvancedSubtensor1(Op):
     """Implement x[ilist] where ilist is a vector of integers."""
 
     def __hash__(self):
         return hash(type(self))
     def __eq__(self, other):
-        type(self) == type(other)
+        return type(self) == type(other)
 
     def make_node(self, x, ilist):
         x_ = as_tensor_variable(x)
@@ -4004,18 +3996,60 @@ class AdvancedSubtensor1(Op):
     def perform(self, node, inp, out_):
         x, i = inp
         out, = out_
+        # Copy always implied by numpy advanced indexing semantic.
         out[0] = x[i]
 
     def grad(self, inputs, grads):
         gz, = grads
-        class NotImplementedOp(Op):
-            # This op should be pruned from the graph.
-            # This Op can be created in a graph,
-            # but it will cause problems if one of your parameters actually depends on it!
-            def make_node(self, *args):
-                return Apply(self, args, [inputs[0].type()])
-        return [NotImplementedOp()(gz)]+[None]*(len(inputs)-1)
+        assert len(inputs)==2
+        return [advanced_inc_subtensor(zeros_like(inputs[0]),gz,inputs[1])]+[None]*(len(inputs)-1)
+
+    def infer_shape(self, node, ishapes):
+        x, ilist = ishapes
+        return [ilist+x[1:]]
+
+advanced_subtensor1 = AdvancedSubtensor1()
+
+class AdvancedIncSubtensor1(Op):
+    """Increments a subtensor using advanced slicing (list of index)"""
+
+    def __hash__(self):
+        return hash(type(self))
+    def __eq__(self, other):
+        return type(self) == type(other)
+
+    def make_node(self, x, y, ilist):
+        x_ = as_tensor_variable(x)
+        y_ = as_tensor_variable(y)
+        ilist_ = as_tensor_variable(ilist)
+
+        assert x_.type.dtype == y_.type.dtype
+        assert x_.type.ndim == y_.type.ndim
+
+        if ilist_.type.dtype[:3] not in ('int', 'uin'):
+            raise TypeError('index must be integers')
+        if ilist_.type.broadcastable != (False,):
+            raise TypeError('index must be vector')
+        if x_.type.ndim == 0:
+            raise TypeError('cannot index into a scalar')
+        if x_.type.broadcastable[0]:
+            # the caller should have made a copy of x len(ilist) times
+            raise TypeError('cannot index into a broadcastable dimension')
+
+        return Apply(self, [x_, y_, ilist_], [x_.type()])
+
+    def perform(self, node, inp, out_):
+        # TODO opt to make this inplace
+        x, y, idx = inp
+        out, = out_
+        x = x.copy()
+        # x[idx] += y don't work if the same index is present many times.
+        # It do it only once
+        for (j,i) in enumerate(idx):
+            x[i] += y[j]
+        out[0] = x
 
+advanced_inc_subtensor = AdvancedIncSubtensor1()
 
 class AdvancedSubtensor(Op):
     """Return a subtensor copy, using advanced indexing.

theano/tensor/tests/test_basic.py

@@ -1576,6 +1576,92 @@ class T_subtensor(unittest.TestCase):
         good[1,0] = numpy.exp(data[1,0])
         self.failUnless(numpy.allclose(gval, good), (gval, good))
 
+    def test_ok_list(self):
+        for data, idx in [(numpy.random.rand(4), [1,0]),
+                          (numpy.random.rand(4,5), [2,3]),
+                          (numpy.random.rand(4,2,3), [0,3]),
+                          (numpy.random.rand(4,2,3), [3,3]),
+                          ]:
+            n = shared(data)
+            t = n[idx]
+            f = function([], t, mode=None)
+            topo = f.maker.env.toposort()
+            assert len(topo) == 1
+            assert isinstance(topo[0].op, theano.tensor.basic.AdvancedSubtensor1)
+            val = f()
+            good = data[idx]
+            self.failUnless(numpy.allclose(val, good), (val, good))
+
+    def test_err_invalid_list(self):
+        n = shared(numpy.asarray(5))
+        self.assertRaises(TypeError, n.__getitem__, [0,0])
+
+    def test_err_invalid_2list(self):
+        # TODO the error message is not clear
+        n = shared(numpy.ones((3,3))*5)
+        self.assertRaises(TypeError, n.__getitem__, ([0,0],[1,1]))
+
+    def test_err_bound_list(self):
+        n = shared(numpy.ones((2,3))*5)
+        t = n[0,4]
+        self.failUnless(isinstance(t.owner.op, Subtensor))
+        self.assertRaises(IndexError, eval_outputs, [t])
+
+    def grad_list_(self, idxs, data):
+        n = shared(data)
+
+        for idx in idxs:
+            idx_ = shared(numpy.asarray(idx))
+            t = n[idx_]
+            gn = grad(sum(exp(t)), n)
+            f = function([], gn, mode=None)
+            topo = f.maker.env.toposort()
+            assert any([isinstance(node.op, AdvancedIncSubtensor1) for node in topo])
+            assert any([isinstance(node.op, AdvancedSubtensor1) for node in topo])
+            gval = f()
+            good = numpy.zeros_like(data)
+            # good[idx] += numpy.exp(data[idx]) don't work when the same index is used many time
+            for i in idx:
+                good[i] += numpy.exp(data[i])
+            self.failUnless(numpy.allclose(gval, good), (gval, good))
+            def fct(t):
+                return sum(exp(t[idx_]))
+            utt.verify_grad(fct, [data])
+
+    def test_grad_list(self):
+        data = numpy.random.rand(3)
+        idxs = [[i] for i in range(data.shape[0])]
+        debug_mode = isinstance(theano.compile.mode.get_default_mode(),
+                                theano.compile.DebugMode)
+        for i in range(data.shape[0]):
+            for j in range(data.shape[0]):
+                idxs.append([i,j])
+        for i in range(data.shape[0]):
+            for j in range(data.shape[0]):
+                for k in range(data.shape[0]):
+                    idxs.append([i,j,k])
+        self.grad_list_(idxs, data)
+
+        data = numpy.random.rand(4,3)
+        self.grad_list_(idxs, data)
+
+        data = numpy.random.rand(5,3,7)
+
+    def test_shape_list(self):
+        #TODO for all type of subtensor shape
+        for data, idx in [(numpy.random.rand(4), [1,0]),
+                          (numpy.random.rand(4,2), [2,3]),
+                          (numpy.random.rand(4,2,3), [0,3]),
+                          (numpy.random.rand(4,2,3), [3,3,1,2,2,]),
+                          ]:
+            n = shared(data)
+            t = n[idx]
+            f = function([], t.shape, mode=None)
+            topo = f.maker.env.toposort()
+            #assert len(topo) == 1
+            #assert isinstance(topo[0].op, theano.tensor.basic.AdvancedSubtensor1)
+            val = f()
+            self.failUnless(numpy.allclose(val, data[idx].shape))
 
 class T_Join_and_Split(unittest.TestCase):
     """