r - 绘制函数和导数函数

Question

我想绘制数据框 (X,Y) data 以及拟合函数和拟合函数的导数。

fit <- lm(data$Y ~ poly(data$X,32,raw=TRUE))
data$fitted_values <- predict(fit, data.frame(x=data$X))

据我所知，这给了我一个 32 次多项式函数 fit，我用它来计算函数值并将它们存储在 data$fitted。使用 ggplot2 绘制这些系列就像一个魅力。

ggplot(data, aes(x=X)) + 
    geom_line(aes(y = Y), colour="red") + 
    geom_line(aes(y = predict), colour="blue")

到目前为止一切顺利。但我也想绘制的是拟合函数 fit 的一阶导数 data$Y'。我感兴趣的是拟合函数的梯度。

我的问题:如何得到fit的导函数？我假设我可以“预测”之后绘制的绝对值。正确吗？

Answer 1

首先，我将创建一些“有点”像你的测试数据

set.seed(15)
rr<-density(faithful$eruptions)
dd<-data.frame(x=rr$x)
dd$y=rr$y+ runif(8,0,.05)

fit <- lm(y ~ poly(x,32,raw=TRUE), dd)
dd$fitted <- fitted(fit)

ggplot(dd, aes(x=x)) + 
    geom_line(aes(y = y), colour="red") + 
    geom_line(aes(y = fitted), colour="blue")

然后，因为你有一个特殊形式的多项式，我们可以通过将每个系数乘以幂并将所有项向下移动来轻松计算导数。这是计算新系数的辅助函数

deriv_coef<-function(x) {
    x <- coef(x)
    stopifnot(names(x)[1]=="(Intercept)")
    y <- x[-1]
    stopifnot(all(grepl("^poly", names(y))))
    px <- as.numeric(gsub("poly\\(.*\\)","",names(y)))
    rr <- setNames(c(y * px, 0), names(x))
    rr[is.na(rr)] <- 0
    rr
}

我们可以像...一样使用它

dd$slope <- model.matrix(fit) %*% matrix(deriv_coef(fit), ncol=1)

现在我可以画画了

ggplot(dd, aes(x=x)) + 
    geom_line(aes(y = y), colour="red") + 
    geom_line(aes(y = fitted), colour="blue") + 
    geom_line(aes(y = slope), colour="green")

我们可以看到拐点对应于导数为零的地方。