Remark: the M estimator of location must satisfy the following equation
If \rho is differentiable, with
\psi(u) = \rho'(\psi) = \text{d}\psi(u)/\text{d}u,
\hat{\mu} is the solution of the equation (estimating equation)
\sum_{i=1}^n \psi \left( \frac{y_i - \hat{\mu}}{\sigma} \right) = 0.
This estimating equation shows the importance of \psi(y). The equation
can be rewritten to provide an algorithm for estimation of \mu as
\sum_{i=1}^n w_i(y_i - \hat{\mu}) = 0 \quad \text{where} \quad w_i = \psi(y_i - \hat{\mu})/(y_i - \hat{\mu}).
This routine computes the value of \mu which satisfies the above
equation. \sigma corresponds to optional input parameter scaleest.
Note that the value of \sigma is kept fixed in each iteration.
The iterative procedure starts with some \hat{\mu}_0, for example the sample median.
In this routine \hat{\mu}_0 corresponds to optional input parameter
initialmu.
Given the estimate \hat{\mu}_k at stage k, compute
w_{i,k} = w(y_i - \hat{\mu}_k), \quad i = 1,\ldots,n
\hat{\mu}_{k+1} = \sum_{i=1}^n w_{i,k}y_i/\sum_{i=1}^n w_{i,k}.