In this paper the normal-skew-normal distribution, proposed by Gomez et al. (Statistics 47:411-421, 2013), is extended to the multivariate case. It is also a special case of SU N n,2 -distribution, recently studied by Arellano-Valle and Genton (Chil J Stat 2:17-34, 2010). We show that the proposed distribution can be expressed as a shape mixture of the multivariate extended skew-normal distribution. Applying this property leads to deriving stochastic representations for the proposed distribution. Also we give some basic properties for this new family. Computational techniques using EM-type algorithms are employed for iteratively computing maximum likelihood estimates. Finally, an application of the new distribution is illustrated using some real data sets.