【问题标题】:Computing the Jacobian matrix for x^2*y, 5*x+sin(y)计算 x^2*y, 5*x+sin(y) 的雅可比矩阵
【发布时间】:2013-12-03 13:06:28
【问题描述】:

如果我有一个 R^n --> R^m 函数,我如何创建它的雅可比矩阵?

例如:

expression( x^2*y, 5*x+sin(y) ) # f : R^2 --> R^2

我想要一个表达式矩阵,例如:

expression(2*x*y) expression(x^2)
expression(5) expression(cos(y))

有什么解决办法吗?

【问题讨论】:

    标签: r function matrix differentiation


    【解决方案1】:

    一开始你的问题似乎有点含糊(尽管我们后来确实提供了更多细节)。有一个deriv 函数允许简单的符号区分。由于您一开始没有提供任何数据,因此在第二条评论之前,我认为您想要象征性的结果。

    ex <- expression( x^2*y, 5*x+sin(y) )
    sapply(ex, deriv, c("x", "y") )
    #-------result--------------------
    expression({
        .expr1 <- x^2
        .value <- .expr1 * y
        .grad <- array(0, c(length(.value), 2L), list(NULL, c("x", 
            "y")))
        .grad[, "x"] <- 2 * x * y
        .grad[, "y"] <- .expr1
        attr(.value, "gradient") <- .grad
        .value
    }, {
        .value <- 5 * x + sin(y)
        .grad <- array(0, c(length(.value), 2L), list(NULL, c("x", 
            "y")))
        .grad[, "x"] <- 5
        .grad[, "y"] <- cos(y)
        attr(.value, "gradient") <- .grad
        .value
    })
    

    您也可以使用它,它会返回一个函数列表:

    sapply(ex, deriv, c("x", "y") , func=TRUE)
    

    [编辑 1:回答第一条评论] 如果您希望对象允许单个矩阵函数,请尝试以下操作:

    res <- sapply(ex, function(ex1) lapply(c("x","y") , 
                                function(arg) deriv(ex1, arg, func=TRUE) ) )
    
    #--------------
    > res[1,1]
    [[1]]
    function (x) 
    {
        .value <- x^2 * y
        .grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
        .grad[, "x"] <- 2 * x * y
        attr(.value, "gradient") <- .grad
        .value
    }
    

    或者对于表达式版本:

    res <- sapply(ex, function(ex1) lapply(c("x","y") , 
                                       function(arg) deriv(ex1, arg) ) )
    #----------------------
    > res[1,1]
    [[1]]
    expression({
        .value <- x^2 * y
        .grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
        .grad[, "x"] <- 2 * x * y
        attr(.value, "gradient") <- .grad
        .value
    })
    

    [编辑 2:回答第二条评论] 要评估表达式,请将其传递给具有适当环境的 eval

    > eval( res[1,1][[1]], envir=list(x=4,y=5) )
    [1] 80
    attr(,"gradient")
          x
    [1,] 40
    

    需要先使用“[”然后使用“[[”来提取矩阵中列表项的第一个(也是唯一一个)元素,这是一个琐碎的问题。

    【讨论】:

    • 我如何访问差异化表达式?我需要多次评估它们。
    • 当我需要评估具体的 x、y 值时,如何评估 res[1,1] 与 x=4 和 y=5?
    • 您现在了解您的问题是如何不完整的吗?见上文。
    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2019-01-07
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多