Before I begin, I have no axe to grind here with anyone. I just want to comment on some stuff I found and leave everyone else to do what they like with that information. I have no doubt Dr. Zak has done a lot of other work on which he bases his theory, I'm not trying to undermine that work, I'm more wanting to make a point about the state of the top journals in science. All my code here is for R.
The data I was looking at relates to an experiment in which participants were asked to invest money in a trustee (a stranger). If they invested, then the total funds for the investor and trustee increased. If the trustee shares the proceeds then both players end up better off, but if the trustee does not repay the investors’ trust by splitting the fund then the trustee ends up better off and the investor worse off. The question is, will investors trust the trustees to split the funds? If they do then they will invest, if not they won’t. Critically, one group of investors were exposed to exogenous oxytocin (N = 29), whereas a placebo group were not (N = 29). The oxytocin group invested significantly more than the placebo control group suggesting that oxytocin had causally influenced their levels of trust in the trustee. This is the sort of finding that the media loves.
The paper reports a few experiments, I want to focus on the data specifically related to trust shown in Figure 2a (reproduced below):
placebo<-c(3, 3, 4, 4, 4, 5, 6, 6, 6, 6, 6, 7, 7, 8, 8, 9, 9, 11, 11, 11, 12, 12, 12, 12, 12, 12)
oxy<-c(3, 4, 4, 6, 6, 7, 8, 8, 8, 8, 9, 9, 10, 10, 10, 11, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12)
The problem being that this gives us only N = 26 in the placebo group. Bollocks!
Ok, well perhaps they reported the wrong N. Here's their table of descriptives:
Let's have a look at the descriptives I get (you'll need the pastecs package):
> round(stat.desc(placebo, basic = F), 1)
median mean SE.mean CI.mean.0.95 var std.dev coef.var
7.5 7.9 0.6 1.3 10.2 3.2 0.4
> round(stat.desc(oxy, basic = F), 1)
median mean SE.mean CI.mean.0.95 var std.dev coef.var
10.0 9.6 0.5 1.1 8.1 2.8 0.3
For the oxytocin group the mean, median and SD match, but for the placebo group they don't. Hmmm. So, there must be missing cases. Based on where the median is, I guessed that the three missing cases might be values of 10. In other words, Figure 2a in the paper, ought to look like this:
So, the placebo group data now looks like this:
placebo<-c(3, 3, 4, 4, 4, 5, 6, 6, 6, 6, 6, 7, 7, 8, 8, 9, 9, 10, 10, 10, 11, 11, 11, 12, 12, 12, 12, 12, 12)
Let's look at the descriptives now:
> round(stat.desc(placebo, basic = F), 1)
median mean SE.mean CI.mean.0.95 var std.dev coef.var
8.0 8.1 0.6 1.2 9.5 3.1 0.4
median mean SE.mean CI.mean.0.95 var std.dev coef.var
8.0 8.1 0.6 1.2 9.5 3.1 0.4
They now match the table in the paper. So, this gives me confidence that I have probably correctly reconstructed the raw data despite Figure 2's best efforts to throw me off of the scent. Of course, I could be wrong. Let's see if my figures match their analysis. They did a Mann-Whitney U that yielded a p-value of 0.05887, which they halved to report one-tailed as .029. Notwithstanding the fact that you probably shouldn't do one-tailed tests, ever, they conclude a significant between group difference. I replicate their p-value, again convincing me that my reconstructed data matches their raw data.
Put data into a data frame and do wilcoxon rank sum test (which is equivalent to Mann-Whitney):
zak<-data.frame(gl(2, 29), c(placebo, oxy))
names(zak)<-c("Group", "MU")
zak$Group<-factor(zak$Group, labels = c("Placebo", "Oxytocin"))
> wilcox.test(MU~Group, zak)
Wilcoxon rank sum test with continuity correction
data: MU by Group
W = 301, p-value = 0.05887
alternative hypothesis: true location shift is not equal to 0
So, all well and good. There are better ways to look at this than a Mann-Whitney test though, so let's use some of Rand Wilcox's robust tests. First off, a trimmed mean and bootstrap:
> yuenbt(placebo, oxy, tr=.2, alpha=.05, nboot = 2000, side = T)
$ci
[1] -3.84344384 0.05397015
$test.stat
[1] -1.773126
$p.value
[1] 0.056
$est.1
[1] 8.315789
$est.2
[1] 10.21053
$est.dif
[1] -1.894737
$n1
[1] 29
$n2
[1] 29
Gives us a similar p-value to the Mann-Whitney and a confidence interval for the difference that contains zero. However, the data are very skewed indeed:
Perhaps we're better off comparing the medians of the two groups. (Medians are no-where near as biased by skew as means):
> pb2gen(placebo, oxy, alpha=.05, nboot=2000, est=median)
[1] "Taking bootstrap samples. Please wait."
$est.1
[1] 8
$est.2
[1] 10
$est.dif
[1] -2
$ci
[1] -6 1
$p.value
[1] 0.179
$sq.se
[1] 2.688902
$n1
[1] 29
$n2
[1] 29
The p-value is now a larger .179 (which even if you decide to halve it won't get you below the, not very magic, .05 threshold for significance). So, the means might be significantly different depending on your view of one-tailed tests, but the medians certainly are not. Of course null hypothesis significance testing is bollocks anyway (see my book, or this blog) so maybe we should just look at the effect size. Or maybe we shouldn't because the group means and SDs will be biased by the skew (SDs especially because you square the initial bias to compute it). Cohen's d is based on means and SDs so if these statistics are biased then d will be too. I did it anyway though, just for a laugh:
> d<-(mean(oxy)-mean(placebo))/sd(placebo)
> d
[1] 0.4591715
I couldn't be bothered to do the pooled estimate variant, but because there is a control group in this study it makes sense to use the control group SD to standardise the mean difference (because the SD in this condition shouldn't be affected my the experimental manipulation). The result? About half a standard deviation difference (notwithstanding the bias). That's not too shabby, although it'd be even more un-shabby if the estimate wasn't biased. Finally, and notwithstanding various objections to using uninformed priors in Bayesian analysis) we can use the BayesFactor packaged to work out a Bayes Factor.
> ttestBF(oxy, y = placebo)
Bayes factor analysis
--------------
[1] Alt., r=0.707 : 1.037284 ±0%
Against denominator:
Null, mu1-mu2 = 0
---
Bayes factor type: BFindepSample, JZS
The Bayes factor is 1.03, which means that there is almost exactly equal evidence for the null hypothesis as the alternative hypothesis. In other words, there is no support for the hypothesis that oxytocin affected trust.
So, Nature, one of the top journals in science, published a paper where the evidence for oxytocin affecting trust was debatable at best. Let's not even get into the fact that this was based on N = 58. Like I said, I'm sure Dr. Zak has done many other studies and I have no interest in trying to undermine what he or his colleagues do, I just want to comment on this isolated piece of research and offer an alternative interpretation of the data. There's lots of stuff I've done myself that with the benefit of experience I'd rather forget about, so, hey, I'm not standing on any kind of scientific high ground here. What I would say is that in this particular study, based on the things least affected by the shape of the data (medians) and a (admittedly using uninformed priors) Bayesian analysis there is not really a lot of evidence that oxytocin affected trust.
