We utilize a novel numerical technique to model star formation in cosmological simulations of galaxy formation - called SUPERSTARS - to simulate a Milky Way-like galaxy with ≳108 star particles to study the formation and evolution of out-of-equilibrium stellar disc structures in a full cosmological setting. In the plane defined by the coordinate and velocity perpendicular to the mid-plane [vertical phase space, (Z, VZ)], stars in solar-like volumes at late times exhibit clear spirals qualitatively similar in shape and amplitude to the Gaia 'snail shell' phase spiral. We show that the phase spiral forms at a lookback time of ~6 Gyr during the pericentric passage of an ~1010$\rm {\rm M}_{\odot }$ satellite on a polar orbit. This satellite stimulates the formation of a resonant wake in the dark matter halo while losing mass at a rate of ~0.5-1 dex per orbit loop. The peak magnitude of the wake-induced gravitational torque at the solar radius is ~8 times that from the satellite, and triggers the formation of a disc warp that wraps up into a vertical phase spiral over time. As the wake decays, the phase spiral propagates several gigayears to present day and can be described as 'ever-present' once stable disc evolution is established. These results suggest an alternative scenario to explain the Gaia phase spiral, which does not rely on a perturbation from bar buckling or a recent direct hit from a satellite.