In this contribution we derive a computational Bayesian approach to NARMAX model identification. The identification algorithm exploits continuing advances in computational processing power to numerically obtain posterior distributions for both model structure and parameters via sampling methods. The main advantage of this approach over other NARMAX identification algorithms is that for the first time model uncertainty is characterised as a byproduct of the identification procedure. The algorithm is based on the reversible jump Markov chain Monte Carlo (RJMCMC) procedure. Key features of the approach are (i) sampling of unselected model terms for testing for inclusion in the model (the birth move), which encourages global searching of the model term space, (ii) sampling of previously selected model terms for testing for exclusion from the model-a naturally incorporated pruning step (the death move), which leads to model parsimony, and (iii) estimation of model and parameter distributions, which are naturally generated in the Bayesian framework. We present a numerical example to demonstrate the algorithm and a comparison with a forward regression method: the results show that the RJMCMC approach is competitive and gives useful additional information regarding uncertainty in both model parameters and structure. (C) 2013 Elsevier Ltd. All rights reserved.