求解最佳拟合多项式和绘图下拉线

我在Windows 10上使用R 3.3.1(64位).我有一个xy数据集,我适合二阶多项式.我想在y = 4处求解x的最佳拟合多项式,并绘制从y = 4到x轴的下拉线.

这将在数​​据帧v1中生成数据:

v1 <- structure(list(x = c(-5.2549, -3.4893, -3.5909, -2.5546, -3.7247, 
-5.1733, -3.3451, -2.8993, -2.6835, -3.9495, -4.9649, -2.8438, 
-4.6926, -3.4768, -3.1221, -4.8175, -4.5641, -3.549, -3.08, -2.4153, 
-2.9882, -3.4045, -4.6394, -3.3404, -2.6728, -3.3517, -2.6098, 
-3.7733, -4.051, -2.9385, -4.5024, -4.59, -4.5617, -4.0658, -2.4986, 
-3.7559, -4.245, -4.8045, -4.6615, -4.0696, -4.6638, -4.6505, 
-3.7978, -4.5649, -5.7669, -4.519, -3.8561, -3.779, -3.0549, 
-3.1241, -2.1423, -3.2759, -4.224, -4.028, -3.3412, -2.8832, 
-3.3866, -0.1852, -3.3763, -4.317, -5.3607, -3.3398, -1.9087, 
-4.431, -3.7535, -3.2545, -0.806, -3.1419, -3.7269, -3.4853, 
-4.3129, -2.8891, -3.0572, -5.3309, -2.5837, -4.1128, -4.6631, 
-3.4695, -4.1045, -7.064, -5.1681, -6.4866, -2.7522, -4.6305, 
-4.2957, -3.7552, -4.9482, -5.6452, -6.0302, -5.3244, -3.9819, 
-3.8123, -5.3085, -5.6096, -6.4557), y = c(0.99, 0.56, 0.43, 
2.31, 0.31, 0.59, 0.62, 1.65, 2.12, 0.1, 0.24, 1.68, 0.09, 0.59, 
1.23, 0.4, 0.36, 0.49, 1.41, 3.29, 1.22, 0.56, 0.1, 0.67, 2.38, 
0.43, 1.56, 0.07, 0.08, 1.53, -0.01, 0.12, 0.1, 0.04, 3.42, 0.23, 
0, 0.34, 0.15, 0.03, 0.19, 0.17, 0.2, 0.09, 2.3, 0.07, 0.15, 
0.18, 1.07, 1.21, 3.4, 0.8, -0.04, 0.02, 0.74, 1.59, 0.71, 10.64, 
0.64, -0.01, 1.06, 0.81, 4.58, 0.01, 0.14, 0.59, 7.35, 0.63, 
0.17, 0.38, -0.08, 1.1, 0.89, 0.94, 1.52, 0.01, 0.1, 0.38, 0.02, 
7.76, 0.72, 4.1, 1.36, 0.13, -0.02, 0.13, 0.42, 1.49, 2.64, 1.01, 
0.08, 0.22, 1.01, 1.53, 4.39)), .Names = c("x", "y"), class = "data.frame", row.names = c(NA, 
-95L))

这是绘制y vs x的代码,绘制最佳拟合多项式,并在y = 4处绘制一条线.

> attach(v1)
> # simple x-y plot of the data
> plot(x,y, pch=16)
> # 2nd order polynomial fit
> fit2 <- lm(y~poly(x,2,raw=TRUE))
> summary(fit2)
> # generate range of numbers for plotting polynomial
> xx <- seq(-8,0, length=50)
> # overlay best fit polynomial
>lines(xx, predict(fit2, data.frame(x=xx)), col="blue")
> # add horizontal line at y=4
> abline(h=4, col="red")
>

从图中可以明显看出y = 4时x约为-2和-6.5,但我想实际求解这些值的回归多项式.

理想情况下,我喜欢从红蓝线交叉点下降到x轴的线(即绘制垂直上升,终止于两个y = 4解).如果那是不可能的话,我会很高兴看到好的旧垂直上升,一直到情节,只要他们在正确的x解决方案值.

此图表示当y> 4时将超出规范的部分,因此我想使用下拉线突出显示将生成符合规范的部分的x值范围.

您可以使用二次公式计算值:

betas <- coef(fit2)    # get coefficients
betas[1] <- betas[1] - 4    # adjust intercept to look for values where y = 4

# note degree increases, so betas[1] is c, etc.
betas
##             (Intercept) poly(x, 2, raw = TRUE)1 poly(x, 2, raw = TRUE)2 
##               8.7555833               6.0807302               0.7319848 

solns <- c((-betas[2] + sqrt(betas[2]^2 - 4 * betas[3] * betas[1])) / (2 * betas[3]), 
           (-betas[2] - sqrt(betas[2]^2 - 4 * betas[3] * betas[1])) / (2 * betas[3]))

solns
## poly(x, 2, raw = TRUE)1 poly(x, 2, raw = TRUE)1 
##               -1.853398               -6.453783 

segments(solns, -1, solns, 4, col = 'green')    # add segments to graph

更简单(如果你能找到)是polyroot:

polyroot(betas)
## [1] -1.853398+0i -6.453783+0i

由于它返回一个复杂的向量,因此as.numeric如果要将其传递给它,则需要将其包装起来segments.

我完全理解这个简单的二次多项式有一个解析解.我向您展示数值解的原因是您在回归设置中提出这个问题.当您有更复杂的回归曲线时,数值解可能总是您的解决方案.

在下面我将使用uniroot函数.如果您不熟悉它,请首先阅读以下简短答案:在R中解决问题.

这是使用您的代码生成的图.你快到了.这是一个根查找问题,您可以在数字上使用uniroot.让我们定义一个函数:

f <- function (x) {
  ## subtract 4
  predict(fit2, newdata = data.frame(x = x)) - 4
  }

从图中可以看出,有两个根,一个在里面[-7, -6],另一个在里面[-3, -1].我们uniroot用来找到两个:

x1 <- uniroot(f, c(-7, -6))$root
#[1] -6.453769

x2 <- uniroot(f, c(-3, -1))$root
#[1] -1.853406

现在,您可以将这些点的垂直线向下拖动到x轴:

y1 <- f(x1) + 4  ## add 4 back
y2 <- f(x2) + 4  

abline(h = 0, col = 4)  ## x-axis
segments(x1, 0, x1, y1, lty = 2)
segments(x2, 0, x2, y2, lty = 2)

你有一个二次方程

0.73198 * x^2 + 6.08073 * x + 12.75558 = 4
OR
0.73198 * x^2 + 6.08073 * x + 8.75558 = 0

您可以使用二次公式来解析分析.R给出了两个根源:

(-6.08073 + sqrt(6.08073^2 -4*0.73198 * 8.75558)) / (2 * 0.73198)
[1] -1.853392
(-6.08073 - sqrt(6.08073^2 -4*0.73198 * 8.75558)) / (2 * 0.73198)
[1] -6.453843

abline(v = c(-1.853392,-6.453843))