Non-parametric stats, bootstrap, permutation testing

Some notes related to presentation

Talk

James Read-Tannock gave talk with general background ideas on non-parametric stats, more classical tests still taught in Psych programmes, and then introduced the ideas of the bootstrap, jackknife, and permutation testing. [Slides]

• Wikipedia article on the general idea 😀 of bootstrapping a la Baron von Münchhausen.
• In the statistics context, the bootrap refers to a method of obtaining the distribution of a statistic / estimator by resampling data.
• The Alpaca shampoo website w/ scrollytelling - the cool visual explainer shown by James during the talk

alpaca-shampoo

• There is a key paper by Efron that introduces the idea Efron (1979). Some bit of this paper are very technical, so might be beyond the scope of what we need here, but a useful resource to dip into for the nitty gritty, if you are interested.
• The plug-in principle came up and we also talked about population vs sample statistics and ideas around
• using the bootstrap procedure to get confidence intervals for a statistic, whereas a related idea, permutation tests, is more suited to hypothesis testing. There is a free online book chapter, which includes a discussion of the accuracy of the boostrap distribution, dependence on sample size, etc.

Efron, B. 1979. “Bootstrap Methods: Another Look at the Jackknife.” The Annals of Statistics 7 (1). https://doi.org/10.1214/aos/1176344552.

Worked examples with data from the Efron & Tibshirani book

Some thoughts on computationally expensive methods. To get some hands on the basics, let’s have a look at a basic example (from Chapter 2 in the book. Table 2.1), in which survival data from mouse in a treatment / control experiment are presented.

I decided to use julia in this notebook, but it shouldn’t be hard to reproduce in R, python, matlab, or your own preferred method of working with and plotting data.

The julia packages required are listed in the Project.toml and Manifest.toml. You should be able to download the folder from the github repo and import Pkg follwed by Pkg.activate(".") and Pkg.instantiate().

using Bootstrap
using Statistics
using Plots
gr(); # activate gr backend

# data from Table 2.1 in Efron & Tibshirani book
treatment = [94,197,16,38,99,141,23];
control = [52,104,146,10,51,30,40,27,46];

plot(treatment,ones(size(treatment)), line=:stem, color=:red, lw=2, label="Treatment")
plot!(control,ones(size(treatment)), line=:stem, 
        color=:blue, lw=2, label="Control",
        size=[600,100], ylim=[0, 1.2], 
        legend = :outertopright, yticks=[0,1],
        title = "Efron & Tibshirani murine example", 
        titlefont = font(10,"Arial"),
        titlefonthaligns=:right,
        tickfont = font(10,"Arial"),
        legendfont = font(10,"Arial"))
xaxis!(label="Survival time (days)")

Bootstrap of the median

Now we can create 10000 bootstrap samples from the treatment and control data - making sure to respect the size of each list/vector, and to use sampling with replacement*.

## basic bootstrap of the "median" treatment data
@time bs1 = bootstrap(median, treatment, BasicSampling(n_boot))

  0.332936 seconds (802.73 k allocations: 40.556 MiB, 99.63% compilation time)

Bootstrap Sampling
  Estimates:
     Var │ Estimate  Bias      StdError
         │ Float64   Float64   Float64
    ─────┼──────────────────────────────
       1 │     94.0  -14.0293   37.3796
  Sampling: BasicSampling
  Samples:  10000
  Data:     Vector{Int64}: { 7 }

The Bootstrap.jl package (see documentation) also has tools for getting confidence intervals, and extracting other info like # of observations, etc

## calculate 95% confidence intervals
cil = 0.95;

## basic CI
bci1 = confint(bs1, BasicConfInt(cil));

## percentile CI
bci2 = confint(bs1, PercentileConfInt(cil));

## BCa CI
bci3 = confint(bs1, BCaConfInt(cil));

## Normal CI
bci4 = confint(bs1, NormalConfInt(cil));

which lets you compare, e.g. the normal confidence interval (94.0, 34.7667175864441, 181.29188241355592), with the percentile confidence interval (94.0, 23.0, 141.0).

To visualise the distribution of bootstrapped medians, we can get them with straps() and use the histogram() function from Plots.jl. The bootstrap estimate of the standard error of the statistic of interest can be obtained from the standard deviation of the bootstrap replications. (The library will do this under the hood).

# look at the distribution (called "straps")
histogram( straps(bs1)[1], title="Distribution of bootstrapped medians [treatment]", 
                label="",
                size=[600,300] )

Bootstrap of other, more unusual statistics

We can also try something a bit more complicated, like boostrap sample treatment and sample and use those to look at the distribution of difference in means.

## basic bootstrap of the "mean"  data
bs_t = bootstrap(x->mean(x), treatment, BasicSampling(n_boot))
bs_c = bootstrap(y->mean(y), control, BasicSampling(n_boot))

# the difference estimator
# calculated from the
t = straps(bs_t)[1];
c = straps(bs_c)[1];
z =  t - c;
nbins = 50;

bins = range(minimum([t;c;z]),maximum([t;c;z]), nbins )
# percentile
pc = quantile(z, [0.025, 0.975])

histogram(t,  bins = bins, alpha=0.8, color=:gray, label="treatment",size=[600,400])
histogram!(c,  bins = bins, alpha=0.8, color=:white, label="control")

… and from those bootstrapped means, we can calculate the difference and look at the distribution, and percentiles of that:

# look at the distribution (called "straps")
histogram( z, bins = bins, 
              title="Distribution of bootstrapped difference", 
              color=:red, 
              label="treatment - control",
              size=[600,400])
vline!(pc, lw=2, label="2.5 and 97.5 percentiles")
vline!([mean(z)], lw=2, color=:black, label="")

Permutation testing

There is a whole ecosystem for doing hypothesis testing in julia which also includes permutation testing (see documentation for HypothesisTests.jl) or dip into R, which will have tons of options, too.

Lots more interesting stuff out there - if you want to add things here, send me a message.

References

There is a pretty comprehensive article on boot strap methods by the inventors of the methods - worth a look if this is something you will use in your work Efron and Tibshirani (1986).

Efron, B., and R. Tibshirani. 1986. “Bootstrap Methods for Standard Errors, Confidence Intervals, and Other Measures of Statistical Accuracy.” Statistical Science 1 (1): 54–75. https://doi.org/10.1214/ss/1177013815.

The classic textbook is quite expensive and a bit harder to get but you can ask around and/or look in the library to check it out.

An Introduction to the Bootstrap

The first couple of chapters are available at Born lab, Harvard Medical School as part of some neurobio class. Check it out. Tibshirani and Efron (1993).

Tibshirani, Robert J, and Bradley Efron. 1993. “An Introduction to the Bootstrap.” Monographs on Statistics and Applied Probability 57 (1): 1–436.