# Some code to illustrate permutation testing for two sample
# comparison. May 8/07

# See data in Table 9.4, Givens and Hoeting

stom <- c( 25, 42, 45, 46, 51, 103, 124, 340, 396, 412, 876, 1112 )
brea <- c( 24, 40, 719, 727, 791, 1166, 1235, 1581, 1804, 3460, 3808 )

# t-test, log scale

ttest <- t.test( log(stom), log(brea), alternative="two.sided", 
                  var.equal=TRUE  )

# wilcoxon test (note, same results for raw or log transformed...)

wtest <- wilcox.test( stom, brea, alternative="two.sided" )

# permutation versions

alldat <- c( log(stom), log(brea) )
sp <- c( rep("stomach", length(stom) ), rep("breast", length(brea) ) )
m <- length(stom); n <- length(brea)

B <- 10^6   ### use smaller B on slowish computer
tstats <- numeric(B)
wstats <- numeric(B)
tstats[1] <- ttest$statistic   # value observed
ranks <- rank(alldat)
wstats[1] <- sum( ranks[1:m] )
for( i in 2:B )
 {
  shuff <- sample(ranks)  # permutes 
  newdat <- alldat[shuff]
  tstats[i] <- ( t.test( newdat[1:m], newdat[(m+1):(m+n)]  , 
                  alternative="two.sided", var.equal=TRUE ) )$statistic
  wstats[i] <- sum( shuff[1:m] )
 }

# take a look
pdf(file="perm2.pdf",height=10,width=10)
par(mfrow=c(3,1), oma=c(0,0,4,0) )

stripchart( split(exp(alldat),sp), cex=2, pch=19, xlab="survival time (days)")
 

hist( tstats, prob=T, 100, 
    main="t statistic (log): permutation & t distribution",
		xlab="t", ylab="" )
abline( v=tstats[1], lwd=2, col="red" )
pv1 <- mean( abs(tstats)  >= abs(tstats[1]) )
text( -4, .3, format( pv1, digits=2 ),cex=1.5 )
text( -4, .2, "(0.017)" ,cex=1.5 )
supp <- seq(-4,4, length=50) 
lines( supp, dt(supp,df=(m+n-2)), lwd=2 )
hist( wstats, 100, prob=T, main="sum of ranks (stomach)" )
abline( v=wstats[1], lwd=2, col="red" )
pv2 <- 2*mean( wstats  <= wstats[1] )
text( 90, .02, format( pv2, digits=2 ),cex=1.5 )

title( main="Two sample comparison: Cameron & Pauling, 1978 PNAS, 75:4538-52",
  outer=T , cex.axis=1.5, cex.main=1.5)
dev.off()

#These data are used in the book "Computational Statistics" 
#by G.H. Givens and J.A. Hoeting.
#
#Name: cancersurvival.dat
#Chapter: 9
#Discussed in: Problem 9.5
#
#Source: Excerpts from Table 1 in E. Cameron and L. Pauling,
#"Supplemental ascorbate in the supportive treatment of cancer:
#Re-evaluation of prolongation of survival times in terminal human
#cancer", Proceedings of the National Academy of Science, 1978, 75,
#4538-4542.
#
#Description: Survival times (days) for patients with terminal stomach
#(disease=1) and breast (disease=2) cancer.  Both groups were treated
#with ascorbate.
#
#Variable labels: survivaltime, disease