楼主: ywh19860616
64608 65

带置信区间的画图问题 [推广有奖]

学术权威

32%

还不是VIP/贵宾

-

0%
威望
0 级
论坛币
870 个
通用积分
3477.1298
学术水平
854 点
热心指数
991 点
信用等级
637 点
经验
116346 点
帖子
3976
精华
0
在线时间
7729 小时
注册时间
2009-9-3
最后登录
2023-12-9
楼主
ywh19860616 发表于 2011-1-10 08:54:22 |只看作者 |坛友微信交流群|倒序 |AI写论文
相似文件 换一批

+2 论坛币
k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群
赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币
我想使用R画出带置信区间的图像，类似于了分位数回归图，见图例。
    文字描述：
    假如我有两个向量X和Y，每个X对应一个Y
   X=[0.1, 0.2, 0.3, 0.4, 0.6, 0.7];
   Y=[0.2,0.3,0.4,0.5,0.6,0.7];
   如果只画X和Y的对应图，那很简单，直接plot(X,Y)就能解决。
  但是现在问题是Y中每个数值都对应一个置信区间，
如[0.15,0.22],[0.25,0.32],[0.35,0.42],[0.45,0.52],[0.55,0.62],[0.65,0.7]
  怎么画如下图形？请大家指点，谢谢
二维码

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0
关键词：画图问题 置信区间 分位数回归 plot 分位数 画图 置信区间

相关帖子

未命名.jpg (15.64 KB)

未命名.jpg

使用道具

沙发
ltx5151 发表于 2011-1-10 10:03:42 |只看作者 |坛友微信交流群
三个选择
1. lines
2. abline
如果真的想要途中的那种填充效果，用contour。
这些函数help里面都介绍的比较详细，lz可以自己查阅一下。
已有 1 人评分学术水平 热心指数 收起 理由
ywh19860616 + 1 + 1 谢谢您的建议

总评分: 学术水平 + 1  热心指数 + 1   查看全部评分

使用道具

藤椅
shenbaiseshatan 在职认证  发表于 2011-1-10 15:02:30 |只看作者 |坛友微信交流群
简单点可以实现类似功能的程序
  1. n<-100
  2. set.seed(123);x<-rnorm(n,5,1);y<-x+rnorm(n,3,.5)
  3. lm.xy<-lm(y~x)
  4. lm.for.plot<-as.data.frame(predict.lm(lm.xy,interval = "confidence"))
  5. attach(lm.for.plot)
  6. index<-order(fit)
  7. plot(sort(fit),pch=20,col="blue")
  8. points(1:n,lwr[index],pch=20,col="red",type='l')
  9. points(1:n,upr[index],pch=20,col="red",type='l')
  10. segments(1:n,lwr[index],1:n,upr[index],col="brown")
复制代码

  1. n<-100
  2. set.seed(123);x<-rnorm(n,5,1);y<-x+rnorm(n,3,.5)
  3. lm.xy<-lm(y~x)
  4. lm.for.plot<-as.data.frame(predict.lm(lm.xy,interval = "confidence"))
  5. attach(lm.for.plot)
  6. index<-order(fit)
  7. index.fill<-seq(1,n,by=2)
  8. points.fill<-c(upr[index],rev(lwr[index]))
  9. plot(c(1:n,n:1),points.fill,type='n')
  10. polygon(c(1:n,n:1),points.fill,border=NA,col="gray")
  11. points(sort(fit),pch=20,col="blue",add=T)
  12. points(1:n,lwr[index],pch=20,col="red",type='l')
  13. points(1:n,upr[index],pch=20,col="red",type='l')
  14. #segments(1:n,lwr[index],1:n,upr[index],col="brown")
复制代码
已有 1 人评分学术水平 热心指数 信用等级 收起 理由
ywh19860616 + 1 + 1 + 1 谢谢您

总评分: 学术水平 + 1  热心指数 + 1  信用等级 + 1   查看全部评分

胜人者有力，自胜者强！

使用道具

板凳
楚韵荆风 学生认证  发表于 2011-1-11 13:35:49 |只看作者 |坛友微信交流群
把置信带下面的数据和上面的数据分别用两个数组y1,y2存储，再用lines命令就可以了，
plot(x,y,type="l")
lines(x,y1)
lines(x,y2)
已有 1 人评分学术水平 热心指数 收起 理由
ywh19860616 + 1 + 1 谢谢

总评分: 学术水平 + 1  热心指数 + 1   查看全部评分

共享是一种彼此的快乐

使用道具

报纸
trier2006 发表于 2011-1-13 08:58:47 |只看作者 |坛友微信交流群
嗯，大家的建议本质好像都是计算出区间值，然后画上去。似乎没有现成的函数用。
最好的医生是自己，最好的药物是时间……

使用道具

地板
shenbaiseshatan 在职认证  发表于 2011-1-14 14:12:08 |只看作者 |坛友微信交流群
5# trier2006
现成的函数有类似的，比如Hmisc包的Dotplot
胜人者有力，自胜者强！

使用道具

7楼
epoh 发表于 2011-3-26 11:18:23 |只看作者 |坛友微信交流群
function rq.fit.panel()
is now superseded by the use of the function rqss() in the R quantreg package,
or the use of method = "lasso" for the function rq() in that package.

Package‘quantreg’Version 4.57  Date 2011-2-26  

#############
library(quantreg)
data(engel)
attach(engel)
summary(rq(foodexp~income,tau = 1:49/50,data=engel))

Call: rq(formula = foodexp ~ income, tau = 1:49/50, data = engel)
tau: [1] 0.02
Coefficients:
            coefficients lower bd  upper bd
(Intercept) 108.60499    103.39624 141.26345
income        0.34664      0.28526   0.34664

Call: rq(formula = foodexp ~ income, tau = 1:49/50, data = engel)
tau: [1] 0.04
............

#Example of plotting of coefficients and their confidence bands
plot(summary(rq(foodexp~income,tau = 1:49/50,data=engel)))
已有 1 人评分学术水平 热心指数 信用等级 收起 理由
ywh19860616 + 1 + 1 + 1 谢谢哈，不过我看过数据集engel，并不是面板结构的，而只是截面的

总评分: 学术水平 + 1  热心指数 + 1  信用等级 + 1   查看全部评分

使用道具

8楼
ywh19860616 发表于 2011-3-26 13:02:39 |只看作者 |坛友微信交流群
##epoh ,哈哈，谢谢啊
我看过rqss（）函数帮助，但是找不到识别面板结构的程序，还不知道怎么使用。
一份耕耘，一份收获。

使用道具

9楼
epoh 发表于 2011-3-28 19:58:36 |只看作者 |坛友微信交流群
panel data:
  package "Ecdat"
  Produc Us States Production
  a panel of 48 observations from 1970 to 1986
  number of observations : 816

##############
library(Ecdat)
data("Produc")  
#n=48, T=17, N=816
x1=log(Produc$pcap)
x2=log(Produc$pc)
x3=log(Produc$emp)
x4=Produc$unemp
X=cbind(x1,x2,x3,x4)
y=log(Produc$gsp)
s= rep(1:48,rep(17,48))
w=1

####case 1
####use function rq.fit.panel()
require(SparseM)
require(quantreg)

X <- as.matrix(X)
p <- ncol(X)
n <- length(levels(as.factor(s)))
N <- length(y)
if(N != length(s) || N != nrow(X))
    stop("dimensions of y,X,s must match")
Z <- as.matrix.csr(model.matrix(~as.factor(s)-1))
xx<- cbind(as(w,"matrix.diag.csr") %x% X,w %x% Z)
rq.fit.sfn(xx,y,tau=0.5)
>$coef
[1] -0.001856897  0.227955745  0.806905899 -0.003254027  2.348518722
[6]  2.531511125  2.416402785  2.646310826  2.554288408  2.598468250
[11]  2.540988532  2.473361510  2.401692284  2.497764134  2.537963266
[16]  2.425204261  2.482897448  2.510841856  2.553142722  2.703008416
[21]  2.404225489  2.537962766  2.486090983  2.567535469  2.455404207
[26]  2.403301315  2.463135907  2.545757570  2.474922409  2.515871136
[31]  2.443668078  2.580372182  2.643673666  2.605560174  2.384878447
[36]  2.541337214  2.468143678  2.587964518  2.509402455  2.432777812
[41]  2.495204924  2.258353874  2.465798076  2.372996183  2.600392620
[46]  2.471721790  2.457974490  2.519673773  2.590197854  2.472063213
[51]  2.470425047  2.862672642


####case 2
####use function rqss()
tt=as.matrix(xx)
fit=rqss(y~tt-1)
summary(fit)
>
Formula:
y ~ tt - 1

Parametric coefficients:
      Estimate    Std. Error   t value    Pr(>|t|)   
tt1  -0.001857   0.040474  -0.046   0.9634   
tt2   0.227956   0.033949   6.715 3.68e-11 ***
tt3   0.806906   0.041542  19.424  < 2e-16 ***
tt4  -0.003254   0.001402  -2.321   0.0205 *  
tt5   2.348519   0.227838  10.308  < 2e-16 ***
tt6   2.531511   0.226443  11.179  < 2e-16 ***
tt7   2.416403   0.216806  11.145  < 2e-16 ***
tt8   2.646311   0.254619  10.393  < 2e-16 ***
tt9   2.554288   0.222018  11.505  < 2e-16 ***
tt10  2.598468   0.220106  11.806  < 2e-16 ***
tt11  2.540989   0.210428  12.075  < 2e-16 ***
tt12  2.473362   0.234957  10.527  < 2e-16 ***
tt13  2.401692   0.226512  10.603  < 2e-16 ***
tt14  2.497764   0.208589  11.975  < 2e-16 ***
tt15  2.537963   0.241893  10.492  < 2e-16 ***
tt16  2.425204   0.227976  10.638  < 2e-16 ***
tt17  2.482897   0.232136  10.696  < 2e-16 ***
tt18  2.510842   0.231957  10.825  < 2e-16 ***
tt19  2.553143   0.231057  11.050  < 2e-16 ***
tt20  2.703008   0.244407  11.059  < 2e-16 ***
tt21  2.404225   0.201550  11.929  < 2e-16 ***
tt22  2.537963   0.228998  11.083  < 2e-16 ***
tt23  2.486091   0.220021  11.299  < 2e-16 ***
tt24  2.567535   0.239121  10.737  < 2e-16 ***
tt25  2.455404   0.235207  10.439  < 2e-16 ***
tt26  2.403301   0.223782  10.739  < 2e-16 ***
tt27  2.463136   0.226028  10.897  < 2e-16 ***
tt28  2.545758   0.227966  11.167  < 2e-16 ***
tt29  2.474922   0.232918  10.626  < 2e-16 ***
tt30  2.515871   0.211395  11.901  < 2e-16 ***
tt31  2.443668   0.201243  12.143  < 2e-16 ***
tt32  2.580372   0.224698  11.484  < 2e-16 ***
tt33  2.643674   0.222489  11.882  < 2e-16 ***
tt34  2.605560   0.253096  10.295  < 2e-16 ***
tt35  2.384878   0.221027  10.790  < 2e-16 ***
tt36  2.541337   0.225034  11.293  < 2e-16 ***
tt37  2.468144   0.238798  10.336  < 2e-16 ***
tt38  2.587965   0.227848  11.358  < 2e-16 ***
tt39  2.509402   0.221705  11.319  < 2e-16 ***
tt40  2.432778   0.239579  10.154  < 2e-16 ***
tt41  2.495205   0.191899  13.003  < 2e-16 ***
tt42  2.258354   0.211827  10.661  < 2e-16 ***
tt43  2.465798   0.222536  11.080  < 2e-16 ***
tt44  2.372996   0.232450  10.209  < 2e-16 ***
tt45  2.600393   0.251702  10.331  < 2e-16 ***
tt46  2.471722   0.214506  11.523  < 2e-16 ***
tt47  2.457974   0.199837  12.300  < 2e-16 ***
tt48  2.519674   0.226890  11.105  < 2e-16 ***
tt49  2.590198   0.241976  10.704  < 2e-16 ***
tt50  2.472063   0.224931  10.990  < 2e-16 ***
tt51  2.470425   0.229340  10.772  < 2e-16 ***
tt52  2.862673   0.230734  12.407  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


  Quantile Fidelity at tau = 0.5  is      10.8697
  Effective Degrees of Freedom = 52        Sample Size = 816
已有 1 人评分学术水平 热心指数 信用等级 收起 理由
ywh19860616 + 1 + 1 + 1 谢谢，我以为rqss函数没有识别面板数据结构的哦，看结果，两个基本一样

总评分: 学术水平 + 1  热心指数 + 1  信用等级 + 1   查看全部评分

使用道具

10楼
king64 发表于 2011-3-28 21:36:19 |只看作者 |坛友微信交流群
#（1）先回归
   lmfit <- lm( dat$Sunday ~ dat$Daily, data=dat )
# （2）  计算置信区间（95% ）---------------------------
     Ypred <- predict( lmfit, interval = "confidence", newdata= dat, level = 0.95 )
#（3） 画图
        plot(dat$Daily, dat$Sunday,
              xlab="Daily", ylab="Sunday",
             main="Figure:  Yhat Predict confidence")
        lines ( Ypred[ ,1 ] ~ dat$Daily, col="red", lwd=1 )
        lines ( Ypred[ ,2 ] ~ dat$Daily, col="gray", lwd=1 )
        lines ( Ypred[ ,3 ] ~ dat$Daily, col="gray", lwd=1 )
已有 1 人评分学术水平 热心指数 信用等级 收起 理由
ywh19860616 + 1 + 1 + 1 谢谢您

总评分: 学术水平 + 1  热心指数 + 1  信用等级 + 1   查看全部评分

使用道具

返回列表
1234567下一页
高级模式
B Color Image Link Quote Code Smilies
您需要登录后才可以回帖 登录 | 我要注册

本版微信群
加好友,备注cda
拉您进交流群
  • 手机版  |
  • 意见反馈  |
  • 帮助  |
  • 新手入门  |
  • 用户手册  |
  • 友情链接  |

    • 如有投资本站、合作意向或投放广告，请联系：13661292478（刘老师）

    联系客服

    邮箱：service@pinggu.org 投诉或不良信息处理：（010-68466864）

    010-68466864

    京ICP备16021002-2号 京B2-20170662号 京公网安备 11010802022788号 论坛法律顾问：王进律师 知识产权保护声明   免责及隐私声明

    GMT+8, 2024-4-25 06:20