In line commutated converter based high-voltage direct current (LCC-HVDC) transmission systems, the transformer saturation can induce harmonic instability, which poses a serious threat to the safe operation of the power system. However, the nonlinear characteristics of the power grids introduced by the transformer saturation considerably limit the application of the conventional analysis methods. To address the issue, this paper derives a linear model for the transformer saturation caused by the DC current due to the converter modulation. Afterwards, the nonlinear characteristics of power grids with the transformer saturation is described by a complex valued impedance matrix. Based on the derived impedance matrix, the system harmonic stability is analyzed and the mechanism of the transformer saturation induced harmonic instability is revealed. Finally, the sensitivity analysis is conducted to find the key factors that influence the system core saturation instability. The proposed impedance model is verified by the electromagnetic transient simulation, and the simulation results corroborate the effectiveness of the proposed impedance model.