Finite mixture models, typically Gaussian mixtures, are well known and widely used as model-based clustering. In practical situations, there are many non-Gaussian data that are heavy-tailed and/or asymmetric. Normal inverse Gaussian (NIG) distributions are normal-variance mean which mixing densities are inverse Gaussian distributions and can be used for both haavy-tail and asymmetry. For NIG mixture models, both expectation-maximization method and variational Bayesian (VB) algorithms have been proposed. However, the existing VB algorithm for NIG mixture have a disadvantage that the shape of the mixing density is limited. In this paper, we propose another VB algorithm for NIG mixture that improves on the shortcomings. We also propose an extension of Dirichlet process mixture models to overcome the difficulty in determining the number of clusters in finite mixture models. We evaluated the performance with artificial data and found that it outperformed Gaussian mixtures and existing implementations for NIG mixtures, especially for highly non-normative data.