Abstract. A theory for nonlinear evolution of the mirror modes near the instability threshold is developed. It is shown that during initial stage the major instability saturation is provided by the flattening of the velocity distribution function in the vicinity of small parallel ion velocities. The relaxation scenario in this case is accompanied by rapid attenuation of resonant particle interaction which is replaced by a weaker adiabatic interaction with mirror modes. The saturated plasma state can be considered as a magnetic counterpart to electrostatic BGK modes. After quasi-linear saturation a further nonlinear scenario is controlled by the mode coupling effects and nonlinear variation of the ion Larmor radius. Our analytical model is verified by relevant numerical simulations. Test particle and PIC simulations indeed show that it is a modification of distribution function at small parallel velocities that results in fading away of free energy driving the mirror mode. The similarity with resonant Weibel instability is discussed.