The problem of concern here is parameter estimation of harmonics in the presence of multiplicative and additive noise. Cyclic statistics are employed to estimate the frequencies and phases, after which the time series is demodulated and cumulants of the noise processes are estimated. The latter are then supplied to linear or nonlinear cumulant-based algorithms to identify ARMA model parameters for the noises. Cyclic statistics and higher order spectra-based approaches are shown to yield the same frequency estimates. Simulation examples illustrate the algorithms.