I got some leads from this StackOverflow thread and Gavin Simpson’s blog.
library(RCurl)
## Warning: package 'RCurl' was built under R version 3.5.1
## Loading required package: bitops
script <- getURL("https://gist.githubusercontent.com/gavinsimpson/ca18c9c789ef5237dbc6/raw/095c9be4d3654b5a8c05aaa6b9037ad1bdab53b3/derivSimulCI.R", ssl.verifypeer = FALSE)
eval(parse(text = script))
Let us simulate some data:
x <- seq(0, 1, length.out = 100)
# setting arbitrary a and b values
a <- 30
b <- 15
# simulate Y with noise
logY <- ((a*x) / (1 + b*x)) + rnorm(length(x), 0, 0.1)
library(mgcv)
m <- gam(logY ~ s(x))
logYpred <- predict(m)
plot(x, logY)
lines(x, logYpred)
Courtesy of Gavin Simpson, we can easily assess the slope by calculating the first derivative of the GAM spline. When the first derivative (i.e. slope / gradient) is not significantly different from zero, we have reached the flat plateau and can (somehow) conclude that it’s the inflection point.
fd <- derivSimulCI(m) # m is the gam model object
## Loading required package: MASS
plot(fd, sizer = TRUE)
Now we locate the first time (minimum x value) that the slope/gradient becomes zero.
CI <- lapply(fd[1], function(x) t(apply(x$simulations, 1, quantile, probs = c(0.025, 0.975))))
first.zero.slope.index <- min(which(sign(CI$x[, "2.5%"]) != sign(CI$x[, "97.5%"])))
fd$eval[first.zero.slope.index]
## [1] 0.1859295
The inflection point, as per our definition, is at x = 0.186!