This paper proposes a two-phase identification approach to Mamdani fuzzy neural networks. The first phase is the system identification which includes a novel forward recursive input-output clustering method for the structure initialization and the gradient descent algorithm for the parameter initialization. The main advantage of the proposed method is that it fits perfectly the special clustering requirement for system identification: coarser clustering in the regions where the identified system is smoother and finer clustering in the regions where the system is more variable or nonlinear. The second phase is the system simplification which includes the accurate similarity analysis and merging method for similar fuzzy rules and the gradient descent algorithm for the parameter finalization. The accurate similarity analysis developed solves the long standing open problem how to compute the exact (rather than approximate) similarity between fuzzy sets and rules with Gaussian membership functions. Numerical experiments based on well-known benchmark data sets are used to verify the effectiveness and accuracy of the proposed approach. (C) 2016 Elsevier B.V. All rights reserved.