问题描述
我正在使用nlme(版本3.1-149)的gls函数使用Cowtan-Way数据集绘制自1970年以来的全球表面温度。我想为趋势线添加95%的置信区间,但遇到了障碍。这是到目前为止我正在运行的代码:
cowtan <- read.table("https://www-users.york.ac.uk/~kdc3/papers/coverage2013/had4sst4_krig_v2_0_0.txt",header=F)
cowtan.sub<-subset(cowtan,V1>=1970)
x<-cowtan.sub$V1
y<-cowtan.sub$V2
library(nlme)
library(forecast)
fit<-lm(y~x)
auto.arima(resid(fit),ic=c("bic"))
cowtan.gls<-gls(y~x,correlation=corARMA(p=2,q=0))
summary(cowtan.gls)
plot(y~x,type="l",xlab="Year",ylab="Temperature anomalies (ºC)",main="Global surface air temperatures since 1970")
curve(cowtan.gls$coefficients[1]+cowtan.gls$coefficients[2]*x,add=TRUE,col="red",lwd=2)
到目前为止,太好了。障碍在于从gls模型中提取趋势线的95%置信区间。我尝试使用
newx<-seq(min(x),max(x),by=0.05)
conf_interval <- predict(cowtan.gls,newdata=data.frame(x=newx),interval="confidence",level=0.95)
lines(newx,conf_interval[,2],lty=2)
,但是收到错误消息“ conf_interval [,2]中的错误:维数不正确”。在检查conf_interval对象时,它具有合适的大小,但缺少我要遵循的上下置信区间。
我已经使用Google和Stacks Overflow搜索了解决方案,但到目前为止却空手而归。使用nlme的interval函数,我获得了每个参数的95%置信区间,但无法将其转换为趋势线的置信区间。
谢谢您的帮助。
解决方法
暂无找到可以解决该程序问题的有效方法,小编努力寻找整理中!
如果你已经找到好的解决方法,欢迎将解决方案带上本链接一起发送给小编。
小编邮箱:dio#foxmail.com (将#修改为@)