# why does R rlm {MASS} return different coefficients almost each time it is called?

I'm noticed that the rlm {MASS} returns almost every time different coefficients, even though I'm using the same parameters and the same data set

I'm calling:

model <- rlm(price ~ ., data = data[,-1], weights = weights,
maxit = 1000000, init = "lts", method = "MM",
psi=psi.huber, acc=0.00001, scale.est="proposal 2", cor = T)


results:

Intercept)  livingArea        area    discrete       dummy
-570.621795   17.169323    2.275109   46.002527  143.812900

(Intercept)  livingArea        area    discrete       dummy
-581.893552   16.828956    3.912192   48.253955  180.875439

(Intercept)  livingArea        area    discrete       dummy
-303.488284   16.747009    2.928579   26.951809  -14.795652


This are the three most frequent results.

Could somebody explain me why does the rlm behave like this?

NOTE: I'm am a beginner in the field of statistics.

• Please edit to add name of the library that rlm function comes from, we can guess it is a MASS library, but it is better to state this explicitly. – Tim Oct 23 '15 at 13:16
• @Tim yes, it is MASS. I edited the post and added this information. – Paul Oct 23 '15 at 13:18

## 1 Answer

Because of the algorithm you choose to compute the initial values (init) (from which the second stage of the algorithm, i.e. the MM steps, start).

Setting init=lts uses the coefficients fitted by the FastLTS algorithm as starting points of the MM iterations. The FastLTS algorithm in turn uses many random starting points so is itself random. Unless you fix the seed argument, you will get different solutions (which is as it should be!).

library(MASS)
model1 <- rlm(stack.loss ~ ., data = stackloss,maxit = 1000000, init = "lts",seed=1, method = "MM", psi=psi.huber, acc=0.00001, scale.est="proposal 2", cor = T)
model2 <- rlm(stack.loss ~ ., data = stackloss,maxit = 1000000, init = "lts",seed=1, method = "MM", psi=psi.huber, acc=0.00001, scale.est="proposal 2", cor = T)
model1$coef-model2$coef


Setting method = "MM" without setting init uses the coefficient fitted by the FastS algorithm as starting points of the MM iterations (source). The FastS algorithm in turn uses many random starting points so is itself random. Unless you fix the seed argument, you will get different solutions (which is as it should be!).

library(MASS)
model1 <- rlm(stack.loss ~ ., data = stackloss,maxit = 1000000, seed=1, method = "MM", psi=psi.huber, acc=0.00001, scale.est="proposal 2", cor = T)
model2 <- rlm(stack.loss ~ ., data = stackloss,maxit = 1000000, seed=1, method = "MM", psi=psi.huber, acc=0.00001, scale.est="proposal 2", cor = T)
model1$coef-model2$coef

• also letting the init=lts away, generates different results – Paul Oct 23 '15 at 14:08
• Did my edit address your comment? – user603 Oct 23 '15 at 18:44