Model Thinking - Segregation Indices for Schelling Output

SegregationIndex.R

# this is an internal function just to make up block-like neighborhoods on a matrix with no natural borders.
# no need to tinker with it unless you're curious
DecideBLocks <- function(A,hoodsize){
	nhood_side <- floor(sqrt(hoodsize))
	Neighborhoods <- matrix(paste(ceiling(col(A)/nhood_side), ceiling(row(A)/nhood_side), sep="-"), nc=ncol(A))	
	nhoods.out <- length(unique(c(Neighborhoods)))
	return(Neighborhoods)
}

# this calculates a segregation index given a matrix of either 1s and 0s or TRUE and FALSE entries, and an equally dimensioned matrix of neighborhoods (character or factor), such as that made by the above function
# this function is also internal, so no need to use it directly unless you're modifying it
SegregationIndexi <- function(M,Neighborhoods){
	
	# grand sums
	sum1 <- sum(M,na.rm=TRUE)
	sum2 <- sum(abs(M-1),na.rm=TRUE)
	
	# neighborhood sums
	sb1 <- tapply(c(M),c(Neighborhoods), sum, na.rm=TRUE)
	sb2 <- tapply(c(M),c(Neighborhoods), function(x){sum(abs(x-1),na.rm=TRUE)})
	
	sum(abs(sb1/sum1 - sb2/sum2))/2
}

# this function is the one that you can actually use (calling the other 2). It takes the standard output from SchellingProcess(), here the argument M (no need to modify the list output), and a given neighborhood size measured in cells, where 16 cells means a 4x4 neighborhood, for instance. If you don't give a square number, the function will round down to a square number. All neighborhoods will be square, except possibly some edges, but it doesn't make much difference.

SchellingSegregation <- function(M,hoodsize=16){
	N <- DecideBLocks(M[[1]][[1]],hoodsize)
	l <- length(M$M)-1
	unlist(lapply(M[[1]][1:l],SegregationIndexi,Neighborhoods=N))
}