In this paper, we develop a new statistical procedure for the comparison of frequency distributions on systems of ordinal indicators, based on a multidimensional fuzzy extension of the first order dominance (FOD) criterion. The procedure, named fuzzy-first order dominance (F-FOD), employs concepts and tools from partially ordered set theory and from fuzzy relational calculus and is designed to overcome the main limitations of previously developed algorithms for FOD analysis. In particular, F-FOD produces full pairwise comparison matrices, allows for partial orderings and rankings of the statistical units to be derived from the input data, is computationally sufficiently light to be applied in most cases of practical interest and is freely available in the R package PARSEC. To illustrate its effectiveness, we also show F-FOD in action on two real datasets concerning health in Denmark and child well-being in the Democratic Republic of Congo.