Вопрос: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.
Я читаю о jacobian Matrix, пытаясь построить один и из того, что я прочитал до сих пор, этот код python следует рассматривать как jacobian. Я понимаю это право?
Ответ №1
Вы можете использовать autograd библиотеку Гарварда (ссылка), где grad и jacobian принимают функцию в качестве аргумента:
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]))
В противном случае вы можете использовать метод jacobian доступный для матриц в 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)
Также вам может быть интересно посмотреть этот низкоуровневый вариант (ссылка).
Ответ №2
Якобиан определяется только для вектор- функций. Вы не можете работать с массивами, заполненными константами, для вычисления якобиана; вы должны знать основную функцию и ее частные производные или численное их приближение. Это очевидно, если учесть, что (частичная) производная от константы (по отношению к чему-то) равна 0.
В Python вы можете работать с символическими математическими модулями, такими как SymPy или SymEngine для вычисления якобиан функций. Здесь простая демонстрация примера из Википедии:
Использование модуля SymEngine:
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 Ответ №3
В python 3 вы можете попробовать пакет sympy:
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’])
который выдает:
Matrix([[2, 3],[2, -3]]) Ответ №4
Вот Python-реализация математического якобиана вектор-функции f(x), которая, как предполагается, возвращает одномерный массив с нулевыми значениями.
import numpy as np def J(f, x, dx=10^-8): n = len(x) func = f(x) jac = np.zeros((n, n)) for j in range(n): #through columns to allow for vector addition Dxj = (abs(x[j])*dx if x[j] != 0 else dx) x_plus = [(xi if k != j else xi+Dxj) for k, xi in enumerate(x)] jac[:, j] = (f(x_plus)-func)/Dxj return jac
Рекомендуется сделать dx ~ 10 ^ -8.
