swh:1:snp:7293efff0a0e3c53bcbd00ffc839903d0e184a6a
Raw File
Tip revision: 3ebe7b76043299bc4b6a541b5be0d7012895ca1c authored by Joao Sollari Lopes on 21 November 2017, 15:04:40 UTC
remove binaries folder
Tip revision: 3ebe7b7
model_choice.r
# @author:    joao lopes
# @workplace: University of Reading
# @date:      20th April 2010

# Function to perform model-choice on a file with rejection step results.
# It also plots a bar plor of the obtained posterior distribution.
# Creates .eps image files in the directory from where it is run

# @arg data_file - file with summary of 'real' data (.trg)
# @arg rej_file  - file with rejection step results (.rej)
# @arg pri_file  - file with sample of priors (.pri)
model_choice <- function(data_file , rej_file ){

	#import VGAM library
	library(VGAM)

	#import Mark Beaumont's scripts
	source("calmod.r")

	#import the .rej file
	abc.rej <- data.matrix(read.table(rej_file))

	#import the .trg files
	target <- data.matrix(read.table(data_file))

	#Bar plot for model-choice (rejection)
	top1<-length(which(abc.rej[,2]==1))/length(abc.rej[,2])
	top2<-length(which(abc.rej[,2]==2))/length(abc.rej[,2])
	print(c(top1,top2))
	write(c(top1,top2),"modelprob_rej.txt")
	barplot(c(top1,top2),names=c("Model1","Model2"))
	abline(h=1.0/2,col="red")

	#Bar plot for model-choice (regression)
	res.topol<-calmod(target,abc.rej[,2],abc.rej[,17:34],1,rej=F)
	print(res.topol$x2)
	write(res.topol$x2,"modelprob_reg.txt")
	barplot(res.topol$x2,names=c("Model1","Model2"))
	abline(h=1.0/2,col="red")
	print("model-choice done.")

}
back to top