The OP links to an abstract by Tableman, Nguyen, and Ernst (2014), who mentioned a method for permutation CIs but didn't actually describe it in their abstract.
Most likely, it was the method from this MS thesis paper by one of the authors:
Nguyen, M.D. (2009). "Nonparametric Inference using Randomization and Permutation Reference Distribution and their Monte-Carlo Approximation" [unpublished MS thesis; Mara Tableman, advisor], Portland State University. Dissertations and Theses. Paper 5927. https://archives.pdx.edu/ds/psu/37406.
Section 2.2, Lemma 1, and Theorem 1 of Nguyen (2009) define a method that matches the description in Tableman, Nguyen, and Ernst (2014)'s abstract. It really does require only a single set of permutations, and it is different from Garthwaite (1996).
I haven't found Nguyen's method implemented in code anywhere, so Emily Tupaj and I have coded it up as an R package, CIPerm:
https://github.com/ColbyStatSvyRsch/CIPerm
https://cran.r-project.org/web/packages/CIPerm/index.html
I'd welcome any feedback.
We also included a brief summary of Nguyen (2009)'s method in our package vignette.
Here's an even briefer overview:
Let $Y$ be a vector of $n$ observations from one group, and let $X$ be $m$ observations from the other group. The difference in sample means is
$$t_0 = \frac{\sum_{j=1}^n Y_j}{n} - \frac{\sum_{i=1}^m X_i}{m}$$
In the standard permutation or randomization test, at each permutation let $k$ denote the number of swapped labels (i.e., $k$ of the $X_i$s are assigned to the $Y$ group and vice versa).
Then the permuted test statistics have the form
$$t_{k,d} = \frac{\sum_{i=1}^k X_i + \sum_{j=k+1}^n Y_j}{n} - \frac{\sum_{j=1}^k Y_j + \sum_{i=k+1}^m X_i}{m}$$
(where $d$ indexes over different permutations with the same value of $k$ --- I should really show the indices for summation changing too but the notation would get messy --- hopefully the idea is clear).
- Naive approach: For a given mean difference $\Delta$, you could test $H_0: \mu_Y-\mu_X = \Delta$ vs $H_A: \mu_Y-\mu_X \neq \Delta$ by replacing all $Y_j$ values with $Y_{j,\Delta}=Y_j-\Delta$ and running the usual test for a null difference of 0. If you do this many times for many different $\Delta$ values, then the failed-to-reject values form a confidence interval.
- Nguyen's method: His Theorem 1 and its Corollaries show that the quantiles of $w_{k,d}$ give you the CI endpoints directly, where
$$w_{k,d} = \frac{t_0-t_{k,d}}{k\left(\frac{1}{n}+\frac{1}{m}\right)}$$
Thus, you don't have to try out different $\Delta$ values or run new permutations. One set of permutations is enough, as long as you track $k$ at each permutation.
For small problems, the naive approach might be fine. But I consulted once for someone who had enough data that a single permutation test run took hours, and they wanted to run it daily. A naive grid search over many $\Delta$ values would have been impractical, but Nguyen's method worked just fine for them.
