加权最小二乘(Weighted Least Squares,WLS)回归及R语言计算
与OLS回归相比,加权最小二乘(Weighted Least Squares,WLS)回归是当残差中方差不变的最小二乘假设被违背(异方差性)时可以考虑的一种方法,即different observations have different residuals。
> #部分示例数据
> head(dat_mod)time pace diameter
1 168.5000 40 5.28
2 179.4762 40 5.28
3 155.8077 50 5.28
4 166.5278 50 5.28
5 148.2500 50 5.28
6 143.5714 60 5.28
#做图查看数据分布
dat_mod$diameter=as.factor(dat_mod$diameter)
ggplot(data=dat_mod,aes(x=pace,y=time,color=diameter,shape=diameter))+geom_point()+scale_y_continuous(limits = c(60,190),n.breaks = 14)+scale_x_continuous(breaks = c(30,40,50,60,90,120,150,180,210,240,270,300,330,360,390))+scale_color_manual(values=c('#228B22','#FFFF00','#FF0000'))
可以看出数据比较符合幂函数模型,所以思路是对y和x都取对数,将非线性模型线性化,然后再用最小二乘法拟合。
#幂函数模型
y3=log(dat_mod$time)
mod3=lm(y3~log(dat_mod$pace))
summary(mod3)
> summary(mod3)Call:
lm(formula = y3 ~ log(dat_mod$pace))Residuals:Min 1Q Median 3Q Max
-0.186191 -0.024459 0.007359 0.031775 0.147288 Coefficients:Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.37989 0.05151 123.86 <2e-16 ***
log(dat_mod$pace) -0.36094 0.01032 -34.97 <2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1Residual standard error: 0.06502 on 63 degrees of freedom
Multiple R-squared: 0.951, Adjusted R-squared: 0.9502
F-statistic: 1223 on 1 and 63 DF, p-value: < 2.2e-16pace3=seq(20,400,length.out=1000)
time_pred3=(pace3^(-0.36094))*exp(6.37989)
#或者mod3$fitted.values
dat_pred3=data.frame(pace3,time_pred3)ggplot()+geom_point(data=dat_mod,aes(x=pace,y=time,color=diameter,shape=diameter),size=3)+scale_y_continuous(limits = c(60,200),n.breaks = 15)+scale_x_continuous(breaks = c(20,30,40,50,60,90,120,150,180,210,240,270,300,330,360,390,400))+scale_color_manual(values=c('#228B22','#1E90FF','#FF0000'))+theme(panel.grid=element_blank(),panel.background=element_rect(fill='white', color='black',size=0.5),axis.text.x=element_text(size=18),axis.title.x=element_text(size=22,face='bold'),axis.text.y=element_text(size=18),axis.title.y=element_text(size=22,face='bold'),legend.title = element_text(size=22,face='bold'),legend.text = element_text(size=18))+labs(x='pace(ug/min)',y='time(s)',color='diameter(cm)',shape='diameter(cm)')+geom_line(data=dat_pred3,aes(x=pace3,y=time_pred3),color='#000000',size=1)
可以看出流速慢的时候残差大,流速越快残差越小,符合WLS使用条件异方差性
OLS是minimize sum(residuals^2),而WLS是minimize sum(w*residuals^2),即将权数与残差平方相乘后再求和,所以要先定义权重。
#定义权重
wt=1/(mod3$residuals)^2#加权最小二乘(WLS)回归,通过参数 weight 指定加权的变量
fit.wls=lm(y3~log(dat_mod$pace),weights = wt)
summary(fit.wls)
###################################################Call:lm(formula = y3 ~ log(dat_mod$pace), weights = wt)Weighted Residuals:Min 1Q Median 3Q Max
-0.9949 -0.9613 1.0009 1.0294 1.4474 Coefficients:Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.372919 0.006742 945.3 <2e-16 ***log(dat_mod$pace) -0.359768 0.001201 -299.7 <2e-16 ***---Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1Residual standard error: 0.9959 on 63 degrees of freedom
Multiple R-squared: 0.9993, Adjusted R-squared: 0.9993
F-statistic: 8.98e+04 on 1 and 63 DF, p-value: < 2.2e-16###################################################pace4=seq(20,400,length.out=1000)
time_pred4=(pace4^(-0.359768))*exp(6.372919)
#或者predict()
dat_pred4=data.frame(pace4,time_pred4)ggplot()+geom_point(data=dat_mod,aes(x=pace,y=time,color=diameter,shape=diameter),size=3)+scale_y_continuous(limits = c(60,200),n.breaks = 15)+scale_x_continuous(breaks = c(20,30,40,50,60,90,120,150,180,210,240,270,300,330,360,390,400))+scale_color_manual(values=c('#228B22','#1E90FF','#FF0000'))+theme(panel.grid=element_blank(),panel.background=element_rect(fill='white', color='black',size=0.5),axis.text.x=element_text(size=18),axis.title.x=element_text(size=22,face='bold'),axis.text.y=element_text(size=18),axis.title.y=element_text(size=22,face='bold'),legend.title = element_text(size=22,face='bold'),legend.text = element_text(size=18))+labs(x='pace(ug/min)',y='time(s)',color='diameter(cm)',shape='diameter(cm)')+geom_line(data=dat_pred4,aes(x=pace4,y=time_pred4),color='#8B008B',size=1)+geom_line(data=dat_pred3,aes(x=pace3,y=time_pred3),color='#000000',size=1)
这是使用了WLS后得到的拟合曲线,虽然r方=0.9993较之前r方=0.9502增加了很多,但是参数和OLS得到的参数基本一样,得到的拟合曲线基本重合,为什么会这样呢?
我们先比较下两种算法得到的残差平方和:
> sum(mod3$residuals^2)
[1] 0.2663267
> sum(fit.wls$residuals^2)
[1] 0.2664723
可以看到由WLS计算得到的残差平方和要大于由OLS计算得到的残差平方和,如果按照皮尔森相关系数的计算公式1-(RSS/TSS)那么WLS计算得到的R方肯定会更小,所以我猜测WLS计算R方的公式应该是
1-(sum(w×residuals∧2) / sum(w×(y-mean(y))∧2)),下面来验证一下
> tss=(wt*(y3-mean(y3))^2)%>%sum()
> rss=(wt*fit.wls$residuals^2)%>%sum()
> 1-rss/tss
[1] 0.9999152
总结
在加权回归分析中, 相关系数不能作为评判权重选择是否最优的指标,因为计算得到的是加权相关系数而不是皮尔森相关系数,但是采用加权最小二乘法可以使拟合曲线更加靠近残差小的数据点。
加权最小二乘(Weighted Least Squares,WLS)回归及R语言计算相关推荐
- R语言计算回归模型残差平方和(Residual Sum of Squares)实战,并基于残差平方和比较模型优劣
R语言计算回归模型残差平方和(Residual Sum of Squares)实战,并基于残差平方和比较模型优劣 目录
- R语言计算线性回归的最小二乘估计
R语言计算线性回归的最小二乘估计 全称:线性回归的最小二乘法(OLS回归),ordinary least square,字面翻译:普通最小平方: 内容:包括三个部分:简单线性回归.多项式回归.多元线性 ...
- R语言计算回归模型每个样本(观察、observation、sample)的DFFITS度量实战:忽略单个观察(样本)时,回归模型所做的预测会发生多大的变化
R语言计算回归模型每个样本(观察.observation.sample)的DFFITS度量实战:忽略单个观察(样本)时,回归模型所做的预测会发生多大的变化 目录
- R语言计算回归模型每个样本(观察、observation、sample)的杠杆值(leverage)实战:如果一个样本的预测变量比其他样本的预测变量值更极端,那么被认为具有很高的杠杆作用
R语言计算回归模型每个样本(观察.observation.sample)的杠杆值(leverage)实战:如果一个样本的预测变量比其他样本的预测变量值更极端,那么被认为具有很高的杠杆作用 目录
- R语言计算F1评估指标实战:F1 score、使用R中caret包中的confusionMatrix()函数为给定的logistic回归模型计算F1得分(和其他指标)
R语言计算F1评估指标实战:F1 score.使用R中caret包中的confusionMatrix()函数为给定的logistic回归模型计算F1得分(和其他指标) 目录
- R语言计算回归模型标准化残差实战(Standardized Residuals):识别回归模型中离群点
R语言计算回归模型标准化残差实战(Standardized Residuals):识别回归模型中离群点 目录
- R语言计算回归模型的SST、SSR以及SSE指标实战
R语言计算回归模型的SST.SSR以及SSE指标实战 目录 R语言计算回归模型的SST.SSR以及SSE指标实战 #仿真数据
- R语言计算回归模型学生化残差(Studentized Residuals)实战:如果样本学生化残差(Studentized Residuals)绝对值大于3则是离群值
R语言计算回归模型学生化残差(Studentized Residuals)实战:如果样本学生化残差(Studentized Residuals)绝对值大于3则是离群值 目录
- R语言计算回归模型每个样本(观察、observation、sample)的DFBETAS值实战:每一个样本对给定系数的估计有多大的影响
R语言计算回归模型每个样本(观察.observation.sample)的DFBETAS值实战:每一个样本对给定系数的估计有多大的影响 目录
最新文章
- Python面向对象编程:深度认识类class
- 字符串与数组的常用方法
- linux服务器配置https访问
- 【Python入门】列表的常用操作,这十张图把它说的明明白白!
- redis入门demo
- SessionListener失败,退出
- linux chrome 安装过程记录
- jvm(2)-OutOfMemoryError 异常(内存溢出异常)
- 推荐两款 GTD 工具
- android java标准时间_Android 时间 日期 相关
- OpenSSL命令---pkcs7
- 关于第五届全国高校新一代信息技术暑假教师培训班的通知
- 最新emoji表情代码大全_2019七夕节最新撩妹句子大全,浪漫的七夕节表情包集锦...
- 【Android】实现应用简单的用户登录界面
- if...else 语句双分支结构 计算分段函数
- Python 批量发送邮件脚本
- UEFI服务器PXE网络安装CentOS7.5
- 试图加载格式不正确的程序 0x8007000b
- java写excel
- Android开发之获取当前展示的activity的包名,类名