在R中计算不规则形状的面积

Question

在R中计算不规则形状的面积

4

我希望计算R中以下图形的封闭区域：

代码如下：

a <- c(0,1,2,2,1,0,-1,-2,-2,-1,0)
b <- c(0,0,0,1,0,0,0,0,1,0,0)
id <- order(a)
AUC <- sum(abs(diff(a[id])*rollmean(b[id],2)))

AUC的结果为0.5。为什么不是1？

当我对不同的向量执行相同操作时，结果如下：代码如下：

a <- c(0,1,1,0,0,0,1,1,0,0,0)
b <- c(1,1,2,1,0,-1,-1,-2,-1,0,1) 
id <- order(a) 
AUC <- sum(abs(diff(a[id])*rollmean(b[id],2)))

结果是1，这正是我所预期的。

最后，也是最重要的一个问题是：我有一个数据框，其中有两列x和y（100行，2列），当将x与y绘制在一起时，形成了一个不规则形状的封闭曲线（循环图），就像这张图片中的一样：我想使用这两列x和y计算该形状的面积，有没有办法做到这点？

更新： 第三个例子（循环图）的数据如下：

structure(c(-0.317857301365341, -0.27254852000745, -0.239116190750992, 
-0.212617899743053, -0.188665601384051, -0.164652457277714, -0.516811363687365, 
-0.704142203990645, -0.925833008271225, -1.17644987082148, -1.44698060870818, 
-1.72501192276032, -1.99567842934217, -2.24329812421979, -2.45336219338653, 
-2.61434899714043, -2.71894439436842, -2.76453987881648, -2.75307392975534, 
-2.69034804478649, -2.58498022688992, -2.44717942836848, -2.28751639917794, 
-2.11583498227994, -1.94040889261091, -1.76740411983085, -1.37034128566268, 
-0.945317248194438, -0.603583991912112, -0.245876237842544, 0.12405281597095, 
0.500705070183803, 1.07928003293753, 1.3786192588209, 1.93625569819096, 
2.32102987210859, 2.67667485884349, 2.9817217810439, 3.21508203608072, 
3.35805065942255, 3.39629843099713, 3.32158207702952, 3.13297595252244, 
2.83749468084523, 2.44998880428773, 1.99219777869042, 1.34133852767483, 
0.905097245217795, 0.332495432968877, 0.000565676359279427, -0.304679355313777, 
-0.676576930379378, -0.799504493858939, -0.877037748519715, -0.907917643885644, 
-0.895574022744618, -0.847333419315191, -0.773304922267926, -0.68506594358595, 
-0.594261154658402, -0.511242949335073, -0.443908244418344, -0.396899995824444, 
-0.371311027504522, -2.95318107124334, -2.91494448796198, -2.80017869044755, 
-2.61894793180667, -2.38576761259564, -2.11805651326137, -1.41368736155233, 
-0.896608829738209, -0.363137342157785, 0.165324041578864, 0.668413502524084, 
1.12765173596156, 1.52716559439287, 1.85443533430032, 2.10104197466128, 
2.26318484127474, 2.34177274880676, 2.34207035998987, 2.2729829455437, 
2.14607786936849, 1.97445291610407, 1.77157880225516, 1.5502403566333, 
1.32167995293426, 1.09501806084454, 0.876985405337858, 0.780951474855999, 
0.715653776169357, 0.631076472349792, 0.572023400386015, 0.534441560957589, 
0.513328184556074, 0.275997813967411, 0.34478353180991, 0.504249489908016, 
0.587434743269948, 0.654235885315243, 0.701450060739606, 0.72746455454756, 
0.732567090012129, 0.718940924224628, 0.690264488986728, 0.650927827300914, 
0.604983707715769, 0.555050419060634, 0.501440460475964, 0.370959009632643, 
0.337726723524854, 0.216270960220213, 0.162458598954556, 0.060364386887637, 
-0.136702597984821, -0.309488796931729, -0.535968167790145, -0.807754519418657, 
-1.11285625385781, -1.43643730574278, -1.76172554792592, -2.0710195962479, 
-2.34678141009815, -2.57280751486803, -2.73543188753935, -2.82464965255479, 
-2.83501170706251), .Dim = c(64L, 2L), .Dimnames = list(c("309", 
"310", "311", "312", "313", "314", "315", "316", "317", "318", 
"319", "320", "321", "322", "323", "324", "325", "326", "327", 
"328", "329", "330", "331", "332", "333", "334", "335", "336", 
"337", "338", "339", "340", "341", "342", "343", "344", "345", 
"346", "347", "348", "349", "350", "351", "352", "353", "354", 
"355", "356", "357", "358", "359", "360", "361", "362", "363", 
"364", "365", "366", "367", "368", "369", "370", "371", "372"
), c("PC1", "PC2")))

- Dori

1

如果您提供圆图的数据点，我们就可以进行测试！ - Waldi

嗨Waldi，谢谢你的回答！我已经在帖子的<updated>部分提供了循环图的点。 - Dori

1

请查看我的编辑，其中包含cyclogram的近似值。 - Waldi

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Waldi · Accepted Answer

不确定您使用的方法和为什么在这种情况下无法正常工作。
另一种选择是使用 sf（simple features），它可以处理面积计算：

library(sf)  

a <- c(0,1,2,2,1,0,-1,-2,-2,-1,0)
b <- c(0,0,0,1,0,0,0,0,1,0,0)

poly <- st_polygon(list(cbind(a,b)))
poly <- st_make_valid(poly)
plot(poly)

st_area(poly)
#> [1] 1

注意使用st_make_valid因为多边形应该是非交叉的，而这在您提供的三个示例中并不是这种情况。
在第一个示例中，多边形被转换为2个多边形和2个线串。

st_make_valid(poly)
GEOMETRYCOLLECTION (MULTIPOLYGON (((-2 1, -1 0, -2 0, -2 1)), ((2 0, 1 0, 2 1, 2 0))), 
                    MULTILINESTRING ((0 0, -1 0), (0 0, 1 0)))

对于圆形图，如果您想要精确的面积，则需要轨迹的交点，而这在提供的数据中并非如此。
作为近似值，您可以使用concaveman：

library(concaveman)
st_polygon(list(concaveman(data)))

plot(st_polygon(list(concaveman(data))), col='grey')
plot(st_linestring(data),add=T,col='red')

st_area(st_polygon(list(concaveman(data))))
[1] 2.693168