Does it make sense to do PCA before a Tree-Boosting model?

Question

All I could find about this was this answer, which verifies my initial intuition that

Decision trees, by virtue of doing recursive splitting of your samples, with splits being based on a single variable, can only generate decision boundaries parallel to the axes of your co-ordinate system. So by rotating the data to directions of maximum variance/diagonalizing your covariance matrix as best you can, it might be easier to put decision boundaries between your class distributions

but starts with

Disclaimer: I'm usually wrong at things.

I want some more validated answer regarding this:

1. Will this improve accuracy?
2. Will this improve time or steps to convergence?
3. Does it make sense in general?

---

Edit:

Please notice I said nothing about dimension reduction using PCA.
If there is a difference in answer depending on weather PCA reduces dimension or not, please also explain why that is.

Answer 1

In general, good features will improve the performance of any model, and should require fewer steps / result in faster convergence. One nice example of this is whether you want to use the distance from the hole for modeling the golf putting probability of success, or whether you design a new feature based on the geometry (hole size, ball size, tolerance for deviation from the optimal angle).

Whether PCA helps depends on whether resulting representation is useful for the modeling problem. That is not in general clear - e.g. if you are predicting some decision that happens in real-life and people use a threshold on a single one of your features, then dimension reducing with PCA can only hurt performance (probably not a situation in which you'd want to use Machine Learning, really). It might also be that some dimension reduction is desirable, but the "ideal" representation would use some kind of very non-linear transformation (as e.g. UMAP might be able to achieve).

My intuition would be that PCA might help the most when there's a huge number of initial features (or a really complex signal), especially if that includes ones that are slightly different, but correlated versions of the same types of features similar

Update: Without dimensionality reduction, here is an example, where we can immediately predict that PCA will not help (first column of figures; orange = outcome 1, green = 0) and where PCA will help (second column; blue = 1, red = 0). We should expect the scenario where xgboost (and other similar tree based methods like random forrest, LightGBM etc.) simply needs to learn a cut-off along a single dimension to be the best scenario for it, while the cases where it needs to learn a diagonal (or even more complex) boundary should be harder for it. In the first column of plots PCA actually takes us away from the right feature representation, while in the second column of plots it takes us to a better representation.

The R code for this is below, but in the first case we get

• Column 1 without PCA: train-logloss:0.011828+0.000477 test-logloss:0.016520+0.008104 (best iteration after 185 iterations)
• Column 1 with PCA: train-logloss:0.013675+0.000435 test-logloss:0.025791+0.011341 (best iteration after 313 iterations)
• Column 2 without PCA: train-logloss:0.022771+0.000451 test-logloss:0.049508+0.009839 (best iteration after 467 iterations)
• Column 2 with PCA: train-logloss:0.019837+0.000593 test-logloss:0.026960+0.009282 (best iteration after 131 iterations)

So, in one case we need fewer iterations and get better accuracy without PCA, while in the other case we need fewer iterations and get better accuracy when using the PCA features.

Of course, you could argue that you should just use both sets of features, but - again - there's no way of knowing that either leaving the features as is or using PCA will be the best way to transform the features. And sure, perhaps there's better hyperparameter choices that would minimize these differences. And sure, the CV-performance after early stopping might be a tiny bit biased, but the general pattern is clear enough that we do not need to do anything more complicated to see it, I think

library(tidyverse)
library(patchwork)
library(xgboost)

set.seed(1234)
simulated = tibble(`original x`=rnorm(5000, mean=0, sd=5),
                   `original y`=rnorm(5000, mean=0, sd=20),
                   x = `original x` * cos(2*pi/8) - `original y` * sin(2*pi/8),
                   y = `original x` * sin(2*pi/8) + `original y` * cos(2*pi/8),
                   outcome1 = (x+rnorm(5000,0,0.25)>5)*1L,
                   outcome2 = (`original x`+rnorm(5000,0,0.25)>5)*1L)

pca = prcomp(simulated %>% dplyr::select(x,y), center = TRUE,scale. = TRUE)

simulated = simulated %>%
  mutate(`PCA x` = pca$x[,1], 
         `PCA y` = pca$x[,2])

p1 = simulated %>%
  ggplot(aes(x=x, y=y, col=factor(outcome1))) + 
  geom_point(alpha=0.2) +
  theme_bw(base_size=18) +
  theme(legend.position="none") +
  scale_color_brewer(palette="Dark2")

p2 = simulated %>%
  ggplot(aes(x=`original x`, y=`original y`, col=factor(outcome1))) + 
  geom_point(alpha=0.2) +
  theme_bw(base_size=18) +
  theme(legend.position="none") +
  scale_color_brewer(palette="Dark2")

p3 =simulated %>%
  ggplot(aes(x=x, y=y, col=factor(outcome2))) + 
  geom_point(alpha=0.2) +
  theme_bw(base_size=18) +
  theme(legend.position="none") +
  scale_color_brewer(palette="Set1")

p4 = simulated %>%
  ggplot(aes(x=`original x`, y=`original y`, col=factor(outcome2))) + 
  geom_point(alpha=0.2) +
  theme_bw(base_size=18) +
  theme(legend.position="none") +
  scale_color_brewer(palette="Set1")

p5 = simulated %>%
  ggplot(aes(x=`PCA x`, y=`PCA y`, col=factor(outcome1))) + 
  geom_point(alpha=0.2) +
  theme_bw(base_size=18) +
  theme(legend.position="none") +
  scale_color_brewer(palette="Dark2")

p6 =simulated %>%
  ggplot(aes(x=`PCA x`, y=`PCA y`, col=factor(outcome2))) + 
  geom_point(alpha=0.2) +
  theme_bw(base_size=18) +
  theme(legend.position="none") +
  scale_color_brewer(palette="Set1")

#(p2 + p4) / 
(p1 + p3) / (p5 + p6)

xgb.cv(
  params = list(booster="gbtree", 
                objective="binary:logistic",
                eta=0.05, 
                max_depth=128,
                min_child_weight=4, 
                subsample=0.65,
                colsample_bytree=1),
  data = data.matrix(simulated %>% dplyr::select(x,y)),
  nrounds=1000,
  nfold=10,
  label = simulated$outcome1,
  showsd = TRUE,
  metrics = "logloss",
  early_stopping_rounds = 20,
  print_every_n=100)


xgb.cv(
  params = list(booster="gbtree", 
                objective="binary:logistic",
                eta=0.05, 
                max_depth=128,
                min_child_weight=4, 
                subsample=0.65,
                colsample_bytree=1),
  data = data.matrix(simulated %>% dplyr::select(x,y)),
  nrounds=1000,
  nfold=10,
  label = simulated$outcome1,
  showsd = TRUE,
  metrics = "logloss",
  early_stopping_rounds = 20,
  print_every_n=100)

xgb.cv(
  params = list(booster="gbtree", 
                objective="binary:logistic",
                eta=0.05, 
                max_depth=128,
                min_child_weight=4, 
                subsample=0.65,
                colsample_bytree=1),
  data = data.matrix(simulated %>% dplyr::select(`PCA x`,`PCA y`)),
  nrounds=1000,
  nfold=10,
  label = simulated$outcome1,
  showsd = TRUE,
  metrics = "logloss",
  early_stopping_rounds = 20,
  print_every_n=100)


xgb.cv(
  params = list(booster="gbtree", 
                objective="binary:logistic",
                eta=0.05, 
                max_depth=128,
                min_child_weight=4, 
                subsample=0.65,
                colsample_bytree=1),
  data = data.matrix(simulated %>% dplyr::select(x,y)),
  nrounds=1000,
  nfold=10,
  label = simulated$outcome2,
  showsd = TRUE,
  metrics = "logloss",
  early_stopping_rounds = 20,
  print_every_n=100)

xgb.cv(
  params = list(booster="gbtree", 
                objective="binary:logistic",
                eta=0.05, 
                max_depth=128,
                min_child_weight=4, 
                subsample=0.65,
                colsample_bytree=1),
  data = data.matrix(simulated %>% dplyr::select(`PCA x`,`PCA y`)),
  nrounds=1000,
  nfold=10,
  label = simulated$outcome2,
  showsd = TRUE,
  metrics = "logloss",
  early_stopping_rounds = 20,
  print_every_n=100)