3D表面插值

Question

我有一个包含纬度、经度和地下测量值的三列数据框。我正在尝试找出如何在我拥有的点之间进行数据插值（这些点是不规则分布的），然后创建整个区域的平滑曲面图。我尝试使用'rgl'包中的'surface3d'函数，但我的结果看起来像一个单独的巨大尖峰。我已经能够用'plot3d'绘制出数据，但我需要更进一步地使用插值填补空白处。有什么想法或建议吗？我也可以考虑使用其他包，只是当时'rgl'似乎是最合适的选择。

编辑：以下是我的数据摘录（含含水层深度的测量值）：

 lat_dd_NAD83 long_dd_NAD83 lev_va_ft
1     37.01030     -101.5006    288.49
2     37.03977     -101.6633    191.68
3     37.05201     -100.4994    159.34
4     37.06567     -101.3292    174.07
5     37.06947     -101.4561    285.08
6     37.10098     -102.0134    128.94

- ffarg

---

您可以在asbio包中使用loess.surf函数。如果您能提供数据，我们可以更好地帮助您。 - adaien

使用gstat包进行Kriging估计看起来是一个显然要尝试的技术。 - Richard Telford

“3D表面”的各种解释都是可能的。您是否试图根据“lat”和“long”的基础插值测量的“lev”值？（我会称之为嵌入三维空间中的2D流形。如果是这样，那么在RHS上具有2个变量的loess建议（尽管在统计包中）似乎是正确的，并且已经有了一个说明性答案：http://stackoverflow.com/questions/15019725/r-3d-surface-plot） - IRTFM

也许你会在 demo(persp) 中找到一些有用的东西。 - RHertel

Answer 1

关于插值，我想补充一点小但可能很重要的说明。

使用非常好用的"akima"软件包，您可以轻松地对数据进行插值：

library(akima)
library(rgl)
# library(deldir)

# Create some fake data
x <- rnorm(100)
y <- rnorm(100)
z <- x^2 + y^2

# # Triangulate it in x and y
# del <- deldir(x, y, z = z)
# triangs <- do.call(rbind, triang.list(del))
# 
# # Plot the resulting surface
# plot3d(x, y, z, type = "n")
# triangles3d(triangs[, c("x", "y", "z")], col = "gray")

n_interpolation <- 200

spline_interpolated <- interp(x, y, z,
                              xo=seq(min(x), max(x), length = n_interpolation),
                              yo=seq(min(y), max(y), length = n_interpolation),
                              linear = FALSE, extrap = TRUE)

x.si <- spline_interpolated$x
y.si <- spline_interpolated$y
z.si <- spline_interpolated$z

persp3d(x.si, y.si, z.si, col = "gray")

样条插值图片（200个步骤）

使用此软件包，您可以轻松更改插值步数等参数。使用此软件包至少需要10个点（点越多越好），才能获得合理的样条插值。线性版本无论点的数量如何都可以正常工作。

P.S. 感谢2554330用户 - 对于deldir不了解，这在某些情况下真的是非常有用的东西。

- Dmitriy

Answer 2

您可以使用 deldir 包来获得您的点的 Delaunay 三角剖分，然后将其转换为绘图所需的数据形式，以供 triangles3d 使用。我不知道这在一个非常大的数据集上会有多有效，但它似乎对于100个点是有效的。

library(deldir)
library(rgl)
# Create some fake data
x <- rnorm(100)
y <- rnorm(100)
z <- x^2 + y^2

# Triangulate it in x and y
del <- deldir(x, y, z = z)
triangs <- do.call(rbind, triang.list(del))

# Plot the resulting surface
plot3d(x, y, z, type = "n")
triangles3d(triangs[, c("x", "y", "z")], col = "gray")

编辑后添加：

R-forge上的rgl版本现在有一个函数可以轻松实现这一点。您现在可以使用此函数生成类似于上面的图形。

library(deldir)
library(rgl)
plot3d(deldir(x, y, z = z))

还有一个函数可以从 deldir() 的输出构建 mesh3d 对象。

- user2554330