2

I am currently working with R and I'm trying to write a function that derives the partial residuals for a multiple linear model. I know that there are existing functions in R but I want to write a function myself.

However, my problem is that executing the function always returns numeric(0):

y <- c(2, 3, 6, 8, 12, 15)
x1 <- c(1, 1.5, 2.5, 4, 7, 9)
x2 <- c(10, 22, 30, 31, 38, 42)

fitlm <- lm(y ~ x1+x2)

PR <- function(fit, predictor)
{
PRes <- as.numeric(fit$residuals) + (fit$coefficients["predictor"]*fit$model$predictor)
return(PRes)
}
# with fit = name of linear model and predictor = name of predictor

PR_x1 <- PR(fitlm, x1)
PR_x2 <- PR(fitlm, x2)

but when I only execute the code outside the function everything works and I get the partial residuals:

as.numeric(fitlm$residuals)+fitlm$coefficients["x1"]*fitlm$model$x1

Does anybody know how to fix this? What am I doing wrong?

Thanks for your help!

Philipp
  • 57
  • 1
  • 8

1 Answers1

4

You are hardcoding the predictor in

fit$coefficients["predictor"]

but i'm guessing you want to use whatever predictor you pass in. This should do the trick:

PR <- function(fit, predictor)
{
  PRes <- fit$residuals + fit$coefficients[predictor] * fit$model[,predictor]
  return(as.numeric(PRes))
}

x <- 1:1000
y <- rnorm(1000)
fit <- lm(y~x)

head(PR(fit,"x"))

prints

> head(PR(fit,"x"))
[1]  0.5813375 -0.2879395  0.1891699  0.8803139  0.9769820  0.5359668

EDIT OP's just provided an example.

If you want to specify predictor itself (i.e. using variable x instead of "x" as string) then

PR2 <- function(fit, predictor)
{
  name <- deparse(substitute(predictor))
  PRes <- fit$residuals + fit$coefficients[name] * fit$model[,name]
  return(as.numeric(PRes))
}

x <- 1:1000
y <- rnorm(1000)
fit <- lm(y~x)

head(PR2(fit,x))

which still prints

> head(PR2(fit,x))
[1]  1.16292881 -1.03540054  1.86687897 -0.02262028 -1.46575723  1.79580590
rbm
  • 3,243
  • 2
  • 17
  • 28
  • 1
    Thanks a lot! This works perfectly fine. So my mistake was just the hardcoding of the predictor? – Philipp May 19 '16 at 12:26
  • No, the other issue was the second part `fit$model$predictor)`. Since `fit$model` is a data frame, then you need to access column stored in variable `predictor` using `fit$model[,predictor]` – rbm May 19 '16 at 14:04