| Library: | robtech |
| See also: | median lms linreg |
| Quantlet: | robmest | |
| Description: | calculates M-estimators in linear model |
| Usage: | {mest,iter}=robmest(y,x,type{,pars,start,prec,maxit}) | |
| Input: | ||
| y | n x 1 vector; dependent variable | |
| x | n x p matrix of independent variables | |
| type | string, type of score function: "huber" or "hampel" | |
| pars | optional scalar(s), parameter(s) of the score function: scalar h: 0 < h < Inf for "huber"; default is 1.5 two comma separated scalars h1,h2 such that 0 < h1 < h2 < Inf for "hampel"; default is h1=1, h2=2*h1 | |
| start | optional p x 1 vector representing the starting estimate, or string "ls" for least squares estimate (default) | |
| prec | optional scalar, precision of the estimate, iteration stop if the change in sum of absolute values of components of the estimate is smaller than prec; default is 0.0001 | |
| maxit | optional scalar, maximum number of iterations; default is 100 | |
| Output: | ||
| mest | p x 1 vector containing the M-estimate | |
| iter | scalar, number of iterations | |
library("robtech")
n=100
randomize(0)
x=matrix(n)~normal(n,2)
tmp=(uniform(n,1)>0.80)*100
true=(1|2|4)
y=x*true+normal(n,1)+tmp.*normal(n,1)
ls=inv(x'*x)*x'*y
{mest1,i1}=robmest(y,x,"huber",1,0|0|0)
{mest2,i2}=robmest(y,x,"huber",1.5,0|0|0,0.000001)
{mest3,i3}=robmest(y,x,"hampel",1,3,"ls",0.000001,15)
true~ls~mest1~mest2~mest3
Contents of _tmp
[1,] 1 4.9191 1.1768 1.2396 1.0012
[2,] 2 3.9026 2.192 2.2628 1.947
[3,] 4 7.5184 4.0706 4.1097 3.9515
Matrix containing in columns
1) the true value of the parameters,
2) least squares estimate
3) Huber's type M-estimate with parameter 1 and starting
estimate (0,0,0),
4) Huber's type M-estimate with parameter 1.5, starting
estimate (0,0,0), and precision 0.000001
5) Hampel's type M-estimate with parameters 1 and 3,
least squares starting estimate,
precision 0.000001 and a maximum of 15 iterations