Directed bond percolation process is an important model in statistical physics. It provides a paramount example of non-equilibrium phase transitions. By now its universal properties are known only up to the second-order of the perturbation theory. Here, our aim is to put forward a numerical technique with anomalous dimensions of directed percolation to the higher orders of perturbation theory. It is based on the perturbative renormalization scheme in deviation from the upper critical dimension, ϵ = 4-d. Universal quantities are expressed in terms of irreducible renormalized Feynman diagrams and there is no need for calculation of renormalization constants. A numerical evaluation of integrals has been performed using the Vegas algorithm from the Cuba library. Whitin the framework, the anomalous dimensions are computed up to three-loop order in ϵ.