BCPA/PathPlot

From QERM Wiki
Jump to: navigation, search

This function (an optional component of the Behavioral Change Point Analysis) takes path "z" at times "t" with windowsweep output "ws" and partition parameter output "pp" and generates a plot of the path, color coded by  \hat\rho , with thickness proportional to  \hat\mu .

PathPlot <- function(t, z, ws, pp, threshold = 10)
{
	n<-names(pp)
	for(i in 1:length(n))
		assign(n[i],pp[,i])

	# Prep rho legend
	rho.int <- round((rho.hat-min(rho.hat,na.rm=1))*1000)+1
	rho.col <- topo.colors(max(rho.int))[rho.int]
	rho.index<-quantile(0:max(rho.int))
	rho.index[1]<-1

	legend.cols<-topo.colors(max(rho.int))[round(rho.index)]	
	legend.rhos<-round(quantile((min(rho.hat)*1000):(max(rho.hat)*1000))/1000,1)

	# Prep mu legend
	sizes <- round(quantile(seq(0,max(mu.hat,na.rm=1),.01)),2)[2:5]

	# get "notable" changepoints and corresponding locations
	x.breaks <- ws$Break
	x.model <- ws$Model

	mids <- hist(x.breaks, breaks = t, plot = F)$mid
	freq <- hist(x.breaks, breaks = t, plot = F)$count

	t.breaks <- mids[freq > threshold]
	z.mid <- (z[-(1:2)]+z[-( (length(z)-1):(length(z)))])/2
	t.mid <- (t[-1]+t[-length(t)])/2
	z.breaks <- z.mid[match(t.breaks, t.mid)+1]


	# plot everything
	plot(z, asp=1, type="l", xlab="X", ylab="Y")
	segments(Re(z[-1]), Im(z[-1]), Re(z[-length(z)]), Im(Z[-length(Z)]), 
		 lwd=mu.hat/max(mu.hat)*12, col="darkgrey")
	segments(Re(z[-1]), Im(z[-1]), Re(z[-length(z)]), Im(Z[-length(Z)]), 
		 lwd=mu.hat/max(mu.hat)*10, col=rho.col)
	points(z.breaks, col=rgb(1,.2,.2,.8), cex=sqrt(freq/max(freq))*5, pch=4, lwd=2)
	points(z,pch=19,cex=0.5)
	lines(z, lwd=0.75) 

	# add a legend
	legend("topright",fill=legend.cols,legend=legend.rhos,
		title=expression(hat(rho)),bty="n",cex=1.2)
	legend("bottomleft",col="darkgrey",legend=round(sizes, -1),
		lty=1, lwd=sizes/max(mu.hat)*10,title=expression(hat(mu)),cex=1.2,bty="n")
}
Personal tools