Nonparemetric tests: how to support the null hypothesis you claim to be testing

Let us assume that we have taken an unbalanced number of independent random samples from 5 different populations, which will be analogous to 5 different locations in this example. Each observation belongs do a unique individual. We have measured some continuous variable -say the concentration of some chemical- in each individual that we sampled. For this example, we will assume that it is perfectly logical to directly compare this variable (i.e., the chemical) across our samples purely as a function of which location (population) that they were sampled from.

I will simulate this data by drawing samples from normal distributions with somewhat similar means and standard deviations:

set.seed(123)
data <- data.frame(group = factor(rep(c(paste0("G",1:5)), c(10,24,10,12,9))),
                   val = c(rnorm(10, mean=1.34,sd=0.17), 
                           rnorm(24, mean = 1.14, sd=0.11),
                           rnorm(10, mean=1.19, sd=0.15),
                           rnorm(12, mean=1.06, sd=0.11),
                           rnorm(9, mean=1.09, sd = 0.10)))

Here, group is the population/location from which observations were sampled, and val is the value of the continuous variable.

Now let us check some sample statistics, calculate standard errors for each group, and plot the distribution of samples, and run a test for normality

library(tidyverse)
se <- function(x) sd(x) / sqrt(length(x))
data%>%
  group_by(group)%>%
  summarise_at(., "val", list(mean=mean,med=median,sd=sd,se=se))%>%
  mutate(across(is.numeric, round, 2))
group  mean   med    sd    se
  <fct> <dbl> <dbl> <dbl> <dbl>
1 G1     1.35  1.33  0.16  0.05
2 G2     1.14  1.15  0.11  0.02
3 G3     1.21  1.17  0.14  0.05
4 G4     1.09  1.06  0.09  0.03
5 G5     1.05  1.06  0.07  0.02

#note we fail this though we "know" these were sampled from normal distributions, but lets go along with it
shapiro.test(data$val)

    Shapiro-Wilk normality test

data:  data$val
W = 0.9394, p-value = 0.003258
#make density plots
data%>%
  group_by(group)%>%
  ggplot(., aes(x=val))+
  geom_density(aes(color=group))

Now from here, we want to know if individuals that were sampled from different locations have different concentrations of this "chemical". We dont meet the assumptions of normality so we have decided to use an omnibus Kruskal-Wallis test:

kruskal.test(data$val,data$group)

    Kruskal-Wallis rank sum test

data:  data$val and data$group
Kruskal-Wallis chi-squared = 23.95, df = 4,
p-value = 8.174e-05

It suggests at least one of the locations is different, so we want to know which ones they are. We will approach this question with Dunn’s test:

#let us ignore the issue of multiple comparisons for the moment, this is just a conceptual example
dunn.test(data$val,data$group)

  Kruskal-Wallis rank sum test

data: x and group
Kruskal-Wallis chi-squared = 23.9499, df = 4, p-value = 0


                           Comparison of x by group                            
                                (No adjustment)                                
Col Mean-|
Row Mean |         G1         G2         G3         G4
---------+--------------------------------------------
      G2 |   3.189730
         |    0.0007*
         |
      G3 |   1.762110  -1.096030
         |     0.0390     0.1365
         |
      G4 |   3.956793   1.396187   2.116328
         |    0.0000*     0.0813    0.0172*
         |
      G5 |   4.250052   1.924417   2.534939   0.586373
         |    0.0000*     0.0272    0.0056*     0.2788

alpha = 0.05
Reject Ho if p <= alpha/2

It appears that we indeed have some "significant differences", but what exactly are there significant differences in? For each of these comparisons, exactly what null hypothesis did we just accept or reject? Of course in practice, we should have a clear answer to this question before conducting an experiment, but again this is just an example.

My understanding is that Dunn’s test compares the average rank for each group using the rank sums from the Kruskal-Wallis test to test the null hypothesis that the average rank of each group is the same, and the alternative hypothesis is that one group stochastically dominates the other. Depending on the specific situation, a significant result can be interpreted as having one group that stochastically dominates the other, meaning that you have a higher probability of randomly selecting a larger observation from one group than the other, or if you can assume that both groups were generated from the same distribution, a significant result would be interpreted as two groups that have different medians. Just about every document I have found states this much with a fair amount of clarity, but they don’t talk about how to tell which case applies to a given situation.

According to the R documentation:
"dunn.test computes Dunn’s test (1964) for stochastic dominance and reports the results among multiple pairwise comparisons after a Kruskal-Wallis test for stochastic dominance among k groups (Kruskal and Wallis, 1952). The interpretation of stochastic dominance requires an assumption that the CDF of one group does not cross the CDF of the other. dunn.test makes m = k(k-1)/2 multiple pairwise comparisons based on Dunn’s z-test-statistic approximations to the actual rank statistics. The null hypothesis for each pairwise comparison is that the probability of observing a randomly selected value from the first group that is larger than a randomly selected value from the second group equals one half"

If I understand this correctly, along with the other information I have provided, in no case does Dunn’s test make inferences about the distributions from which the data were drawn. In fact, to interpret Dunn’s test, we require another approach to estimate if the data for each group was generated from the same distribution in the first place.
So my question is how do we know, or how do we support, our claim to the specific null hypothesis we have tested in each case for the data above?