# RobCov

robCov computes covariance matrix of robust regression coefficients

## Syntax

• out=RobCov(X,scaledres,scaleest)example
• out=RobCov(X,scaledres,scaleest,Name,Value)example

## Description

Under some regularity conditions, robust (S and MM) estimates are asymptotically normal, thereby allowing for Wald-type tests and confidence intervals. The covariance matrix of the estimated parameters $cov(\hat \beta)= q^2 \times \sigma^2 \times v \times V_X^{-1}$

consists of four parts: 1) $q$ a correction factor for the scale estimate;

2) $\sigma$ the scale parameter.

3) $v$ a correction factor depending on the $\psi$ function which is used;

4) $V_X$= a matrix part. For OLS $V_X=X'X$. Given that in robust regression we give a weight to each observation, the matrix $X'X$ should be replaced by something like $X'WX$, where $W$ is a diagonal matrix containing the weights assigned to each observation.

The purpose of this function is to provide the user with different options for the estimate of $cov(\hat \beta)$ where $\hat \beta$ is a vector of regression coefficients obtained using S or MM estimation and a particular $\rho$ function.

 out =RobCov(X, scaledres, scaleest) Compare the 5 estimates of cov matrix.

 out =RobCov(X, scaledres, scaleest, Name, Value) Compare t stat from S and MM estimator.

## Examples

expand all

### Compare the 5 estimates of cov matrix.

n=200;
p=3;
state1=123456;
randn('state', state1);
X=randn(n,p);
y=randn(n,1);
kk=10;
ycont = y;
ycont(1:kk)=ycont(1:kk)+7;
[outS]=Sreg(ycont,X);
rhofunc='optimal';
bdp=0.5;
out=RobCov(X,outS.residuals,outS.scale);
disp('Compare 5 estimates of cov(\hat beta)')
disp('covrob'); disp(out.covrob)
disp('--------')
disp('covrob1'); disp(out.covrob1)
disp('--------')
disp('covrob2'); disp(out.covrob2)
disp('--------')
disp('covrob3'); disp(out.covrob3)
disp('--------')
disp('covrob4'); disp(out.covrob4)
disp('--------')
disp('covrobc'); disp(out.covrobc)

### Compare t stat from S and MM estimator.

n=200;
p=3;
state1=123456;
randn('state', state1);
X=randn(n,p);
y=randn(n,1);
kk=10;
ycont = y;
ycont(1:kk)=ycont(1:kk)+7;
[outS]=Sreg(ycont,X);
rhofunc='optimal';
bdp=0.5;
out=RobCov(X,outS.residuals,outS.scale,'rhofunc',rhofunc,'bdp',0.5);
covrobS=out.covrob;
covrobS1=out.covrob1;
covrobS2=out.covrob2;
covrobS3=out.covrob3;
covrobS4=out.covrob4;
% Compute robust S t-statistics
tstatS=outS.beta./sqrt(diag(covrobS));
tstatS1=outS.beta./sqrt(diag(covrobS1));
tstatS2=outS.beta./sqrt(diag(covrobS2));
tstatS3=outS.beta./sqrt(diag(covrobS3));
tstatS4=outS.beta./sqrt(diag(covrobS4));
eff=0.95;
outMM=MMregcore(ycont,X,outS.beta,outS.scale);
out=RobCov(X,outMM.residuals,outS.scale,'rhofunc',rhofunc,'eff',eff);
covrobMM=out.covrob;
covrobMM1=out.covrob1;
covrobMM2=out.covrob2;
covrobMM3=out.covrob3;
covrobMM4=out.covrob4;
tstatMM=outMM.beta./sqrt(diag(covrobMM));
tstatMM1=outMM.beta./sqrt(diag(covrobMM1));
tstatMM2=outMM.beta./sqrt(diag(covrobMM2));
tstatMM3=outMM.beta./sqrt(diag(covrobMM3));
tstatMM4=outMM.beta./sqrt(diag(covrobMM4));
disp('tstat from S')
disp([tstatS tstatS1 tstatS2 tstatS3 tstatS4])
disp('--------')
disp('tstat from MM')
disp([tstatMM tstatMM1 tstatMM2 tstatMM3 tstatMM4])
qhat=out.q;
disp('tstat from MM after correction for sigma')
disp([tstatMM/qhat tstatMM1/qhat tstatMM2/qhat tstatMM3/qhat tstatMM4/qhat])

## Related Examples

expand all

### An example of the need of using covrobc.

XX is a 50-by-2 matrix generated by standard random normal.

% Estimate beta using Sreg and compare the four estimators of
% cov (\hat \beta).
XX=[0.8115    1.3104
-0.4477    0.8110
0.9438    1.9121
-0.2346    0.7058
1.8364    3.1318
0.2911    0.7837
0.3669    1.7571
0.4560    1.0876
-1.0831    0.0679
0.7018    1.0923
1.6976    2.7744
-0.0624    1.2015
0.2016    0.7747
-0.1542    1.3095
-0.9416   -0.5251
0.4873    0.8116
-0.3561    0.9261
0.5185    2.1219
-0.7674   -0.0514
0.3360    1.8620
-0.2083    1.3796
-0.1700    0.9053
1.9152    2.4369
-0.9912   -0.6885
0.5939    2.0449
-1.1842   -0.8506
1.1375    1.6184
-0.2482    1.2104
0.7936    2.2313
-0.2744    1.5259
1.1455    2.2739
-0.0637    1.2618
-0.5054    0.0536
0.7282    1.9411
-0.0703    0.8090
-1.8186   -0.8725
-2.0621   -0.5299
-1.7750   -0.7765
-0.8553   -0.3311
1.8966    2.9885
-0.1076    0.5106
0.9055    1.1283
0.3086    1.5491
-0.0487    0.4964
-1.6901   -0.7348
-0.7253    0.9436
0.2124    0.6549
-0.2748   -0.2009
-1.2804   -0.7036
0.3627    1.3027];
% Define X and y;
y=XX(:,2);
X=XX(:,1);
outS=Sreg(y,X,'rhofunc','bisquare','bdp',0.5);
out=RobCov(X,outS.residuals,1,'rhofunc','bisquare','bdp',0.5);
disp('covrob'); disp(out.covrob)
disp('--------')
disp('covrob1'); disp(out.covrob1)
disp('--------')
disp('covrob2'); disp(out.covrob2)
disp('--------')
disp('covrob3'); disp(out.covrob3)
disp('--------')
disp('covrob4'); disp(out.covrob4)
disp('--------')
disp('covrobc'); disp(out.covrobc)
% plot X and y and add LS and robust regression line
plot(X,y,'o')
xlabel('X')
ylabel('y')
hLS=lsline;
hLS.Color='r';
hLS.LineStyle = ':';
hLS.LineWidth = 2;
hROB=refline([outS.beta(2) outS.beta(1)]);
hROB.LineWidth = 2;
legend({'' 'LS','Robust'},'Location','Best')
Total estimated time to complete S estimate:  0.69 seconds
covrob
0.6426   -0.0036
-0.0036    0.7248

--------
covrob1
29.6355   -0.4573
-0.4573   29.4373

--------
covrob2
29.6292   -0.1654
-0.1654   33.4196

--------
covrob3
0.2391   -0.0037
-0.0037    0.2375

--------
covrob4
1.0e-03 *

0.2843   -0.0069
-0.0069    0.2487

--------
covrobc
0.6426   -0.0036
-0.0036    0.7248



## Input Arguments

### X — Data matrix of explanatory variables (also called 'regressors') of dimension (n x p-1). Rows of X represent observations, and columns represent variables.

Data Types: single | double

### scaledres — Scaled residuals.Vector. n-times-1 vector containing scaled residuals $r_i/\hat \sigma$.

Data Types: single | double

### scaleest — robust estimate of the scale. Scalar.

Robust estimate of sigma ($\hat \sigma$).

Data Types: single | double

### Name-Value Pair Arguments

Specify optional comma-separated pairs of Name,Value arguments. Name is the argument name and Value is the corresponding value. Name must appear inside single quotes (' '). You can specify several name and value pair arguments in any order as  Name1,Value1,...,NameN,ValueN.

Example:  'intercept',false , 'eff',0.99 , 'bdp',0.4 , 'rhofunc','optimal' , 'rhofuncparam',5 

### intercept —Indicator for constant term.true (default) | false.

Indicator for the constant term (intercept) in the fit, specified as the comma-separated pair consisting of 'Intercept' and either true to include or false to remove the constant term from the model.

Example:  'intercept',false 

Data Types: boolean

### eff —nominal efficiency.scalar.

Scalar defining nominal efficiency (i.e. a number between 0.5 and 0.99). The default value is 0.95 Asymptotic nominal efficiency is: $(\int \psi' d\Phi)^2 / (\psi^2 d\Phi)$

Example:  'eff',0.99 

Data Types: double

### bdp —breakdown point.scalar.

It measures the fraction of outliers the algorithm should resist. In this case any value greater than 0 but smaller or equal than 0.5 will do fine.

Note that given bdp nominal efficiency is automatically determined.

REMARK: just one between bdp and eff must be specified. If both of them are specified an error is produced. If both of them are not specified the defulat is tu use the tuning constant associated to a nominal efficiency of 0.95.

Example:  'bdp',0.4 

Data Types: double

### rhofunc —rho function.string.

String which specifies the rho function which must be used to weight the residuals.

Possible values are 'bisquare' 'optimal' 'hyperbolic' 'hampel'.

'bisquare' uses Tukey's $\rho$ and $\psi$ functions.

See TBrho and TBpsi.

'optimal' uses optimal $\rho$ and $\psi$ functions.

See OPTrho and OPTpsi.

'hyperbolic' uses hyperbolic $\rho$ and $\psi$ functions.

See HYPrho and HYPpsi.

'hampel' uses Hampel $\rho$ and $\psi$ functions.

See HArho and HApsi.

'mdpd' uses Minimum Density Power Divergence $\rho$ and $\psi$ functions.

See PDrho and PDpsi.

'AS' uses Andrews' sine $\rho$ and $\psi$ functions.

See ASrho and ASpsi.

The default is bisquare

Example:  'rhofunc','optimal' 

Data Types: double

### rhofuncparam —Additional parameters for the specified rho function.scalar | vector.

For hyperbolic rho function it is possible to set up the value of k = sup CVC (the default value of k is 4.5).

For Hampel rho function it is possible to define parameters a, b and c (the default values are a=2, b=4, c=8)

Example:  'rhofuncparam',5 

Data Types: single | double

## Output Arguments

### out — description Structure

A structure containing the following fields

Value Description
covrob

p-times-p (if intercept is true else is (p-1)-by-(p-1)) matrix containing asymptotic variance covariance matrix of regression coefficients. covrob uses equation (4.49) of p. 101 of Maronna et al. (2006) namely: $\mbox{covrob} = cov( \hat \beta) = \hat \sigma^2 \hat v (X'X)^{-1}$ where $\hat v = \frac{n}{n-p} n\frac{\sum_{i=1}^n \psi(r_i/\hat \sigma)^2}{[\sum_{i=1}^n \psi'(r_i/\hat \sigma)]^2}$

covrob1

p-times-p (if intercept is true else is (p-1)-by-(p-1)) matrix containing asymptotic variance covariance matrix of regression coefficients. covrob1 uses equation (7.81) of p. 171 of Huber and Ronchetti (2009) with $(X'X)^{-1}$ replaced by $(X' W X)^{-1}$ namely: $\mbox{covrob1} = K^2 \hat v (X' W X)^{-1};$ where $K=1+\frac{p}{n} \frac{var(\psi' (r/\hat \sigma))}{ \left[ \sum_{i=1}^n \psi'(r_i/\hat \sigma)/n \right]^2}= 1+[CV(\psi')]^2 p/n$. The notation $CV(\psi')$ stands for the coefficient of variation of $\psi'$, where $\psi$ and $\psi'$ are, respectively, the first and second derivatives of the $\rho$ function.

covrob2

p-times-p (if intercept is true else is (p-1)-by-(p-1)) matrix containing asymptotic variance covariance matrix of regression coefficients. covrob2 uses equation (7.81) of p. 171 of Huber and Ronchetti (2009) with $X'X$ and $K^2$ namely: $\mbox{covrob2} = K^2 \hat v (X' X)^{-1};$

covrob3

p-times-p (if intercept is true else is (p-1)-by-(p-1)) matrix containing asymptotic variance covariance matrix of regression coefficients. covrob uses equation (7.82) of p. 171 of Huber and Ronchetti (2009).

namely:

$\mbox{covrob3} = K \frac{n}{n-p} \frac{\sum_{i=1}^n \psi (r_i/\hat \sigma)^2}{\sum_{i=1}^n \psi'(r_i/\hat \sigma)} (X' W X)^{-1};$

covrob4

p-times-p (if intercept is true else is (p-1)-by-(p-1)) matrix containing asymptotic variance covariance matrix of regression coefficients. covrob4 uses equation (7.83) of p. 171 of of Huber and Ronchetti (2009).

namely:

$\mbox{covrob4} = \frac{1}{n-p} K^{-1} \sum_{i=1}^n \psi(r_i/\hat \sigma)^2 (X' W X)^{-1} X'X (X' W X)^{-1};$

covrobc

p-times-p (if intercept is true else is (p-1)-by-(p-1)) matrix containing best variance covariance matrix of regression coefficients. covrobc checks the coefficient of variation of the derivative of function $\psi(r_i/\hat \sigma)$. If it is greater than the threshold proposed in Salini, et al. (2020), it uses covrob else covrob2.

q

scalar. Correction for scale estimate (see Maronna and Yohai CSDA 2010). It is defined as $q=1+\frac{p}{2n} \frac{a}{b \times c}$ where $a = \frac{1}{n} \sum_{i=1}^n \psi(r_i/\hat \sigma)^2$ $b = \frac{1}{n} \sum_{i=1}^n \psi'(r_i/\hat \sigma)$ $c = \frac{1}{n} \sum_{i=1}^n \psi(r_i/\hat \sigma) r_i/\hat \sigma$

## References

Maronna, R.A., Martin D. and Yohai V.J. (2006), "Robust Statistics, Theory and Methods", Wiley, New York.

Huber, P.J. and Ronchetti, E.M. (2009), "Robust Statistics, 2nd Edition", Wiley.

Maronna, R.A., and Yohai V.J. (2010), Correcting MM estimates for fat data sets, "Computational Statistics and Data Analysis", Vol. 54, pp. 3168-3173.

Koller, M. and W. A. Stahel (2011), Sharpening wald-type inference in robust regression for small samples, "Computational Statistics & Data Analysis", Vol. 55, pp. 2504-2515.

Croux, C., Dhaene G., and Hoorelbeke D. (2003), Robust standard errors for robust estimators. Technical report, Dept. of Applied Economics, KU Leuven.

Salini, S., Laurini, F., Morelli, G., Riani M. and Cerioli A. (2022), Covariance matrices of S robust regression estimators, Journal of Statistical Computation and Simulation, Vol. 92, pp. 724-747, https://doi.org/10.1080/00949655.2021.1972300