Nearly 15 years of radial velocity (RV) monitoring and direct imaging enabled the detection of two giant planets orbiting the young, nearby star $beta$ Pictoris. The $delta$ Scuti pulsations of the star, overwhelming planetary signals, need to be carefully suppressed. In this work, we independently revisit the analysis of the RV data following a different approach than in the literature to model the activity of the star. We show that a Gaussian Process (GP) with a stochastically driven damped harmonic oscillator kernel can model the $delta$ Scuti pulsations. It provides similar results as parametric models but with a simpler framework, using only 3 hyperparameters. It also enables to model poorly sampled RV data, that were excluded from previous analysis, hence extending the RV baseline by nearly five years. Altogether, the orbit and the mass of both planets can be constrained from RV only, which was not possible with the parametric modelling. To characterize the system more accurately, we also perform a joint fit of all available relative astrometry and RV data. Our orbital solutions for $beta$ Pic b favour a low eccentricity of $0.029^{+0.061}_{-0.024}$ and a relatively short period of $21.1^{+2.0}_{-0.8}$ yr. The orbit of $beta$ Pic c is eccentric with $0.206^{+0.074}_{-0.063}$ with a period of $3.36pm0.03$ yr. We find model-independent masses of $11.7pm1.4$ and $8.5pm0.5$ M$_{Jup}$ for $beta$ Pic b and c, respectively, assuming coplanarity. The mass of $beta$ Pic b is consistent with the hottest start evolutionary models, at an age of $25pm3$ Myr. A direct direction of $beta$ Pic c would provide a second calibration measurement in a coeval system.