r - 使用 `D` 获取回归多项式的导数时出错

Question

我在二维图形上有点。我想找到适合这个模型的最佳第三多项式并得到它的一阶导数。但我无法得到D功能工作。这是一个简单的例子:

a <- 0:10
b <- c(2, 4, 5, 8, 9, 12, 15, 16, 18, 19, 20)
plot(a, b)
m1 <- lm(b ~ a + I(a ^ 2) + I(a ^ 3))
s <- coef(m1)

## try to get 1st derivative of the regression polynomial
D(expression(s[1] + s[2] * a + (a ^ 2) * s[3] + (a ^ 3) * s[4]), "a")

Error in D(expression(s[1] + s[2] * a + (a^2) * s[3] + (a^3) * s[4]), :

Function '[' is not in the derivatives table

我想避免通过差分计算数值导数。感谢帮助!

Answer 1

您看到的错误消息“ Function ' [ ' 不在导数表 中是因为 D只能识别一组特定的符号操作函数。您可以在 ?D 中找到它们:

The internal code knows about the arithmetic operators ‘+’, ‘-’,
‘*’, ‘/’ and ‘^’, and the single-variable functions ‘exp’, ‘log’,
‘sin’, ‘cos’, ‘tan’, ‘sinh’, ‘cosh’, ‘sqrt’, ‘pnorm’, ‘dnorm’,
‘asin’, ‘acos’, ‘atan’, ‘gamma’, ‘lgamma’, ‘digamma’ and
‘trigamma’, as well as ‘psigamma’ for one or two arguments (but
derivative only with respect to the first).  (Note that only the
standard normal distribution is considered.)

而 "["实际上是 R 中的一个函数(阅读 ?Extract 或 ?"[" )。

要演示类似的行为，请考虑:

s <- function (x) x

D(expression(s(x) + x ^ 2), name = "x")
# Error in D(expression(s(x) + x^2), name = "x") : 
#  Function 's' is not in the derivatives table

在这里，即使 s已定义为简单函数， D对它无能为力。

我最近对 ​​ Function for derivatives of polynomials of arbitrary order (symbolic method preferred) 的回答已经解决了您的问题.我的三个答案中提供了三种方法，其中没有一个基于数值导数。我个人更喜欢 the one using outer (LaTeX 数学公式的唯一答案)，至于多项式，一切都是精确的。

要使用该解决方案，请使用函数 g在那里，并指定参数 x通过要评估导数的值(比如 0:10 )和 pc由您的多项式回归系数 s .默认情况下， nderiv = 0L所以多项式本身就像 predict.lm(m1, newdata = list(a = 0:10)) 一样返回被称为。但曾经 nderiv指定，您将获得回归曲线的精确导数。

a <- 0:10
b <- c(2, 4, 5, 8, 9, 12, 15, 16, 18, 19, 20)
plot(a, b)
m1 <- lm(b ~ a + I(a ^ 2) + I(a ^ 3))
s <- coef(m1)
#(Intercept)           a      I(a^2)      I(a^3) 
# 2.16083916  1.17055167  0.26223776 -0.02020202 

## first derivative at your data points
g(0:10, s, nderiv = 1)
# [1] 1.1705517 1.6344211 1.9770785 2.1985237 2.2987568 2.2777778 2.1355866
# [8] 1.8721834 1.4875680 0.9817405 0.3547009

其他备注:建议您使用 poly(a, degree = 3, raw = TRUE)而不是 I() .他们在这里做同样的事情，但 poly更简洁，如果你想要交互，让它更容易，比如 How to write interactions in regressions in R?