Basic Statistics in R

Statistics let you summarize a set of measurements and judge whether a pattern in your data is real or just noise. In this lesson you will summarize gene expression numbers, compare two groups of samples, and measure how two genes move together. Every function used here is built into R, so there is nothing to install.

Descriptive Statistics

# <- stores a value under a name so you can reuse it later.
# c() combines several numbers into a vector, which is one ordered list of values.
# Each number below is one sample's expression level for a gene.
control <- c(4.2, 5.1, 4.8, 6.0, 5.5, 4.9)

mean(control)     # average: the sum of the values divided by how many there are
median(control)   # the middle value once the numbers are sorted; with an even count
                  # like this one (6 values) it is the average of the two middle
                  # values, here 4.9 and 5.1, which gives 5
sd(control)       # standard deviation: how spread out the values are

[1] 5.083333
[1] 5
[1] 0.6177918

Comparing Two Groups

# A second group of samples, where the same gene was measured after treatment.
treated <- c(6.5, 7.1, 6.8, 7.9, 7.2, 6.9)

# t.test() asks whether the two groups' averages differ by more than random noise.
# Read the p-value and the 95 percent confidence interval in the printout below.
t.test(control, treated)

	Welch Two Sample t-test

data:  control and treated
t = -6.2287, df = 9.3904, p-value = 0.0001281
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -2.699102 -1.267565
sample estimates:
mean of x mean of y 
 5.083333  7.066667 

Measuring Correlation

# Two genes measured across the same six samples, in matching order.
geneA <- c(4.2, 5.1, 4.8, 6.0, 5.5, 4.9)
geneB <- c(3.9, 4.7, 4.5, 5.8, 5.2, 4.6)

# cor() returns one number from -1 to 1 describing how the two move together.
# +1 is a perfect upward line, 0 is no linear link, -1 is a perfect downward line.
cor(geneA, geneB)

# cor.test() adds a p-value and a confidence interval for that correlation.
cor.test(geneA, geneB)

[1] 0.9962565

	Pearson's product-moment correlation

data:  geneA and geneB
t = 23.049, df = 4, p-value = 2.099e-05
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.9645830 0.9996099
sample estimates:
      cor 
0.9962565 

Both t.test() and cor.test() report a p-value and a 95 percent confidence interval. The p-value is the probability of seeing a result at least this extreme if there were truly no difference between groups (or no correlation), so a small value, conventionally below 0.05, is evidence against that no-effect assumption, while the confidence interval gives the plausible range for the true value. A small p-value does not tell you the effect is large or biologically important, and it cannot rescue biased sampling or a tiny, unrepresentative dataset, so always report the interval and the size of the effect alongside it.