我已按照以下模型参数计算/绘制了线性和95%CI的图
lm <- lm(cars$speed~cars$dist)
conf <- predict(lm, interval='confidence')
conf <- cbind(cars,conf)
CI <- as.data.frame(confint(lm))
library(ggplot2)
plot<-ggplot(conf,aes(dist,speed)) +
geom_line(aes(y=fit),color='black') +
geom_line(aes(y=lwr),color='red',linetype='dashed') +
geom_line(aes(y=upr),color='red',linetype='dashed')
plot
我想知道计算图上上下限(红线)的方程是什么?我以为可以使用confint()
函数中的值来计算这些值?我试图像这样计算lwr和upr值,但没有得到相同的结果。
lower <- CI[1,1] + CI[2,1]*cars$dist
upper <- CI[1,2] + CI[2,2]*cars$dist
这是使用以下公式在lm.predict
中计算置信区间的方式:
可以如下实现:
my.lm <- lm(cars$speed~cars$dist)
intercept <- model.matrix(delete.response(terms(my.lm)), cars)
fit.values <- c(Xp %*% coef(my.lm))
data.fit <- data.frame(dist=cars$dist, fit=fit.values)
# compute t-value
tval <- qt((1-0.95)/2, df=nrow(data.fit)-2)
# compute Sxx
Sxx <- sum((data.fit$dist - mean(data.fit$dist))^2)
# compute MSres
MSres <- sum(lm$residuals^2)/(nrow(data.fit)-2)
# calculate confidence interval
CI <- data.frame(t(apply(data.fit, 1, FUN = function(row){
sqrt(MSres * (1/nrow(data.fit) + (as.numeric(row[1]) - mean(data.fit$dist))^2/Sxx)) * tval * c(1, -1) + as.numeric(row[2])
})))
names(CI) <- c("lwr","upr")
head(CI)
# lwr upr
#1 6.917090 10.31299
#2 8.472965 11.40620
#3 7.307526 10.58483
#4 10.764584 13.08820
#5 9.626909 12.23906
#6 8.472965 11.40620
您可以将结果与从predict
获得的结果进行比较。
希望有帮助。