想象一下，我有一个 3 列矩阵
x、y、z，其中 z 是 x 和 y 的函数。
我知道如何绘制这些点的“散点图”
plot3d(x,y,z)
但是如果我想要一个表面，我必须使用其他命令，例如 surface3d 问题是它不接受与 plot3d 相同的输入，它似乎需要一个矩阵
(nº elements of z) = (n of elements of x) * (n of elements of x)
我怎样才能得到这个矩阵？我已经尝试过使用命令 interp，就像我在需要使用等高线图时所做的那样。
如何在不计算此矩阵的情况下直接从 x,y,z 绘制曲面？如果我有太多点，这个矩阵就太大了。
干杯
如果您的 x 和 y 坐标不在网格上，那么您需要将 x、y、z 表面插值到一个网格上。您可以使用任何地统计软件包（geoR、gstat 等）或更简单的技术（例如反距离加权）通过克里金法来完成此操作。
我猜你提到的“interp”函数来自 akima 包。请注意，输出矩阵与输入点的大小无关。如果需要，您可以在输入中有 10000 个点并将其插入到 10x10 网格中。默认情况下 akima::interp 在 40x40 网格上执行此操作：
require(akima)
require(rgl)
x = runif(1000)
y = runif(1000)
z = rnorm(1000)
s = interp(x,y,z)
> dim(s$z)
[1] 40 40
surface3d(s$x,s$y,s$z)
因为它的随机数据，这看起来很尖刻和垃圾。希望你的数据不是！
您可以使用该函数outer()来生成它。
看一下函数的演示persp()，它是一个基本的图形函数，用于绘制曲面的透视图。
这是他们的第一个例子：
x <- seq(-10, 10, length.out = 50)
y <- x
rotsinc <- function(x,y) {
sinc <- function(x) { y <- sin(x)/x ; y[is.na(y)] <- 1; y }
10 * sinc( sqrt(x^2+y^2) )
}
z <- outer(x, y, rotsinc)
persp(x, y, z)
这同样适用于surface3d()：
require(rgl)
surface3d(x, y, z)
你可以看看使用 Lattice。在这个例子中，我定义了一个网格，我想在上面绘制 z~x,y。它看起来像这样。请注意，大部分代码只是构建我使用线框函数绘制的 3D 形状。
变量“b”和“s”可以是 x 或 y。
require(lattice)
# begin generating my 3D shape
b <- seq(from=0, to=20,by=0.5)
s <- seq(from=0, to=20,by=0.5)
payoff <- expand.grid(b=b,s=s)
payoff$payoff <- payoff$b - payoff$s
payoff$payoff[payoff$payoff < -1] <- -1
# end generating my 3D shape
wireframe(payoff ~ s * b, payoff, shade = TRUE, aspect = c(1, 1),
light.source = c(10,10,10), main = "Study 1",
scales = list(z.ticks=5,arrows=FALSE, col="black", font=10, tck=0.5),
screen = list(z = 40, x = -75, y = 0))
rgl很棒，但是需要进行一些实验才能使轴正确。
如果你有很多点，为什么不从它们中随机抽样，然后绘制结果表面。您可以根据来自相同数据的样本添加多个表面，以查看采样过程是否严重影响您的数据。
所以，这是一个非常可怕的功能，但它做了我认为你想要它做的事情（但没有采样）。给定一个矩阵 (x, y, z)，其中 z 是高度，它将绘制点和表面。限制是每个 (x,y) 对只能有一个 z。因此，在自身上环回的飞机会引起问题。
将plot_points = T绘制制作表面的各个点 - 这对于检查表面和点是否实际相遇很有用。将plot_contour = T在 3d 可视化下方绘制 2d 等高线图。将颜色设置rainbow为提供漂亮的颜色，其他任何东西都会将其设置为灰色，但是您可以更改该功能以提供自定义调色板。无论如何，这对我有用，但我确信它可以被整理和优化。打印出verbose = T很多输出，我用来在函数中断时调试函数。
plot_rgl_model_a <- function(fdata, plot_contour = T, plot_points = T,
verbose = F, colour = "rainbow", smoother = F){
## takes a model in long form, in the format
## 1st column x
## 2nd is y,
## 3rd is z (height)
## and draws an rgl model
## includes a contour plot below and plots the points in blue
## if these are set to TRUE
# note that x has to be ascending, followed by y
if (verbose) print(head(fdata))
fdata <- fdata[order(fdata[, 1], fdata[, 2]), ]
if (verbose) print(head(fdata))
##
require(reshape2)
require(rgl)
orig_names <- colnames(fdata)
colnames(fdata) <- c("x", "y", "z")
fdata <- as.data.frame(fdata)
## work out the min and max of x,y,z
xlimits <- c(min(fdata$x, na.rm = T), max(fdata$x, na.rm = T))
ylimits <- c(min(fdata$y, na.rm = T), max(fdata$y, na.rm = T))
zlimits <- c(min(fdata$z, na.rm = T), max(fdata$z, na.rm = T))
l <- list (x = xlimits, y = ylimits, z = zlimits)
xyz <- do.call(expand.grid, l)
if (verbose) print(xyz)
x_boundaries <- xyz$x
if (verbose) print(class(xyz$x))
y_boundaries <- xyz$y
if (verbose) print(class(xyz$y))
z_boundaries <- xyz$z
if (verbose) print(class(xyz$z))
if (verbose) print(paste(x_boundaries, y_boundaries, z_boundaries, sep = ";"))
# now turn fdata into a wide format for use with the rgl.surface
fdata[, 2] <- as.character(fdata[, 2])
fdata[, 3] <- as.character(fdata[, 3])
#if (verbose) print(class(fdata[, 2]))
wide_form <- dcast(fdata, y ~ x, value_var = "z")
if (verbose) print(head(wide_form))
wide_form_values <- as.matrix(wide_form[, 2:ncol(wide_form)])
if (verbose) print(wide_form_values)
x_values <- as.numeric(colnames(wide_form[2:ncol(wide_form)]))
y_values <- as.numeric(wide_form[, 1])
if (verbose) print(x_values)
if (verbose) print(y_values)
wide_form_values <- wide_form_values[order(y_values), order(x_values)]
wide_form_values <- as.numeric(wide_form_values)
x_values <- x_values[order(x_values)]
y_values <- y_values[order(y_values)]
if (verbose) print(x_values)
if (verbose) print(y_values)
if (verbose) print(dim(wide_form_values))
if (verbose) print(length(x_values))
if (verbose) print(length(y_values))
zlim <- range(wide_form_values)
if (verbose) print(zlim)
zlen <- zlim[2] - zlim[1] + 1
if (verbose) print(zlen)
if (colour == "rainbow"){
colourut <- rainbow(zlen, alpha = 0)
if (verbose) print(colourut)
col <- colourut[ wide_form_values - zlim[1] + 1]
# if (verbose) print(col)
} else {
col <- "grey"
if (verbose) print(table(col2))
}
open3d()
plot3d(x_boundaries, y_boundaries, z_boundaries,
box = T, col = "black", xlab = orig_names[1],
ylab = orig_names[2], zlab = orig_names[3])
rgl.surface(z = x_values, ## these are all different because
x = y_values, ## of the confusing way that
y = wide_form_values, ## rgl.surface works! - y is the height!
coords = c(2,3,1),
color = col,
alpha = 1.0,
lit = F,
smooth = smoother)
if (plot_points){
# plot points in red just to be on the safe side!
points3d(fdata, col = "blue")
}
if (plot_contour){
# plot the plane underneath
flat_matrix <- wide_form_values
if (verbose) print(flat_matrix)
y_intercept <- (zlim[2] - zlim[1]) * (-2/3) # put the flat matrix 1/2 the distance below the lower height
flat_matrix[which(flat_matrix != y_intercept)] <- y_intercept
if (verbose) print(flat_matrix)
rgl.surface(z = x_values, ## these are all different because
x = y_values, ## of the confusing way that
y = flat_matrix, ## rgl.surface works! - y is the height!
coords = c(2,3,1),
color = col,
alpha = 1.0,
smooth = smoother)
}
}
在add_rgl_model没有选项的情况下执行相同的工作，但将表面覆盖到现有的 3dplot 上。
add_rgl_model <- function(fdata){
## takes a model in long form, in the format
## 1st column x
## 2nd is y,
## 3rd is z (height)
## and draws an rgl model
##
# note that x has to be ascending, followed by y
print(head(fdata))
fdata <- fdata[order(fdata[, 1], fdata[, 2]), ]
print(head(fdata))
##
require(reshape2)
require(rgl)
orig_names <- colnames(fdata)
#print(head(fdata))
colnames(fdata) <- c("x", "y", "z")
fdata <- as.data.frame(fdata)
## work out the min and max of x,y,z
xlimits <- c(min(fdata$x, na.rm = T), max(fdata$x, na.rm = T))
ylimits <- c(min(fdata$y, na.rm = T), max(fdata$y, na.rm = T))
zlimits <- c(min(fdata$z, na.rm = T), max(fdata$z, na.rm = T))
l <- list (x = xlimits, y = ylimits, z = zlimits)
xyz <- do.call(expand.grid, l)
#print(xyz)
x_boundaries <- xyz$x
#print(class(xyz$x))
y_boundaries <- xyz$y
#print(class(xyz$y))
z_boundaries <- xyz$z
#print(class(xyz$z))
# now turn fdata into a wide format for use with the rgl.surface
fdata[, 2] <- as.character(fdata[, 2])
fdata[, 3] <- as.character(fdata[, 3])
#print(class(fdata[, 2]))
wide_form <- dcast(fdata, y ~ x, value_var = "z")
print(head(wide_form))
wide_form_values <- as.matrix(wide_form[, 2:ncol(wide_form)])
x_values <- as.numeric(colnames(wide_form[2:ncol(wide_form)]))
y_values <- as.numeric(wide_form[, 1])
print(x_values)
print(y_values)
wide_form_values <- wide_form_values[order(y_values), order(x_values)]
x_values <- x_values[order(x_values)]
y_values <- y_values[order(y_values)]
print(x_values)
print(y_values)
print(dim(wide_form_values))
print(length(x_values))
print(length(y_values))
rgl.surface(z = x_values, ## these are all different because
x = y_values, ## of the confusing way that
y = wide_form_values, ## rgl.surface works!
coords = c(2,3,1),
alpha = .8)
# plot points in red just to be on the safe side!
points3d(fdata, col = "red")
}
所以我的方法是，尝试用你所有的数据来做（我很容易绘制从~15k点生成的表面）。如果这不起作用，请取几个较小的样本并使用这些函数一次将它们全部绘制出来。
也许现在已经晚了，但是跟随 Spacedman，您是否尝试过 duplicate="strip" 或任何其他选项？
x=runif(1000)
y=runif(1000)
z=rnorm(1000)
s=interp(x,y,z,duplicate="strip")
surface3d(s$x,s$y,s$z,color="blue")
points3d(s)