The analysis and design process in large scale practical systems are computationally tedious because more number of variables are involved. Model reduction techniques are simplification methods based on physical considerations or by mathematical approaches employed to realize reduced models for the original high order systems. There are many reduction techniques available in the literature for single variable systems, but there are only few methods for reduction of multi variable systems. However, the methods related to single variable can be extended to reduction of linear Multi Input Multi Output (MIMO) systems. In this paper, a concept of interlacing property based on Hermite – Biehler theorem is discussed which includes the theory of interlacing of odd and even polynomials. In this method the denominator of MIMO system is reduced by interlacing property and reduced numerator is obtained by making use of matching of coefficients of higher order system. This method is applied for model reduction of 10^th order multi variable linear time invariant model of a power system. To prove the validity of the proposed algorithm, the performance indices and stability responses like settling time, overshoot, etc., are compared to that of the other existing methods.