hiatus bug triggered by level

[tsdye]

Describe the bug
A clear and concise description of what the bug is.

hiatus returns an error at a one level, but exits correctly with a slightly different level.

To Reproduce
A minimal reproducible example (AKA a reprex).

hiatus(eve, level=0.950)

Error in if (p[2, i] == p[1, i]) return(no_hiatus) : 
  the condition has length > 1

hiatus(eve, level=0.949)

|                                         |   start |     end |         duration |
|-----------------------------------------+---------+---------+------------------|
| Austral Islands settlement-Paleotsunami | 1177.91 | 1423.34 | 246.427260273972 |

hiatus(eve, level=0.951)

|                                         |      start |              end |         duration |
|-----------------------------------------+------------+------------------+------------------|
| Austral Islands settlement-Paleotsunami | 1178.58041 | 1423.33991657977 | 245.756766853743 |

Expected behavior
A clear and concise description of what you expected to happen.

All levels should exit correctly.

Test output
Output of tinytest::test_package("ArchaeoPhases")
tinytest::test_package("ArchaeoPhases")
test_activity.R............... 2 tests OK 0.3s
test_allen.R.................. 11 tests OK 15ms
test_boundaries.R............. 1 tests OK 0.2s
test_events.R................. 1 tests OK 30ms
test_hiatus.R................. 2 tests OK 0.2s
test_interval.R............... 4 tests OK 93ms
test_occurrence.R............. 1 tests OK 58ms
test_older.R.................. 2 tests OK 26ms
test_phases.R................. 2 tests OK 0.2s
test_summary.R................ 2 tests OK 0.2s
test_tempo.R.................. 1 tests OK 0.2s
test_transition.R............. 1 tests OK 68ms
All ok, 30 results (1.7s)

Screenshots
If applicable, add screenshots to help explain your problem.

Additional context
Add any other context about the problem here.

This bug is also present in ArchaeoPhases 1.8.
I can send the eve data, if need be.

[nfrerebeau]

I can't reproduce it on my side, could you send me your eve data?

[tsdye]

Here they are, saved with write.csv. Column names didn't make it into the csv file, but V1 is Austral Islands settlement and V2 is Paleotsunami.

eve.csv

FYI, here is the procedure used to create eve, with a slightly altered return value from the output sent yesterday:

  load_archaeophases(ver)
  library(ArchaeoPhases)
  events <- read_oxcal(mcmc)
  ind <- which(names(events) %in% c(early, late))
  eve <- events[, ind]
  ret <- if (ver == "1.8") dates_hiatus(eve[,1], eve[,2], level = lev) else hiatus(eve, level=lev)
  if (ver == "1.8") ret$duration else round(as.data.frame(ret), 0)

[nfrerebeau]

It seems to be a corner case. To determine the hiatus, the duration of all possible intervals is calculated to determine the interval with the maximum length (if any). In your case, two intervals have the same duration with very slightly different end points (less than a year).

An easy (ugly) solution is to add very low random noise (of the order of a thousandth of a day, which is insignificant for archaeological data):

eve <- read.csv(file = "eve.csv")
eve <- eve + rnorm(nrow(eve)) / 1000
eve <- as_events(eve, calendar = NULL)
plot(eve)

hia950 <- hiatus(eve, level=0.950)
as.data.frame(hia950)

hia949 <- hiatus(eve, level=0.949)
as.data.frame(hia949)

hia951 <- hiatus(eve, level=0.951)
as.data.frame(hia951)

Clearly not a viable solution...

Maybe we should keep only the first interval? Or compute the mean endpoints of all possible intervals?

[tsdye]

Thanks for the analysis. In my experience the bug bites frequently. I haven't worked through the code and can't recommend a solution. Perhaps check with Anne?TomSent from my iPhoneOn Nov 6, 2024, at 6:31 AM, N. Frerebeau ***@***.***> wrote: It seems to be a corner case. To determine the hiatus, the duration of all possible intervals is calculated to determine the interval with the maximum length (if any). In your case, two intervals have the same duration with very slightly different end points (less than a year). An easy (ugly) solution is to add very low random noise (of the order of a thousandth of a day, which is insignificant for archaeological data): eve <- read.csv(file = "eve.csv") eve <- eve + rnorm(nrow(eve)) / 1000 eve <- as_events(eve, calendar = NULL) plot(eve) hia950 <- hiatus(eve, level=0.950) as.data.frame(hia950) hia949 <- hiatus(eve, level=0.949) as.data.frame(hia949) hia951 <- hiatus(eve, level=0.951) as.data.frame(hia951) Clearly not a viable solution... Maybe we should keep only the first interval? Or compute the mean endpoints of all possible intervals? —Reply to this email directly, view it on GitHub, or unsubscribe.You are receiving this because you authored the thread.Message ID: ***@***.***>

[tsdye]

After line 31 of hiatus.R, make sure i is length 1 by keeping the first interval.

i <- if(length(i) > 1) i[1] else i

This corrects the error and seems to yield a reasonable result.

Any thoughts?

[nfrerebeau]

This seems to be a possible solution. I've implemented it as you suggest. You can install the development version from GitHub, to see if it works.

[tsdye]

I'm still able to push changes :/

Please review the change I pushed, which solves the error.

If the change is OK from a statistical point of view, then please feel free to close this issue.

Thanks for your help, and for keeping me as a contributor.

[nfrerebeau]

It seems OK, thanks! We can keep this issue open for a while longer, just to make sure it's the best approach.