What is the best way to identify outliers in multivariate data? Suppose I have a large set of multivariate data with at least three variables. How can I find the outliers? Pairwise scatterplots won't work as it is possible for an outlier to exist in 3 dimensions that is not an outlier in any of the 2 dimensional subspaces.
I am not thinking of a regression problem, but of true multivariate data. So answers involving robust regression or computing leverage are not helpful.
One possibility would be to compute the principal component scores and look for an outlier in the bivariate scatterplot of the first two scores. Would that be guaranteed to work? Are there better approaches?
 A: It may be an overshoot, but you may train an unsupervised Random Forest on the data and use the object proximity measure to detect outliers. More details here.
A: For moderate dimensions, like 3, then some sort of kernel cross-validation technique as suggested elsewhere seems reasonable and is the best I can come up with.
For higher dimensions, I'm not sure that the problem is solvable; it lands pretty squarely into 'curse-of-dimensionality' territory.  The issue is that distance functions tend to converge to very large values very quickly as you increase dimensionality, including distances derived from distributions.  If you're defining an outlier as "a point with a comparatively large distance function relative to the others", and all your distance functions are beginning to converge because you're in a high-dimensional space, well, you're in trouble.
Without some sort of distributional assumption that will let you turn it into a probabilistic classification problem, or at least some rotation that lets you separate your space into "noise dimensions" and "informative dimensions", I think that the geometry of high-dimensional spaces is going to prohibit any easy -- or at least robust -- identification of outliers.
A: I'm not sure what you mean when you say you aren't thinking of a regression problem but of "true multivariate data".  My initial response would be to calculate the Mahalanobis distance since it doesn't require that you specify a particular IV or DV, but at its core (as far as I understand it) it is related to a leverage statistic.
A: I'm not aware that anyone is doing this, but I generally like to try dimensionality reduction when I have a problem like this.  You might look into a method from manifold learning or non-linear dimensionality reduction.
An example would be a Kohonen map.  A good reference for R is "Self- and Super-organizing Maps in R: The kohonen Package".
A: Have a look at the mvoutlier package which relies on ordered robust mahalanobis distances, as suggested by @drknexus.
A: I think Robin Girard's answer would work pretty well for 3 and possibly 4 dimensions, but the curse of dimensionality would prevent it working beyond that. However, his suggestion led me to a related approach which is to apply the cross-validated kernel density estimate to the first three principal component scores. Then a very high-dimensional data set can still be handled ok.
In summary, for i=1 to n


*

*Compute a density estimate of the first three principal component scores obtained from the data set without Xi. 

*Calculate the likelihood of Xi for the density estimated in step 1.
call it Li. 


end for
Sort the Li (for i=1,..,n) and the outliers are those with likelihood below some threshold. I'm not sure what would be a good threshold -- I'll leave that for whoever writes the paper on this! One possibility is to do a boxplot of the log(Li) values and see what outliers are detected at the negative end.
A: You can find a pedagogical summary of the various methods available in (1)
For some --recent-- numerical comparisons of the various methods listed there, you can check 
(2) and (3).
there are many older (and less exhaustive) numerical comparisons, typically found in books. You will find one on pages 142-143 of (4), for example.
Note that all the methods discussed here have an open source R implementation, mainly through the rrcov
 package.


*

*(1) P. Rousseeuw and M. Hubert (2013)  High-Breakdown Estimators of
Multivariate Location and Scatter.

*(2) M. Hubert, P. Rousseeuw, K. Vakili (2013).
Shape bias of robust covariance estimators: an empirical study.
Statistical Papers.

*(3) K. Vakili and E. Schmitt (2014). Finding multivariate outliers with FastPCS. Computational Statistics & Data Analysis.

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

A: I would do some sort of "leave one out testing algorithm" (n is the number of data):
for i=1 to n


*

*compute a density estimation of the data set obtained by throwing $X_i$ away. (This density estimate should be done with some assumption if the dimension is high, for example, a gaussian assumption for which the density estimate is easy: mean and covariance)

*Calculate the likelihood of $X_i$ for the density estimated in step 1. call it $L_i$. 


end for 
sort the $L_i$ (for i=1,..,n) and use a multiple hypothesis testing procedure to say which are not good ... 
This will work if n is sufficiently large... you can also use "leave k out strategy" which can be more relevent when you have "groups" of outliers ...
A: I didn't see anybody mention influence functions.  I first saw this idea in Gnanadesikan's multivariate book.
In one dimension an outlier is either an extremely large or an extremely small value.  In multivariate analysis it is an observation removed from the bulk of the data.  But what metric should we use to define extreme for the outlier?  There are many choices. The Mahalanobis distance is just one.  I think that looking for every type of outlier is futile and counterproductive.  I would ask why do you care about the outlier?  In estimating a mean they can have a great deal of influence on that estimate.  Robust estimators downweight and accommodate outliers but they do not formally test for them.  Now in regression, the outliers--like leverage points--could have large effects on the slope parameters in the model.  With bivariate data they can unduly influence the estimated correlation coefficient and in three or more dimensions the multiple correlation coefficient.
Influence functions were introduced by Hampel as a tool in robust estimation and Mallows wrote a nice unpublished paper advocating their use.  The influence function is a function of the point you are at in n-dimensional space and the parameter.  It essentially measures the difference between the parameter estimate with the point in the calculation and with the point left out. Rather than go to the trouble of doing the calculation of the two estimates and taking the difference, often you can derive a formula for it.  Then the contours of constant influence tell you the direction that is extreme with respect to the estimate of this parameter and hence tell you where in the n-dimensional space to look for the outlier.
For more you can look at my 1983 paper in the American Journal of Mathematical and Management Sciences titled "The influence function and its application to data validation."   In data validation we wanted to look for outliers that affected the intended use of the data.  My feeling is that you should direct your attention to outliers that greatly affect the parameters you are interested in estimating and not care so much about others that don't. 
A: You can find candidates for "outliers" among the support points of the minimum volume bounding ellipsoid.  (Efficient algorithms to find these points in fairly high dimensions, both exactly and approximately, were invented in a spate of papers in the 1970's because this problem is intimately connected with a question in experimental design.)
A: I novel approach I saw was by IT Jolliffe Principal Components Analysis.  You run a PCA on your data (Note: PCA can be quite a useful data exploration tool in its own right), but instead of looking at the first few Principal Components (PCs), you plot the last few PCs.  These PCs are the linear relationships between your variables with the smallest variance possible.  Thus they detect "exact" or close to exact multivariate relationships in your data.
A plot of the PC scores for the last PC will show outliers not easily detectable by looking individually at each variable.  One example is for height and weight - some who has "above average" height and "below average" weight would be detected by the last PC of height and weight (assuming these are positively correlated), even if their height and weight were not "extreme" individually (e.g. someone who was 180cm and 60kg).
A: One of the above answers touched in mahalanobis distances.... perhaps anpther step further and calculating simultaneous confidence intervals would help detect outliers!
