Okay, I'm a bit late to this party, but while I agree with what dsaxton says in the first paragraph, I think the second paragraph gets lost.
Re-randomization works very well to specify the null distribution for a large variety of statistics. However, you've managed to cause a problem by combining two pathological distributions (point distributions centered on 9 and 10 respectively) with the median -- a statistic which is perhaps at its least useful where there are only two possible values because it can become very unstable.
I'm going to try to walk through comparisons for several sample sizes to show what is happening here. It should help explain dsaxton's insight that the consistency is where the real statistical power lies.
Imagine we took one ride on each bus. We get one 9 and one 10.
We randomize 10,000 times to conduct inference. In half of them, the positions switch, in half they don't. Thus if we measured medians, half the time the difference in medians will be -1, and half the time it will be 1. Similarly for means, half the time the difference in means will be -1 and half the time it will be 1.
Now imagine we took 10 rides on each bus, resulting in ten 10s and ten 9s.
We re-randomize.
This time, most of the randomizations result in having about five of each 10 and 9 in each sample. The means will form normal (really a shifted binomial) distributions around 9.5 for each sample, giving a difference centered on 0.
The difference in medians can occasionally be 0 -- if we actually get five of each time in each sample -- giving medians in each sample of 9.5, but its more likely to have a slight imbalance. That slight imbalance makes the medians 9 and 10 or 10 and 9. Thus most of the time the difference of medians will be either -1 or 1, which is similar to our real result, giving the extra high p-value.
It may seem like continuing to raise the number of bus rides should fix this problem, but while that makes the mean more stable -- and fixes the null firmly around 0, it actually destabilizes the median. It becomes less and less likely to get that exact match, and so the middle ground disappears.
Okay. Maybe that made sense. I'm going to include some R code to make this concrete.
n = 10
a = rep(10,n) #initial samples
b = rep(9,n)
joint.sample = c(a,b) #Combining samples for ease
bootstraps = 10000 #Number of replications
est.mean = mean(a) - mean(b) #Estimate of treatment
boot.mean = replicate(bootstraps, {
new.sample = sample(joint.sample)
mean(new.sample[1:n]) - mean(new.sample[1:n+n])
}) #Simply resamples and takes means of the two groups
CI.mean = quantile(boot.mean,prob=c(0.025,0.975) #Calculates a CI
pval.mean = mean(boot.mean >= est.mean)*2 #Two-sided p-value
#Same things but with median
est.median = median(a)-median(b)
boot.median = replicate(bootstraps, {
new.sample = sample(joint.sample)
median(new.sample[1:n]) - median(new.sample[1:n+n])
})
CI.median = quantile(boot.median,prob=c(0.025,0.975)
pval.median = mean(boot.median >= est.median)*2
That should give results for you that show that randomization with a mean would strongly reject that these were the same. Feel free to fiddle with the sample size n to see how that affects things, but mostly, for such a clear cut case as this, it doesn't take a large sample to spot the difference. You should also be able to reject using a median -- but you would need a different pair of distributions such that the medians moved around a bit more. Anything continuous should do I think, and then its a matter of sample size.
One note of caution. I used the defaults for the sample function here to dictate whether I was going with or without replacement. In general you want to think really hard about which sampling type you're using because that can and will affect results.