r - 如何计算 R 中圆拟合的预测区间

Question

我希望用公式 > r² = (x-h)²+(y-k)² 拟合圆来计算半径的预测区间。 r-圆的半径，x，y，是高斯坐标，h，k，标记拟合圆的中心。

# data
x <- c(1,2.2,1,2.5,1.5,0.5,1.7)
y <- c(1,1,3,2.5,4,1.7,0.8)
# using nls.lm from minpack.lm (minimising the sum of squared residuals)
library(minpack.lm)

residFun <- function(par,x,y) {
  res <- sqrt((x-par$h)^2+(y-par$k)^2)-par$r
  return(res)
}
parStart <- list("h" = 1.5, "k" = 2.5, "r" = 1.7)
out <- nls.lm(par = parStart, x = x, y = y, lower =NULL, upper = NULL, residFun)

问题是，predict() 不适用于 nls.lm，因此我尝试使用 nlsLM 计算圆拟合。 (我可以手动计算它，但在创建我的 Designmatrix 时遇到困难)。`

这就是我接下来尝试的:

dat = list("x" = x,"y" = y)
out1 <- nlsLM(y ~ sqrt(-(x-h)^2+r^2)+k, start = parStart )

结果是:

Error in stats:::nlsModel(formula, mf, start, wts) : 
  singular gradient matrix at initial parameter estimates

问题 1a:nlsLM() 如何处理圆拟合？ (优点是通用的 predict() 可用。 问题 1b:如何获得圆拟合的预测区间？

线性回归示例(这是我想要的圆形回归)

attach(faithful)     
eruption.lm = lm(eruptions ~ waiting) 
newdata = data.frame(waiting=seq(45,90, length = 272)) 
# confidence interval
conf <- predict(eruption.lm, newdata, interval="confidence") 
# prediction interval
pred <- predict(eruption.lm, newdata, interval="predict")
# plot of the data [1], the regression line [1], confidence interval [2], and prediction interval [3]
plot(eruptions ~ waiting)
lines(conf[,1] ~ newdata$waiting, col = "black") # [1]
lines(conf[,2] ~ newdata$waiting, col = "red") # [2]
lines(conf[,3] ~ newdata$waiting, col = "red") # [2]
lines(pred[,2] ~ newdata$waiting, col = "blue") # [3]
lines(pred[,3] ~ newdata$waiting, col = "blue") # [3]

亲切的问候

编辑摘要:

编辑1:重新排列nlsLM中的公式，但参数(h，k，r)结果现在在out和out1中不同...

Edit2:添加了 2 个维基百科链接，用于澄清所使用的术语:(参见下文)

confidence interval

prediction interval

Edit3:对问题的一些改写

Edit4:添加了线性回归的工作示例

Answer 1

我很难弄清楚你想做什么。让我说明一下数据是什么样的以及有关“预测”的一些内容。

plot(x,y, xlim=range(x)*c(0, 1.5), ylim=range(y)*c(0, 1.5))
lines(out$par$h+c(-1,-1,1,1,-1)*out$par$r, # extremes of x-coord
      out$par$k+c(-1,1,1,-1 ,-1)*out$par$r, # extremes of y-coord
      col="red")

那么我们所说的“预测区间”是什么？ (我确实意识到您正在考虑一个圆圈，如果您只想在这个背景上绘制一个圆圈，那也很容易。)

lines(out$par$h+cos(seq(-pi,pi, by=0.1))*out$par$r, #center + r*cos(theta)
      out$par$k+sin(seq(-pi,pi, by=0.1))*out$par$r, #center + r*sin(theta)
      col="red")