Semi-empirical models are commonly used to predict vortex-induced vibrations (VIV). Empirical parameters are often based on forced harmonic motion tests of short rigid cylinders in one degree of freedom, which is a significant simplification of the three dimensional VIV experienced by flexible cylinders. Still, the models may provide satisfactory estimates of the response of slender marine structures exposed to current. At least in terms of root mean square (r.m.s.) of displacements and strains, dominating frequency and dominating mode. However, VIV of long slender beams, especially in sheared current, show large response variability in time and space, which appears to be of a random nature. Semi-empirical models that predict steady-state VIV are hence unable to reflect the non-stationary response observed in experiments, even though the time and space averaged result might be well represented.In this work, a previously published semi-empirical time domain model for cross-flow vortex-induced vibrations is modified to describe the stochastic nature of the response of long slender beams subjected to stationary current. The mean non-dimensional frequency of the synchronization model (previously constant) is taken to be a simplified Gaussian process, where standard deviation and spectral frequencies are input. The stochastic synchronization model allows the response to jump between different eigenfrequencies and corresponding modes, and still predict mean and r.m.s. values close to the deterministic model. Its performance is verified through simulation of an experiment with a long riser in sheared flow. The response sensitivity of standard deviation and spectral frequencies is investigated. It indicates that for a proper choice of empirical coefficients, the chaotic response of the riser can be quite realistically simulated in terms of frequency variation and, to some extent, amplitude modulation.