如何拟合知道R中最大值的二次模型？

Question

我有一个数据框x和y,我知道y的最大值.我想将这些数据拟合到二次模型中.我怎么能在R知道最大值？如果我不知道最大值,我会用lm(y~x + I(x ^ 2))拟合.任何人都可以对此有所了解吗？提前致谢!

Answer 1

你必须尽量减少受约束条件影响的平方和; lm不允许这样的约束,所以你必须使用泛型优化函数,例如optim.这是一种可行的方式.

弥补一些数据.在这里,我会说已知的最大值是50.

set.seed(5)
d <- data.frame(x=seq(-5, 5, len=51))
d$y <- 50 - 0.3*d$x^2 + rnorm(nrow(d))
M <- 50

使用给定的二次和线性系数得到x点的二次曲线并给出最大值M的函数.微积分是直截了当的; 请参阅duffymo的答案以获取详细信息.

qM <- function(a, b, x, M) {
  c <- M - (3*b^2)/(4*a)
  a*x^2 + b*x + c
}

创建一个函数,得到给定二次和线性系数的二次曲线与d中的数据之间的平方和.

ff <- function(ab, d, M) {
  p <- qM(ab[1], ab[2], d$x, M)
  y <- d$y
  sum((p-y)^2)
}

得到普通的lm适合用作起始值.

m0 <- lm(y ~ I(x^2) + x, data=d)
start <- coef(m0)[2:3]

优化ff函数中的系数.

o <- optim(start, ff, d=d, M=M)
o$par

制作一个图表,显示拟合的最大值为50; 原来的lm适合不.

plot(d)
xs <- seq(-5, 5, len=101)
lines(xs, predict(m0, newdata=data.frame(x=xs)), col="gray")
lines(xs, qM(o$par[1], o$par[2], xs, M))
abline(h=50, lty=3)