In this paper, we construct and analyze a new second-order exponential method for simulating genetic regulatory systems. For such systems involving highly nonlinearities, a novel splitting technique is proposed to improve the stability and thus further increase accuracy and efficiency of the method. This idea can be also used with high-order exponential Runge–Kutta methods. Another possibility is to use exponential Rosenbrock methods which are derived based on dynamic linearization of the system. In particular, we additionally suggest a second- and fourth-order exponential Rosenbrock methods for genetic regulatory systems depending on their structure and properties. The linear stability analysis of the proposed methods is also presented. Our numerical experiments on a 2-gene model, 3-gene repressilator and a 13-gene budding yeast cell cycle model show the effectiveness of the proposed methods and technique over several widely used methods for genetic regulatory systems in the literature.