用“疯狂”计算偏导数
Compute partial derivatives with `madness`
madness
包,如前所述 here,非常适合 R 中的 autodiff。
我现在想计算导数 wrt y 的导数 wrt x。
$\frac{\partial}{\partial x}\frac{\partial}{\partial y}xy$
如何使用 madness
完成此操作?
更新:实际上在这里我猜它是因数..也许这可以通过将两个导数相乘来解决?也许只有当 x 是 y 的函数时这才会困难。
这是在madness
中使用numderiv
函数的方法:
library(madness)
dxdy <- function(x, y, f) {
dy <- function(x, y) {
dvdx(f(x, y))
}
numderiv(dy, x, y = y)
}
x <- matrix(1:3, nrow = 1)
y <- matrix(3:1, ncol = 1)
# identity matrix, as expected
dvdx(dxdy(madness(x), madness(y), function(x, y) x%*%y))
#> [,1] [,2] [,3]
#> [1,] 1 0 0
#> [2,] 0 1 0
#> [3,] 0 0 1
x <- matrix(2, ncol = 1)
y <- matrix(3, ncol = 1)
dvdx(dxdy(madness(x), madness(y), function(x, y) y^x))
#> [,1]
#> [1,] 9.591674
# compare to analytical solution
y^(x-1)*(x*log(y) + 1)
#> [,1]
#> [1,] 9.591674
x <- matrix(1:3, ncol = 1)
y <- matrix(3:1, ncol = 1)
dvdx(dxdy(madness(x), madness(y), function(x, y) sum(y^x)))
#> [,1] [,2] [,3]
#> [1,] 2.098612 0.000000 0
#> [2,] 0.000000 4.772589 0
#> [3,] 0.000000 0.000000 1
# compare to analytical solution
y^(x-1)*(x*log(y) + 1)
#> [,1]
#> [1,] 2.098612
#> [2,] 4.772589
#> [3,] 1.000000
x <- matrix(1:3, ncol = 1)
y <- matrix(3:1, ncol = 1)
dvdx(dxdy(madness(x), madness(y), function(x, y) sum(sin(x*y))))
#> [,1] [,2] [,3]
#> [1,] -1.413352 0.000000 0.000000
#> [2,] 0.000000 2.373566 0.000000
#> [3,] 0.000000 0.000000 -1.413353
# compare to analytical solution
cos(x*y) - x*y*sin(x*y)
#> [,1]
#> [1,] -1.413353
#> [2,] 2.373566
#> [3,] -1.413353
请注意,madness
导数是每个 x
(列)相对于每个 xy
(行)的 y
部分的部分。
madness
包,如前所述 here,非常适合 R 中的 autodiff。
我现在想计算导数 wrt y 的导数 wrt x。
$\frac{\partial}{\partial x}\frac{\partial}{\partial y}xy$
如何使用 madness
完成此操作?
更新:实际上在这里我猜它是因数..也许这可以通过将两个导数相乘来解决?也许只有当 x 是 y 的函数时这才会困难。
这是在madness
中使用numderiv
函数的方法:
library(madness)
dxdy <- function(x, y, f) {
dy <- function(x, y) {
dvdx(f(x, y))
}
numderiv(dy, x, y = y)
}
x <- matrix(1:3, nrow = 1)
y <- matrix(3:1, ncol = 1)
# identity matrix, as expected
dvdx(dxdy(madness(x), madness(y), function(x, y) x%*%y))
#> [,1] [,2] [,3]
#> [1,] 1 0 0
#> [2,] 0 1 0
#> [3,] 0 0 1
x <- matrix(2, ncol = 1)
y <- matrix(3, ncol = 1)
dvdx(dxdy(madness(x), madness(y), function(x, y) y^x))
#> [,1]
#> [1,] 9.591674
# compare to analytical solution
y^(x-1)*(x*log(y) + 1)
#> [,1]
#> [1,] 9.591674
x <- matrix(1:3, ncol = 1)
y <- matrix(3:1, ncol = 1)
dvdx(dxdy(madness(x), madness(y), function(x, y) sum(y^x)))
#> [,1] [,2] [,3]
#> [1,] 2.098612 0.000000 0
#> [2,] 0.000000 4.772589 0
#> [3,] 0.000000 0.000000 1
# compare to analytical solution
y^(x-1)*(x*log(y) + 1)
#> [,1]
#> [1,] 2.098612
#> [2,] 4.772589
#> [3,] 1.000000
x <- matrix(1:3, ncol = 1)
y <- matrix(3:1, ncol = 1)
dvdx(dxdy(madness(x), madness(y), function(x, y) sum(sin(x*y))))
#> [,1] [,2] [,3]
#> [1,] -1.413352 0.000000 0.000000
#> [2,] 0.000000 2.373566 0.000000
#> [3,] 0.000000 0.000000 -1.413353
# compare to analytical solution
cos(x*y) - x*y*sin(x*y)
#> [,1]
#> [1,] -1.413353
#> [2,] 2.373566
#> [3,] -1.413353
请注意,madness
导数是每个 x
(列)相对于每个 xy
(行)的 y
部分的部分。