set.seed(1)

s <- 1e5
K <- seq (5, 100, 5)

original <- augmented <- matrix(0, length(K), 1)
rownames(original) <- rownames(augmented) <- K

row <- 0
for (k in K) {
	row <- row + 1
	for (sim in 1:s) {
		#Generate data from a pupil on k items
		theta <- rnorm (1)
		b <- runif (k, -1,1)
		x <- sum (rlogis (k) <= theta - b)
		
		#Generate a proposal distribution
		z = rlogis (k + 1, location = c(b, 0), scale = 1)
		j = x + 1
		o = order(z)
		theta.ast = z[o[j]] #choose (x + 1)-largest as proposed point
		
		#Compute log of acceptance probability
		alpha = 0
		if(o[j] <= k) { 
			alpha = alpha + log ( 1 + exp (theta.ast - b[o[j]]))
			alpha = alpha + log ( 1 + exp (theta))
			alpha = alpha - log ( 1 + exp (theta - b[o[j]]))
			alpha = alpha - log ( 1 + exp (theta.ast))
		}
		
		if (log (runif (1)) <= alpha) 
			augmented [row] <- augmented [row] + 1/s
		
		#Original SVE algorithm
		theta.ast <- rnorm (1)
		x.ast <- sum (rlogis (k) <= theta.ast - b)
		if (log (runif (1)) <= (theta - theta.ast)*(x.ast - x)) 
			original [row] <- original [row] + 1/s
	}
}

#Figure 6
par (mfrow=c(1,2), cex.main = 1.5, mar = c(4, 5, 0, 0) + 0.1, mgp = c(3.5, 1, 0),
cex.lab = 1.5, font.lab = 2, cex.axis = 1.3, bty = "n", las = 1)

plot (K, original, ylim = c(0,1) , xlim = c(0,100), lwd = 2, lty = 1, 
  xlab = " ", ylab = " ", pch = 20)
mtext ("Number of items",1, line = 3, cex = 1.5, font = 2)
mtext ("Acceptance rate", 2, line = 3.5, cex = 1.5, font = 2, las = 0)
segments (0,0.2,100,0.2,lty = 2, col = "gray", lwd = 0.5)
segments (0,0.4,100,0.4,lty = 2, col = "gray", lwd = 0.5)
segments (0,0.6,100,0.6,lty = 2, col = "gray", lwd = 0.5)
segments (0,0.8,100,0.8,lty = 2, col = "gray", lwd = 0.5)
segments (0,1,100,1,lty = 2, col = "gray", lwd = 0.5)

plot(K, augmented, ylim = c(0,1) , xlim = c(0,100), lwd = 2, lty = 1, xlab = " ",
  ylab = " ", pch = 20)
mtext ("Number of items",1, line = 3, cex = 1.5, font = 2)
mtext ("Acceptance rate", 2, line = 3.5, cex = 1.5, font = 2, las = 0)
segments (0,0.2,100,0.2,lty = 2, col = "gray", lwd = 0.5)
segments (0,0.4,100,0.4,lty = 2, col = "gray", lwd = 0.5)
segments (0,0.6,100,0.6,lty = 2, col = "gray", lwd = 0.5)
segments (0,0.8,100,0.8,lty = 2, col = "gray", lwd = 0.5)
segments (0,1,100,1,lty = 2, col = "gray", lwd = 0.5)