# Mscale

Mscale finds the M estimator of the scale

## Syntax

• sc=Mscale(u, psifunc)example
• sc=Mscale(u, psifunc, initialsc)example
• sc=Mscale(u, psifunc, initialsc, tol)example
• sc=Mscale(u, psifunc, initialsc, tol, maxiter)example

## Description

 sc =Mscale(u, psifunc) Example of M estimate of scale.

 sc =Mscale(u, psifunc, initialsc) Estimate of scale using Hampel rho function.

 sc =Mscale(u, psifunc, initialsc, tol) Use of options initialsc, tol and maxiter.

 sc =Mscale(u, psifunc, initialsc, tol, maxiter) Compare scale estimate using two differnt link functions.

## Examples

expand all

### Example of M estimate of scale.

M estimate of the scale using Tukey biweight rho function with a value of c associated to a breakdown point of 0.5.

psifunc=struct;
psifunc.class='TB';
bdp=0.5;
c=TBbdp(bdp,1);
% kc = E(rho) = sup(rho)*bdp
kc=c^2/6*bdp;
psifunc.c1=c;
psifunc.kc1=kc;
n=10000;
shift=5;
u=2*randn(n,1);
u(1:10)=u(1:10)+shift;
s=Mscale(u,psifunc)

### Estimate of scale using Hampel rho function.

M estimate of the scale using Hampel rho function with a value of c associated to a breakdown point of 0.5

psifunc=struct;
psifunc.class='HA'
abc=[1.5 3.5 8];
bdp=0.5;
c=HAbdp(bdp,1,abc);
% kc = E(rho) = sup(rho)*bdp
kc=HArho(c*abc(3),[c, abc])*bdp;
psifunc.c1=[c abc];
psifunc.kc1=kc;
n=10000;
shift=5;
u=3*randn(n,1);
u(1:10)=u(1:10)+shift;
s=Mscale(u,psifunc)

### Use of options initialsc, tol and maxiter.

M estimate of the scale using Tukey biweight rho function with a value of c associated to a breakdown point of 0.5.

psifunc=struct;
psifunc.class='TB';
bdp=0.5;
c=TBbdp(bdp,1);
% kc = E(rho) = sup(rho)*bdp
kc=c^2/6*bdp;
psifunc.c1=c;
psifunc.kc1=kc;
n=10000;
shift=5;
u=2*randn(n,1);
u(1:10)=u(1:10)+shift;
s=Mscale(u,psifunc,3,1e-7,20)

### Compare scale estimate using two differnt link functions.

psifunc=struct;
psifunc.class='HA'
abc=[1.5 3.5 8];
bdp=0.5;
c=HAbdp(bdp,1,abc);
% kc = E(rho) = sup(rho)*bdp
kc=HArho(c*abc(3),[c, abc])*bdp;
psifunc.c1=[c, abc];
psifunc.kc1=kc;
n=10000;
shift=100;
u=3*randn(n,1);
u(100:200)=u(100:200)+shift;
sHA=Mscale(u,psifunc)
% kc = E(rho) = sup(rho)*bdp
psifunc=struct;
psifunc.class='TB';
c=TBbdp(bdp,1);
psifunc.c1=c;
kc=c^2/6*bdp;
psifunc.kc1=kc;
sTB=Mscale(u,psifunc,3,1e-7,20);
sMLE=std(u);
cate=categorical({'Robust scale (Hampel)' 'Robust scale (TB)' 'Non robust scale'});
bar(cate,[sHA sTB sMLE])

## Input Arguments

### u — residuals or Mahalanobis distances. Vector.

n x 1 vector which contains the (unscaled) residuals or Mahalanobis distances

Data Types: single | double

### psifunc — rho (psi) function. Structure.

A structure specifying the class of rho (psi) function to use, the consistency factor, and the value associated with the Expectation of rho in correspondence of the consistency factor.

psifunc must contain the following fields:

Value Description
class

string identyfing the rho (psi) function to use.

Admissible values for class are 'bisquare' (TB), 'optimal', (OPT) 'hyperbolic' (HYP), 'hampel' (HA) 'power divergence' (PD)

c1

consistency factor (and other parameters) associated to required breakdown point or nominal efficiency.

More precisely, psifunc.c1(1) contains consistency factor associated to required breakdown point or nominal efficiency psifunc.c1(2:end) contain other parameters associated with the rho (psi) function.

For example, if psifunc.class='hampel', c1(2:4) must contain parameters (a, b and c) of Hampel rho function.

kc1

Expectation of rho associated with c1 to get a consistent estimator at the model distribution kc1=E(rho)

Example: Example - psifunc.class='TB';psifunc.c1=1.5476;psifunc.kc1=0.1996

Data Types: struct

### initialsc — scalar. The initial estimate of the scale.

If not defined, scaled MAD of vector |u| is used.

Example: 'initialsc',0.34 

Data Types: double

### tol — scalar. The tolerance for controlling convergence.

If not defined, tol is fixed to 1e-7.

Example: 'tol',1e-10 

Data Types: double

### maxiter — scalar. Maximum number of iterations to find the scale.

If not defined, maxiter is fixed to 200.

Example: 'maxiter',100 

Data Types: double

## Output Arguments

### sc —M-estimate of the scale. Scalar

Robust M estimate of scale.

This routine is called by Taureg and Sreg.m

u = residuals or Mahalanobis distances (note that u is kept fixed in each iteration).

Remark: the M estimator of scale must satisfy the following equation

$(1/n) \sum_{i=1}^n \rho((u_i/c)/s) = kc$

This routine computes the value of $s$ which satisfies the above equation.

Huber, P.J. and Ronchetti, E.M. (2009), "Robust Statistics, 2nd Edition", Wiley. [equation 7.119, p. 176].