Listings 8.9 to 8.12

May 12, 2013

the code in Listings 8.9 to 8.12 is correct, but it can have numerical problems when sample size is large. It is better to do the code in logs. Here is Listing 8.9 with that alternative, and it would be straightforward to do same adjustment to the other 3.

bern_test=function(p,v){
# calculation with numerical problems when sample size is large
	a=p^(sum(v))*(1-p)^(length(v)-sum(v))
	b=(sum(v)/length(v))^(sum(v))*(1-(sum(v)/length(v)))^(length(v)-sum(v))
#	return(-2*log(a/b))

# more stable calculation when sample size is large
	lv=length(v)
	sv=sum(v)

	al=log(p)*sv+log(1-p)*(lv-sv)
	bl=log(sv/lv)*sv +log(1-sv/lv)*(lv-sv)
	return(-2*(al-bl))
}

---