In this study, a recursive algorithm is presented for identification of both linear and non-linear blocks of the Hammerstein–Wiener (H–W) models. An iterative sampling schema is used, along with non-parametric describing of non-linear functions by Gaussian processes, to consider the Gaussian distributed noise and disturbances. Different aspects of solving the H–W problem are discussed, covering the formulation of the non-linear functions, uniqueness assumptions and the effect of stochastic disturbances on the solution. The proposed method is used to identify an H–W model for the systems encounter all possible disturbances consists of input and output measurement noise, disturbances in the linear block and the additive disturbances to the interconnecting signals of the linear and non-linear blocks. Effectiveness and simplicity of the method are shown by an example. The results show that the accuracy of the proposed method is acceptable in presence of noise and disturbances.