#bootstrap analysis on Genetic Distances p.dist<-function(obj,nperm){ obs<-pairwise.fst(obj) test<-sapply(1:nperm,function(e)pairwise.fst(df2genind(X = swans[,4:16], sep = "/", ind.names = swans[, 1], pop = sample(swans[, 2]), NA.char = "NA", ploidy = 2, type = "codom"))) pvalues<-sapply(1:nrow(test),function(e)as.randtest(test[e,],obs[e])$pvalue) return(pvalues) } #Black_Swans_SpaceMixAnalysis run.spacemix.analysis( n.fast.reps = 0, fast.MCMC.ngen = 1e5, fast.model.option = "source_and_target", long.model.option = "source_and_target", data.type = "counts", fast.likelihood.option = "wishart", long.likelihood.option = "wishart", proj.mat.option=NULL, sample.frequencies=NULL, mean.sample.sizes=NULL, counts = swanpop@tab, sample.sizes = size, sample.covariance=NULL, target.spatial.prior.scale=NULL, source.spatial.prior.scale=NULL, spatial.prior.X.coordinates = coorpop[,2], spatial.prior.Y.coordinates = coorpop[,1], round.earth = FALSE, long.run.initial.parameters=NULL, k = nrow(swanpop@tab), #number of samples loci = 13, #number of loci ngen = 1e7, #number of MCMC gnereations for the long MCMC printfreq = 1e3, #frequency with which updates are printed samplefreq = 1e4, #frequency with which samples are logged from the MCMC (basically the thinning) mixing.diagn.freq = 50, savefreq=5e5, #frequency with which MCMC_output object is saved directory=NULL, #directory into which you want output to be saved prefix="Swan_example") #prefix to be attached to all output files plot(last.params$population.coordinates[1:last.params$k,],type='n', ylab="LAT",xlab="LONG", xlim=c(min(coorpop[,2]),max(coorpop[,2])), ylim=c(min(coorpop[,1]),max(coorpop[,1]))) text(coorpop[,c(2,1)],labels=c("BA","NZ","QL","TS","WA","WL","WE")) text(last.params$population.coordinates[1:last.params$k,], col=pop.cols,labels=c("BA","NZ","QL","TS","WA","WL","WE")) arrows(x0 = last.params$population.coordinates[(last.params$k+1):(2*last.params$k),1], y0 = last.params$population.coordinates[(last.params$k+1):(2*last.params$k),2], x1 = last.params$population.coordinates[1:last.params$k,1], y1 = last.params$population.coordinates[1:last.params$k,2], lwd = 1, col = pop.cols)