## Thursday, March 22, 2012

### Comparing Banzhaf and Shapley-Shubik power indices

Last week I analyzed Shapley-Shubik power index in R. I got several requests to write a code calculating Banzhaf power index. Here is the proposed code.
Again I use data from Warsaw School of Economics rector elections (the details are in my last post). I give the code for calculation of Shapley-Shubik and Banzhaf power indices below.

# Constituency list
const <- c(30, 22, 27, 27, 41, 2 + 11, 38 + 5, 1 + 9)

# Shapley-Shubik power index
library(gtools)
perms <- permutations(8, 8)
outcome <- apply(perms, 1, function(x) {
x[sum(cumsum(const[x]) < 107) + 1] })
sspi <- prop.table(table(outcome))

# Banzhaf power index
subsets <- length(const) - 1
subs <- matrix(FALSE, subsets, subsets)
for (i in 1:subsets) {
subs[,i] <- rep(c(rep(FALSE, 2 ^ (i - 1)),
rep(TRUE, 2 ^ (i - 1))),
2 ^ (subsets - i))
}
banzhaf <- function(i) {
other <- const[-i]
part.sum <- apply(subs, 1, function(x) { sum(other[x]) } )
sum((part.sum < 106.5) & ((part.sum + const[i]) > 106.5))
}
bpi <- prop.table(sapply(1:8, banzhaf))

# power index comparison
names(bpi) <- c("C_1","C_2", "C_3", "C_4", "C_5",
"C_67", "C_89", "C_1011")
barplot(rbind(bpi, sspi, const / sum(const))[,order(const)],
col=c(2,4,1), beside = TRUE,