Joint stiffness is often represented by a parallel cascade model. The present study proposes a new approach to identify the parameters of such model structures from nonlinear frequency response functions. At first, a harmonic probing technique is used to derive the linear and higher-order frequency response functions (called the generalized frequency response functions (GFRFs)) of systems represented by parallel cascade models. The computation of the GFRFs is a recursive procedure where each lower order GFRF contains no effects from higher order terms. Thus the parameter estimation problem can be formulated in a linear least squares framework where the parameters corresponding to nonlinearities of different orders can be estimated independently, beginning with first order and then building up to include the nonlinear terms using the weighted complex orthogonal estimator, which is a modified version of the standard orthogonal least squares, that accommodates complex data. Simulation results are included to demonstrate that the proposed method can successfully estimate the parameters of the system under the effects of significant levels of noise.