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