<p>This paper proposed a new methodology to address the key problem in model order reduction methods of setting random values of lower and upper bounds and initial values of the parameters in optimization-based approaches. The moth flame optimization (MFO) method is utilized for the model order reduction of linear time-invariant systems, wherein the search space boundaries are found using a novel strategy with the classical balanced truncation technique. Both the numerator and denominator coefficients of the desired reduced-order system are found using the proposed optimization approach. The integral square error (ISE) is employed as the objective function in the optimization of single-input–single-output (SISO) systems while a novel objective function is framed using ISE for the multi-input–multi-output (MIMO) systems. The key advantage of using balanced truncation-based search space boundaries ensures targeted search with potential solutions and stability of the reduced-order model. Further, the disadvantage of steady-state error of the balanced truncation is overcome using a gain adjustment factor. The overall methodology takes very less simulation time while keeping all the necessary parameters of the reduced system close to those of the original system. To test the efficacy of the proposed methodology, five real-world high-order transfer function systems with two SISO systems, two MIMO systems and one discrete-time system are considered and compared with existing methods through several error metrics and time- and frequency-domain specifications. Further, a high-order benchmark state-space model is considered to demonstrate the efficacy of proposed methodology. It has been found that the proposed methodology results in significant reduction of ISE and improvement in matching of step responses, preserving stability of the reduced-order models.</p>