In this study, we apply our MACE potential to predict bond dissociation energy (BDE), similarly to how the AM1 method is used in the CYP metabolism prediction described by Tyzack et al. We focus on the hydrogen abstraction reactions facilitated by CYPs, whose rate-limiting step is characterised by the aliphatic C-H bond dissociation energy. It is possible to choose BDE to be the target and predicted directly by the ML model, for example, similarly to how st. John et al. were modelling bond dissociation enthalpies. In such a case, a database of BDEs is computed with an ab initio method and an ML model is created to predict those BDEs based on, for example, the structures' 2D representations (SMILES strings). Our approach instead chooses to retain the mechanistic aspect of computing the BDEs and instead uses an ML force field as a computationally efficient DFT alternative to drive the geometry optimisation and evaluate the total energies from which the BDE is computed. Therefore, no BDE data has been included in the MLIP training set and, because BDE is an equilibrium property, the training set only includes structures sampled from 500 K MD trajectories and no structures from the full hydrogen abstraction paths. By fitting an ML interatomic potential, rather than predicting the property of interest directly, the model can be additionally applied to compute other (thermodynamic) quantities, such as bond dissociation enthalpies.