# NOT RUN {
library(HelpersMG)
# _______________________________________________________________
# right censored distribution with gamma distribution
# _______________________________________________________________
# Detection limit
DL <- 100
# Generate 100 random data from a gamma distribution
obc <- rgamma(100, scale=20, shape=2)
# remove the data below the detection limit
obc[obc>DL] <- +Inf
# search for the parameters the best fit these censored data
result <- cutter(observations=obc, upper_detection_limit=DL,
cut_method="censored")
result
plot(result, xlim=c(0, 150), breaks=seq(from=0, to=150, by=10))
# _______________________________________________________________
# The same data seen as truncated data with gamma distribution
# _______________________________________________________________
obc <- obc[is.finite(obc)]
# search for the parameters the best fit these truncated data
result <- cutter(observations=obc, upper_detection_limit=DL,
cut_method="truncated")
result
plot(result, xlim=c(0, 150), breaks=seq(from=0, to=150, by=10))
# _______________________________________________________________
# left censored distribution with gamma distribution
# _______________________________________________________________
# Detection limit
DL <- 10
# Generate 100 random data from a gamma distribution
obc <- rgamma(100, scale=20, shape=2)
# remove the data below the detection limit
obc[obc<DL] <- -Inf
# search for the parameters the best fit these truncated data
result <- cutter(observations=obc, lower_detection_limit=DL,
cut_method="censored")
result
plot(result, xlim=c(0, 200), breaks=seq(from=0, to=200, by=10))
# _______________________________________________________________
# left and right censored distribution
# _______________________________________________________________
# Generate 100 random data from a gamma distribution
obc <- rgamma(100, scale=20, shape=2)
# Detection limit
LDL <- 10
# remove the data below the detection limit
obc[obc<LDL] <- -Inf
# Detection limit
UDL <- 100
# remove the data below the detection limit
obc[obc>UDL] <- +Inf
# search for the parameters the best fit these censored data
result <- cutter(observations=obc, lower_detection_limit=LDL,
upper_detection_limit=UDL,
cut_method="censored")
result
plot(result, xlim=c(0, 150), col.DL=c("black", "grey"),
col.unobserved=c("green", "blue"),
breaks=seq(from=0, to=150, by=10))
# _______________________________________________________________
# Example with two values for lower detection limits
# corresponding at two different methods of detection for example
# with gamma distribution
# _______________________________________________________________
obc <- rgamma(50, scale=20, shape=2)
# Detection limit for sample 1 to 50
LDL1 <- 10
# remove the data below the detection limit
obc[obc<LDL1] <- -Inf
obc2 <- rgamma(50, scale=20, shape=2)
# Detection limit for sample 1 to 50
LDL2 <- 20
# remove the data below the detection limit
obc2[obc2<LDL2] <- -Inf
obc <- c(obc, obc2)
# search for the parameters the best fit these censored data
result <- cutter(observations=obc,
lower_detection_limit=c(rep(LDL1, 50), rep(LDL2, 50)),
cut_method="censored")
result
# It is difficult to choose the best set of colors
plot(result, xlim=c(0, 150), col.dist="red",
col.unobserved=c(rgb(red=1, green=0, blue=0, alpha=0.1),
rgb(red=1, green=0, blue=0, alpha=0.2)),
col.DL=c(rgb(red=0, green=0, blue=1, alpha=0.5),
rgb(red=0, green=0, blue=1, alpha=0.9)),
breaks=seq(from=0, to=200, by=10))
# ___________________________________________________________________
# left censored distribution comparison of normal, lognormal and gamma
# ___________________________________________________________________
# Detection limit
DL <- 10
# Generate 100 random data from a gamma distribution
obc <- rgamma(100, scale=20, shape=2)
# remove the data below the detection limit
obc[obc<DL] <- -Inf
# search for the parameters the best fit these truncated data
result_gamma <- cutter(observations=obc, lower_detection_limit=DL,
cut_method="censored", distribution="gamma")
result_gamma
plot(result_gamma, xlim=c(0, 250), breaks=seq(from=0, to=250, by=10))
result_lognormal <- cutter(observations=obc, lower_detection_limit=DL,
cut_method="censored", distribution="lognormal")
result_lognormal
plot(result_lognormal, xlim=c(0, 250), breaks=seq(from=0, to=250, by=10))
result_normal <- cutter(observations=obc, lower_detection_limit=DL,
cut_method="censored", distribution="normal")
result_normal
plot(result_normal, xlim=c(0, 250), breaks=seq(from=0, to=250, by=10))
compare_AICc(gamma=result_gamma,
lognormal=result_lognormal,
normal=result_normal)
# ___________________________________________________________________
# Test for similarity in gamma left censored distribution between two
# datasets
# ___________________________________________________________________
obc1 <- rgamma(100, scale=20, shape=2)
# Detection limit for sample 1 to 50
LDL <- 10
# remove the data below the detection limit
obc1[obc1<LDL] <- -Inf
obc2 <- rgamma(100, scale=10, shape=2)
# remove the data below the detection limit
obc2[obc2<LDL] <- -Inf
# search for the parameters the best fit these censored data
result1 <- cutter(observations=obc1,
distribution="gamma",
lower_detection_limit=LDL,
cut_method="censored")
logLik(result1)
plot(result1, xlim=c(0, 200),
breaks=seq(from=0, to=200, by=10))
result2 <- cutter(observations=obc2,
distribution="gamma",
lower_detection_limit=LDL,
cut_method="censored")
logLik(result2)
plot(result2, xlim=c(0, 200),
breaks=seq(from=0, to=200, by=10))
result_totl <- cutter(observations=c(obc1, obc2),
distribution="gamma",
lower_detection_limit=LDL,
cut_method="censored")
logLik(result_totl)
plot(result_totl, xlim=c(0, 200),
breaks=seq(from=0, to=200, by=10))
compare_AICc(Separate=list(result1, result2),
Common=result_totl, factor.value=1)
compare_BIC(Separate=list(result1, result2),
Common=result_totl, factor.value=1)
# }
Run the code above in your browser using DataLab