01 曲线拟合法

曲线拟合法只适用于单一变量与目标函数之间的关系分析，而曲面拟合则多用于二维变量与目标函数之间关系的分析。所谓曲面拟合，就是根据实验测试数据，求取函数f(x,y)与变量x及y之间的解析式，使其通过或近似通过所有的实验测试点。也就是说，使所有实验数据点能近似地分布在函数f(xy)所表示的空间曲面上。
曲面拟合通常采用两种方式，即插值方式和逼近方式来实现。两者的共同点是均利用曲面上或接近曲面的一组离散点寻求良好的曲面方程。两者主要的区别是:插值方式得到的方程所表示的曲面全部通过这组数据点，比如LOESS曲面拟合;而逼近方式，只要求在某种准则下其方程表示的曲面与这组数据点接近即可，比如多项式曲面拟合。逼近方式一般使用最小二乘法实现。最小二乘法是一种逼近理论，也是采样数据进行拟合时最常用的一种方法。曲面一般不通过已知数据点，而是根据拟合的曲面在取样处的数值与实际值之差的平均和达到最小时求得，它的主旨思想就是使拟合数值与实际数值之间的偏平方和达到最小。
在R中，rgl 包的surface3d0函数、scatterplot3d 包的plane3d 0函数、lttie 包的wireframe()函数，或者plot3D包的persp3DQ函数都可以绘制曲面拟合。

02 多项式拟合


library(plot3D)
library(gridExtra)
library(reshape2)
library(RColorBrewer)

colormap<-colorRampPalette(rev(brewer.pal(11,'RdYlGn')))(100)
mydata <- read.csv("Surface_Data.csv", sep= ",", header=T)

#多项式拟合z=f(x,y)=a+bx+cy+dxx+eyy
x <- mydata$x
y <- mydata$y
z <- mydata$variable
x2<-x*x
y2<-y*y
poly_z <- lm(z ~ x + y +x2+y2)
print(poly_z)

#设定为30X30的网格数据(x, y)，并根据拟合方程求其数值
N<-30
xmar <- seq(min(x),max(x),(max(x)-min(x))/N)
ymar <- seq(min(y),max(y),(max(y)-min(y))/N)
Grid_xy<-expand.grid(list(x=xmar,y=ymar))
Grid_xy$x2<-Grid_xy$x*Grid_xy$x
Grid_xy$y2<-Grid_xy$y*Grid_xy$y
Grid_z <- predict.lm(poly_z, newdata=Grid_xy)  


pred_z<-matrix(Grid_z, length(xmar), length(ymar))


persp3D (xmar, ymar, pred_z,
      theta = 150, phi = 40, d=3, 
      col = colormap, 
      scale = TRUE, border = "black", 
       bty = "f",box = TRUE,ticktype = "detailed",#nticks=5,
      ylab = "0-60 mph (sec)", 
      xlab = "Gax Mileage (mpg)",
      zlab="Power (KW)",
      clab="Power (KW)",
      zlim=c(20,180),
      add=TRUE,
      colkey = list(length = 0.5, width = 1))

效果如下：

03 loess回归式拟合

mydata <- read.csv("Surface_Data.csv", sep= ",", header=T)

x <- mydata$x
y <- mydata$y
z <- mydata$z


elev.loess <- loess(z ~ x * y, mydata,span=0.95)
print(elev.loess)


xmar <- seq(min(x),max(x),(max(x)-min(x))/30)
ymar <- seq(min(y),max(y),(max(y)-min(y))/30)
# get fitted (interpolated) values
elev.interp <- predict(elev.loess, expand.grid(list(x=xmar,y=ymar)))

pred_z<-matrix(elev.interp, length(xmar),length(ymar))

# 显示曲面网格，网格边线颜色为洋红，显示box框线
persp3D (xmar, ymar, pred_z,
         theta = 150, phi = 40, d=3, 
         col = colormap, 
         scale = TRUE, border = "black", 
         bty = "f",box = TRUE,ticktype = "detailed",#nticks=5,
         xlab = "0-60 mph (sec)", 
         ylab = "Gax Mileage (mpg)",
         zlab="Power (KW)",
         zlim=c(20,180))