在R中绘制3D曲面

问题描述

我有一个lm(log(z) ~ x*y)形式的lm.

我可以使用s3d <- scatterplot3d(x,y,log(z))制作数据的3D散点图.

如何绘制回归函数的表面?

我知道我可以为简单的线性模型做s3d$plane3d(lm(log(z) ~ x + y)，但这不适用于非平面表面.例如s3d$plane3d(lm(log(z) ~ x * y)无法正常工作.

解决方案

因此，首先，如果您提供示例数据集(或实际数据，可能会更好)和以下示例的可重现性，您很有可能会获得帮助.您的代码.否则，您会强迫我们这样做-我认为这就是为什么您的问题在6个小时内被忽略的原因.有关更多详细信息，请参见此链接.

第二，像这样使用lm(...)和formula=log(z)~x*y是一个非常糟糕的主意.线性建模基于以下假设:响应中的误差(在您的情况下为log(z))通常以恒定方差分布.如果您的z-data具有具有恒定方差的正态分布误差，那么log(z)肯定不会.这是一个经典的错误；正确的方法是使用带有family=poisson的广义线性建模(请参见glm程序包等).

最后，您的问题.下面的代码创建了一个3D散点图，该散点图覆盖在响应面上.它使用rgl程序包，该程序包生成可旋转的3D图.

在这里，我已经为这些点着色，以使表面以下的点为红色，而上方的部分为绿色，并在每个点到表面之间添加了垂线.

棘手的是，在surface3d(...)中，x和y参数是对应于网格的向量，而z是矩阵，每个x值一行，每个y值一行.

对于实际数据，您可能需要调整open3d(...)中的scale=...参数.

# create sample dataset - you have this already,,,
set.seed(1)   # for reproducible example
df <- data.frame(x=sample(1:50,50)/50,y=sample(1:50,50)/50)
df$z <- with(df,exp(4*x + 2*y - 6*x*y + 6)+rnorm(50,sd=500))
fit <- lm(log(z) ~ x*y,df)   # bad, bad, bad - don't do this!!!

# you start here...
df$pred <- predict(fit)
# set up matrix of z-values
x <- seq(min(df$x),max(df$x),len=100)
y <- seq(min(df$y),max(df$y),len=100)
plot.df <- expand.grid(x=x,y=y)
plot.df$z <- predict(fit,newdata=plot.df)
library(reshape2)
z <- dcast(plot.df,x~y,value.var="z")[-1]

# plot the points, the fitted surface, and droplines
library(rgl)
colors <- 2.5+0.5*sign(residuals(fit))
open3d(scale=c(1,1,0.2))
points3d(df$x,df$y,log(df$z),col=colors)
surface3d(x,y,as.matrix(z),col="blue",alpha=.2)
apply(df,1,function(row)lines3d(rep(row[1],2),rep(row[2],2),c(log(row[3]),row[4]),col=colors))
axes3d()
title3d(xlab="X",ylab="Y",zlab="log(Z)")

I have a lm of the form lm(log(z) ~ x*y).

I can make a 3D scatterplot of the data using s3d <- scatterplot3d(x,y,log(z)).

How do I plot the surface of my regression function?

I know I can do s3d$plane3d(lm(log(z) ~ x + y) for a simple linear model, but this doesn't work for nonplanar surfaces; e.g. s3d$plane3d(lm(log(z) ~ x * y) won't work.

解决方案

So first, you are much more likely to get help if you provide a sample dataset (or your actual data, which would be better) and a reproducible example of your code. Otherwise you are forcing us to do it - which I think is why your question was ignored for 6 hours. See this link for more details.

Second, using lm(...) like that, with formula=log(z)~x*y is an exceptionally bad idea. Linear modeling is based on the assumption that the errors in the response (log(z) in your case) are normally distributed with constant variance. If your z-data has normally distributed error with constant variance, then log(z) certainly will not. This is a classic mistake; the right way to do this is using generalized linear modeling (see the glm package, among others), with family=poisson.

Finally, to your question. The code below creates a 3D scatter plot overlaid onto a response surface. It uses the rgl package, which generates rotatable 3D plots.

Here I've colored the points so that those below the surface are red and those above are green, and added dropline from each point to the surface.

The tricky bit is that in surface3d(...), the x and y arguments are vectors corresponding to a grid, and z is a matrix with one row for every x value and one column for every y value.

For your real data you may need to tweak the scale=... argument in open3d(...).

# create sample dataset - you have this already,,,
set.seed(1)   # for reproducible example
df <- data.frame(x=sample(1:50,50)/50,y=sample(1:50,50)/50)
df$z <- with(df,exp(4*x + 2*y - 6*x*y + 6)+rnorm(50,sd=500))
fit <- lm(log(z) ~ x*y,df)   # bad, bad, bad - don't do this!!!

# you start here...
df$pred <- predict(fit)
# set up matrix of z-values
x <- seq(min(df$x),max(df$x),len=100)
y <- seq(min(df$y),max(df$y),len=100)
plot.df <- expand.grid(x=x,y=y)
plot.df$z <- predict(fit,newdata=plot.df)
library(reshape2)
z <- dcast(plot.df,x~y,value.var="z")[-1]

# plot the points, the fitted surface, and droplines
library(rgl)
colors <- 2.5+0.5*sign(residuals(fit))
open3d(scale=c(1,1,0.2))
points3d(df$x,df$y,log(df$z),col=colors)
surface3d(x,y,as.matrix(z),col="blue",alpha=.2)
apply(df,1,function(row)lines3d(rep(row[1],2),rep(row[2],2),c(log(row[3]),row[4]),col=colors))
axes3d()
title3d(xlab="X",ylab="Y",zlab="log(Z)")

这篇关于在R中绘制3D曲面的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持！