The statistics of unsteady turbulence with uniform stratification $N$ (Brunt–Väisälä frequency) and shear $\alpha ({=}\,{\rm d}U_{1}/{\rm d}x_{3})$ are analysed over the entire time range ($0 \,{<}\, \alpha t \,{<}\, \infty$) using rapid distortion theory (RDT) over a wide range of Richardson number ${\hbox{\it Ri}} ({=}\,N^{2}/\alpha^{2})$, and initial conditions. The solutions are found to be described by the Legendre functions of complex degree with pure-imaginary argument and are compared with previously published results of both direct numerical simulations (DNS) and experiments. In the initial stage of development many of the characteristics are similar to those in stratified flow with no shear, since the turbulence is determined by $\hbox{\it Nt}$ at the leading order, and the effects of vertical shear $\alpha$ generally appear at higher order. It is shown how in developing turbulence for ${\hbox{\it Ri}} \,{>}\, 0$ and ${\hbox{\it Ri}} \,{>}\, 0.25$ respectively, oscillatory momentum and positive and negative density fluxes develop. Above a critical value of ${\hbox{\it Ri}}_{\hbox{\scriptsize\it crit}} ({\sim}\,0.3)$, their average values are persistently countergradient. This structural change in the turbulence is the primary mechanism whereby stable stratification reduces the fluxes and the production of variances. It is quite universal and differs from the energy and stability mechanisms of Richardson (1926) and Taylor (1931). The long-time asymptotics of the energy ratio ER$({=}\,\hbox{\it PE/VKE}$) of the potential energy to the vertical kinetic energy generally decreases with ${\hbox{\it Ri}}({\geq}\,0.25)$, reaching the smallest value of $3/2$ when there is no shear (${\hbox{\it Ri}} \,{\rightarrow}\, \infty$). For strong mean shear (${\hbox{\it Ri}}\,{<}\,0.25$), RDT significantly overestimates ER since (as in unstratified shear flow) it underestimates the vertical kinetic energy VKE. The RDT results show that the asymptotic values of the energy ratio ER and the normalized vertical density flux are independent of the initial value of ER, in agreement with DNS. This independence of the initial condition occurs because the ratios of the contributions from the initial values $P\!E_{0}$ and $K\!E_{0}$ are the same for PE and VKE and can be explained by the linear processes. Stable stratification generates buoyancy oscillations in the direction of the energy propagation of the internal gravity wave and suppresses the generation of turbulence by mean shear. Because the shear distorts the wavenumber fluctuations, the low-wavenumber spectrum of the vertical kinetic energy has the general form $E_{33}(k)\,{\propto}\, (\alpha tk)^{-1}$, where $(L_{X} \alpha t)^{-1} \,{\ll}\,k \,{\ll}\,L_{X}^{-1}$ ($L_{X}$: integral scale). The viscous decay is controlled by the shear, so that the components of larger streamwise wavenumber $k_{1}$ decay faster. Then, combined with the spectrum distortion by the shear, the energy and the flux are increasingly dominated by the small-$k_{1}$ components as time elapses. They oscillate at the buoyancy period $\pi/N$ because even in a shear flow the components as $k_{1} \,{\rightarrow}\, 0$ are weakly affected by the shear. The effects of stratification $N$ and shear $\alpha$ at small scales are to reduce both VKE and PE. Even for the same ${\hbox{\it Ri}}$, larger $N$ and $\alpha$ reduce the high-wavenumber components of VKE and PE. This supports the applicability of the linear assumption for large $N$ and $\alpha$. At large scales, the stratification and shear effects oppose each other, i.e. both VKE and PE decrease due to the stratification but they increase due to the shear. We conclude that certain of these unsteady results can be applied directly to estimate the properties of sheared turbulence in a statistically steady state, but others can only be applied qualitatively.