import numpy as np
a = np.array([[1,2,3],
[4,5,6],
[7,8,9]])
b = np.array([[1,2,3]]).T
c = a.dot(b) #function
jacobian = a # as partial derivative of c w.r.t to b is a.
I am reading about jacobian Matrix, trying to build one and from what I have read so far, this python code should be considered as jacobian. Am I understanding this right?
The Jacobian matrix collects all first-order partial derivatives of a multivariate function that can be used for backpropagation. The Jacobian determinant is useful in changing between variables, where it acts as a scaling factor between one coordinate space and another.
You can use the Harvard autograd
library (link), where grad
and jacobian
take a function as their argument:
import autograd.numpy as np
from autograd import grad, jacobian
x = np.array([5,3], dtype=float)
def cost(x):
return x[0]**2 / x[1] - np.log(x[1])
gradient_cost = grad(cost)
jacobian_cost = jacobian(cost)
gradient_cost(x)
jacobian_cost(np.array([x,x,x]))
Otherwise, you could use the jacobian
method available for matrices in sympy
:
from sympy import sin, cos, Matrix
from sympy.abc import rho, phi
X = Matrix([rho*cos(phi), rho*sin(phi), rho**2])
Y = Matrix([rho, phi])
X.jacobian(Y)
Also, you may also be interested to see this low-level variant (link). MATLAB provides nice documentation on its jacobian
function here.
UPDATE: Note that the autograd
library has since been rolled into jax
, which provides functions for computing forward and inverse Jacobian matrices (link).
The Jacobian is only defined for vector-valued functions. You cannot work with arrays filled with constants to calculate the Jacobian; you must know the underlying function and its partial derivatives, or the numerical approximation of these. This is obvious when you consider that the (partial) derivative of a constant (with respect to something) is 0.
In Python, you can work with symbolic math modules such as SymPy
or SymEngine
to calculate Jacobians of functions. Here's a simple demonstration of an example from Wikipedia:
Using the SymEngine
module:
Python 2.7.11 (v2.7.11:6d1b6a68f775, Dec 5 2015, 20:40:30) [MSC v.1500 64 bit (AMD64)] on win32
Type "help", "copyright", "credits" or "license" for more information.
>>>
>>> import symengine
>>>
>>>
>>> vars = symengine.symbols('x y') # Define x and y variables
>>> f = symengine.sympify(['y*x**2', '5*x + sin(y)']) # Define function
>>> J = symengine.zeros(len(f),len(vars)) # Initialise Jacobian matrix
>>>
>>> # Fill Jacobian matrix with entries
... for i, fi in enumerate(f):
... for j, s in enumerate(vars):
... J[i,j] = symengine.diff(fi, s)
...
>>> print J
[2*x*y, x**2]
[5, cos(y)]
>>>
>>> print symengine.Matrix.det(J)
2*x*y*cos(y) - 5*x**2
In python 3, you can try sympy package:
import sympy as sym
def Jacobian(v_str, f_list):
vars = sym.symbols(v_str)
f = sym.sympify(f_list)
J = sym.zeros(len(f),len(vars))
for i, fi in enumerate(f):
for j, s in enumerate(vars):
J[i,j] = sym.diff(fi, s)
return J
Jacobian('u1 u2', ['2*u1 + 3*u2','2*u1 - 3*u2'])
which gives out:
Matrix([[2, 3],[2, -3]])
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With