An unsupervised machine learning algorithm based on nonnegative matrix factor Two-dimensional deconvolution (NMF2D) with approximated optimum model order is proposed. The proposed algorithm adapted under the hybrid framework that combines the generalized EM algorithm with multiplicative update. As the number of parameters in the NMF2D grows exponentially the number of frequency basis increases linearly, the issues of model-order fitness, initialization, and parameters estimation become ever more critical. This paper proposes a variational Bayesian method to optimize the number of components in the NMF2D by using the Gamma-Exponential process as the observation-latent model. In addition, it is shown that the proposed Gamma-Exponential process can be used to initialize the NMF2D parameters. Finally, the paper investigates the issue and advantages of using different window length. Experimental results for the synthetic convolutive mixtures and live recordings verify the competence of the proposed algorithm.