Dimension reduction regression in R

Regression is the study of the dependence of a response variable y on a collection predictors p collected in x. In dimension reduction regression, we seek to find a few linear combinations β1x,...,βdx, such that all the information about the regression is contained in these linear combinations. If d is very small, perhaps one or two, then the regression problem can be summarized using simple graphics; for example, for d=1, the plot of y versus β1x contains all the regression information. When d=2, a 3D plot contains all the information. Several methods for estimating d and relevant functions of β1,..., βdhave been suggested in the literature. In this paper, we describe an R package for three important dimension reduction methods: sliced inverse regression or sir, sliced average variance estimates, or save, and principal Hessian directions, or phd. The package is very general and flexible, and can be easily extended to include other methods of dimension reduction. It includes tests and estimates of the dimension , estimates of the relevant information including β1,..., βd, and some useful graphical summaries as well.


Introduction
In the general regression problem, we have a response # of dimension $ (usually, $ &% (' ) and ) - dimensional predictor 0 , and the goal is to learn about how the conditional distributions 1 32 4# 65 0 87 as 0 varies through its sample space.In parametric regression, we specify a functional form for the condi- tional distributions that is known up to a few parameters.In nonparametric regression, no assumptions are made about 1 , but progress is really only possible if the dimensions ) and $ are small.
Dimension reduction regression is one intermediate possibility between the parametric and nonparametric extremes.In this setup, we assume without loss of information that the conditional distributions can be indexed by 9 linear combinations, or for some probably unknown ) A@ B9 matrix C 1 32 D# 65 0 ©7 E% GF H2 D# 65 C PI Q0 ©7 (1) or in words, all the information in 0 about # is contained in the 9 linear combinations C I 0 .This representation always holds trivially, by setting C R% TS , the ) U@ V) identity matrix, and so the usual goal is to find the C of lowest possible dimension for which this representation holds.If (1) holds for a particular C , then it also holds for C 3W X% YC a`, where `is any full rank matrix, and hence the unique part of the regression summary is the subspace that is spanned by C , which we denote b c2 4C a7 .Cook (1998)  In this paper, we discuss software for estimating the subspace b 2 4C 7 spanned by C , and tests concerning the dimension 9 based on dimension reduction methods.This software was written using R, but can also be used with Splus.Most, but not all, of the methods available here are also included in the Xlisp-Stat program Arc, Cook and Weisberg (1999).The R platform allows use of dimension reduction methods with existing statistical methods that are not readily available in Xlisp-Stat, and hence in Arc.For example, the R functions are more suitable for Monte Carlo experimentation than Arc.In addition, R includes a much wider array of options for smoothing, including multidimensional smoothers.On the other hand, Arc takes full advantage of the dynamic graphical capabilities of Xlisp-Stat, and at least for now the graphical summaries of dimension reduction regression are clearly superior in Arc.Thus, there appears to be good reason to have these methods available using both platforms.Cook (1998) provides the most complete introduction to this area, and we leave the details of all the methods described here to that book.See also Cook and Weisberg (1994) for a more gentle introduction to dimension reduction.In this paper we give only the barest outline of dimension reduction methodology, concentrating on the software.Suppose we have data 2 40 ¡ £¢ # ¡ 7 , for ¤ % R' ¢ ¦¥ ¦¥ §¥ ¨¢ © that are independent and collected into a matrix and a vector if $ % ' and a matrix if $ ' .In addition, suppose we have nonnegative weights ¢ ¦¥ ¦¥ §¥ ¨¢ whose sum is © ; if unspecified, we take all the ¡ % &' .Generally following Yin (2000), a procedure for estimating b c2 4C a7 and for obtaining tests concerning 9 is: 1. Scale and center as where (0 % 43 ¡ 0 ¡ 65 3 ¡ is the vector of weighted column means, !% "7 98 A@ CB ©2 ¡ 7 , and is any square root of the inverse of the sample covariance matrix for (for example, using a singular value decomposition) and (0 is a vector of weighted sample means.In this scaling, the rows of have zero mean and identity sample covariance matrix. 2. Use the scaled and centered data to find a dimension ) @ U) symmetric matrix , and . In some applications (like phd) the eigenvalues may be negative.

4.
A test that the dimension 9 % 9 XW against the alternative that 9 Y 9 XW is based on a partial sum of eigenvalues of the form: where 0 e is a method-specific term, and q is generally equal to 1, but it is equal to 2 for phd.The distribution of these partial sums depends on assumptions and on the method of obtaining 0 H . 5. Given 9 , the estimate of b c2 4C a7 is the span of the first 9 eigenvectors.When viewed as a sub- space of r , the basis for this estimated subspace is . These directions can then be back-transformed to the -scale.Given the estimate of b 2 4C 7 , graphical methods can be used to recover information about 1 , or about particular aspects of the conditional distributions, such as the conditional mean function.

Usage
The R package for dimension reduction called dr can be obtained from the Comprehensive R Archive Network, CRAN, at http://www.R-project.org.The primary function in this package is called dr.Its usage is modeled after the R function lm that is used for fitting linear models, and so many of the arguments that can be used in lm can also be used in dr.The following statements will examine dimension reduction using sir applied to the Australian Institute of Sport data: > library(dr) > data(ais) > attach(ais) > i1 <-dr(LBM~Ht + Wt + log(RCC) + WCC, method="sir", nslices=8) First, provide a formula giving the response on the left side, and the predictors on the right side.
The response can be a vector of © elements or for sir and save it can be a matrix with © rows and $ columns.The model statement on the right side can include factors and interactions, or transformations, so the syntax is very general.The remaining arguments are all optional, and the values shown are all defaults.method can equal either "sir", "save", "phd", "phdres", "phdq" or "phdy".Adding other methods is not difficult, and is described in Section 3.6.The number of slices is relevant only for the methods that use slicing, sir and save, and has default equal to the larger of 8 and the number of predictors plus 3. If the response has $ columns, then the argument nslices should be a vector of $ elements.If it is specified as a number rather than a vector, then that number will give the total number of cells, approximately.For example, if $ U% ¡ and nslices=8, the program will slice ¢ £ ¥¤ §¦ slices along each of the two response variables for a total of ¦ @ ¦ % ©¨cells.The argument numdir specifies the number of directions to display in printed output and the number of tests to compute.The default is 4. Keywords inherited from lm include weights, which if set should be a vector of the same length of the response of positive numbers; contrasts, which specifies how to turn factors into variables; na.action, which specifies how to "handle" missing values.The default for this argument is na.omit, which will reduce the dataset by eliminating all cases with one or more missing value.The argument subset can be set to a list of case indices to be used in the computations; all other cases are ignored.
Brief printed output from the dr method is obtained by simply typing the name of the object created, i1 in the above example: This output repeats most of the basic output, plus it also gives information on the slices and on tests of dimension, if available.The class of the object created by dr depends on the value of the argument method.For example, if method="sir", the object is of class sir if the response is univariate and msir if the response is multivariate.If method="save", the object is of class save or msave.All these objects inherit from the class dr.

Sliced inverse regression
Sliced inverse regression, or sir, was proposed by Li (1991); see Cook (1998, Chapter 11).In sir, we make use of the fact that given certain assumptions on the marginal distribution of 1 , the inverse regression problem 1 32 ¡ 85 # 7 I b 2 4C 7 .The general computational outline for sir is as follows: 1. Examine 1 32 ¡ 85 # 7 by dividing the range of into ¢ slices, each with approximately the same number of observations.With a multivariate response ( has $ columns), divide the range of @ Y $ ¥ §¥ ¦¥ @ Y ¤£ into ¢ cells.For example, when $ % ¦ , and we slice into 3 slices, $ into 2 slices, and ¤¥ into 4 slices, we will have ¢ % © @ ¦ @ §¦ a% © ¨¦ cells.The number of slices or cells ¢ is a tuning parameter of the procedure.
2. Assume that within each slice or cell 1 32 ¡ 85 # 7 is approximately constant.Then the expected value of the within-slice vector of sample means will be a vector in b 2 C a7 .

Form the
¢ @ ) matrix whose ¤ -th row is the vector of weighted sample means in the ¤ -th slice.
The matrix 0 H is the ) A@ a) sample covariance matrix of these sample mean vectors.
sir thus concentrates on the mean function E2 5 7 , and ignores any other dependence.
The output given in the last section is an example of typical output for sir.First is given the eigenvectors and eigenvalues of 0 H ; the eigenvectors have been back-transformed to the original -scale.
Assuming that the dimension is 9 , the estimate of b c2 4C a7 is given by the first 9 eigenvectors.Also given along with the eigenvectors is the square of the correlation between the ols fitted values and the first 9 principal directions.The first direction selected by sir is almost always about the same as the first direction selected by ols, as is the case in the example above.
For sir, Li (1991) provided asymptotic tests of dimension based on partial sums of eigenvalues, and these tests are given in the summary.The tests have asymptotic Chi-square distributions, with the number of degrees of freedom shown in the output.
Examining the tests shown in the final output, we see that the test of 9 % © versus 9 " © has a very small ) -value, so we would reject 9 Y% © .The test for 9 % ' versus 9 ' has ) -value near © ¥ © © ' , suggesting that 9 is at least 2. The test for 9 % versus 9 © has ) -value of about © ¥ ¦ ' , so we suspect that 9 a% © for this problem.This suggests that further analysis of this regression problem can be done based on the 3D graph of the response versus the linear combinations of the predictors determined by the first two eigenvectors, and the dimension of the problem can be reduced from 4 to 2 without loss of information.See Cook (1998), and Cook and Weisberg (1994, 1999), for further examples and interpretation.
When the response is multivariate, the format of the call is: The summary for a multivariate response is similar: 1 For the methods discussed here, we generally need the linearity condition, the E ¥ !H¥ " $# is a linear function for !¥ % , where is any matrix, and This condition holds, for example, if the marginal distribution of 1 is normal, although that assumption is much stronger than is needed.See Cook (1998), Cook andWeisberg (1994, 1999) for more discussion.The test statistics are the same as in the univariate response case, as is the interpretation of the eigenvalues and vectors.The output gives the squared correlation of each of the responses with the eigenvectors.

Sliced average variance estimation
Sliced average variance estimation, or save, was proposed by Cook and Weisberg (1991).As with sir, we slice the range of into ¢ slices, but rather than compute the within-slice mean we compute withinslice covariance matrices.If ¡ is the weighted within slice sample covariance matrix in slice ¤ , then the matrix 0 H is given by where ¡ ¡ is the sum of the weights in the slice; if all weights are equal, then the ¡ ¡ are just the number of observations in each slice.save looks at second moment information and may miss first-moment information, particularly it may miss linear trends.
Output for save is similar to sir, except that no asymptotic tests have been developed.However, tests of dimension based on a permutation test are available; see Section 4. Li (1992) proposed the method called principal Hessian directions, or pHd.This method examines the matrix 0 H

Principal Hessian direction
given by where F ¡ is either equal to # ¡ (method phdy) or F ¡ is an ols residual (method phd or phdres), and ¡ is the weight for the ¤ -th observation (recall again that we assume that 3 ¡ % © , and if this is not satisfied the program rescales the weights to meet this condition).While all methods produce The column of tests called "normal theory" were proposed by Li (1992) and require that the predictors are normally distributed.These statistics are asymptotically distributed as Chi-square, with the degrees of freedom shown.When the method is phdres additional tests are provided.Since this method is based on residuals, it gives tests concerning the central subspace for the regression of the residuals on rather than the response on .The subspace for this residual regression may be, but need not be, smaller than the subspace for the original regression.For example, the column marked "Indep.test" is essentially a test of 9 U% © versus 9 © described by Cook (1998) for the residual regression.Should the significance level for this test be large, we might conclude that the residual regression subspace is of dimension zero.From this we have two possible conclusions: (1) the dimension of the response regression may be 1 if using the residuals removed a linear trend, or (2) the dimension may be 0 if the residuals did not remove a linear trend.
Similarly, if the significance level for the independence test is small, then we can conclude that the dimension is at least 1.It could be one if the method is picking up a nonlinear trend in the OLS direction, but it will be 2 if the nonlinearity is in some other direction.
The independence test and the final column, also from Cook (1998), use the same test statistic, but different distributions based on different assumptions.Significance levels are obtained by comparing the statistic to the distribution of a random linear combination of Chi-square statistics, each with one df.These statistics do not require normality of the predictors.The way the significance levels in this column are approximated using the method of Wood (1989).Li (1992) proposed an alternative method of estimating the Hessian matrix based on quadratic regression, as follows: (1) fit the ols regression of the response on the full second-order model based on all ) predictors; (2) set 0 H to be the ) @ ) matrix whose 2 ¤ ¢ ¡ 7 element is the estimated coefficient of 0 ¡ 0 g .

H
, proceed as with other phd methods.

Multivariate phd
No multivariate extensions of phd have yet been proposed, and so the response must be univariate for any phd method.

Adding other methods to dr
You may skip this section unless you are interested in adding additional dimension reduction methods to dr.The dr function is designed to be flexible and to make adding additional methods to the package as easy as possible.Here are the steps you need to follow: 1. Select a name for your method.For example suppose you select the name "poly" to mean that you will be estimating H by fitting polynomials to the inverse plots each of the predictors versus the response.When you call dr with method="poly", an object of class poly will be created.If your response is a matrix, the object will be of type mpoly, which inherits from poly.
2. You will need to write a function called dr.fit.M.poly that estimates the H matrix for your method.You can model this function after the function dr.fit.M.sir shown in Table 1.The important arguments that are passed to this function include z, which is an © @ ) rotated and cen- tered data matrix with no missing values; y, which is the response vector or matrix, and w, which is the vector of weights, if any.If your method requires other parameters, for example setting a degree of a polynomial, simply add the argument degree=2 to the list of function arguments.This sets the default value of degree equal to 2. The ". . ." argument in the functions allow you to add the degree argument when you call the function dr.Your function must return a list, including the argument M, which is the matrix of interest.
H can be either a square matrix leading to an analysis of eigenvalues and eigenvectors, as in the example for sir, or it can be a rectangular matrix, leading to use of singular values and vectors.All entries in the list will become attributes of the resulting object.For example, if your list is list(z=z,degree=degree), when you create an object like > i1 <-dr(LBM~Ht+Wt+RCC+WCC,method="poly",degree=3) Table 1: The dr.fit.M method for sir.
3. If your method works differently when y is a matrix, write a method called dr.fit.M.mpoly to do the computations for this case.For sir, the only difference between univariate and multivariate responses is in the way the slices are obtained, and the method dr.slices works for either case.As a result, a separate multivariate fit method is not required for sir. .For example, the function dr.fit.y.phdy returns the left-hand side variable from the formula specified in the call to dr, while dr.fit.y.phdres returns the residuals from the regression of the response on the predictors.You may also want to write a summary.dr.poly method if the default summary is not adequate for your needs.Cook (1998) and Yin (2000) discuss permutation tests of dimension that can be used with a dimension reduction method.These are implemented in the function dr.permutation.test.The function requires the name of the object.The number of permutations defaults to 50 and the number of directions defaults to 3. Increasing either can increase the computation time required to obtain the solution.The permutation test results for the example are very similar to the asymptotic results given earlier.

Graphical methods
The function call plot(i1) returns a scatterplot matrix of the response and the principal directions.
The call plot(i1,mark.by.y=TRUE) produces a scatterplot matrix of the principal directions, but with points marked (colored) according to the value of the response; the point marking does not work with Splus.
The function call dr.coplot(i1) returns a plot of the response versus the first principal direction, with a separate panel conditioning on the value of the second principal direction (a coplot).The call dr.coplot(i1,mark.by.y=T), gives a coplot of the first direction versus the second direction conditioning on the third direction and using color to encode the information about the response.
The function rotplot is a generic function that allows looking at a number of static views of a 3D plot.The call > rotplot(dr.directions(m1,1:2),dr.y(m1),number=16) gives 16 views of the 3D plot of the response versus linear combinations of the first two principal directions.

Weights
Weights are generally used in dimension reduction methods to make the resulting weighted sample closer to a normal distribution than the original sample.Cook (1998, Section 8.4) discusses the method that is implemented here.When weights are present, they are used in centering the data and computing the covariance matrix, and they are used in computing the objective matrix H for phd.Weights may be provided by the user with the weights argument.If weights=NULL, the default, no weighting is used.
The function dr.weights is used to estimate weights using the algorithm described by Cook (1998, Sec. 8.4).There are several other arguments that control how the weights are computed, as described below, and on the help page for the function dr.weights.The algorithm works as follows: 1.For an © @ B) data matrix , find estimates and ¡ of the mean and covariance matrix.For this purpose, in R the function cov.rob in the lqs package is used, while in Splus the function covRob in the robust package is used; in either case the needed package will be loaded automatically.If you do not want to use one of these routines, you must rewrite the function robust.center.scale to use your preferred code.In R, the method of computing and ¡ is determined by the argument covmethod.If covmethod="classical", the usual estimator is used for ¡ , but is estimated by medians.If method="mve", the default, or method="mcd", the covariance matrix is estimated by the minimum volume ellipsoid method and the minimum determinant method, respectively.These latter two also return a robust estimate of center.Any tuning parameters for the method to compute the robust estimate of and ¡ can be passed from the call to dr.See the documentation for cov.rob for a description of these additional parameters.All the defaults are sensible, so most users will not need to use these additional parameters.

Compute the matrix
. If the data were normally distributed N2 ¢ ¡ 7 , the rows of would be like a sample from N2 ¡© ¢ S 7 .
3. Obtain a random vector ¡ from the N2 © ¢ £¢ $ S 7 distribution.The parameter sigma=1 is a tuning parameter that can be set in the call to dr, and values near 1 or slightly smaller seem appropriate.Find the row of that is closest to ¡ (the code uses Euclidean distance), and increase a counter for that row by 1.
4. The argument nsamples determines the number of times this last step is repeated; the default is nsamples=10*dim(x) [1] where is the © @ ) data matrix; this number may be too small.5.Return a vector of weights given by the value of the counter divided by nsamples and multiplied by © , so the sum of the weights will be © .

Miscellaneous functions included in dr
The function dr.direction takes two arguments, the object name, and which directions are wanted (for, example, 1:3 returns the first three directions).It returns the matrix ¥¤ , scaled to have unit column length unless the argument norm is false, where ¤ gives the specified eigenvectors.
The function calls dr.x(i1) and dr.y(i1) return the model matrix and the response, respectively.dr.z(i1,center=T,rotate=T) returns the centered and rotated matrix from .You can also use dr.z by explicity providing a matrix in place of the first argument, and if necessary a vector of weights as a second argument.
The routine used for slicing is called dr.slices(y,h) to slice y into h slices.If y has ) columns and h has ) elements, then slicing is done recursively.The first column of y is sliced into h[1] slices.
Within each of these slices, the second column of y is sliced into h[2] slices, giving h[1]*h [2] slices.This process is then repeated for any additional columns.If h is a scalar, then each dimension is sliced into the smallest integer larger than h £# T slices.For example, if ) % § and h% £ , then each dimension has 3 slices for a total of 9.

Bug Reports
Please send bug reports to sandy@stat.umn.edu.9 Splus dr has been tested with Splus versions 5 and 6 under Linux.If you want to use dr with Splus, you must either: (1) set the variable whichengine <-"s6" before you load the file dr.R, or else you must set this value on line 21 of the file.Otherwise, the same code is used for R and Splus.All options available in R seem to work in Splus, except that the argument mark.by.y on plots and coplots and the function markby do not work with Splus.

Description
The function dr implements dimension reduction methods, including SIR, SAVE and pHd.

Arguments
formula a symbolic description of the model to be fit.The details of the model are the same as for lm.
data an optional data frame containing the variables in the model.By default the variables are taken from the environment from which 'dr' is called.
subset an optional vector specifying a subset of observations to be used in the fitting process.
weights an optional vector of weights to be used where appropriate.
na.action a function which indicates what should happen when the data contain 'NA's.The default is 'na.omit,' which will force calculations on a complete subset of cases.
method This character string specifies the method of fitting."sir" specifies sliced inverse regression and "save" specifies sliced average variance estimation."phdy" uses principal hessian directions using the response as suggested by Li, and "phdres" uses the LS residuals as suggested by Cook.Other methods may be added contrasts an optional list.See the 'contrasts.arg' of 'model.matrix.default'.
numdir Maximum number of directions to consider ... additional items that may be required or permitted by some methods.nslices is the number of slices used by sir and save.

Details
The general regression problem studies ¢¡ ¡ ¤£ £ ¦¥ , the conditional distribution of a response ¡ given a set of predictors £ .This function provides methods for estimating the dimension and central subspace of a general regression problem.That is, we want to find a ¢ ¨ § a matrix © such that ¢¡ ¡ ¤£ £ ¦¥ ¢¡ ¡ ¤£ © ¥ £ ¦¥ Both the dimension and the subspace ¡ © ¥ are unknown.These methods make few assumptions.All the methods available in this function estimate the unknowns by study of the inverse problem, ¢¡ £ £ ¡ ¥ .In each, a kernel matrix is estimated such that the column space of should be close to the central subspace.
Eigenanalysis of is then used to estimate the central subspace.Objects created using this function have appropriate print, summary and plot methods.
Weights can be used, essentially to specify the relative frequency of each case in the data.Empirical weights that make the contours of the weighted sample closer to elliptical can be computed using dr.weights.This will usually result in zero weight for some cases.The function will set zero estimated weights to missing.
Several functions are provided that require a dr object as input.dr.permutation.testsuses a permutation test to obtain significance levels for tests of dimension.dr.coplot allows visualizing the results using a coplot of either two selected directions conditioning on a third and using color to mark the response, or the resonse versus one direction, conditioning on a second direction.plot.drprovides the default plot method for dr objects, based on a scatterplot matrix.

Value
Returns an object that inherits from dr (the name of the type is the value of the method argument), with attributes: M A matrix that depends on the method of computing.The column space of M should be close to the central subspace. evalues The eigenvalues of M (or squared singular values if M is not symmetric). evectors The eigenvectors of M (or of M'M if M is not square and symmetric) ordered according to the eigenvalues. numdir The maximum number of directions to be found.The output value of numdir may be smaller than the input value.

ols.coef
Estimated regression coefficients, excluding the intercept, for the (weighted) LS fit.
ols.fit LS fitted values.
slice.info output from sir.slice, used by sir and save.
method the dimension reduction method used.
Other returned values repeat quantities from input.

Details
The basic outline is: (1) Estimate a mean m and covariance matrix S using a possibly robust method; (2) For each iteration, obtain a random vector from N(m,sigma*S).Add 1 to a counter for observation i if the i-th row of the data matrix is closest to the random vector; (3) return as weights the sample faction allocated to each observation.If you set the keyword weights.only to T on the call to dr, then only the list of weights will be returned.

Value
Returns a list of n weights, some of which may be zero.

Description
These routines provide default plotting methods for dimension reduction regression.... arguments passed to plot or coplot.In particular, if the car library is available, the argument panel=panel.carwill add smoothers to a coplot.

Value
Produces a scatterplot matrix (plot) or coplot (dr.coplot) of the specified directions in an dimension reduction regression

Description
This function draws several 2D views of a 3D plot, sort of like a spinning plot.

Arguments
x a matrix with 2 columns giving the horizontal axes of the full 3D plot.
y the vertical axis of the 3D plot.
theta a list of rotation angles ... additional arguments passed to coplot

Value
returns a graph object.

Value
dr.
x returns an § ¢ matrix of terms, excluding the intercept, constructed from the dimension reduction regression object.dr.y returns the response.dr.z returns a possibly centered and scaled version of x.

See Also
dr Examples data(ais) attach(ais) m1 <-dr(LBM~Ht+Wt+RCC+WCC) dr.x(m1) dr.y(m1) The following documentation is from the package lqs by Brian Ripley, and is included here for convenience:

Description
Compute a multivariate location and scale estimate with a high breakdown point -this can be thought of as estimating the mean and covariance of the good part of the data.cov.mve and cov.mcd are compatibility wrappers.

Arguments
x a matrix or data frame.
cor should the returned result include a correlation matrix?
quantile.used the minimum number of the data points regarded as good points.
method the method to be used -minimum volume ellipsoid, minimum covariance determinant or classical product-moment.Using cov.mve or cov.mcd forces mve or mcd respectively.
nsamp the number of samples or "best" or "exact" or "sample".If "sample" the number chosen is min(5*p, 3000), taken from Rousseeuw and Hubert (1997).If "best" exhaustive enumeration is done up to 5000 samples: if "exact" exhaustive enumeration will be attempted however many samples are needed.
seed the seed to be used for random sampling: see RNGkind.The current value of .Random.seed will be preserved if it is set.

Details
For method "mve", an approximate search is made of a subset of size quantile.usedwith an enclosing ellipsoid of smallest volume; in method "mcd" it is the volume of the Gaussian confidence ellipsoid, equivalently the determinant of the classical covariance matrix, that is minimized.The mean of the subset provides a first estimate of the location, and the rescaled covariance matrix a first estimate of scatter.The Mahalanobis distances of all the points from the location estimate for this covariance matrix are calculated, and those points within the 97.5% point under Gaussian assumptions are declared to be good.The final estimates are the mean and rescaled covariance of the good points.
The rescaling is by the appropriate percentile under Gaussian data; in addition the first covariance matrix has an ad hoc finite-sample correction given by Marazzi.
For method "mve" the search is made over ellipsoids determined by the covariance matrix of p of the data points.For method "mcd" an additional improvement step suggested by Rousseeuw and van Driessen (1997) is used, in which once a subset of size quantile.usedis selected, an ellipsoid based on its covariance is tested (as this will have no larger a determinant, and may be smaller).

Value
estimates of vectors in b 2 C a7 , the residual methods are more suitable for tests of dimension.SeeCook (1998, Chapter 12)for details.Output for phd is again similar to sir, except for the tests.Here is the output for the same setup as before, but for method phdres:
These functions estimate weights to apply to the rows of a data matrix to make the resulting weighted matrix as close to multivariate normality as possible.This probably will not work with Splus.If classical is selected, the usual estimate of the covariance matrix is used, but the center is the medians, not the means.
Arguments objecta dimension reduction regression object name, or an n by p matrix sigma A tuning parameter, with default 1, usually in the range .2 to 1.0 covmethod covmethod is passed as the argument method to the function cov.rob in the required package lqs.The choices are "classical", "mve" and "mcd".
(only is cor = TRUE) the estimate of the correlation matrix.singmessagegiving number of singular samples out of total crit the value of the criterion on log scale.For MCD this is the determinant, and for MVE it is proportional to the volume.bestthe subset used.For MVE the best sample, for MCD the best set of size quantile.used. cor