In ball-type robust universal hypothesis testing (UHT), the null hypothesis is a set of probability distributions constrained by a ball of radius \(r > 0\) denoted \(B(\mathcal {P}_0,r)\) based on the cumulative density function of the nominal distribution \(\mathcal {P}_0\) . A major limitation is that this method is originally designed only for one-dimensional distributions. To overcome this limitation, this paper proposes a new method to deal with multidimensional samples. For this purpose, first of all, new bounds are defined in the multidimensional domain. Later, a new mathematical programming model based on the transformed region of \(B(\mathcal {P}_0,r)\) , namely empirical multidimensional robust UHT problem based on Kullback–Leibler divergence is proposed for ball-type robust UHT. The power of the new testing method combined with different types of bounds was then demonstrated by a computational study. This method fills the research gap by enabling ball-type robust UHT for multidimensional samples and is flexible in that it can be used with different type of bounds.