I am dealing with fairly highly dimensional data (spectroscopic data), which I preprocess with PCA, then apply the chosen learning algorithm. When I use linear regression, using the Caret package, I get an $R^2=0.83$ on my test set of 300 PCA components however, when I use random forests on the same set:

modelFit <- train(train$age~.,data=trainPC[,1:300],method = "rf")

I get worse performance, with $R^2$ much lower on the test set. RMSE is also not as good for random forests vs linear regression.

For clarity - I am using Caret's in-built cross validation, and am evaluating the performance on an independent test set.

I get similar results when I use xgboost, neural networks etc. I know that different machine learning algorithms can perform better in different circumstances, but I am nonetheless surprised that linear regression appears to be performing best, and makes me suspect that I am approaching things incorrectly, or making a basic error with using the Caret package.

I have also tried using PLS (the main feature selection method used for my particular example) and get an $R^2=0.83$. The literature suggests that neural nets/non-linear-methods should outperform the linear methods I am using, but at the moment, I cannot reproduce this.

Does anyone have any suggests/tips on this? I also attach the data set I am using (I am trying to predict age) - data link

I also attach my current R script that I am using to do the learning - any help would be much appreciated! - R script

  • The purpose of principal components regression is that you use a smaller number of components to avoid overfitting. A plot of $k$ vs. AIC from using the first $k$ components usually helps in choosing $k$. – Frank Harrell Sep 30 '15 at 12:29
  • Another thing, if you're going to be using random forests then you might think about using the ranger package--its implementation is much faster than randomForest. – Steve S Sep 30 '15 at 14:41

Definitely check on how $R^2$ is being evaluated and whether it is uniform across the different algorithms. Beyond that, a few thoughts occur to me:

  • If your features have a smooth, nearly linear dependence on the covariates, then linear regression will model the dependence better than random forests, which will basically approximate a linear curve with an ugly irregular step function. If the dependence is multivariate linear and smooth, with $v$ significant covariates producing the dependence, the fit performance of random forests can be expected to get worse and worse for larger and larger $v$. RF has a much greater ability than a single decision tree to model linearity, since we are adding tree predictions together - but still, it's just not very efficient to approximate a high-dimensional linear relationship with a series of step functions. I think this is the most likely theoretical explanation for RF underperforming linear regression.

  • ...that said, the step functions will get nicer as you add more trees to the random forest. So you may want to consider increasing the number of trees a lot - as high as you are willing to wait for. Maybe see if caret will allow you to track the $R^2$ as the number of trees increases? I would also expect boosted trees to potentially work well here, but you may have to fit a ton of trees and fiddle a lot with the parameters. There's no guarantee though - if this is a really high-dimensional smooth relationship, it could take an unfathomably large number of trees to approximate it well.

  • Preprocessing with PCA is not ideal for high-dimensional learners like random forests. For best results, they should generally get direct access to the full dataset. If the raw dataset has sparsity of main effects, this nice property could easily be destroyed by PCA. This would make even worse the problem of my first bullet, where RF is trying to approximate a high-dimensional linear relationship by step functions.

  • Advanced machine learning algorithms are favored for their generalization performance, not their ability to closely fit training data and produce a high $R^2$. You are performing an interesting diagnostic and that is worth doing, but be aware that this diagnostic is not a good measure of final performance.

  • Random Forests are not the only useful algorithm out there for high-dimensional data. You could also try elastic net regression either on your PCA processed data or on the raw dataset. Elastic net gives you the ability to model smooth relationships, as well as the ability to perform high-dimensional feature selection using the $\ell_1$ regression coefficient penalty.

  • 2
    +1 for the 4th point, I would rather spend time on measuring the prediction accuracy of the model than measuring the $R^2$ – Metariat Sep 30 '15 at 12:08
  • @Paul Thanks very much for your well-considered answer. I agree that $R^2$ isn't necessarily a good measure of performance - I could have replaced all my references to it with RMSE, and all that is said in my points above would still stand. I am now trying using a random forest on my data without doing any pre-processing. Would you also recommend using neural nets without doing PCA first? Again, many thanks! – ben18785 Sep 30 '15 at 12:58
  • I would check out Steve S's answer before doing anything else. Just to clarify, are you measuring performance by predicting on the training data, or an independent set of test data? Are you using cross-validation? – Paul - MM DataSci Sep 30 '15 at 13:07
  • Hi, I will do. I am measuring performance on an independent test set. I am using Caret's in-built cross validation - sorry for the lack of clarity here! – ben18785 Sep 30 '15 at 13:18

Check out caret's function findLinearCombos and run it on your data. The object returned is a list--the second element is a vector of indices which can safely be removed from your data set since they are linear combinations of other columns (which, now that I think about it makes total sense because you only had 500 observations and, as a result, your model was totally over-specified).

When you pass in your data frame the resulting vector has a thousand elements. Basically, you can go right ahead and delete columns 501 to 1500.

Do that from the get-go and try re-building those models.

A couple quick side notes:

  • If you include the 'data' argument, then only refer to the columns by their names in the formula (i.e. write age ~. instead of data$age ~.).
  • If you plan on doing PCA then you might want to consider setting it up via the preProcess argument in the train function (you can tweak the threshold in trainControl function). (Since you're already using carat anyway you might as well have it handle your pre-processing, too--it gives you one less thing to worry about).

By the way, try checking out this code.

prepped <- preProcess(x=df[, -1], method=c('center', 'scale', 'pca'), 
pca <- predict(prepped, df[, -1])
pca$age <- df$age
fit <- ln(age ~ ., data=pca)

num.unique <- length(unique(df$age))

cols <- rainbow(num.unique)[factor(df$age)]

plot(fit, col=cols)
  • Very good sleuthing - multicollinearity that bad could really hurt RF in training. This should be higher than my answer! – Paul - MM DataSci Sep 30 '15 at 12:40
  • @Paul: Ha! I barely even finished training a model when I saw your answer--I couldn't believe you had written so much so quickly! I guess it's sort of like point #4: Your advice generalizes well; I just happened to notice a quirk of the data set. – Steve S Sep 30 '15 at 13:40
  • Thanks very much for your thorough answer. Whilst that plot looks like it was generated through some weird transformation of the data, it is actually due to the fact that the age variable only takes on a discrete set of values (13 unique values). If the actual age is 1 day, then the predicted value - 1 day will always lie on a straight line with a slope of minus 1. Does that make sense? – ben18785 Sep 30 '15 at 13:54
  • @SteveS - I have now added a larger data set (the actual one rather than the smaller test set I was using before), to the problem above. I think that findLinearCombos was only finding perfect collinearity due to the number of columns being greater than the number of predictors. – ben18785 Sep 30 '15 at 14:11
  • @ben18785: Ok, that makes sense. I don't really follow that last comment but I'm also really tired... If you run the Breusch-Pagan test (bptest from lmtest) then the resulting p-value is pretty minuscule. – Steve S Sep 30 '15 at 14:28

Your Answer


By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Not the answer you're looking for? Browse other questions tagged or ask your own question.