[R] Possible bug in rlm

[Francis Bursa]

Dear all,

I believe I have found a bug in rlm in the MASS package. Specifically, the scale estimate can be wrong when there are no outliers. The following code snippet is an example:

dose <- c(0,1,2,0,1,2)
response <- c(0.659,1.633,3.621,1.803,3.093,4.424)
line <- c(1,1,1,2,2,2)
k2 <- seq(1.5,5,by=0.01)
repNA <- rep(NA,length(k2))
scale <- repNA
niter <- repNA
for (i in 1:length(k2)){
  rlm.fit <- rlm(response~dose+factor(line), psi=psi.huber,k=1.345,
                         scale.est="proposal 2",k2=k2[i])
  scale[i] <- rlm.fit$s
  niter[i] <- length(rlm.fit$conv)
}
plot(k2,scale,type="b",col=niter)

For this dataset there are no outliers, so I would expect the scale to be a smooth function of k2 once k2 is reasonably large, certainly for k2 > 2. However, there is a funny jump in the scale estimate around k2 = 2.4, just at the point where the number of iterations to convergence falls from 3 to 1.

Looking at the source code, it appears that on each iteration, the scale is updated, then the parameters, and then a check for convergence is carried out just for the parameters, not the scale. So I would guess that in the range around k2=2.5 convergence is being reached when in fact the scale estimate hasn't converged.

I am using MASS version 7.3-33 and R version 3.1.0 on Windows.

I am not sure how common this issue is but there does not seem to be anything special about my dataset so it could be quite generic. Am I right that this is a bug?

Many thanks,
Francis Bursa

------------------ 
Francis Bursa
Statistician

Quantics Consulting Ltd
28 Drumsheugh Gardens
Edinburgh
EH3 7RN
 
Telephone: +44 (0) 131 440 2781 ext 207
 
Quantics - complex data into clear results                         
Quantics is an ISO 9001 registered company                                   
www.quantics.co.uk

Please note that the contents of this e-mail (including any attachments) are confidential and may be legally privileged. If you are not the intended recipient you may not read, copy, distribute or make any other use of this email or its contents. If received in error, please tell us immediately by telephone on +44 (0) 131 440 2781 quoting the name of the sender and the intended recipient, then delete it from your system. Thank you.