Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Compute the Jacobian matrix in Python

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?

like image 218
filtertips Avatar asked Mar 29 '18 09:03

filtertips


People also ask

What is Jacobian matrix in machine learning?

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.


3 Answers

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).

like image 115
Adam Erickson Avatar answered Sep 19 '22 23:09

Adam Erickson


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:

enter image description here

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
like image 37
dROOOze Avatar answered Sep 18 '22 23:09

dROOOze


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]])
like image 38
OliverQ Avatar answered Sep 19 '22 23:09

OliverQ