# In PCA, is there a systematic way of dropping variables to maximise the segregation of two populations?

I am trying to investigate using principal component analysis whether it is possible to guess with good confidence from which population ("Aurignacian" or "Gravettian") a new datapoint came from. A datapoint is described by 28 variables, most of which are relative frequencies of archaeological artefacts. The remaining variables are computed as ratios of other variables.

Using all variables, the populations partially segregate (subplot (a)), but there is still some overlap in their distribution (90% t-distribution prediction ellipses, though I am not sure I can assume normal distribution of populations). I thought it was therefore not possible to predict with good confidence the origin of a new datapoint:

Removing one variable (r-BEs), the overlap becomes much more important, (subplots (d), (e), and (f)), as the populations do not segregate in any paired PCA plots: 1-2, 3-4, ..., 25-26, and 1-27. I took this to mean r-BEs was essential for separating the two populations, because I thought that taken together, these PCA plots represent 100% of the "information" (variance) in the dataset.

I was hence extremely surprised to notice that the populations actually did segregate almost completely if I dropped all but a handful of variables:

Why is this pattern not visible when I perform a PCA on all variables? With 28 variables, there are 268,435,427 ways of dropping a bunch of them. How may one find those that will maximise population segregation and best allow guessing the population of origin of new datapoints? More generally, is there a systematic way of finding "hidden" patterns like these?

EDIT: Per amoeba's request, here are the plots when the PCs are scaled. The pattern is clearer. (I realise I'm being naughty by continuing to knock out variables, but the pattern this time resists to the knock-out of r-BEs, implying the "hidden" pattern is picked up by the scaling):

• +1 for very pretty figures. Are you doing your PCA on the covariance matrix or on the correlation matrix, i.e. are all the variables normalized or not? Can it be that you did not normalize your variables, and G1/2/7/8/9 variables have much less variance than the variables that the "full" PCA is primarily picking up? May 21, 2017 at 20:35
• @amoeba Thanks for the compliment :) The percentages are not normalised. I tried scaling just now and it does seem to help. I had decided against scaling at first because I thought it might add noise: very rare artefacts with low percentages might be amplified, and with them the effect of chance. It's still surprising, however, that the G1/2/7/8/9's effect cannot be seen even in the least important variables. May 21, 2017 at 21:04
• I mean in the least important PCs May 21, 2017 at 21:10
• I don't find it particularly surprising. What you seem to have here, is a group of variables (Gs) that are all strongly positively correlated between each other (see second row of second figure) and can well predict population identity. However, these variables have less variance than some of the other variables, and are not picked up when you do PCA on all variables. At the same time, it so happens that they do not get into one single component of PCA-28 (they get spread out into multiple PCs), that's why you don't see this separation on any of the pairwise scatter plots. May 21, 2017 at 21:14
• Thanks. Looks nice. The arrows should not be of the same length; after normalization they are of the same length in the full 28-dimensional space, but once you project everything on the 2D, they will have different lengths. The longer the arrow -- the stronger the contribution of the corresponding variable to the PC1/2. May 22, 2017 at 20:59

Principal Components (PCs) are based on the variances of the predictor variables/features. There is no assurance that the most highly variable features will be those that are most highly related to your classification. That is one possible explanation for your results. Also, when you limit yourself to projections onto 2 PCs at a time as you do in your plots, you might be missing better separations that exist in higher-dimensional patterns.

As you are already incorporating your predictors as linear combinations in your PC plots, you might consider setting this up as a logistic or multinomial regression model. With only 2 classes (e.g., "Aurignacian" versus "Gravettian"), a logistic regression describes the probability of class membership as a function of linear combinations of the predictor variables. A multinomial regression generalizes to more than one class.

These approaches provide important flexibility with respect both to the outcome/classification variable and to the predictors. In terms of the classification outcome, you model the probability of class membership rather than making an irrevocable all-or-none choice in the model itself. Thus you can for example allow for different weights for different types of classification errors based on the same logistic/multinomial model.

Particularly when you start removing predictor variables from a model (as you were doing in your examples), there is a danger that the final model will become too dependent on the particular data sample at hand. In terms of predictor variables in logistic or multinomial regression, you can use standard penalization methods like LASSO or ridge regression to potentially improve the performance of your model on new data samples. A ridge-regression logistic or multinomial model is close to what you seem to be trying to accomplish in your examples. It is fundamentally based on principal components of the feature set, but it weights the PCs in terms of their relations to the classifications rather than by the fractions of feature-set variance that they include.

• +1. It is worth emphasizing that using Lasso penalty will provide sparse (and hence more interpretable) solutions, which is what OP seems to be after here. May 21, 2017 at 20:37
• Many thanks for your suggestions. I had indeed noticed that one can miss a PC1-PC3 pattern if looking only at sequential PCs (PC1-PC2, PC3-PC4). What I am looking for is a statement of the type "For every future data point, I can correctly assign membership to class A with XX% confidence if it is from class A, and to class B with YY% confidence if it is from class B". Would logistic regression and ridge-regression allow me to get that? May 21, 2017 at 21:32
• Logistic regression provides statements like: "Given these predictor values, the probability is XX that this case is in Class A." That's often the most useful for prediction, with a class-membership cutoff typically of 1/2. A statement like the one you seek also requires knowledge of the distributions of predictor values among cases of Classes A and B. As your data sample is an estimate of those distributions, you can use cross-validation or bootstrapping to help generate such a statement. ISLR provides helpful details.
– EdM
May 22, 2017 at 14:55
• @TheThunderChimp your situation isn't different from many classification problems, in which the class helps determine the values of the "predictor" variables and you are trying to infer the underlying class. A LASSO or elastic net approach to the logistic regression should help with correlated variables. Perfect separation often is seen in logistic regression, not necessarily related to variable correlations. Follow the hauck-donner-effect tag on this site for advice; this answer may be particularly helpful.
– EdM
May 22, 2017 at 21:42
• @TheThunderChimp I don't know how well that approach would generalize to more than 2 dimensions. This page might be even more useful than what I linked in an earlier comment.
– EdM
May 23, 2017 at 21:27