A novel approach is presented for discovering PDEs that govern the motion of satellites in space. The method is based on SINDy, a data-driven technique capable of identifying the underlying dynamics of complex physical systems from time series data. SINDy is utilized to uncover PDEs that describe the laws of physics in space, which are non-deterministic and influenced by various factors such as drag or the reference area (related to the attitude of the satellite). In contrast to prior works, the physically interpretable coordinate system is maintained, and no dimensionality reduction technique is applied to the data. By training the model with multiple representative trajectories of LEO - encompassing various inclinations, eccentricities, and altitudes - and testing it with unseen orbital motion patterns, a mean error of around 140 km for the positions and 0.12 km/s for the velocities is achieved. The method offers the advantage of delivering interpretable, accurate, and complex models of orbital motion that can be employed for propagation or as inputs to predictive models for other variables of interest, such as atmospheric drag or the probability of collision in an encounter with a spacecraft or space objects. In conclusion, the work demonstrates the promising potential of using SINDy to discover the equations governing the behaviour of satellites in space. The technique has been successfully applied to uncover PDEs describing the motion of satellites in LEO with high accuracy. The method possesses several advantages over traditional models, including the ability to provide physically interpretable, accurate, and complex models of orbital motion derived from high-entropy datasets. These models can be utilised for propagation or as inputs to predictive models for other variables of interest.