Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Plotting a function and a derivative function

Tags:

plot

r

derivative

I would like to plot a dataframe (X,Y) data together with a fitted function and the derivative of the fitted function.

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

As far as I understood, this gives me a polynomial function of the 32nd degree, fit, that I use to calculate the function values and store them in data$fitted. Plotting these series works like a charm with ggplot2.

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

Plot

So far so good. But what I'm would like to plot too is the first derivative, data$Y', of the fitted function fit. What I'm interested in is the gradient of the fitted function.

My Question: How can I get the derivative function of fit? I assume I can then "predict" the absolute values for plotting afterwards. Correct?

like image 507
Oliver Avatar asked Jun 20 '26 05:06

Oliver


1 Answers

First, i'll create some test data that "kind of" looks like yours

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")

enter image description here

Then, because you have a special form of a polynomial we can easily calculate the derivative by multiplying each of the coefficients by the power and shifting all the terms down. Here's a helper function to calcualte the new coefficients

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
}

which we can use like...

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

And now I can plot

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")

enter image description here

and we can see that the inflection points correspond to places where the derivative is zero.

like image 159
MrFlick Avatar answered Jun 23 '26 02:06

MrFlick



Donate For Us

If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!