Elegant way to express residual calculation in R?

   Hi All,

   I am illustrating a simple, two-way ANOVA using the following data and I'm
   having difficulty in expressing the predicted values succinctly in R.

   X<- data.frame(read.table(textConnection("
       Machine.1    Machine.2    Machine.3
       53           61           51
       47           55           51
       46           52           49
       50           58           54
       49           54           50"
   ), header=TRUE))
   rownames(X)<- paste("Operator.", 1:nrow(X), sep="")
   print(X)

   # I'd like to know if there is a more elegant way to calculate the residuals
   # than the following, which seems to be rather a kludge. If you care to read
   # the code you'll see what I mean.

   machine.adjustment<-  colMeans(X) - mean(mean(X))    # length(machine.adjustment)==3
   operator.adjustment<- rowMeans(X) - mean(mean(X))    # length(operator.adjustment)==5
   X.predicted<- numeric(0)
   for (j in 1:ncol(X))
   {
       new.col<- mean(mean(X)) + operator.adjustment + machine.adjustment[j]
       X.predicted<- cbind(X.predicted, new.col)
   }
   print(X.predicted)
   X.residual<- X - X.predicted
   SS.E<- sum( X.residual^2 )

It seems like there ought to be some way of doing that a little bit cleaner ...

Thanks,

Jack.


---------------------------------

Bring photos to life! New PhotoMail  makes sharing a breeze.
       [[alternative HTML version deleted]]

Try this:

Xm <- as.matrix(X)
X.lm <- lm(c(Xm) ~ factor(col(Xm)) + factor(row(Xm)))
sum(resid(X.lm)^2)

or if the idea was to do it without using lm try replacing
your calculation of X.predicted with this:

X.predicted <-
 outer(operator.adjustment, machine.adjustment, "+") + mean(mean(X))


On 2/27/06, John McHenry <john_d_mchenry at yahoo.com> wrote:

   Hi All,

   I am illustrating a simple, two-way ANOVA using the following data and I'm
   having difficulty in expressing the predicted values succinctly in R.

   X<- data.frame(read.table(textConnection("
       Machine.1    Machine.2    Machine.3
       53           61           51
       47           55           51
       46           52           49
       50           58           54
       49           54           50"
   ), header=TRUE))
   rownames(X)<- paste("Operator.", 1:nrow(X), sep="")
   print(X)

   # I'd like to know if there is a more elegant way to calculate the residuals
   # than the following, which seems to be rather a kludge. If you care to read
   # the code you'll see what I mean.

   machine.adjustment<-  colMeans(X) - mean(mean(X))    # length(machine.adjustment)==3
   operator.adjustment<- rowMeans(X) - mean(mean(X))    # length(operator.adjustment)==5
   X.predicted<- numeric(0)
   for (j in 1:ncol(X))
   {
       new.col<- mean(mean(X)) + operator.adjustment + machine.adjustment[j]
       X.predicted<- cbind(X.predicted, new.col)
   }
   print(X.predicted)
   X.residual<- X - X.predicted
   SS.E<- sum( X.residual^2 )

It seems like there ought to be some way of doing that a little bit cleaner ...

Thanks,

Jack.


---------------------------------

Bring photos to life! New PhotoMail  makes sharing a breeze.
       [[alternative HTML version deleted]]