R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型
原文链接:http://tecdat.cn/?p=20531
原文出处:拓端数据部落公众号
===
在标准线性模型中,我们假设
文章图片
。当线性假设无法满足时,可以考虑使用其他方法。
- 多项式回归
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/208dde079c874f9ea5ccca918cfcc845.png)
文章图片
同样,在标准线性模型方法(使用GLM的条件正态分布)中,参数
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/1cf82e349fab481982af9e703f18412c.png)
文章图片
可以使用最小二乘法获得,其中
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/9683d0883ab043378f690b94c068da14.png)
文章图片
在
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/ad1d13bf9b1f4127af276fd510a0af86.png)
文章图片
。
即使此多项式模型不是真正的多项式模型,也可能仍然是一个很好的近似值
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/ec099e3e22dc46679e55d3502ed87d8c.png)
文章图片
。实际上,根据 Stone-Weierstrass定理,如果
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/f943bee4af904ab09bb2ff4633b997b0.png)
文章图片
在某个区间上是连续的,则有一个统一的近似值
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/f9e13211972d447e8d1b7a89959e02d6.png)
文章图片
,通过多项式函数。
仅作说明,请考虑以下数据集
db = data.frame(x=xr,y=yr)
plot(db)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/4b83c7891f244ab7a6ee23f893e0b827.jpg)
文章图片
与标准回归线
reg = lm(y ~ x,data=https://www.it610.com/article/db)
abline(reg,col="red")
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/aa7e63d0f6594c6aa786f02d49112adf.jpg)
文章图片
考虑一些多项式回归。如果多项式函数的次数足够大,则可以获得任何一种模型,
reg=lm(y~poly(x,5),data=https://www.it610.com/article/db)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/02c0ffb65e2042b1b718d3d3d7caa8c8.jpg)
文章图片
但是,如果次数太大,那么会获得太多的“波动”,
reg=lm(y~poly(x,25),data=https://www.it610.com/article/db)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/8906831f4d3848e7a66841edd27437d1.jpg)
文章图片
并且估计值可能不可靠:如果我们更改一个点,则可能会发生(局部)更改
yrm=yr;
yrm\[31\]=yr\[31\]-2
lines(xr,predict(regm),col="red")
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/fac3e4e4707241e6b3fc3d3649dd8e77.jpg)
文章图片
- 局部回归
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/641c1122b5004eb3a49c3d28a1f3f4fa.png)
文章图片
,为什么不使用局部回归?
使用加权回归可以很容易地做到这一点,在最小二乘公式中,我们考虑
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/a6bd1a9fedc444fdb8a46d182c6e8127.png)
文章图片
- 在这里,我考虑了线性模型,但是可以考虑任何多项式模型。在这种情况下,优化问题是
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/f657edb097784792ab342ba283e50f4d.png)
文章图片
可以解决,因为
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/7a297202b84846b88bcfb44a4a3dd08a.png)
文章图片
例如,如果我们想在某个时候进行预测 , 考虑
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/6585220ecc3d4990b45e68ae63b2dadd.png)
文章图片
。使用此模型,我们可以删除太远的观测值,
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/d20dc52f606d4265a6f8bd9c8ac00ebe.jpg)
文章图片
更一般的想法是考虑一些核函数
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/b09f4edde3514d2fa8ea225e07238c31.png)
文章图片
给出权重函数,以及给出邻域长度的一些带宽(通常表示为h),
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/d41cd350e5204aba87e2e24047552911.png)
文章图片
这实际上就是所谓的 Nadaraya-Watson 函数估计器
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/dfdeae5e5fae454aa4f63049ea1455e0.png)
文章图片
。
在前面的案例中,我们考虑了统一核
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/5135828138f04b2081f532a1912d4b0d.png)
文章图片
,
但是使用这种权重函数具有很强的不连续性不是最好的选择,尝试高斯核,
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/e5af1abed3f94b7daf1c87f1f41a4cb5.png)
文章图片
这可以使用
w=dnorm((xr-x0))
reg=lm(y~1,data=https://www.it610.com/article/db,weights=w)
在我们的数据集上,我们可以绘制
w=dnorm((xr-x0))
plot(db,cex=abs(w)*4)
lines(ul,vl0,col="red")
axis(3)
axis(2)
reg=lm(y~1,data=https://www.it610.com/article/db,weights=w)
u=seq(0,10,by=.02)
v=predict(reg,newdata=data.frame(x=u))
lines(u,v,col="red",lwd=2)
在这里,我们需要在点2进行局部回归。下面的水平线是回归(点的大小与宽度成比例)。红色曲线是局部回归的演变
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/fe8ea63392de4168a6a6f9739e586d64.jpg)
文章图片
让我们使用动画来可视化曲线。
但是由于某些原因,我无法在Linux上轻松安装该软件包。我们可以使用循环来生成一些图形
name=paste("local-reg-",100+i,".png",sep="")
png(name,600,400)
for(i in 1:length(vx0)) graph (i)
然后,我使用
当然,可以考虑局部线性模型,
return(predict(reg,newdata=https://www.it610.com/article/data.frame(x=x0)))}
甚至是二次(局部)回归,
lm(y~poly(x,degree=2), weights=w)
当然,我们可以更改带宽
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/a3b274f5073c4f6cbd18ea41c1029e49.jpg)
文章图片
请注意,实际上,我们必须选择权重函数(所谓的核)。但是,有(简单)方法来选择“最佳”带宽h。交叉验证的想法是考虑
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/536351fc41424058992802cd8fb0ddee.png)
文章图片
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/f8b9312ba579459087bfab26256c4939.png)
文章图片
是使用局部回归获得的预测。
我们可以尝试一些真实的数据。
library(XML)
data = https://www.it610.com/article/readHTMLTable(html)
整理数据集,
plot(data$no,data$mu,ylim=c(6,10))
segments(data$no,data$mu-1.96*data$se,
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/82a9ac25019b415d9e3757fd0737b581.jpg)
文章图片
我们计算标准误差,反映不确定性。
for(s in 1:8){reg=lm(mu~no,data=https://www.it610.com/article/db,
lines((s predict(reg)/[1:12/]
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/99232fb6860c4776879f09b536682823.jpg)
文章图片
所有季节都应该被认为是完全独立的,这不是一个很好的假设。
smooth(db$no,db$mu,kernel = "normal",band=5)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/00ffacc5c7ab42ca8bb795e7997309b7.jpg)
文章图片
我们可以尝试查看带宽较大的曲线。
db$mu\[95\]=7plot(data$no,data$mulines(NW,col="red")
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/393507646ac5487582fc84de1e60b58a.jpg)
文章图片
样条平滑 接下来,讨论回归中的平滑方法。假设
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/32781c4f518d4cda91647cf275fa7a66.jpg)
文章图片
,
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/2a86ea0384134e4da8ad24538491bf56.jpg)
文章图片
是一些未知函数,但假定足够平滑。例如,假设
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/42b07b3c33ee4515a4dd4c167eeddf95.jpg)
文章图片
是连续的,
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/e1de8e0f71134ed085ba1201ebacfd48.jpg)
文章图片
存在,并且是连续的,
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/38a57a887ce4407c9f2ffffa9f661ce1.jpg)
文章图片
存在并且也是连续的等等。如果
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/42b07b3c33ee4515a4dd4c167eeddf95.jpg)
文章图片
足够平滑,可以使用泰勒展开式。 因此,对于
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/7b4a43fb75a34a118de714d7b9c1b432.jpg)
文章图片
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/48f9326756d647e2bac772133e250df8.jpg)
文章图片
也可以写成
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/8fd34e2b9e6e4f0f96c36c4bf606d7c1.jpg)
文章图片
第一部分只是一个多项式。
使用 黎曼积分,观察到
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/8e275c56c3e446d8b36f84d011f339bf.jpg)
文章图片
因此,
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/2a8d4618bbfb4dfb822ebec90c5426f5.jpg)
文章图片
我们有线性回归模型。一个自然的想法是考虑回归
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/4ff69759f5064538a3cddda024257a5e.gif)
文章图片
,对于
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/e9d923ee91004598b7b126920e322d02.gif)
文章图片
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/ea8cf245b2e84ee8bc03ade9480e45c4.jpg)
文章图片
给一些节点
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/e098f4c8454a4ea386662030dad9cd24.gif)
文章图片
。
plot(db)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/902e68e9d8534bffb42547b82b030099.jpg)
文章图片
如果我们考虑一个节点,并扩展阶数1,
B=bs(xr,knots=c(3),Boundary.knots=c(0,10),degre=1)
lines(xr\[xr<=3\],predict(reg)\[xr<=3\],col="red")
lines(xr\[xr>=3\],predict(reg)\[xr>=3\],col="blue")
可以将用该样条获得的预测与子集(虚线)上的回归进行比较。
lines(xr\[xr<=3\],predict(reg)\[xr<=3lm(yr~xr,subset=xr>=3)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/72231e4ef6434adf8614e2f263572109.jpg)
文章图片
这是不同的,因为这里我们有三个参数(关于两个子集的回归)。当要求连续模型时,失去了一个自由度。观察到可以等效地写
lm(yr~bs(xr,knots=c(3),Boundary.knots=c(0,10)
回归中出现的函数如下
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/7b541aa6d0194eb5ab05fc57ea4b35f2.jpg)
文章图片
现在,如果我们对这两个分量进行回归,我们得到
matplot(xr,Babline(v=c(0,2,5,10),lty=2)
【R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型】如果加一个节点,我们得到
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/8b1f92f1105f4331b47bfa9916e70c94.jpg)
文章图片
预测是
lines(xr,predict(reg),col="red")
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/18d337366ea24e019541854e81ece784.jpg)
文章图片
我们可以选择更多的节点
lines(xr,predict(reg),col="red")
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/7912b6abc1cc4fe68fe01ec07aa23307.jpg)
文章图片
我们可以得到一个置信区间
polygon(c(xr,rev(xr)),c(P\[,2\],rev(P\[,3\]))points(db)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/ada0846cffec4634801bdf76f3e697cc.jpg)
文章图片
如果我们保持先前选择的两个节点,但考虑泰勒的2阶的展开,我们得到
matplot(xr,B,type="l")
abline(v=c(0,2,5,10),lty=2)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/5f439cb78c0e45cf986788e5fdd47cd0.jpg)
文章图片
如果我们考虑常数和基于样条的第一部分,我们得到
B=cbind(1,B)
lines(xr,B\[,1:k\]%*%coefficients(reg)\[1:k\],col=k-1,lty=k-1)
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/24d3faea773a4967b9f9dec2e9d54bb4.jpg)
文章图片
如果我们将常数项,第一项和第二项相加,则我们得到的部分在第一个节点之前位于左侧,
k=3
lines(xr,B\[,1:k\]%*%coefficients(reg)\[1:k\]
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/c31ed2b14b8249eca464de69513c6d6d.jpg)
文章图片
通过基于样条的矩阵中的三个项,我们可以得到两个节点之间的部分,
lines(xr,B\[,1:k\]%*%coefficients(reg)\[1:k\]
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/aa67134a54e14a859ad6968d591dbca2.jpg)
文章图片
最后,当我们对它们求和时,这次是最后一个节点之后的右侧部分,
k=5
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/44aee2eafc7f49e4867a153822e3a3f4.jpg)
文章图片
这是我们使用带有两个(固定)节点的二次样条回归得到的结果。可以像以前一样获得置信区间
polygon(c(xr,rev(xr)),c(P\[,2\],rev(P\[,3\]))points(db)
lines(xr,P\[,1\],col="red")
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/8e80be2cabe1415baa08b90b0ead657c.jpg)
文章图片
使用函数
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/9a60aae17a0c4e8e8311efcbc1ead17b.gif)
文章图片
,可以确保点的连续性
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/3ef2c9f96dbc46adbc73ff50e16a11ba.gif)
文章图片
。
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/549c11e5040b4dfb9cc77cb5b1662b64.jpg)
文章图片
再一次,使用线性样条函数,可以增加连续性约束,
lm(mu~bs(no,knots=c(12*(1:7)+.5),Boundary.knots=c(0,97),lines(c(1:94,96),predict(reg),col="red")
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/0e9685c2607d4a8682672101688f6b5e.jpg)
文章图片
但是我们也可以考虑二次样条,
abline(v=12*(0:8)+.5,lty=2)lm(mu~bs(no,knots=c(12*(1:7)+.5),Boundary.knots=c(0,97),
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/cba38dfa631e43d7b8f8e2d690921df2.jpg)
文章图片
![R语言中的多项式回归、局部回归、核平滑和平滑样条回归模型](https://img.it610.com/image/info9/6dd90f58dbea4094831e454a2b3efa62.jpg)
文章图片
最受欢迎的见解
1.R语言多元Logistic逻辑回归 应用案例
2.面板平滑转移回归(PSTR)分析案例实现分析案例实现")
3.matlab中的偏最小二乘回归(PLSR)和主成分回归(PCR)
4.R语言泊松Poisson回归模型分析案例
5.R语言回归中的Hosmer-Lemeshow拟合优度检验
6.r语言中对LASSO回归,Ridge岭回归和Elastic Net模型实现
7.在R语言中实现Logistic逻辑回归
8.python用线性回归预测股票价格
9.[R语言如何在生存分析与Cox回归中计算IDI,NRI指标
](http://tecdat.cn/r%E8%AF%AD%E... "R语言如何在生存分析与Cox回归中计算IDI,NRI指标")
推荐阅读
- 热闹中的孤独
- JS中的各种宽高度定义及其应用
- 我眼中的佛系经纪人
- 《魔法科高中的劣等生》第26卷(Invasion篇)发售
- Android中的AES加密-下
- 【生信技能树】R语言练习题|【生信技能树】R语言练习题 - 中级
- 放下心中的偶像包袱吧
- 一起来学习C语言的字符串转换函数
- C语言字符函数中的isalnum()和iscntrl()你都知道吗
- C语言浮点函数中的modf和fmod详解