5

我想知道是否有人知道如何计算我的置信椭圆内的蓝色阴影区域。任何建议或看的地方都非常感谢。另外，我希望找到一个通用公式，因为椭圆不一定必须位于应用中的那个区域（即，椭圆可能更大）。如果有帮助，这是我的图片代码：

library(car)

x = c(7,4,1)
y = c(4,6,7)

plot(x,y,xlim=c(0,10),ylim=c(0,10))
rect(x,y,max(x)+100000,max(y)+100000,col="lightblue",border=NA)
points(x,y,col="red",pch=19)

ellipse(center=c(3.5,5),shape=matrix(c(1,.5,.5,2),nrow=2),radius=3,col="green")

在此处输入图像描述

4

1 回答 1

5

如果您可以将椭圆和蓝色区域转换为SpatialPolygons对象，那么使用rgeos包中的函数来计算它们的面积、它们的交点面积以及所有其他类型的有趣量就很容易了。

不幸的是，以正确的形式获取对象需要一些繁重的工作：

library(car)
library(sp)
library(rgeos)

## Function for creating a SpatialPolygons object from data.frame of coords
xy2SP <- function(xy, ID=NULL) {
    if(is.null(ID)) ID <- sample(1e12, size=1)
    SpatialPolygons(list(Polygons(list(Polygon(xy)), ID=ID)),
                    proj4string=CRS("+proj=merc"))
}

## Ellipse coordinates
plot.new()  # Needed by ellipse() 
ell <- ellipse(center=c(3.5,5),shape=matrix(c(1,.5,.5,2),nrow=2),radius=3)
dev.off()   # Cleaning up after plot.new()

## Three rectangles' coordinates in a length-3 list
x <- c(7,4,1)
y <- c(4,6,7)
mx <- max(x) + 1e6
my <- max(y) + 1e6
rr <- lapply(1:3, function(i) {
    data.frame(x = c(x[i], x[i], mx, mx, x[i]),
               y = c(y[i], my, my, y[i], y[i]))
})


## Make two SpatialPolygons objects from ellipse and merged rectangles
ell <- xy2SP(ell)
rrr <- gUnionCascaded(do.call(rbind, lapply(rr, xy2SP)))

## Find area of their intersection
gArea(gIntersection(ell, rrr))
# [1] 10.36296
于 2013-05-14T23:16:38.163 回答