GPU Acceleration in GBDT

[BoanSy]

Dear authors,

I have to do 4 dimensions reweighting for 10M events using GBDT. With the CPU only, the process takes a full day.

My question is : Is it possible to write or is there already an interface for GBDT GPU acceleration?

Many thanks for your reply.

Best regards,
BaSy

[arogozhnikov]

Hi @BoanSy,

There is no GPU acceleration. That could be an interesting problem to tackle for someone aching for well-defined pieces of code to be optimized.

I should say, there is a lot to be optimized on CPU - code is pure-numpy, thus single-threaded, so as the first attempt
**Update: ** I was confused here. Pure-numpy refers to losses in uniform boosting.

1. profile the code to find out if it is tree building too slow or loss computation (by e.g. training with just 5-10 trees).
2. if trees - significant optimization gains possible by using feature binning
3. if loss is slow - worth replacing numpy in computations with torch/jax and see if those can use multiprocessing for the same operations.

[BoanSy]

Hi @arogozhnikov, Thank you for your careful reply.

I think my previous question in an overly general way, and I'd like to try to specify it here.

Firstly, what I am trying to achieve is to modify GBreweighter.fit while keeping writing of the tree and the creation of the dataset as unchanged as possible.

Also I think although each decision tree requires inputs from the training results of the previous one, resulting in the impossibility of parallel arithmetic. I am wondering if it is possible to apply automatic differentiation when solving for residuals, for example as implemented in the TF or JAX.

These are my preliminary thoughts and I welcome any comments you may have.

[arogozhnikov]

modify GBreweighter.fit while keeping writing of the tree and the creation of the dataset as unchanged as possible.

GBReweighter may only be stuck on tree building. Target computation should be fast.
I am surprised by your runtime, becuase 15M with 4 variables and default setting takes ~18 minutes on my laptop.

I am wondering if it is possible to apply automatic differentiation when solving for residuals, for example as implemented in the TF or JAX.

Automatic differentiation is not helpful in GB (in general, and in GBReweighter in particular).

My recommendation is:

1. Check how long it takes default settings, and compare to some other machine.
2. Tune arguments, e.g. gb_args = {‘subsample’: 0.1, ‘max_features’: 0.75}, n_estimators=100, learning_rate=0.1, min_samples_leaf=200, use small trees.

If you're interested in speeding up GBReweighter (e.g. as a research direction), I'd recommend looking at how HistGradientBoostingRegressor works and deriving that idea, because binning (aka histogram-ming) provides a dramatic speed up on large datasets.

[BoanSy]

Sorry for taking so long to respond to you, for it took me a while to do the tests as you suggestions.

My test run as follows:

1: Dataset size: 800k events for original and target one.

2: Weight: I inject weight factor containing positive and negative values to BOTH dataset.

3; GBDT setting: n_estimators=200; learning_rate=0.1; max_depth=3; min_samples_leaf=1000; gb_args={'subsample': 0.7}. Ohters are set by default.

The time spent on the fit is : 472s.

Meanwhile, if use HistGradientBoostingRegressor, the time is 1-3s approximately, but the stability of the result is far less than default, resulting from that the parameters are more difficult to tune, e.g. max_bins value.

In addtion, I try to combine the lightgbm or xgboost with hep_ml. It has been done since code is pure-numpy. However, none of them support negative weight factor inputs, which doesn't help me in my physical analysis, for which I need to extract my signal events from the data via SPlot.

So I'm asking if you have any other good suggestions.

BaSy

[saur-m]

Hello @arogozhnikov, @BoanSy,

let me jump to this thread as I am working on this topic together with @BoanSy.

At this moment we are indeed mainly interested in enhancing the reweighter library to make it faster when working with larger and larger samples of O(10^6-10^8) events and also with several variables/dimensions to be reweighted. So as the code is a pure numpy so far, my idea was to modify original numpy-only code with TensorFlow or PyTorch and then be able to access their built-in GPU libraries in order to offload problematic computation and speed up the code.
Another alternative I am looking lately, would be to use numba to check if it could get additional performance from the numpy code.

Some time ago I spent quite some time linking GBReweighter logic with Cython and some C++ libraries, but it become a huge unsustainable mess depending on too many 3rd-party code. So, I would say that is not the way how to make this code faster.

Additionally I am thinking if there is a way how to optimise memory usage (a very common problem with python code), as now we are reweighting samples with a size of a few up to tens of gigabytes. And with the with LHC Run 3 data we can only expect that these numbers will go even higher. As hep_ml, and especially GBReweighter, become very important tool for many analyses at LHCb (and let me add it is a very useful tool!) I would be very happy for any ideas how to speed up and modernise GBReweighter (even other classes in hep_ml later on). Ideally, after making it working for our case, we would like to include these changes even in the official repository, if that would be fine with you.

Thanks.

[arogozhnikov]

and let me add it is a very useful tool!

Thanks, I'm very glad to hear that.

original numpy-only code with TensorFlow or PyTorch and then be able to access their built-in GPU libraries in order to offload problematic computation and speed up the code

Sorry, I see where I produced that confusion in the thread.
GBReweighter has pure-numpy code and tree-building code (just invokes sklearn).

• pure-numpy part is easy to optimize with TF/torch/jax, but I doubt it contributes significantly.
• tree building will take > 95% of time, and TF/torch/jax can't speed this up unfortunately.

Meanwhile, if use HistGradientBoostingRegressor, the time is 1-3s approximately, but the stability of the result is far less than default, resulting from that the parameters are more difficult to tune, e.g. max_bins value.

Just leave default 255, that's memory efficient, fast, and I doubt that at any point you'll need to be more precise than 0.5 / 255 ~ 0.001 of quantile. Very few cases when one needs that, and I can't recall any in HEP.

Quite seriously, I fully recommend you folks to implement that Hist Trick - it is simple and dramatically efficient (I used it in a number of projects, and you see speed up in sklearn is dramatic). It will get you farther than e.g. just moving computations to H100.

You'll need a binner:
https://github.com/scikit-learn/scikit-learn/blob/2a2772a87b6c772dc3b8292bcffb990ce27515a8/sklearn/ensemble/_hist_gradient_boosting/binning.py#L70

and a TreeGrower:
https://github.com/scikit-learn/scikit-learn/blob/7f9bad99d6e0a3e8ddf92a7e5561245224dab102/sklearn/ensemble/_hist_gradient_boosting/grower.py#L138

If you want a minimal reference code (sklearn implementation covers many additional cases of no interest to you), see my implementation in python2 (obviously, unmaintained). Will take some time to read as it uses bit twiddling.