In this study, an improved hybrid coupled model (hereafter referred to as “HCM_emd”) is developed by parameterizing empirically the temperature of subsurface water entrained into the mixed layer (Te), in which a tropical Pacific OGCM is coupled with a statistical atmospheric wind stress anomaly model that is constructed from a singular value decomposition (SVD) analysis through a SSTA submodel where Te is estimated empirically by means of the statistical relation with the sea level anomalies (SLA) simulated by the OGCM. As a comparison, a standard hybrid coupled model is also set up (hereafter referred to as “HCM_std”), in which the same OGCM and statistical atmospheric model are coupled directly. The long time free integrations show that there are many biases in the simulated interannual variability in HCM_std, which are observed in many other coupled models, such as the strong SST anomalies confined in a narrow band of the equator around the date line, but too weak anomalies in the eastern equatorial Pacific and along the coast of South America with a higher frequency oscillation (quasi-biennial oscillation). In contrast, HCM_emd performs much better. The simulated ENSO variability distributes much more realistically in both longitude and latitude directions with a frequency of 4 years, which is much closer to the observation. In particular, the seasonal phase-locked ENSO phenomenon is also reproduced in this improved coupled model. The further analysis shows that the simulated ENSO in HCM_emd is consistent to a certain extent with the delayed oscillator, western Pacific oscillator, recharge-discharge oscillator and advective-reflective oscillator from different aspects, implying that multi-mechanisms may operate at the same time during ENSO cycles.