I am trying to determine the population median, 10th, 20th, 30th, and 40th percentiles (+/- 95%CI) based on 71 samples of sediment total DDT(ug/Kg).
Those data that I am using have 71 observations total, including 60 censored values (84.5% censoring). Censored values include 59 values below three different detection levels (<1, <0.05, <0.267-ug/Kg) and one estimated value (0.463-ug/Kg J) ---- between method detection level (MDL) and reporting level (RL). 11 observations were above detection level.
Is it possible to even calculate any defensible summary statistics with this level of censoring?
Rather than trying to calculate summary statistics an alternative would be to describe the population as having 84.5 percent values below detection and 90th and 95th percentile values of X and Y. However using the NADA library in R the the output for the 95%CI values does not seem to make sense. Thank you for your help and feedback.
Below is the a REPREX of my code in R with data:
#Original values:
#<1,<1,<1,<1,<1,<1,<1,<1,<1,<1,<1,<1,<1,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,
#<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.267,<0.267,<0.267,<0.267,<0.267,<0.267,
#<0.267,<0.267,<0.267,<0.267,<0.267,<0.267,<0.267,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,
#<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,<0.05,
#0.463J,0.11,0.2,0.25,0.3,0.39,0.671,0.93,1.72,3.81,5.2,10.3
#Load data into data frame with one column containing values and another indicating censoring
tddt <-c(1,1,1,1,1,1,1,1,1,1,1,1,1,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,
0.05,0.05,0.05,0.05,0.267,0.267,0.267,0.267,0.267,0.267,0.267,0.267,0.267,0.267,0.267,
0.267,0.267,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,0.05,
0.05,0.05,0.05,0.05,0.463,0.11,0.2,0.25,0.3,0.39,0.671,0.93,1.72,3.81,5.2,10.3)
tddtcens <- c(rep(T,60),rep(F,11))
tddt.df <- data.frame(tddt,tddtcens)
#NADA Library
#Kaplan-Meier - nonparametric method
library(NADA)
attach(tddt.df)
#Data are not normally distributed
shapiro.test(tddt)
censummary(tddt,tddtcens)
censtats(tddt,tddtcens)
tddtkm = cenfit(tddt,tddtcens)
tddtkm
plot(tddtkm)
#Robust ROS - semi-parametric method
tddtros = cenros(tddt,tddtcens)
tddtros
plot(tddtros)
#Maximum Likelyhood Estimate - parametric method - not applicable
tddtmle = cenmle(tddt,tddtcens)
tddtmle
#All 3 Methods
#Maximum likelyhood, kaplan meier, robust regression on order statistics
censummary(tddt,tddtcens)
censtats(tddt,tddtcens)
#Confidence Interval
tddtkm.ci = cenfit(tddt,tddtcens,conf.int=0.95)
tddtkm.ci
mean(tddtkm.ci)
median(tddtkm.ci)
#S4 Class object - access to slot w median + 95CI
tddtkm@survfit
#Quantiles
quantile(tddtkm, c(.1,.2,.3,.4,.5,.8,.82,.85,.90,.95), conf.int=T)
tddt.km.df <- data.frame(quantile(tddtkm, c(.1,.2,.3,.4,.5,.82,.85,.90,.95), conf.int=T))
tddt.km.df