r - 关于 R 中椭球体的澄清

Question

不知道是否有人可以帮我解决一些关于R包cluster中的函数ellipsoidhull的问题。我用它来查找包含一系列二维点的最小椭圆。例如

library(cluster)
d <- matrix(c(1,2,3,1,3,2),ncol=2)
e <- ellipsoidhull(d)

该函数计算椭圆 OK，返回一个包含椭圆中心和协方差矩阵的结构。

summary(e)
## 'ellipsoid' in 2 dimensions:
## center = ( 2 2 ); squared ave.radius d^2 =  2 
## and shape matrix =
##        [,1]    [,2]
##   [1,] 0.66667 0.33333
##   [2,] 0.33333 0.66667
##     hence, area  =  3.6276 

问题

a) 如何使用这些数据来检查给定点是否属于椭圆？

b) 如何使用这些数据来计算从给定点到椭圆的距离？

Answer 1

我们可以尝试以下方法:

library(cluster)
d <- matrix(c(1,2,3,1,3,2),ncol=2)
e <- ellipsoidhull(d)
eg <- eigen(e$cov)
axes <- sqrt(eg$values)
angle <- atan(eg$vectors[1,1]/eg$vectors[2,1]) # angle of major axis with x axis

# check if the point (xp, yp) belongs to the ellipse with parameters a,b,... with tolerance eps
belongs.to <- function (xp, yp, a, b, x0, y0, alpha, eps=1e-3) {
  return(abs((cos(alpha)*(xp-x0)+sin(alpha)*(yp-y0))^2/a^2+(sin(alpha)*(xp-x0)-cos(alpha)*(yp-y0))^2/b^2 - 1) <= eps)
} 

# check if the point (xp, yp) is inside the ellipse with parameters a,b,...
is.inside <- function (xp, yp, a, b, x0, y0, alpha) {
  return((cos(alpha)*(xp-x0)+sin(alpha)*(yp-y0))^2/a^2+(sin(alpha)*(xp-x0)-cos(alpha)*(yp-y0))^2/b^2 <= 1)
}

# plot ellipse
plot(e$loc, xlim=c(0,4), ylim=c(0,4), main = "ellipsoidhull", xlab='x', ylab='y')
lines(predict(e), col="blue")
points(rbind(e$loc), col = "red", cex = 3, pch = 13)

x0 <- e$loc[1] # centroid locations
y0 <- e$loc[2]  
a <- sqrt(e$d2) * axes[1]  # major axis length
b <- sqrt(e$d2) * axes[2]  # minor axis length

alpha <- angle
xp <- 3
yp <- 2.9
is.inside(xp, yp, a, b, x0, y0, alpha)
# [1] TRUE
points(xp, yp, pch=19, col='green')
xp <- 3
yp <- 3.1
is.inside(xp, yp, a, b, x0, y0, alpha)
# [1] FALSE
points(xp, yp, pch=19, col='blue')
xp <- 3
yp <- 3
belongs.to(xp, yp, a, b, x0, y0, alpha)
# [1] TRUE
points(xp, yp, pch=19, col='pink')


# distance of a point from the center of the ellipse
sqrt((xp-x0)^2+(yp-y0)^2)