title: Vectorized Partial Derivative of Multivariate Function in Python
body: <p>There was a phenomenal answer posted by alko for computing a partial derivative of a multivariate function numerically in <a href="https://stackoverflow.com/questions/20708038/scipy-misc-derivative-for-mutiple-argument-function">this thread</a>.</p>
  
  <p>I have a follow-up question now about enhancing this function to accept an array of input values. I have some code where I'm looping through a big long list of n-dimensional points, calculating the partial derivatives with respect to each variable, and this becomes quite computationally expensive.</p>
  
  <p>It's easy enough to vectorize the function in question with <code>np.vectorize</code>, but it causes issues with the <code>partial_derivative</code> wrapper:</p>
  
  <pre><code>from scipy.misc import derivative
  import numpy as np
  
  def foo(x, y):
      return(x**2 + y**3)
  
  def partial_derivative(func, var=0, point=[]):
      args = point[:]
      def wraps(x):
          args[var] = x
          return func(*args)
      return derivative(wraps, point[var], dx=1e-6)
  
  vfoo = np.vectorize(foo)
  
  &gt;&gt;&gt;foo(3,1)
  &gt;&gt;&gt;10
  &gt;&gt;&gt;vfoo([3,3], [1,1])
  &gt;&gt;&gt;array([10,10])
  
  &gt;&gt;&gt;partial_derivative(foo,0,[3,1])
  &gt;&gt;&gt;6.0
  &gt;&gt;&gt;partial_derivative(vfoo,0,[[3,3], [1,1]])
  &gt;&gt;&gt;TypeError: can only concatenate list (not "float") to list
  </code></pre>
  
  <p>The last line should ideally return <code>[6.0, 6.0]</code>. In this case the two arrays supplied to the vectorized function <code>vfoo</code> are essentially zipped up pairwise, so <code>([3,3], [1,1])</code> gets transformed into two points, <code>[3,1]</code> and <code>[3,1]</code>. This seems to get mangled when it gets passed to the function <code>wraps</code>. The point that it ends up passing to the function <code>derivative</code> is <code>[3,3]</code>. In addition, there's obviously the <code>TypeError</code> thrown up. </p>
  
  <p>Does anyone have any recommendations or suggestions? Has anyone ever had a need to do something similar?</p>
  
  <h2>Edit</h2>
  
  <p>Sometimes I think posting on SO is just what it takes to break a mental block. I think I've got it working for anyone who might be interested:</p>
  
  <pre><code>vfoo = np.vectorize(foo)
  
  foo(3,1)
  X = np.array([3,3])
  Y = np.array([1,1])
  vfoo(X, Y)
  partial_derivative(foo,0,[3,1])
  partial_derivative(vfoo,0,[X, Y])
  </code></pre>
  
  <p>And the last line now returns <code>array([ 6.,  6.])</code></p>