五年前的帖子里写过
1. hist3d https://yifanyang.wordpress.com/2011/09/18/my-r-hist3d/
2. binplot.3d(你要的应该是这个) modified https://yifanyang.wordpress.com/2011/05/
时过境迁,帖子里面的外部链接可能失效...代码中引号请自行更改
library(rgl)
binplot.3d<-function(x,y,z,alpha=1,topcol=“#ff0000”,sidecol=“#aaaaaa”){
save <- par3d(skipRedraw=TRUE)
on.exit(par3d(save))
x1<-c(rep(c(x[1],x[2],x[2],x[1]),3),rep(x[1],4),rep(x[2],4))
z1<-c(rep(0,4),rep(c(0,0,z,z),4))
y1<-c(y[1],y[1],y[2],y[2],rep(y[1],4),rep(y[2],4),rep(c(y[1],y[2],y[2],y[1]),2))
x2<-c(rep(c(x[1],x[1],x[2],x[2]),2),rep(c(x[1],x[2],rep(x[1],3),rep(x[2],3)),2))
z2<-c(rep(c(0,z),4),rep(0,8),rep(z,8)
)
y2<-c(rep(y[1],4),rep(y[2],4),rep(c(rep(y[1],3),rep(y[2],3),y[1],y[2]),2)
)
rgl.quads(x1,z1,y1,col=rep(sidecol,each=4),alpha=alpha)
rgl.quads(c(x[1],x[2],x[2],x[1]),rep(z,4),c(y[1],y[1],y[2],y[2]),
col=rep(topcol,each=4),alpha=1)
rgl.lines(x2,z2,y2,col=“#000000”)
}
hist3d<-function(x,y=NULL,nclass=“auto”,alpha=1,col=“#ff0000”,scale=10){
save <- par3d(skipRedraw=
TRUE
)
on.exit(par3d(save))
xy <- xy.coords(x,y)
x <- xy$x
y <- xy$y
n<-length(x)
if(nclass == “auto”){
nclass<-ceiling(sqrt(nclass.Sturges(x)))
}
breaks.x <- seq(min(x),max(x),length=(nclass+1))
breaks.y <- seq(min(y),max(y),length=(nclass+1))
z<-matrix(0,(nclass),(nclass))
for(i in 1:nclass){
for(j in 1:nclass){
z[i, j]<-(1/n)*sum(x < breaks.x[i+1]& y < breaks.y[j+1]& x >= breaks.x[i]& y >= breaks.y[j])
binplot.3d(c(breaks.x[i],breaks.x[i+1]),c(breaks.y[j],breaks.y[j+1]),
scale*z[i,j],alpha=alpha,topcol=col)
}
}
}
################################################################################
rgl.open()
rgl.bg(color=“black”)
rgl.light()
mat<-read.csv('demo.csv',header=T)
hist3d(mat$x,mat$y,alpha=0.4,nclass=10,scale=30000)