DiegoIE
DiegoIE

Reputation: 107

R plot confidence interval lines with a robust linear regression model (rlm)

I need to plot a Scatterplot with the confidence interval for a robust linear regression (rlm) model, all the examples I had found only work with LM.

This is my code:

model1 <- rlm(weightsE$brain ~ weightsE$body)
newx <- seq(min(weightsE$body), max(weightsE$body), length.out=70)
newx<-as.data.frame(newx)
colnames(newx)<-"brain"
conf_interval  <- predict(model1, newdata = data.frame(x=newx), interval = 'confidence',
                          level=0.95)

#create scatterplot of values with regression line 
plot(weightsE$body, weightsE$body)
abline(model1)

#add dashed lines (lty=2) for the 95% confidence interval
lines(newx, conf_interval[,2], col="blue", lty=2)
lines(newx, conf_interval[,3], col="blue", lty=2)

but the results of predict don't produce a straight line for the upper and lower level, they are more like random predictions.

Upvotes: 0

Views: 826

Answers (1)

Allan Cameron
Allan Cameron

Reputation: 174468

You have a few problems to fix here.

  1. When you generate a model, don't use rlm(weightsE$brain ~ weightsE$body), instead use rlm(brain ~ body, data = weightsE). Otherwise, the model cannot take new data for predictions. Any predictions you get will be produced from the original weightsE$body values, not from the new data you pass into predict
  2. You are trying to create a prediction data frame with a column called "brain', but you are trying to predict the value of "brain", so you need a column called "body"
  3. newx is already a data frame, but for some reason you are wrapping it inside another data frame when you do newdata = data.frame(x=newx). Just pass newx.
  4. You are plotting with plot(weightsE$body, weightsE$body), when it should be plot(weightsE$body, weightsE$brain)

Putting all this together, and using a dummy data set with the same names as your own (see below), we get:

library(MASS)

model1 <- rlm(brain ~ body, data = weightsE)

newx <- data.frame(body = seq(min(weightsE$body), 
                              max(weightsE$body), length.out=70))

conf_interval  <- predict(model1, newdata = data.frame(x=newx), 
                          interval = 'confidence',
                          level=0.95)

#create scatterplot of values with regression line 
plot(weightsE$body, weightsE$brain)
abline(model1)

#add dashed lines (lty=2) for the 95% confidence interval
lines(newx$body, conf_interval[, 2], col = "blue", lty = 2)
lines(newx$body, conf_interval[, 3], col = "blue", lty = 2)

enter image description here

Incidentally, you could do the whole thing in ggplot in much less code:

library(ggplot2)

ggplot(weightsE, aes(body, brain)) + 
  geom_point() + 
  geom_smooth(method = MASS::rlm)

enter image description here

Reproducible dummy data

data(mtcars)
weightsE <- setNames(mtcars[c(1, 6)], c("brain", "body"))
weightsE$body <- 10 - weightsE$body

Upvotes: 1

Related Questions