71 update the bmd lower boundary censoring code min conc = 0 for lower threshold

[gracezhihuizhao]

Bug: BMD lower censoring threshold should be set based on the lowest experimental dose. However if the data being passed into tcplfit2_core and tcplhit2_core include the untreated control group (conc = 0), the BMD lower bound censoring is not working as intended because the current code min(conc) will return 0.

• Updated the BMD censoring data such that if the data contain untreated control group (i.e. conc = 0), BMD lower bond censoring will based on the lowest experimental concentration that is not 0.
  • Added a unit test to ensure BMD censoring is working appropriately after the update.

[sedavid01]

Code used to check the lower bmd censoring update

Load the current update for BMD lower boundaries

devtools::load_all()

Simulated Data

X <- rep(seq(0,10,length.out = 10),each = 5)

set.seed(111) PARAMS <- data.frame( tp = rep(15,5), ga = c(0.25,1,4,9,20), p = c(0.15,runif(4,min = 0.2,max = 0.7)), er = rep(0.1,5) ) seeds <- c(842,1213,112,922,203) Y <- matrix(nrow = 5,ncol = 5*10) for(i in 1:5){ set.seed(seeds[i]) Y[i,] <- tcplfit2::hillfn(ps = as.vector(t(PARAMS[i,])),x = X) + rt(n = length(X),df = 4) }

Dose-response Fitting

tfit <- vector(length = 5,mode = "list") for(i in 1:5){ tfit[[i]] <- tcplfit2::tcplfit2_core(conc = X,resp = as.vector(t(Y[i,])), cutoff = 2, force.fit = TRUE, bidirectional = TRUE, verbose = FALSE) }

Dose-response Hitting

thit_bound <- vector(length = 5,mode = "list") thit_nobounds <- vector(length = 5,mode = "list") for(i in 1:5){ thit_nobounds[[i]] <- tcplfit2::tcplhit2_core(tfit[[i]],conc = X,resp = as.vector(t(Y[i,])),cutoff = 2,onesd = sd(Y[1:5]),bmr_scale = 1.349) thit_bound[[i]] <- tcplfit2::tcplhit2_core(tfit[[i]],conc = X,resp = as.vector(t(Y[i,])),cutoff = 2,onesd = sd(Y[1:5]),bmr_scale = 1.349,bmd_low_bnd = 0.1,bmd_up_bnd = 10) }

ci_diff <- c() for(i in 1:5){ bound <- thit_bound[[i]]$bmdu - thit_bound[[i]]$bmdl nobound <- thit_nobounds[[i]]$bmdu - thit_nobounds[[i]]$bmdl if(i == 1){ci_diff <- bound - nobound}else{ci_diff <- c(ci_diff,bound-nobound)} }

View(rbind.data.frame(thit_bound[[1]],thit_nobounds[[1]])) View(rbind.data.frame(thit_bound[[2]],thit_nobounds[[2]])) View(rbind.data.frame(thit_bound[[3]],thit_nobounds[[3]])) View(rbind.data.frame(thit_bound[[4]],thit_nobounds[[4]])) View(rbind.data.frame(thit_bound[[5]],thit_nobounds[[5]]))

[sedavid01]

@brown-jason, changes are merged to 'master' please verify when you get the chance and if any branch clean-up is needed please either do so or provide guidance/instructions. Thanks, -Sarah