The Lagrangian L for gravity waves of finite amplitude in an N-layer, stratified shear flow is constructed as a functional of the generalized coordinates qv(t)≡ {qnv(t)} of the N + 1 interfaces, where the qnv are the Fourier coefficients in the expansion of the interfacial displacement ηv(x, t) in a complete, orthogonal set {ψn(x)}. The explicit expansion of L is constructed through fourth-order in the qnν and $\dot{q}_{n}^{\nu}$. Progressive interfacial waves and Kelvin–Helmholtz instability in a two-layer fluid are examined, and the earlier results of Drazin (1970), Nayfeh & Saric (1972) and Weissman (1979) are extended to finite depth. It is found that the pitchfork bifurcation associated with the critical point for Kelvin–Helmholtz instability, which is supercritical for infinitely deep layers, may be subcritical (inverted) for finite depths. The evolution equations that govern Kelvin–Helmholtz waves in the parametric neighbourhood of this critical point are shown to be equivalent to those for a particle in a two-parameter, central force field. The effect of surface tension is examined in an Appendix. Finally, the wave motion forced by flow over a sinusoidal bottom (as in Thorpe's tilting tank) is examined and the corresponding resonance curves and Hopf bifurcations determined. Numerical integrations reveal that stable limit cycles exist in some parametric neighbourhoods of these bifurcations. Period doubling was observed but did not lead to chaotic motion.