Skip to content

Commit 65887a3

Browse files
committed
Fixed issue when only one topology in 95% CI of topo space
Thanks to Arun Prasanna for bringing this to our attention. There was a bug in makeplot.pairs that was causing RWTY to barf when the chain wasn’t exploring enough topologies to have any variation in the 95% CI of topology space. Basically it was trying to color a histogram with categories of zero width. It’s still not producing density plots for correlations between topology and other parameters under these circumstances, but those would basically be a vertical line anyway. Might be worth revisiting, but for now it’s at least functional.
1 parent 0a57d68 commit 65887a3

File tree

1 file changed

+11
-3
lines changed

1 file changed

+11
-3
lines changed

R/makeplot.pairs.R

Lines changed: 11 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -48,7 +48,7 @@ makeplot.pairs <- function(chains, burnin = 0, treedist = 'PD', params = NA, str
4848

4949
problems = setdiff(params, param.names)
5050
if(length(problems)>0){
51-
stop(paste(c("The following names you suppied are not parameters in your parameter table ", "'", problems, "'")))
51+
stop(paste(c("The following names you suppied are not parameters in your parameter table ", "'", paste(problems, " "), "'")))
5252
}
5353
param.names = params
5454
}
@@ -81,7 +81,6 @@ do.pairs.plot <- function(chain, burnin = 0, params, treedist){
8181
# this is so that the distances look nicer in the plots
8282
focal.tree = chains[[1]]$trees[1]
8383
distances = tree.distances.from.first(chains, burnin, focal.tree = focal.tree, treedist = treedist)
84-
8584
ptable$topological.distance = distances$topological.distance
8685

8786
points <- function(data, mapping, ...) {
@@ -96,7 +95,16 @@ do.pairs.plot <- function(chain, burnin = 0, params, treedist){
9695
var = data[,as.character(mapping$x)]
9796
lower = quantile(var, c(0.025))
9897
upper = quantile(var, c(0.975))
99-
fill = as.numeric(cut(var, c(lower, upper)))
98+
if(lower == upper){
99+
ci.width <- round(.025 * length(var))
100+
fill = c(rep(NA, ci.width),
101+
rep(1, length(var) - 2 * ci.width),
102+
rep(NA, ci.width))
103+
}
104+
else{
105+
fill = as.numeric(cut(var, c(lower, upper)))
106+
}
107+
100108
fill[which(is.na(fill))] = 'red'
101109
fill[which(fill == 1)] = 'blue'
102110

0 commit comments

Comments
 (0)