如何计算在R中椭圆的交集体积

3

我想知道如何计算两个椭圆之间的交集,例如在此图表中所示的versicolor和virginca之间的交集的体积: PCA on iris dataset 该图表是基于这个教程的以下mwe绘制的:

data(iris)
log.ir <- log(iris[, 1:4])
ir.species <- iris[, 5]
ir.pca <- prcomp(log.ir, center = TRUE, scale. = TRUE)

library(ggbiplot)
g <- ggbiplot(ir.pca, obs.scale = 1, var.scale = 1, 
          groups = ir.species, ellipse = TRUE,
          circle = TRUE)
g <- g + scale_color_discrete(name = '')
g <- g + theme(legend.direction = 'horizontal', 
           legend.position = 'top')
print(g)

我将按照以下方式获取椭圆的协方差和中心点:

我获取椭圆的协方差和中心点如下:

setosa.cov <- cov(ir.pca$x[ir.species=="setosa",])
versicolor.cov <- cov(ir.pca$x[ir.species=="versicolor",])
virginica.cov <- cov(ir.pca$x[ir.species=="virginica",])
setosa.centre <- colMeans(ir.pca$x[ir.species=="setosa",])
versicolor.centre <- colMeans(ir.pca$x[ir.species=="versicolor",])
virginica.centre <- colMeans(ir.pca$x[ir.species=="virginica",])

但是我已经无计可施了 :-|

编辑: 根据下面@carl-witthoft的指示,这里有一个使用siar::overlap的例子:

library(siar)
setosa <- ir.pca$x[ir.species=="setosa",]
versicolor <- ir.pca$x[ir.species=="versicolor",]
virginica <- ir.pca$x[ir.species=="virginica",]

overlap.fun <- function(data.1, data.2){
   dimensions <- ncol(data.1)
   for(i in 1:(dimensions-1)){
      overlap.out <- overlap(data.1[,i], data.1[,i+1], data.2[,i], data.2[,i+1], steps = 5)
      out$overlap[i] <- overlap.out$overlap
      out$area1[i] <- overlap.out$area1
      out$area2[i] <- overlap.out$area2
   }
   return(out)
}

overlap.fun(versicolor, virginica)
返回值:
$overlap
[1] 0.01587977 0.48477088 0.08375927
$area1
[1]1.020596 1.04614461 0.08758691                 
$area2
[1] 1.028594 1.1535106 0.1208483

奇怪的是,当我进行百分比计算时,结果并不真正对应于ggbiplot PCA中的椭球体:

tmp <- overlap(versicolor[,1], versicolor[,2], virginica[,1], virginica[,2], steps = 5)
virginica.percentage <- round(x=(tmp$overlap/tmp$area2*100), digits = 2)
versicolor.percentage <- round(x=(tmp$overlap/tmp$area1*100), digits = 2)
> virginica.percentage [1] 1.54
> versicolor.percentage[1] 1.56

这比上图1中所示的要少得多。 但是最好在这里开一个新线程讨论。


基本方法是找到交点,计算“上”和“下”曲线的积分,并取距离。您需要将其拆分以确保每个积分都在单值函数范围内。话虽如此,我似乎记得CRAN上有一两个包包括这种交叉区域计算。自然地,我想不起来是哪些了 :-( - Carl Witthoft
1个回答

4

可能的工具:

 spatstat::overlap.owin , geo::geointersect, siar::overlap .

你可能会问,而且你应该问——“他是如何如此快速地得到这些答案的呢?”
获取包sos并输入???overlap即可。

非常感谢您的快速回答!我尝试使用siar::overlap,因为它对我来说是最易懂的。library(siar) versicolor <- ir.pca$x[ir.species=="versicolor",] virginica <- ir.pca$x[ir.species=="virginica",] overlap(versicolor[,1], versicolor[,2], virginica[,1], virginica[,2], steps = 5)返回: '$overlap' [1] 0.01587977
$area1 [1] 1.020596
$area2
[1] 1.028594
- raumkundschafter
@sesselastronaut。很高兴它能够正常工作 - 如果它完全满足您的需求,请将答案标记为“已接受”。 - Carl Witthoft
@carl-wittcroft,就像示例所示,它确实适用于前两个主成分。我想知道如何在鸢尾花数据集上的所有四个主成分中执行此操作?将它们全部比较并相加吗? - raumkundschafter
spatstat::overlap.owin 似乎只接受“Windows(类‘owin’对象)或可被转换为as.owin的数据。”。如何将椭圆转换为这样的对象可能是另一个主题。geo::geointersect同样可以与多边形一起使用。 - raumkundschafter

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接