We propose a hybrid method combining partial differential equation (PDE) and Monte Carlo (MC) techniques to obtain efficient estimates of statistics for plastic deformation related to kinematic hardening models driven by transient coloured noise. Our approach employs a control variate strategy inspired by [CPAM, 75 (3), 455-492, 2022] and relies on a class of PDEs with non-standard boundary conditions, which we derive here. The solutions of those PDEs represent the statistics of models driven by transient white noise and are significantly easier to solve than the coloured noise version. Our approach uses a coupling between the white-noise-driven process and the coloured-noise-driven process, yielding a variance-reduced estimator through control variate techniques. We apply our method to threshold-crossing probabilities, which are used as failure criteria known as ultimate and serviceability limit states under non-stationary excitation. Our contribution provides solid grounds for such calculations and is significantly more computationally efficient in terms of variance reduction compared to standard MC simulations.