# Tasimboko females first birth interval (from marriage to first baby in month)
# 2 generations with interval censoring
# https://doi.org/10.24454/jps.27.0_7
bid <- data.frame(
 gen = factor(c(rep(1, 21), rep(2, 19)), labels=c("old", "new")),
 start = c(25, 37, 22, 30, 44, 13, 90, 24, 49, 25, 37, 25, 8, 101, 37, 22, 
  40, 15, 23, 24, 28, 81, 1, 1, 17, 14, 22, 37, 12, 15, 19, 26, 29, 23,
  20, 18, 9, 9, 3, 2),
 last = c(47, 59, 33, 41, 55, 35, 101, 35, 71, 47, 59, 47, 19, 101, 37, 22,
  40, 15, 23, 24, 28, 92, 23, 23, 17, 14, 22, 37, 12, 15, 19, 26, Inf, Inf,
  Inf, Inf, Inf, Inf, Inf, Inf),
 ef = c(3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1,
  3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0))

library(survival)
rKM <- survfit(Surv(start, last, ef, type="interval") ~ gen, data=bid)
print(rKM)
summary(rKM)
png("BI_TF_KM.png", width=800, height=800)
par(cex=1.5)
plot(rKM, lty=1:2, lwd=2, col=2:3, mark.time=TRUE, conf.int=FALSE,
 main="Birth intervals between 1st and 2nd babies in\n Tasimboko females by 1st baby\n was before 1985 or not")
legend("topright", lty=1:2, lwd=2, col=2:3,
 legend=c("1st baby before 1985", "1st baby in 1986 or later"))
dev.off()
rWB <- survreg(Surv(start, last, ef, type="interval") ~ gen, 
 data=bid, dist="weibull")
print(rWB)
summary(rWB)
# rCPH <- coxph(Surv(start, last, ef, type="interval") ~ gen, data=bid)
# It gives the following error
# > rCPH <- coxph(Surv(start, last, ef, type="interval") ~ gen, data=bid)
# Error in coxph(Surv(start, last, ef, type = "interval") ~ gen, data = bid) : 
#  Cox model doesn't support "interval" survival data
library(icenReg)
rCPH <- ic_sp(Surv(start, last, ef, type="interval") ~ gen, data=bid,
              model = "ph",
              bs_samples = 500)
summary(rCPH) # same as print(rCPH)
#
