agnostic<-function(onoff,total,prior,uniform) {
if (missing(total) || total == '')
total = rep(1,length(onoff))
if (sum(total) == 0)
return(rep(NA,length(onoff)))
if (missing(prior) || prior == '')
prior = 0.5
if (missing(uniform) || uniform != T)
uniform = F
totalBins = length(onoff)
bins = length(onoff)
for (i in 1:bins)
if (total[i] > 0)
lastControl = i
bins = lastControl
onoff = onoff[1:bins]
total = total[1:bins]
if (onoff[bins] > 0)
return(rep(0,totalBins))
withData = 0
for (i in 1:bins)
if (onoff[i] > 0)
withData = withData + 1
if (withData <= 1)
return(rep(NA,totalBins))
first = 0
last = 0
for (i in 1:bins)
if (onoff[i] > 0) {
if (first == 0)
first = i
last = i
}
n = sum(onoff)
# simple likelihood of the data assuming no extinction
# subtract one from the denominator because the last occurrence must be
# held in place
c = choose(sum(total[1:last]) - 1 , n - 1)
aliveLike = c / choose(sum(total) , n)
# key computation: find likelihoods assuming variable-length
# durations before extinction
# the likelihood contribution of the hypothesized extinctions
# before the last appearance is of course zero
# the prior extinction probability is a function of the exponential
# decay coefficient mu
extinctLike = rep(0,totalBins)
mu = log(1 - prior) / bins
for (i in (last+1):bins) {
if (uniform == T)
extinctLike[i] = prior * c / choose(sum(total[1:(i-1)]) , n) / bins
else
extinctLike[i] = c / choose(sum(total[1:(i-1)]) , n) * (exp((i - 1) * mu) - exp(i * mu))
}
# compute interval-by-interval posteriors
post = extinctLike / (sum(extinctLike) + (1 - prior) * aliveLike)
# compute running total
post = cumsum(post)
if (bins < totalBins)
post[(bins+1):totalBins] = post[length(post)]
return(post)
}