I don't think it is going a one sample Z; to me it looks like it is a test against a certain set of priors.
I'm confused, why are you doing a one sample Z using binomial data as your source data? You could simply create a distribution of N successes and see what quantile your actual data was in. However the above method doesn't look like a permutation test per-se to me; as your code doesn't actually permute the observed values between 1 and 0.
That being said, let me comment on your code - to me it looks as though z is defined as
zeta = (mean(mn) - mu) / (sqrt(var/nx))
z[i] <- zeta
Thus, each score in Z is like a Z score of the randomly created binomial vector using the priors you've selected as the null hypothesis. But then you compare that Z to abs(mx); where mx is defined as the mean of your observed binomial vector. At the very least this looks like a problem to me. Z scores should be either compared to some other Z score or means should be compared to means.
As I alluded to above, it is odd that you'd put all of this under a structure of a Z-test. The Z score is nominally a linear transform of the differences between means, as such the result of a test like this should be the same whether you use a Z score or simply look at the differences between means.
Moreover, what you are doing seems a like an attempt to test the observed value against some priors rather than an actual one sample permutation test. What you want to test against is something like permbinom (code provided below) where it could be the case for each observed value that it either was a success or it was not a success. This is in-line with Tukey's classic example with the lady who claimed she could tell whether tea or milk was added first. Critically different from your test, the assumption of this permutation test is that the null hypothesis is fixed at p = .5.
permbinom <- function(x)
{
newx <- x
nx <- length(x)
change <- rbinom(n=nx,size=1,prob=.5)
#This code is readable but inefficient
#Swap the values between 1 and 0 if change == 1
for (i in 1:nx)
{
if ((change[i] == 1) & (x[i] == 1)) {newx[i] <- 0}
if ((change[i] == 1) & (x[i] == 0)) {newx[i] <- 1}
}
return(newx)
}
permtest <- function(x,nsim=2000)
{
permref <- rep(NA,nsim)
obsn <- sum(x)
for (i in 1:nsim)
{
permref[i] <- sum(permbinom(x))
}
pval <- min(c(
(sum(permref > obsn)*2/nsim),
(sum(permref < obsn)*2/nsim)
))
return(pval)
}
I'm not 100% confident regarding how I'm calculating the p-value here; so if someone would kindly correct me if I'm doing it wrong I'll incorporate that as an edit.
For reference, here is a faster permutation function for one-sample tests of binomial data.
permbinomf <- function(x)
{
newx <- x
nx <- length(x)
change <- rbinom(n=nx,size=1,prob=.5)
#This code is readable but inefficient
#Swap the values between 1 and 0 if change == 1
newx <- x + change
newx <- ifelse(newx==2,0,newx)
return(newx)
}
Edit: The question is also put forth, "What about an arbitrary subset with a one-sample z-test?". That would also work, assuming you had a large enough sample to subset. However, it would not be a permutation test, it would be more akin to a bootstrap.
Edit 2: Perhaps the most important answer to your question, is this: You are doing something acceptable (if you fix the Z vs mean computational error noted above), but you aren't doing what you think you are doing. You are comparing your results to results where the null hypothesis is true. This is essentially a Monte-Carlo simulation and if you correct the math (and I suggest you also simplify it) it is an acceptable technique for testing your hypothesis. Also note, my answer above is for a two-tailed test. As noted in the other question, you are ignoring the nesting of binomial observations under participants but independence isn't an assumption in a permutation or monte-carlo test so you should be more or less fine. Though, as also noted there you ignore the possibility that some people are doing better than chance and others are performing at chance.
If you permute the individual values then you are testing the combined null hypothesis that there is no difference between groups and that there is no structure within values from the same individual. So if you reject the null you don't know if it is because the groups differ or because there is structure within an individual.
I would suggest permuting individuals and keeping all values with the individual.
I don't know of a program that does this out of the box, but writing a permutation test in R is often easier/quicker than figuring out a precanned program.
Here is a simple example:
library(nlme)
Odist <- split(Orthodont$distance, Orthodont$Subject)
out <- replicate( 1999, { tmp <- sample(Odist);
mean( unlist( tmp[1:16] ) ) - mean( unlist( tmp[17:27] ) ) } )
out <- c(out, mean(unlist(Odist[1:16])) - mean(unlist(Odist[17:27])) )
hist(out)
abline( v=out[2000] )
mean(out >= out[2000])
Best Answer
There is the possibility of using the
coinpackage for this type of stuff. See its webpage and the accepted answer to this question.An implementation for this type of stuff would be the following.
(Note that you need to use an high number for B in your real example)
Furthermore, you could also use a paired t.test. I don't see any reasons against it:
I think there are is one import thing you haven't discussed so far: With this implementation you would have controlled for the dyads, but not for the individuals (e.g., F1, independent of her male interaction partner). This is a serious threat to your analysis.
I know that there is a bunch of stuff an dyad analyses in psychology and related fields. Unfortunately I cannot give you a real pointer. But you should definitely check this stuff out before finishing your analyses. A quick search on rseek.org returns at least a package called
dyadand this webpage.