6.3 – Combinations and permutations

Introduction

In statistical inference, an important foundational skill is to be able to count the possible outcomes of an experiment — the sample space. For example, when drawing on a sample of individuals for study from specific groups (strata) of a population, combinations may be calculated; where it is unwarranted to assume a particular theoretical probability distribution applies, permutation  tests create a reference distribution by calculating all possible values of a test statistic under rearrangements of the observed data.

The simplest place to begin with probability is to define an event — an event is an outcome — a “Heads” from a toss of coin, or the chance that the next card dealt is an ace in a game of Blackjack. In these simple cases we can enumerate how likely these events are given a number of chances or trials. For example, if the coin is tossed once, then the coin will either land “Heads” or “Tails” (Fig 1).

Heads (left) and Tails (right) of a USA Susan B. Anthony dollar.

Figure 1. Heads (left) and Tails (right) of a USA Susan B. Anthony dollar.

We can then ask what is the likelihood of ten “Heads” in ten tosses of a fair coin?

And this gets us to some basic rules of probability:

  • if successive events are independent, then we multiply the probability
  • if events are not independent, then the probability adds
  • do we sample without replacement, ie, an object can only be used once?
  • do we sample with replacement, ie, an object gets placed “back in the deck” and can be used again and again

Combinations and Permutations

Both combinations and permutations pertain to how samples are selected from a population. Combinations are used if the events or objects of the sample are all of the same kind, whereas for permutations, it is implied that the events or objects are not all of the same kind.

Figure 2. A combination lock

Combinations ignore the order of the events; permutations are ordered combinations. For replacement, the formula for permutations is simply n^r

In the case of the ten “Heads,” in ten successive trials, the probability is 1/2 “ten-times” or \left ( \frac{1}{2} \right )^{10} = 0.0009766 (in R, just type 0.5^10 or  (1/2)^10 at the R prompt).

Examples of permutations in biology include,

Given the four nucleotide bases of DNA, adenine (A), cytosine (C), guanine (G), and thymine (T), how many codons are possible?

    \begin{align*} 4^{3} = 64 \end{align*}

where the three (3) refers to the codon, three nucleotide genetic code. Codons are trinucleotide sequences of DNA or RNA that correspond to a specific amino acid. For transcription initiation of a protein-coding gene, a key sequence must be present in the promoter region of the gene. In Archaea and Eukaryotes, that sequence is called the TATA box. How often do we expect to find by chance the heptamer (ie, seven bases) TATA box sequence in the human genome?

Note 1: The consensus TATA box sequence, TATAWAW (W can be either Adenine or Thymine, per Standard IUB/IUPAC nucleic acid code). The TATA box is the binding site for transcription factors, particularly the TATA binding protein.

    \begin{align*} 4^{7} = 16384 \end{align*}

Thus, we would expect at random to find TATA box sequences every 16,384 bases along the genome. For the human genome of 3.3 billion bases then we would expect at random more than 200,000 TATA box consensus sequences.

Another way to put this is to ask how many combinations of ten tosses gets us ten “Heads” then weigh this against all possible outcomes of tossing a coin ten times — how many times do we get zero Heads? Two “Heads”? And so on. This introduces the combinatorial.

    \begin{align*} C\left ( n,r \right ) = \frac{n!}{\left (n-r \right )!r!} \end{align*}

where n is the number of choices (the list of events) and r is how many times you choose or trials. The exclamation point is called the factorial. The factorial instructs you to multiply the number over and over again, “down to one.” For example, if the number is 10, then

    \begin{align*} 10! = 10 \cdot 9 \cdot 8 \cdot 7 \cdot 6 \cdot 5 \cdot 4 \cdot 3 \cdot 2 \cdot 1 = 3628800 \end{align*}

In R, just type factorial(10) at the R prompt.

factorial(10)
[1] 3628800

An example

Consider an extrasensory perception, ESP, experiment. We place photocopies of two kinds of cards, “Queen” and “King,” into ten envelopes, one card per envelope (Fig 2). Thus, an envelope either contains a “Queen” or a “King.”

Darwin playing cards

Figure 2. Playing cards with images commemorating 150th anniversary of Charles Darwin’s Origin of Species. (Design John R. C. White, Master of the Worshipful Company of Makers of Playing Cards 2008 to 2009.)  

We shuffle the envelopes containing the cards and volunteers guess the contents of the envelopes, one at a time. We score the correct choices after all envelopes had been guessed. (It’s really hard to set these kinds of experiments up to rule out subtle cues and other kinds of experimental error! Diaconis 1978, Pigliucci 2010, pp. 77 – 83.) By chance alone we would expect 50% correct (eg, one way to game the system would have been for the volunteer simply to guess “Queen” every time; since we had placed five “Queen” cards the volunteer would end up being right 50%).

What are the possible combinations?

Let’s start with one correct; the volunteer could be right by guessing the first one correctly, then getting the next nine wrong, or equivalently, the first choice was wrong, but the second was correct followed by eight incorrect choices, and so on.

We need math!

Let n = 10 and r = 1

Using the formula above we have

    \begin{align*} C\left ( 10,1 \right ) = \frac{10!}{\left (10-1 \right )!1!} \end{align*}

or ten ways to get one out of ten correct. This is the number of combinations without replacement.

In R, we can use the function choose(), included in the base install of R. At the R prompt type

choose(10,1)

and R returns

[1] 10

For permutations, use choose(n,k)*factorial(k).

To get permutations for larger sets, you’ll need to write a function or take advantage of functions written by others and available in packages for R. See gtools package, which contains programming tools for R. There are several packages that can be used to expand capabilities of working with permutations and combinations. For example, install a package and library called combinat and then run the combn() function, together with the dim() function, which will return

dim(combn(10,1))[2]
[1] 10

R explained

In the combn() the “10” is the total number of envelopes and the “1” is how many correct guesses. We also used the dim() function to give us the size of the result. The dim() function returns two numbers, the number of rows and the number of columns — combn() returns a matrix and so dim() saves you the trouble of counting the outcomes — the “[2]” tells R to return the total number of columns in the matrix created by combn(). Here’s the output from combn() , but without the dim() function

combn(10,1)
     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,]    1    2    3    4    5    6    7    8    9    10

in order, we see that there is one case where the success came the first time [1,1], another case where success came on the second try [1,2] and so on.

We can write a simple function to return the combinations

for (many in seq(0, 10, by = 1)){
    print(choose(10, many))
}

an inelegant function, but it works well enough, returning

[1] 1
[1] 10
[1] 45
[1] 120
[1] 210
[1] 252
[1] 210
[1] 120
[1] 45
[1] 10
[1] 1

Note that there’s only one way to get zero correct, only one way to get all ten correct.

For completion, what would be the permutations? Use the permn() function (same combinat library) along with length() function to get the number of permutations

length(permn(10))
[1] 3628800

Continue the example

Back to our ESP problem.  If we continue with the formula, how many ways to get two correct? Three correct? And so on, we can show the results in a bar graph (Fig 3).

Figure 4. A simple bar chart

Figure 3. Bar chart of the combinations of correct guesses out of 10 attempts (graph was presented in Chapter 4.1).

The graph in Figure 3 was made in R

Combos <- seq(0,10, by=1)
HowMany <- c(1,10,45,120,210,252,210,120,45,10,1)
barplot(HowMany, names.arg = Combos, xlab = "Number correct", ylab = "Count",col = "darkblue")

If you recall, we had two students score eight out of ten. Is this evidence of “ESP?” It certainly seems pretty good, but how rare is this to score 8 out of 10? The total number of outcomes was

    \begin{align*} 1+10+45+ \cdots +45+10+1 = 1024 \ ways \end{align*}

 

Eight out of ten could be achieved in 45 different ways, so

    \begin{align*} \frac{45}{1024} = 0.043945 \end{align*}

 

So, eight out of ten is pretty unlikely! Is it evidence for ESP powers in our Biostatistics class?  In fact, at the traditional Type I error rate of 5% used in biology and other sciences to evaluate inferences from statistical tests (See Chapter 8), we would say that this observation was statistically significant. Given that this would be an extraordinary claim, this should give you an important clue that statistical significance is not the same thing as evidence for the claim of ESP.

Note 2: This is a good example where frequentist approach may not be the best approach — the Bayesian prior probability for ESP is vanishingly low. See O’Hagan (2006) for a similar ESP problem worked from a Bayesian perspective.

In other words, Is our result biologically significant? Probably not, but I’ll keep my eyes on you, just in case.

Conclusions

Combinations or permutations?

Combinations refers to groups of n things taken k times without repetition. Note that order does not matter, just the combination of things. Permutations, on the other hand, specifically relate the number of ways a particular arrangement can show up.

Questions

  1. Calculate the combinations, from zero correct to ten correct, from our ESP experiment, ie, confirm the numbers reported in Figure 3.
  2. Consider our ESP tests based on guessing cards. Let’s say that one subject repeatedly reports correct guesses at a rate greater than expected by chance. Why or why not should we view this as evidence the person may have extrasensory perception?
  3. Return to our TATA box example. The number of protein-coding genes in the human genome is presently estimated to be 20,620 (NCBI August 2025). Given that each gene has one promoter sequence, what proportion of TATA boxes are not associated with protein-coding genes?
  4. Assuming a standard 52-card deck, how many possible 5-card poker hands — the set of cards a player holds during a round of play — are there?

Quiz Chapter 6.3

Combinations and permutations


Chapter 6 contents

4.1 – Bar (column) charts

Introduction

Bar chart, also called a bar graph, or column chart in many spreadsheet apps (e.g., Google Sheets, Microsoft Excel, Libreoffice Calc — “bar charts” is reserved for horizontal bar plots),  are used to compare counts among two or more categories, i.e., an alternative to pie charts (Fig. 1). 

Number of SNP by DNA element and functional class. A. Pie chart. B. Bar chart.

Figure 1. Single nucleotide variants for human gene ACTB by DNA and functional element (data collected 19 May 2022 from NCBI SNP database with Advanced search query). A. Pie chart. Note that slide for “exon – nonsense” is not visible. B. Bar chart – color coded bars to facilitate comparison with pie chart.

Although bar charts are common in the literature (Cumming et al 2007; Streit and Gehlenborg 2014), bar charts may not be a good choice for comparisons of ratio scale data (Streit and Gehlenbor 2014). Bar charts for ratio data are misleading. Parts of the range implied by the bar may never have been observed: the bars of the chart always start at zero. Box (whisker) plots are better for comparisons of ratio scale data and are presented in the next section of this chapter. That said, I will go ahead and present how to create bar charts for both count — generally considered acceptable — and ratio scale data, for which their use is controversial.

Purpose of the bar chart

Like all graphics, a bar chart should tell a story. The purpose of displaying data is to give your readers a quick impression of the general differences among two or more groups of the data. For counts, that’s where the bar chart comes in. The bar chart is preferred over the pie chart because differences are represented by lengths of the bars in the bar chart. Differences among categories in a pie chart are reflected by angles, and it seems that humans are much better at judging lengths than angles.

Bar chart with counts

Figure 4. A simple bar chart

Figure 2. A simple bar chart

myCombo <- seq(0,10, by=1) 
myCounts <- choose(10, myCombo)    #combinations
barplot(myCounts, names.arg = myCombo, xlab = "Number correct", ylab = "Count",col = "darkblue")

 A stacked bar chart is used to 

myCombo <- seq(0,10, by=1) 
myCounts <- choose(10, myCombo)    #combinations
barplot(myCounts, names.arg = myCombo, xlab = "Number correct", ylab = "Count",col = "darkblue")

The luxury ship RMS Titanic, which sunk 15 April 1912, More than 1500 souls were lost. Public domain image, Wikipedia

Figure 3. The luxury ship RMS Titanic, which sunk 15 April 1912, More than 1500 souls were lost. Public domain image, Wikipedia. Click image to view full sized image.

Stacked bar chart, data set TitanicSurvival in package carData.

igure 6. A stacked bar chart, survival Titanic

Figure 4. A stacked bar chart, survival Titanic

Barplot(passengerClass, by=survived, style="divided", legend.pos="above", xlab="passengerClass", ylab="Frequency")

Bar charts with error bars

Although many data visualization specialists argue against the bar chart — they tend to obscure data — their use is well established. For the familiar bar chart with ratio scale data, the X (horizontal) axis displays the categories of one variable (e.g., location, or treatment group). You plot groups to emphasize comparisons. The Y (vertical) axis then is the mean for each group.

You need error bars with your bar chart to convey variability. A bar chart plus error bars is sometimes called (pejoratively) a dynamite plot. If the mean is displayed, some measure of precision should be (must be?) displayed (Cumming et al 2007). And, as you should recall by now, your choices are standard deviation (Chapter 3.2), standard error of the mean (SEM) (Chapters 3.2, 3.4), or confidence interval (see Chapter 3.4). It is strongly advised, that without a representation of precision, one should not interpret trends or group differences from representations of means (i.e., height of bars) alone.

The bar charts on this page are means plus or minus the standard error of the mean, + SEM. I’m not advocating here for the bar chart — box plots with efforts to show the actual data are better. However, bar charts remain common and may be more familiar to certain audiences. We’ll discuss which choice to make.

Examples

The Copper_rats_PMID3357063 dataset will be used for the next series of graphs (Data set). Refer to Mike’s Workbook for Biostatistics Part07 to review how to import the data.

A portion of the data set is shown below

head(Copper)
         Diet     Body    Heart    Liver
1 Adequate-Cu 320.1381 1.125037 10.259657
2 Adequate-Cu 329.6879 1.158982 9.843295
3 Adequate-Cu 327.9838 1.090374 9.855975
4 Adequate-Cu 334.6669 1.118183 9.942997
5 Adequate-Cu 338.3134 1.172636 9.860971
6 Adequate-Cu 345.4608 1.056183 8.885820

The data set consists of organ weights (heart, liver) from rats fed a diet adequate in copper, deficient in copper, and then a third group who received the adequate diet from perspective of amount of copper, but calorie restricted to match the decreased feeding rates of the rats fed the copper deficient diet. Copper is an essential trace element in our diet. The data set was simulated from descriptive statistics (means, standard deviation, number of subjects) of published data by sampling from a normal distribution. (Table 1, Ovecka et al 1988).

On we go with some graphs.

Note 1: R has many options to create bar charts, and especially ggplot2 can be used to great advantage, but there is a learning curve. One of the great things about R is that folks help each other by sharing code. For example, http://www.cookbook-r.com/Graphs/Plotting_means_and_error_bars_(ggplot2)/

But, I’m still learning about R graphics, and for bar charts with error bars, I find other packages more straight-forward. So, I’ll use this moment to point out that making a good graph is more about the end product then the particular tools. I have been using other tools for years to make my graphics, so I tend to default on these options first. Graphs presented here are mostly from Veusz software program (pronounced “views”). I like it because the software allows me to edit any of the elements of the graph.

Here’s a typical looking bar chart with ratio scale data: means ± SEM. Let’s look at them more critically.

Figure 7. A bar chart with error bars (standard error of the mean)

Figure 5. A bar chart with error bars (standard error of the mean).

Figure 8. Another bar chart with standard errors of mean

Figure 6. Another bar chart. with standard errors of mean

When making comparisons, make sure the axes have the same scale, or consider putting the graphs together.

Figure 9. Bar chart that allows for a comparison among levels of a a factor (organs, liver vs. heart).

Figure 7. Bar chart that allows for a comparison among levels of a a factor (organs, liver vs. heart).

Analysis note: For variables that covary with body size (for example, on average, tall people generally are heavier too), is how to present (and analyze) the data in such a way that body size is accounted for. Here, the solution was to express organ weight as the ratio of organ weight to body weight for the mice.  This may or may not be a good solution, and the answer is too complicated for us now (has to do with a thing called allometric scaling), but I wanted to at least present the issue and show how the graphics can be improved to handle some of these concerns.

Here’s the organ weights again, but taken as ratio of body mass.

 

Figure 10. Same chart as in Figure 9 but on ratios.

Figure 8. Same chart as in Figure 6, but on ratios.

Note 2. Lets be clear about expectations of you for statistics class. Now, R (and Rcmdr) do lots of graphics, pretty much anything you want, but it is not as friendly as it could be. For BI311 homework, the default graphics available via Rcmdr will generally be adequate for assignments. R and Rcmdr has many bar chart options, but there isn’t a straightforward way to get the error bars, unless you are willing to enter some code to the command line or learn a particular package (like gplots or ggplot2).

How to make a bar chart with error bars in R

Option 1. First, let’s try a work-around. Instead of an error bar option for bar chart menu, Rcmdr provides a plot of means that allows you to plot with error bars. These are equivalent graphs, the “bar chart” and the “plot of means”, though you should favor the bar chart format for publishing.

Rcmdr: Graphs → Plot of means…

Figure 11. Rcmdr menu popup for Plots Means

Figure 9. Rcmdr menu popup for Plots Means

Here, Rcmdr takes the data and calculates the mean and your choice of standard errors or deviations, confidence intervals, or no error bars. the resulting graph is below.

Figure 12. Plot of means

Figure 10. Plot of means

That’s an ugly graph (Fig. 10). Functional, good enough for data exploration and preliminary results, and certainly good enough for a Biostatistics homework or report. Additionally, connecting the dots here is a no-no. It implies that if we had measured categories between “adequate copper” and “deficient copper,” then the points would fall on those lines. That would be a complete guess. So, why did I include the connecting lines? That was the default setting for the command, and it makes the point — think before you click. One argument for connecting points in a graph is that it makes it easier for the reader to visualize trends.

This graph (Fig. 10) is fine for exploring data, but you will want to do better for publication.

Let’s make some better graphs with R

Once you are ready to go beyond the default settings available in Rcmdr, there is tremendous functionality in R for graphics. To access R’s potential, you’ll need to get into the commands a bit. I’m going to continue to try and shield you from the programming aspects of R, but from time to time you really need to see what is possible with R. Graphics is one such area. I use the package gplots, with 23 different graphing functions (type at R prompt ?gplots to call up the manual pages).

gplots should be among the packages on your R installation; if not, then install the package and run library(gplots) to complete the installation. We’ll try the barchart2 function.

But first, we need to get means for each of our groups.

At the R command prompt:

hrtWt <- tapply(Dataset$HeartWt, list(Group=Dataset$Group), mean, na.rm=TRUE) 

This code extracts means from our HeartWt variable for each Group, then  stores the three (in this case, because our data set has 3 groups) in the place holder I had called hrtWt. To verify that the three means are there, type “hrtWt” without the quotes then enter.

You should see

hrtWt Group Cu adequate Cu deficient Pair-fed 
               1.200000     1.566667 0.900000

R functions used: tapply, list, mean; na.rm was not needed but would be used to remove all missing values (recall during our import phase we were asked how missing observations were noted in our file; the default is NA).

Next, I want to apply standard error bars

stdDEV <- tapply(Dataset$HeartWt, list(Group=Dataset$Group), sd, na.rm=TRUE)
cil <- hrtWt-(stdDEV/sqrt(3))
ciu <- hrtWt+(stdDEV/sqrt(3))

I used cil and ciu to designate the lower cil and upper ciu values for my ± SEM (standard error of the mean). ciu stands for “confidence interval lower;” ciu stands for “confidence interval upper.”

Finally, here’s the plot command

barplot2(hrtWt, beside = TRUE, main=c("Mice fed different amounts of copper in diet"), col = c("blue", "red", "green"), xlab="Copper in diet", ylab="Heart mass (g)", ylim = c(0, 2), plot.ci=TRUE, ci.l=cil, ci.u=ciu)

Now, draw a box around the graphic

box()

Whew!

What does your new graph look like? My graph is below (Fig. 11).

Figure 13. A bar chart using barplot2

Figure 11. A bar chart using barplot2.

This works and the point is that once the script is written it easy to make small changes as you need in the future to make nice graphs.

If you are impatient like me, I like a GUI option, at least to start crafting the graph. The Rcmdr plugin KMggplot2 provides a good set of tools to make bar charts with error bars. An even better option I think is to use a software package that is designed for graphics, at least simple graphics like a bar chart. I use SciDAVis and Veusz for simple graphs like pie charts and bar charts; much easier to control.

ggplot2 bar charts with error bars

Nevertheless, here’s how to make a bar chart with error bars using ggplot2 (Fig. 11). First, we need to create a statistics summary. The script printed here was modified from scripts at R Graph Cookbook website.

require(plyr)
summarySE <- function(data=NULL, measurevar, groupvars=NULL, na.rm=FALSE,
     conf.interval=.95, .drop=TRUE) {
length2 <- function (x, na.rm=FALSE) {
     if (na.rm) sum(!is.na(x))
     else length(x)
}

#returns a vector with N, mean, and sd
datac <- ddply(data,groupvars, .drop=.drop,
    .fun = function(xx,col){
      c(N=length2(xx[[col]],na.rm=na.rm),
      mean=mean(xx[[col]],na.rm=na.rm),
      sd=sd(xx[[col]],na.rm=na.rm)
     )
    },
    measurevar
   )

#Rename the "mean" column
datac <- rename(datac,c("mean"=measurevar))

#Calculate the standard error of the mean
datac$se <- datac$sd/sqrt(datac$N)

#Get confidence interval
ciMult <- qt(conf.interval/2 + .5, datac$N-1)
datac$ci <- datac$se*ciMult

return(datac)
}

Applying this function to the BMI dataset yields the following output.

sumBMI <- summarySE(BMI, measurevar="BMI", groupvars=c("Sex", "Smoke"));sumBMI

   Sex Smoke   N      BMI       sd        se       ci 
1    F    No  23 26.83567 7.610271 1.5868511 3.290928 
2    F   Yes  14 25.76133 4.658625 1.2450698 2.689810 
3    M    No  10 26.35731 3.363575 1.0636557 2.406156 
4    M   Yes  27 26.71879 4.675631 0.8998256 1.849618

Now we are ready to make the bar chart with error bars

ggplot(tgc,aes(x=Smoke,y=BMI,fill=Sex)) + 
 geom_bar(position=position_dodge(),stat="identity",color="black") + 
 geom_errorbar(aes(ymin=BMI,ymax=BMI+se),width=.2,position=position_dodge(.9))

 

Figure 14. A barchart from ggplot2

Figure 12. A barchart from ggplot2

Questions

1. Why should you use box plots and not bar charts to display comparisons for a ratio scale variable between categories? Obtain a copy of the article by Streit and Gehlenbor 2014 — it’s free! After reading, summarize the pro and cons for box plots over bar charts with error bars.

2. Enter the following data into R. The data are sulfate levels in water, parts per million. 

type = c("Palolo Stream","Chaminade tap water", "Aquafina","Dasani") 
sulfateppm =c(11, 14, 5, 12)
try = data.frame(type,sulfateppm) 
byWater = tapply(try$sulfateppm,list(Group=try$type),mean)

Make a simple bar chart using the boxplot2 function in gplots package. 

3. Change the range of values on the vertical axis to 0, 20

4. Change the color of the bars from gray to blue

5. Add a label to the vertical axis, “Sulfates, ppm” (without the quotes)

6. Add a box around the graph.

Quiz – Chapter 4.1

Bar charts

Data set

DietBodyHeartLiver
Adequate-Cu320.13811.12503710.259657
Adequate-Cu329.68791.1589829.843295
Adequate-Cu327.98381.0903749.855975
Adequate-Cu334.66691.1181839.942997
Adequate-Cu338.31341.1726369.860971
Adequate-Cu345.46081.0561838.88582
Adequate-Cu343.0891.08126110.166647
Adequate-Cu328.34031.11127810.124185
Adequate-Cu324.97231.18919410.158402
Adequate-Cu325.23781.147159.939521
Deficient-Cu195.50521.909737.907565
Deficient-Cu182.78091.8236728.430167
Deficient-Cu184.37011.6322497.619104
Deficient-Cu193.78671.8317658.742489
Deficient-Cu180.04171.7103677.975879
Deficient-Cu208.53492.4956238.652445
Deficient-Cu182.30481.2620537.257726
Deficient-Cu203.04132.1536398.081782
Deficient-Cu193.38291.9860287.807328
Deficient-Cu195.05231.769758.297611
Pair-fed211.08580.69113436.251177
Pair-fed210.40410.69280677.696669
Pair-fed208.59690.69119016.973803
Pair-fed209.33330.70392116.629303
Pair-fed208.88890.70774866.038704
Pair-fed208.29940.70045356.606877
Pair-fed209.45240.69155436.228888
Pair-fed210.26990.69844976.638466
Pair-fed208.81420.72148476.353705
Pair-fed209.29770.68486566.536642

Chapter 4 contents