This script analyses the data set data_DK5208ML.csv. This data set contains the results of an experiment testing whether DK5208 cheats on some of the evolved lines from MyxoEE-III (experimentally evolved descendants of GJV1 and GJV2, hereafter referred to as motility lines).

Step 1

We import the data and make a plot to visualize spore production.

## 'data.frame':    108 obs. of  6 variables:
##  $ replicate : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ strain    : chr  "DK5208" "GJV1" "GJV1" "GJV1" ...
##  $ count     : chr  "alone" "alone" "cheater" "total" ...
##  $ antibiotic: chr  "oxy" "none" "oxy" "none" ...
##  $ dilution  : int  1 5 3 5 5 3 5 5 3 5 ...
##  $ cfus      : num  0 53 532 156 161 345 143 171 55 188 ...
## 'data.frame':    99 obs. of  8 variables:
##  $ replicate : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ strain    : chr  "DK5208" "GJV1" "GJV1" "GJV1" ...
##  $ count     : chr  "alone" "alone" "cheater" "total" ...
##  $ antibiotic: chr  "oxy" "none" "oxy" "none" ...
##  $ dilution  : int  1 5 3 5 5 3 5 5 3 5 ...
##  $ cfus      : num  0 53 532 156 161 345 143 171 55 188 ...
##  $ numspores : num  0 5300000 532000 15600000 16100000 345000 14300000 17100000 55000 18800000 ...
##  $ logspores : num  0 6.72 5.73 7.19 7.21 ...

Step 2

First, we calculate the final frequencies of DK5208 in mixture, then we calculate the relative fitness parameter Wij.

Wij = log(number spores of cheater in mixture/number initial cells of cheater in mixture) - log(number spores of partner in mixture/number initial cells of partner in mixture).

We replace 0 spore counts with 9s. If I saw 0 spores at the lowest dilution plated (10^-1), then 0.9 is the most conservative estimate for the number of spores produced below the detection limit, giving a total of 9 spores produced.

Step 3

Fit a general linear model using ANOVA and then check contrasts. This model excludes the mix with GJV1.

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## strain       8 20.322  2.5402   17.18 5.93e-07 ***
## Residuals   18  2.662  0.1479                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The “strain” variable explains some of the variability of the data set, so we follow up with t-tests.

First, we do a single t-test to see if DK5208 cheats on GJV1, as expected from previous research.

## 
##  One Sample t-test
## 
## data:  mutantdata[which(mutantdata$strain %in% "GJV1"), ]$wij
## t = 3.8, df = 3, p-value = 0.016
## alternative hypothesis: true mean is greater than 0
## 95 percent confidence interval:
##  0.265897      Inf
## sample estimates:
## mean of x 
## 0.6984528

p = 0.016. DK5208 cheats on GJV1, as expected.

Test for a general trend among the motility lines with a sign test (one-tailed exact binomial test). We saw that the natural isolates reduced Wij of DK5208, so we test DK5208’s mean Wij is lower with the motility lines than with GJV1. This is true for 8 out of the 9 motility lines in the experiment.

## 
##  Exact binomial test
## 
## data:  1 and 9
## number of successes = 1, number of trials = 9, p-value = 0.01953
## alternative hypothesis: true probability of success is less than 0.5
## 95 percent confidence interval:
##  0.0000000 0.4291355
## sample estimates:
## probability of success 
##              0.1111111

In this case, “success” means that DK5208 had a similar or greater Wij when mixed with that motility line compared to with GJV1. That happened with 1 motility line out of 9. Our hypothesis, based on what we saw with the natural isolates, was that the probability of “success” should be less than 0.5, since we expect it to be more likely for DK5208 to have a lower Wij when mixed with the motility lines. Based on the p-value of 0.02, we are able to reject the null hypothesis that the probability of “success” is 0.5.

Now we see if Wijs from mixes with the motility lines are different from 0, ie is it cheating like on GJV1 or suppressed like by the natural isolates? (9 t-tests, Bonferroni-Holm corrected for multiple testing).

## [1] "E1" "E2" "E3" "E4" "E5" "E6" "E7" "E8" "E9"
## [1] 0.2021 0.1227 0.9094 0.1227 0.9094 0.1365 0.2258 0.1381 0.1381

Wij is not shown to be different from zero: 0.12 < p < 0.91, 9 2-sided t-tests against 0 corrected using the Bonferroni-Holm method.

Mean Wijs with E2 and E4:

## 
##  One Sample t-test
## 
## data:  mutantdata[which(mutantdata$strain %in% "E2"), ]$wij
## t = -7.2529, df = 2, p-value = 0.01848
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -1.1649965 -0.2974322
## sample estimates:
##  mean of x 
## -0.7312143
## 
##  One Sample t-test
## 
## data:  mutantdata[which(mutantdata$strain %in% "E4"), ]$wij
## t = -5.9305, df = 2, p-value = 0.02727
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -3.6145560 -0.5749826
## sample estimates:
## mean of x 
## -2.094769

E2: -0.73 E4: -2.09

Test whether any of the Wij values are different from that with GJV1 using a Dunnett test.

## 
##   Simultaneous Tests for General Linear Hypotheses
## 
## Multiple Comparisons of Means: Dunnett Contrasts
## 
## 
## Fit: aov(formula = wij ~ strain, data = mutantdata)
## 
## Linear Hypotheses:
##                Estimate Std. Error t value Pr(>|t|)    
## E1 - GJV1 == 0 -0.10755    0.29188  -0.368    1.000    
## E2 - GJV1 == 0 -1.42967    0.29188  -4.898   <0.001 ***
## E3 - GJV1 == 0 -0.73411    0.29188  -2.515    0.125    
## E4 - GJV1 == 0 -2.79322    0.29188  -9.570   <0.001 ***
## E5 - GJV1 == 0 -0.73399    0.29188  -2.515    0.125    
## E6 - GJV1 == 0  0.10327    0.29188   0.354    1.000    
## E7 - GJV1 == 0 -0.33167    0.29188  -1.136    0.854    
## E8 - GJV1 == 0 -0.06609    0.29188  -0.226    1.000    
## E9 - GJV1 == 0 -0.13039    0.29188  -0.447    1.000    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)

Both E2 and E4 are significantly different from GJV1 (p < 0.001 in both cases).

Output figures