Compute the Minimum Covariance Determinant (MCD) estimator, a robust multivariate location and scale estimate with a high breakdown point, via the ‘Fast MCD’ or ‘Deterministic MCD’ (“DetMcd”) algorithm.
covMcd(x, cor = FALSE, raw.only = FALSE,
alpha =, nsamp =, nmini =, kmini =,
scalefn =, maxcsteps =,
initHsets = NULL, save.hsets = FALSE, names = TRUE,
seed =, tolSolve =, trace =,
use.correction =, wgtFUN =, control = rrcov.control())a matrix or data frame.
should the returned result include a correlation matrix?
Default is cor = FALSE.
should only the “raw” estimate be returned, i.e., no (re)weighting step be performed; default is false.
numeric parameter controlling the size of the subsets
over which the determinant is minimized; roughly alpha*n,
(see ‘Details’ below)
observations are used for computing the determinant. Allowed values
are between 0.5 and 1 and the default is 0.5.
number of subsets used for initial estimates or "best",
"exact", or "deterministic". Default is nsamp = 500.
For nsamp = "best" exhaustive enumeration is done, as long as
the number of trials does not exceed 100'000 (= nLarge).
For "exact", exhaustive enumeration will be attempted however
many samples are needed. In this case a warning message may be
displayed saying that the computation can take a very long time.
For "deterministic", the deterministic MCD is computed; as
proposed by Hubert et al. (2012) it starts from the \(h\) most
central observations of six (deterministic) estimators.
for \(n \ge 2 \times n_0\),
\(n_0 := \code{nmini}\), the algorithm splits the data into
maximally kmini (by default 5) subsets, of size
approximately, but at least nmini. When nmini*kmini < n,
the initial search uses only a subsample of size nmini*kmini.
The original algorithm had nmini = 300 and kmini = 5
hard coded.
for the deterministic MCD: function to
compute a robust scale estimate or character string specifying a
rule determining such a function. The default, currently
"hrv2012", uses the recommendation of Hubert, Rousseeuw and
Verdonck (2012) who recommend Qn
for \(n < 1000\) and scaleTau2 for larger n. Alternatively,
scalefn = "v2014", uses that rule with cutoff \(n = 5000\).
maximal number of concentration steps in the deterministic MCD; should not be reached.
NULL or a \(K x h\) integer matrix of initial
subsets of observations of size \(h\) (specified by the indices in
1:n).
(for deterministic MCD) logical indicating if the
initial subsets should be returned as initHsets.
initial seed for random generator, like
.Random.seed, see rrcov.control.
numeric tolerance to be used for inversion
(solve) of the covariance matrix in mahalanobis.
logical (or integer) indicating if intermediate results
should be printed; defaults to FALSE; values \(\ge 2\)
also produce print from the internal (Fortran) code.
whether to use finite sample correction
factors; defaults to TRUE.
a character string or function, specifying
how the weights for the reweighting step should be computed. Up to
April 2013, the only option has been the original proposal in (1999),
now specified by wgtFUN = "01.original" (or via
control). Since robustbase version 0.92-3, Dec.2014,
other predefined string options are available, though experimental,
see the experimental .wgtFUN.covMcd object.
a list with estimation options - this includes those
above provided in the function specification, see
rrcov.control for the defaults. If control is
supplied, the parameters from it will be used. If parameters are
passed also in the invocation statement, they will override the
corresponding elements of the control object.
An object of class "mcd" which is basically a
list with components
the final estimate of location.
the final estimate of scatter.
the (final) estimate of the correlation matrix (only if
cor = TRUE).
the value of the criterion, i.e., the logarithm of the determinant. Previous to Nov.2014, it contained the determinant itself which can under- or overflow relatively easily.
the best subset found and used for computing the raw
estimates, with length(best) == quan =
h.alpha.n(alpha,n,p).
mahalanobis distances of the observations using the final estimate of the location and scatter.
weights of the observations using the final estimate of the location and scatter.
a vector of length two containing the consistency correction factor and the finite sample correction factor of the final estimate of the covariance matrix.
the raw (not reweighted) estimate of location.
the raw (not reweighted) estimate of scatter.
mahalanobis distances of the observations based on the raw estimate of the location and scatter.
weights of the observations based on the raw estimate of the location and scatter.
a vector of length two containing the consistency correction factor and the finite sample correction factor of the raw estimate of the covariance matrix.
the input data as numeric matrix, without NAs.
total number of observations.
the size of the subsets over which the determinant is minimized (the default is \((n+p+1)/2\)).
the number of observations, \(h\), on which the MCD is
based. If quan equals n.obs, the MCD is the classical
covariance matrix.
character string naming the method (Minimum Covariance
Determinant), starting with "Deterministic" when
nsamp="deterministic".
(for the deterministic MCD) contains indices from 1:6 denoting which of the (six) initial subsets lead to the best set found.
(for the deterministic MCD) for each of the initial subsets, the number of C-steps executed till convergence.
the call used (see match.call).
The minimum covariance determinant estimator of location and scatter
implemented in covMcd() is similar to R function
cov.mcd() in MASS. The MCD method looks for
the \(h (> n/2)\) (\(h = h(\alpha,n,p) =\)
h.alpha.n(alpha,n,p)) observations (out of \(n\))
whose classical covariance matrix has the lowest possible determinant.
The raw MCD estimate of location is then the average of these \(h\) points,
whereas the raw MCD estimate of scatter is their covariance matrix,
multiplied by a consistency factor (.MCDcons(p, h/n)) and (if
use.correction is true) a finite sample correction factor
(.MCDcnp2(p, n, alpha)), to make it consistent at the
normal model and unbiased at small samples. Both rescaling factors
(consistency and finite sample) are returned in the length-2 vector
raw.cnp2.
The implementation of covMcd uses the Fast MCD algorithm of
Rousseeuw and Van Driessen (1999) to approximate the minimum
covariance determinant estimator.
Based on these raw MCD estimates, (unless argument raw.only is
true), a reweighting step is performed, i.e., V <- cov.wt(x,w),
where w are weights determined by “outlyingness” with
respect to the scaled raw MCD. Again, a consistency factor and
(if use.correction is true) a finite sample correction factor
(.MCDcnp2.rew(p, n, alpha)) are applied.
The reweighted covariance is typically considerably more efficient
than the raw one, see Pison et al. (2002).
The two rescaling factors for the reweighted estimates are returned in
cnp2. Details for the computation of the finite sample
correction factors can be found in Pison et al. (2002).
Rousseeuw, P. J. and Leroy, A. M. (1987) Robust Regression and Outlier Detection. Wiley.
Rousseeuw, P. J. and van Driessen, K. (1999) A fast algorithm for the minimum covariance determinant estimator. Technometrics 41, 212--223.
Pison, G., Van Aelst, S., and Willems, G. (2002) Small Sample Corrections for LTS and MCD, Metrika 55, 111--123.
Hubert, M., Rousseeuw, P. J. and Verdonck, T. (2012) A deterministic algorithm for robust location and scatter. Journal of Computational and Graphical Statistics 21, 618--637.
cov.mcd from package MASS;
covOGK as cheaper alternative for larger dimensions.
# NOT RUN {
data(hbk)
hbk.x <- data.matrix(hbk[, 1:3])
set.seed(17)
(cH <- covMcd(hbk.x))
cH0 <- covMcd(hbk.x, nsamp = "deterministic")
with(cH0, stopifnot(quan == 39,
iBest == c(1:4,6), # 5 out of 6 gave the same
identical(raw.weights, mcd.wt),
identical(which(mcd.wt == 0), 1:14), all.equal(crit, -1.045500594135)))
## the following three statements are equivalent
c1 <- covMcd(hbk.x, alpha = 0.75)
c2 <- covMcd(hbk.x, control = rrcov.control(alpha = 0.75))
## direct specification overrides control one:
c3 <- covMcd(hbk.x, alpha = 0.75,
control = rrcov.control(alpha=0.95))
c1
## Martin's smooth reweighting:
## List of experimental pre-specified wgtFUN() creators:
## Cutoffs may depend on (n, p, control$beta) :
str(.wgtFUN.covMcd)
cMM <- covMcd(hbk.x, wgtFUN = "sm1.adaptive")
ina <- which(names(cH) == "call")
all.equal(cMM[-ina], cH[-ina]) # *some* differences, not huge (same 'best'):
stopifnot(all.equal(cMM[-ina], cH[-ina], tol = 0.2))
# }
Run the code above in your browser using DataLab