An optimization method is presented for the closed-loop identification of first-order-plus-time-delay (FOPTD) transfer function models of multivariable systems using step responses. A standard least-squares optimization method is used to obtain the parameters of the FOPTD models by matching the closed-loop step responses of the model with those of the actual process. A simple method is proposed to obtain the initial guess values for the transfer function model parameters from the process main and interaction responses. The effects of measurement noise and controller settings on the identified model parameters were also studied. This method was applied to stable FOPTD and higher-order transfer function models of multivariable systems. The proposed method considerably reduces the computational time (by about a factor of 15) for the optimization when compared with the genetic algorithm method reported by Viswanathan et al. (Ind. Eng. Chem. Res. 2001, 40, 2818-2826).