beast2 icon indicating copy to clipboard operation
beast2 copied to clipboard

Unable to plot EBSP in R or use EBSPAnalyser with output of EBSP tutorial

Open TeresaPegan opened this issue 3 years ago • 1 comments

Hello, I would like to use beast's EBSP functionalities and I have been following along the EBSP tutorial available from the website. I tried with both my own data and with the data from the tutorial, and followed the tutorial instructions exactly when using the tutorial data. In both cases, I cannot get either the "plotEBSP.R" script (available in the tutorial) or the EBSPAnalyser to work. Here is what happens in each case with the tutorial data.

R:

> source("/Applications/BEAST 2.6.4/plotEBSP.R")
> plotEBSP("~/Downloads/ebsp2-tut/mystery-mammal/mystery.log", useHPD=FALSE, log="y")
Error in changeTimes[[j]] : subscript out of bounds
Called from: is.unsorted(vec)

EBSPAnalyser:

Processing /Users/teresapegan/Downloads/ebsp2-tut/mystery-mammal/mystery.log
 skipping 1000 line


time	mean	median	95HPD lower	95HPD upper
Exception: Index -1 out-of-bounds for length 0

I dug into the R code a bit more. Here is the traceback:

> traceback()
5: is.unsorted(vec)
4: findInterval(t, changeTimes) at plotEBSP.R#3
3: getPopSize(allTimes[i], changeTimes[[j]], changePops[[j]], isLinear) at plotEBSP.R#105
2: processEBSPdata(df, isLinear) at plotEBSP.R#158
1: plotEBSP("~/Downloads/ebsp2-tut/mystery-mammal/mystery.log", 
       useHPD = FALSE, log = "y")

The problem is in the getPopSize() function, but its root seems to be further back in the processEBSPdata() function. getPopSize() is failing because the object changeTimes is an empty list. Within processEBSPdata(), it looks like changeTimes should be filled in by this snippet of code:

 changeTimes <- list()
    changePops <- list()

    Nmedian <- rep(0, nTimes)
    NupperHPD <- rep(0, nTimes)
    NlowerHPD <- rep(0, nTimes)
    NupperCPD <- rep(0, nTimes)
    NlowerCPD <- rep(0, nTimes)

    for (i in 1:frameLen) {
        theseChangeTimes <- NULL
        theseChangePops <- NULL
        p <- 0
        for (j in 1:nTimes) {
            pair = strsplit(as.character(df[i,1+j]), ":")[[1]]
            time <- as.double(pair[1])
            allTimes[j] <- allTimes[j] + time

            if (length(pair)>1) {
                p <- p + 1
                theseChangeTimes[p] <- time
                theseChangePops[p] <- as.double(pair[2])
            }
        }

        changeTimes[[i]] <- theseChangeTimes
        changePops[[i]] <- theseChangePops
    }

The problem seems to be that length(pair) is never greater than one, so changeTimes[[i]] is always left NULL. In the (j in 1:nTimes) loop, the function seems to be expecting to find values within the log file that can be separated into two objects by a ":", but there aren't any?

I'm not too sure where to go from here or what exactly is triggering the EBSPAnalyser error, although there also it seems like it's trying to interact with an empty object or something. I have attached the xml file I generated from the tutorial.

Let me know if there are any ideas for troubleshooting this, thanks! -Teresa mystery.log

TeresaPegan avatar May 22 '21 15:05 TeresaPegan

Hi Teresa, I cannot replicate your error. Can you please confirm that your mystery.log is the EBSPLogger file (note that by default, unless you have changed it, this is called EBSP.log)? There is an issue with the current code when calling plotEBSPTimesHist() but the error is different (non-character argument).

carlopacioni avatar Jun 15 '21 23:06 carlopacioni