This research paper uses the Bayesian approach to establish an approximate method to specify the four orders of multivariate seasonal autoregressive moving average (SARMA) models. The proposed methodology consists of four coherent consecutive steps. The first step is the approximation of the likelihood function of the model's parameters by a matrix Normal–Wishart form. The second step is to use a proposed semi-intermediate Bayesian procedure to have initial estimates for the four model orders. The third step is to combine the approximate likelihood function and the initial orders prior with one of the matrix Normal–Wishart prior density or Jeffreys' vague prior to developing an approximate joint posterior probability mass function of the orders of the model in a simple form. The last step is to evaluate the posterior probabilities over the range of the four orders and pick out the values of the orders at which the joint probability mass function reaches the highest probability to be the identified orders of the multivariate seasonal time series being analyzed. To test the adequacy of the proposed methodology, two simulation studies with three different prior orders have been conducted. The numerical results showed that the proposed Bayesian methodology is adequate to identify the orders of multivariate SARMA models for medium and large time series lengths.