From d5b42cf8944f16a4b42ae4ea42396cf6291eeb6c Mon Sep 17 00:00:00 2001 From: PharmCat <13901158+PharmCat@users.noreply.github.com> Date: Mon, 28 Nov 2022 16:28:49 +0300 Subject: [PATCH 1/2] remove nca --- Project.toml | 6 +- docs/make.jl | 5 +- docs/src/nca.md | 132 +-- docs/src/ncatable.jl | 1 - docs/src/pd.md | 9 +- docs/src/pk.md | 11 - docs/src/pkplot.md | 30 - src/ClinicalTrialUtilities.jl | 5 +- src/pk.jl | 1505 --------------------------------- src/plots.jl | 426 ---------- test/csv/nanpk.csv | 22 - test/csv/pddata.csv | 14 - test/csv/pkdata.csv | 22 - test/csv/pkdata2.csv | 161 ---- test/csv/upkdata.csv | 6 - test/internal.jl | 22 - test/pktest.jl | 1264 --------------------------- test/test.jl | 50 -- test/testdata.jl | 13 +- test/validation_report.jmd | 142 ---- test/weave.jl | 5 - 21 files changed, 14 insertions(+), 3837 deletions(-) delete mode 100644 docs/src/ncatable.jl delete mode 100644 docs/src/pk.md delete mode 100644 docs/src/pkplot.md delete mode 100644 src/pk.jl delete mode 100644 src/plots.jl delete mode 100644 test/csv/nanpk.csv delete mode 100644 test/csv/pddata.csv delete mode 100644 test/csv/pkdata.csv delete mode 100644 test/csv/pkdata2.csv delete mode 100644 test/csv/upkdata.csv delete mode 100644 test/pktest.jl delete mode 100644 test/validation_report.jmd delete mode 100644 test/weave.jl diff --git a/Project.toml b/Project.toml index 3e81c6f..eab9520 100644 --- a/Project.toml +++ b/Project.toml @@ -1,14 +1,13 @@ name = "ClinicalTrialUtilities" uuid = "535c2557-d7d0-564d-8ff9-4ae146c18cfe" authors = ["Vladimir Arnautov (mail@pharmcat.net)"] -version = "0.6.6" +version = "0.7.0" [deps] DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0" Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" -RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" Roots = "f2b01f46-fcfa-551c-844a-d8ac1e96c665" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" @@ -17,8 +16,7 @@ StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" DataFrames = "1" Distributions = "0.20, 0.21, 0.22, 0.23, 0.24, 0.25" QuadGK = "1, 2" -RecipesBase = "1" -Roots = "1.0, 2" +Roots = "1, 2" SpecialFunctions = "1, 2" StatsBase = "0.30, 0.31, 0.32, 0.33" julia = "1" diff --git a/docs/make.jl b/docs/make.jl index 31b60be..9af94fa 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -12,10 +12,7 @@ makedocs( "Power" => "power.md", "Confidence intervals" => "ci.md", "Descriptive statistics" => "ds.md", - "NCA" => ["NCA" => "nca.md", - "Pharmakokinetics" => "pk.md", - "Pharmacodynamics" => "pd.md", - "PK Plots" => "pkplot.md"], + "NCA" => "nca.md", "Randomization" => "random.md", "Simulations" => "sim.md", "Utilities" => "utils.md", diff --git a/docs/src/nca.md b/docs/src/nca.md index 8aa93f5..963c351 100644 --- a/docs/src/nca.md +++ b/docs/src/nca.md @@ -1,133 +1,7 @@ # Non-compartment analysis -*See [MetidaNCA.jl](https://github.com/PharmCat/MetidaNCA.jl). MetidaNCA.jl will replace CTU in plans.* +Further development of NCA PK/PD will be based on [MetidaNCA.jl](https://github.com/PharmCat/MetidaNCA.jl) package. -NCA analysis based on following steps: +All NCA PK/PD functions moved to [MetidaNCA](https://github.com/PharmCat/MetidaNCA.jl). -1. Loadind data to DataFrame; -2. Constructing subjects list with pkimport/pdimport; -3. Run NCA; -4. Exporting to DataFrame; -5. Descriptive statistics / HTML export. - -### Description - -```julia -include("ncatable.jl") -``` -#### AUC - -```math -AUC = \sum_{n=1}^N AUC_{n} -``` - -Where `AUCn` - partial AUC. - -Linear trapezoidal rule - -```math -AUC\mid_{t_1}^{t_2} = \delta t \times \frac{C_1 + C_2}{2} - -AUMC\mid_{t_1}^{t_2} = \delta t \times \frac{t_1 \times C_1 + t_2 \times C_2}{2} -``` - -Logarithmic trapezoidal rule - -```math -AUC\mid_{t_1}^{t_2} = \delta t \times \frac{ C_2 - C_1}{ln(C_2/C_1)} - -AUMC\mid_{t_1}^{t_2} = \delta t \times \frac{t_2 \times C_2 - t_1 \times C_1}{ln(C_2/C_1)} - \delta t^2 \times \frac{ C_2 - C_1}{ln(C_2/C_1)^2} -``` - -Linear interpolation rule - -```math -C_x = C_1 + \frac{(t_x-t_1)\times(C_2 - C_1)}{t_2 - t_1} -``` - -Logarithmic interpolation rule - -```math -C_x = exp\left(ln(C_1) + \frac{(t_x-t_1)\times(ln(C_2) - ln(C_1))}{t_2 - t_1}\right) -``` - -#### MRTlast - -```math -MRTlast = AUMClast / AUClast -``` - -#### ```math \lambda_z ```- elimination constant - -#### HL - -```math -HL = ln(2) / \lambda_z -``` - -#### AUCinf - -```math -AUCinf = AUClast + \frac{Clast}{\lambda_z} -``` - -#### AUMCinf - -```math -AUMCinf = AUMClast + \frac{tlast\times Clast}{\lambda_z} + \frac{Clast}{\lambda_z^2} -``` - -#### AUCpct - -```math -AUCpct = (AUCinf - AUClast) / AUCinf * 100.0% -``` - -#### Accumulation index - -```math -Accind = \frac{1}{1 - exp(-\lambda_z \tau)} -``` - -### nca! -```@docs -ClinicalTrialUtilities.nca! -``` - -### Scenario - -1 Loading DataFrame - -```julia -using CSV, DataFrames, ClinicalTrialUtilities -pkdata2 = CSV.File("pkdata.csv") |> DataFrame -``` -2 Subject list - -```julia - pkds = pkimport(pkdata2, [:Subject, :Formulation]; time = :Time, conc = :Concentration) -``` - -3 NCA analysis with default settings - -```julia - pk = nca!(pkds) -``` - -4 Exporting - -```julia - ncadf = DataFrame(pk; unst = true) -``` - -5 Descriptive statistics - -```julia - ds = ClinicalTrialUtilities.descriptive(ncadf, stats = [:n, :mean, :sd], sort = [:Formulation]) -``` - -6 Exporting - -```julia - dsdf = ClinicalTrialUtilities.DataFrame(ds; unst = true) -``` +Validation included. diff --git a/docs/src/ncatable.jl b/docs/src/ncatable.jl deleted file mode 100644 index 8b13789..0000000 --- a/docs/src/ncatable.jl +++ /dev/null @@ -1 +0,0 @@ - diff --git a/docs/src/pd.md b/docs/src/pd.md index 9795a6a..d0a5432 100644 --- a/docs/src/pd.md +++ b/docs/src/pd.md @@ -1,6 +1,7 @@ ### Pharmacodynamics -#### pdimport -```@docs -ClinicalTrialUtilities.pdimport -``` +Further development of NCA PK/PD will be based on [MetidaNCA.jl](https://github.com/PharmCat/MetidaNCA.jl) package. + +All NCA PK/PD functions moved to [MetidaNCA](https://github.com/PharmCat/MetidaNCA.jl). + +Validation included. diff --git a/docs/src/pk.md b/docs/src/pk.md deleted file mode 100644 index 8dc1934..0000000 --- a/docs/src/pk.md +++ /dev/null @@ -1,11 +0,0 @@ -### Pharmacokinetics - -#### pkimport -```@docs -ClinicalTrialUtilities.pkimport -``` - -#### auc_sparse -```@docs -ClinicalTrialUtilities.auc_sparse -``` diff --git a/docs/src/pkplot.md b/docs/src/pkplot.md deleted file mode 100644 index 19bf9b3..0000000 --- a/docs/src/pkplot.md +++ /dev/null @@ -1,30 +0,0 @@ -### PK plots - -#### pkplot -```@docs -ClinicalTrialUtilities.pkplot -``` - -#### Examples - - -```@example pkplots -using ClinicalTrialUtilities, DataFrames, CSV, Plots -pkdatapath = joinpath(dirname(pathof(ClinicalTrialUtilities)))*"/../test/csv/pkdata2.csv" -pkdata = CSV.File(pkdatapath) |> DataFrame -pkds = ClinicalTrialUtilities.pkimport(pkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) -plot1 = pkplot(pkds[1], legend = false); -plot2 = pkplot(pkds; pagesort = [:Formulation], typesort = [:Subject])[1]; -``` - -Plot for subject: - -```@example pkplots -plot1 -``` - -First plot for DataSet: - -```@example pkplots -plot2 -``` diff --git a/src/ClinicalTrialUtilities.jl b/src/ClinicalTrialUtilities.jl index 7dbcbd9..42bd94e 100644 --- a/src/ClinicalTrialUtilities.jl +++ b/src/ClinicalTrialUtilities.jl @@ -63,8 +63,7 @@ include("info.jl") #Descriptive statistics include("descriptives.jl") -#PK -include("pk.jl") + #Export include("export.jl") #Randomization @@ -73,8 +72,6 @@ include("randomization.jl") include("utilities.jl") #Show include("show.jl") -#Plots -include("plots.jl") #Deprecated include("deprecated.jl") diff --git a/src/pk.jl b/src/pk.jl deleted file mode 100644 index 6d26083..0000000 --- a/src/pk.jl +++ /dev/null @@ -1,1505 +0,0 @@ -#Pharmacokinetics -#Makoid C, Vuchetich J, Banakar V. 1996-1999. Basic Pharmacokinetics. -#!!! - -isnanormissing(x::Number) = isnan(x) -isnanormissing(x::Missing) = true - - -struct KelData - s::Array{Int, 1} - e::Array{Int, 1} - a::Array{Float64, 1} - b::Array{Float64, 1} - r::Array{Float64, 1} - ar::Array{Float64, 1} - function KelData(s, e, a, b, r, ar)::KelData - new(s, e, a, b, r, ar)::KelData - end - function KelData()::KelData - new(Int[], Int[], Float64[], Float64[], Float64[], Float64[])::KelData - end -end -function Base.push!(keldata::KelData, s, e, a, b, r, ar) - push!(keldata.s, s) - push!(keldata.e, e) - push!(keldata.a, a) - push!(keldata.b, b) - push!(keldata.r, r) - push!(keldata.ar, ar) -end - -mutable struct ElimRange - kelstart::Int - kelend::Int - kelexcl::Vector{Int} - function ElimRange(kelstart::Int, kelend::Int, kelexcl::Vector{Int})::ElimRange - if kelstart < 0 throw(ArgumentError("Kel start point < 0")) end - if kelend < 0 throw(ArgumentError("Kel endpoint < 0")) end - if any(x -> x < 0, kelexcl) throw(ArgumentError("Exclude point < 0")) end - new(kelstart, kelend, kelexcl)::ElimRange - end - function ElimRange(kelstart::Int, kelend::Int) - ElimRange(kelstart, kelend, Vector{Int}(undef, 0)) - end - function ElimRange() - ElimRange(0, 0, Vector{Int}(undef, 0)) - end -end - -struct DoseTime - dose::Number - time::Number - tau::Number - function DoseTime(;dose::Number = NaN, time::Number = 0, tau::Number = NaN) - new(dose, time, tau)::DoseTime - end - function DoseTime(dose) - new(dose, 0, NaN)::DoseTime - end - function DoseTime(dose, time) - new(dose, time, NaN)::DoseTime - end - function DoseTime(dose, time, tau) - new(dose, time, tau)::DoseTime - end -end - - -#Plasma PK subject -mutable struct PKSubject <: AbstractSubject - time::Vector{Float64} - obs::Vector{Float64} - kelauto::Bool - kelrange::ElimRange - dosetime::DoseTime - keldata::KelData - sort::Dict - function PKSubject(time::Vector, conc::Vector, kelauto::Bool, kelrange::ElimRange, dosetime::DoseTime, keldata::KelData; sort = Dict()) - new(time, conc, kelauto, kelrange, dosetime, keldata, sort)::PKSubject - end - function PKSubject(time::Vector, conc::Vector, kelauto::Bool, kelrange, dosetime; sort = Dict()) - PKSubject(time, conc, kelauto, kelrange, dosetime, KelData(); sort = sort) - end - function PKSubject(time::Vector, conc::Vector; sort = Dict()) - PKSubject(time, conc, true, ElimRange(), DoseTime(NaN, 0 * time[1]), KelData(); sort = sort) - end - function PKSubject(data; time::Symbol, conc::Symbol, timesort::Bool = false, sort = Dict()) - if timesort sort!(data, time) end - #Check double time - PKSubject(copy(data[!,time]), copy(data[!,conc]), true, ElimRange(), DoseTime(), KelData(); sort = sort) - end -end - -#Urine PK subject -mutable struct UPKSubject <: AbstractSubject - stime::Vector - etime::Vector - conc::Vector - vol::Vector - time::Vector - obs::Vector - kelauto::Bool - kelrange::ElimRange - dosetime::DoseTime - keldata::KelData - sort::Dict - function UPKSubject(stime::Vector, etime::Vector, conc::Vector, vol::Vector; sort = Dict()) - new(stime, etime, conc, vol, (stime .+ etime) / 2, (conc .* vol) ./ (etime .- stime), true, ElimRange(), DoseTime(NaN, 0 * stime[1]), KelData(), sort)::UPKSubject - end -end - -#PD subject -mutable struct PDSubject <: AbstractSubject - time::Vector - obs::Vector - bl::Float64 - th::Float64 - sort::Dict - function PDSubject(time, resp, bl, th; sort = Dict()) - new(time, resp, bl, th, sort)::PDSubject - end - function PDSubject(time, resp; sort = Dict()) - new(time, resp, 0, NaN, sort)::PDSubject - end -end - -mutable struct PKPDProfile{T <: AbstractSubject} <: AbstractData - subject::T - method - result::Dict - sort::Dict - function PKPDProfile(subject::T, result; method = nothing) where T <: AbstractSubject - new{T}(subject, method, result, subject.sort) - end -end -""" - LimitRule(;lloq::Real = NaN, btmax::Real = NaN, atmax::Real = NaN, nan::Real = NaN, rm::Bool = false) - - - LimitRule(lloq::Real, btmax::Real, atmax::Real, nan::Real, rm::Bool) - -Rule for PK subject. - -STEP 1 (NaN step): replace all NaN values with nan reyword value (if nan !== NaN); -STEP 2 (LLOQ step): replace values below lloq with btmax value if this value befor Tmax or with atmax if this value after Tmax (if lloq !== NaN); -STEP 3 (remove NaN): rm == true, then remove all NaN values; -""" -struct LimitRule{T <: Real} - lloq::T - btmax::Float64 - atmax::Float64 - nan::Float64 - rm::Bool - function LimitRule(lloq::T, btmax, atmax, nan, rm::Bool) where T <: Real - new{T}(lloq, btmax, atmax, nan, rm)::LimitRule - end - function LimitRule(lloq, btmax, atmax, nan) - LimitRule(lloq, btmax, atmax, nan, false) - end - function LimitRule(lloq, btmax, atmax) - LimitRule(lloq, btmax, atmax, NaN, true) - end - function LimitRule(;lloq = NaN, btmax = NaN, atmax = NaN, nan = NaN, rm::Bool = false) - LimitRule(lloq, btmax, atmax, nan, rm) - end -end - - -function Base.getindex(a::PKPDProfile{T}, s::Symbol) where T <: AbstractSubject - return a.result[s] -end - -#= -function Base.getindex(a::DataSet{PKPDProfile}, i::Int) - return a.data[i] -end -function Base.getindex(a::DataSet{PKPDProfile}, i::Int, s::Symbol)::Real - return a.data[i].result[s] -end -function Base.getindex(a::DataSet{PKPDProfile}, d::Pair) - for i = 1:length(a) - if d ∈ a[i].subject.sort return a[i] end - end -end -function Base.getindex(a::DataSet{T}, i::Int) where T <: AbstractSubject - return a.data[i] -end -=# - -#------------------------------------------------------------------------------- -#= -function obsnum(data::T) where T <: AbstractSubject - return length(data.time) -end -function obsnum(keldata::KelData) - return length(keldata.a) -end -=# -function length(data::T) where T <: AbstractSubject - return length(data.time) -end -function length(keldata::KelData) - return length(keldata.a) -end -#------------------------------------------------------------------------------- -function Base.show(io::IO, obj::DataSet{PKSubject}) - println(io, "DataSet: Pharmacokinetic subject") - println(io, "Length: $(length(obj))") - for i = 1:length(obj) - println(io, "Subject $(i): ", obj[i].sort) - end -end -function Base.show(io::IO, obj::DataSet{PDSubject}) - println(io, "DataSet: Pharmacodynamic subject") - println(io, "Length: $(length(obj))") - for i = 1:length(obj) - println(io, "Subject $(i): ", obj[i].sort) - end -end -function Base.show(io::IO, obj::DataSet{PKPDProfile}) - println(io, "DataSet: Pharmacokinetic/Pharmacodynamic profile") - println(io, "Length: $(length(obj))") - for i = 1:length(obj) - println(io, "Subject $(i): ", obj[i].subject.sort) - end -end -function Base.show(io::IO, obj::PKSubject) - println(io, " Pharmacokinetic subject") - println(io, "Observations: $(length(obj)); ", obj.dosetime) - println(io, obj.kelrange) - println(io, "Time Concentration") - for i = 1:length(obj) - println(io, obj.time[i], " => ", obj.obs[i]) - end -end -function Base.show(io::IO, obj::PKPDProfile) - println(io, "PK/PD Profile") - println(io, "Type: $(obj.subject)") - println(io, "Method: $(obj.method)") - print(io, "Results:") - for k in keys(obj.result) - print(io, " $(k)") - end - println(io, ";") -end -function Base.show(io::IO, obj::PDSubject) - println(io, "Pharmacodynamic subject") - println(io, "Observations: $(length(obj))") - println(io, "Time Responce") - for i = 1:length(obj) - println(io, obj.time[i], " => ", obj.obs[i]) - end -end -function Base.show(io::IO, obj::ElimRange) - print(io, "Elimination range: $(obj.kelstart) - $(obj.kelend) ") - if length(obj.kelexcl) > 0 - print(io, "Exclusions: $(obj.kelexcl[1])") - if length(obj.kelexcl) > 1 for i = 2:length(obj.kelexcl) print(io, ", $(obj.kelexcl[i])") end end - print(io, ".") - else - print(io, "No exclusion.") - end -end -function Base.show(io::IO, obj::KelData) - m = copy(obj.s) - m = hcat(m, obj.e) - m = hcat(m, obj.a) - m = hcat(m, obj.b) - m = hcat(m, obj.r) - m = hcat(m, obj.ar) - println(io, "Elimination table:") - printmatrix(io, m) -end -#------------------------------------------------------------------------------- -function notnanormissing(x) - if x === NaN || x === missing return false else true end -end - """ - cmax - tmax - tmaxn - """ - function ctmax(time::Vector, obs::Vector) - if length(time) != length(obs) throw(ArgumentError("Unequal vector length!")) end - if length(obs[notnanormissing.(obs)]) > 0 cmax = maximum(obs[notnanormissing.(obs)]) else cmax = NaN end - tmax = NaN - tmaxn = 0 - for i = 1:length(time) - if obs[i] !== missing - if obs[i] == cmax - tmax = time[i] - tmaxn = i - break - end - end - end - return cmax, tmax, tmaxn - end - #= - function ctmax(data, conc::Symbol, time::Symbol) - return ctmax(data[!, time], data[!, conc]) - end - =# - function ctmax(data::PKSubject) - return ctmax(data.time, data.obs) - end -#= - function ctmax(data::PKSubject, dosetime) - s = 0 - for i = 1:length(data) - 1 - if dosetime >= data.time[i] && dosetime < data.time[i+1] s = i; break end - end - if dosetime == data.time[1] return ctmax(data.time, data.obs) end - mask = trues(length(data)) - cpredict = linpredict(data.time[s], data.time[s+1], data.dosetime.time, data.obs[s], data.obs[s+1]) - mask[1:s] .= false - cmax, tmax, tmaxn = ctmax(data.time[mask], data.obs[mask]) - if cmax > cpredict return cmax, tmax, tmaxn + s else return cpredict, data.dosetime.time, s end - end -=# -function dropbeforedosetime(time::Vector, obs::Vector, dosetime::DoseTime) - s = 0 - for i = 1:length(time) - @inbounds if time[i] < dosetime.time || isnan(obs[i]) s = i else break end - end - return time[s+1:end], obs[s+1:end] -end - """ - Range for AUC - return start point number, end point number - """ - function ncarange(time::Vector, dosetime, tau) - tautime = dosetime + tau - s = 1 - e = length(time) - for i = 1:length(time) - 1 - if dosetime > time[i] && dosetime <= time[i+1] s = i+1; - break - end - end - if tautime >= time[end] - e = length(time) - else - for i = s:length(time) - 1 - if tautime >= time[i] && tautime < time[i+1] e = i; break end - end - end - return s, e - end - - #--------------------------------------------------------------------------- - function slope(x::AbstractVector{T1}, y::AbstractVector{T2}) where T1 where T2 - if length(x) != length(y) throw(ArgumentError("Unequal vector length!")) end - n = length(x) - if n < 2 throw(ArgumentError("n < 2!")) end - T = promote_type(T1, T2) - Σxy = Σx = Σy = Σx2 = Σy2 = zero(T) - @inbounds for i = 1:n - Σxy += x[i] * y[i] - Σx += x[i] - Σy += y[i] - Σx2 += x[i]^2 - Σy2 += y[i]^2 - end - a = (n * Σxy - Σx * Σy)/(n * Σx2 - Σx^2) - b = (Σy * Σx2 - Σx * Σxy)/(n * Σx2 - Σx^2) - r2 = (n * Σxy - Σx * Σy)^2/((n * Σx2 - Σx^2)*(n * Σy2 - Σy^2)) - if n > 2 - ar = 1 - (1 - r2)*(n - 1)/(n - 2) - else - ar = NaN - end - return a, b, r2, ar - end #end slope - #Linear trapezoidal auc - function linauc(t₁, t₂, c₁, c₂) - return (t₂-t₁)*(c₁+c₂)/2 - end - #Linear trapezoidal aumc - function linaumc(t₁, t₂, c₁, c₂) - return (t₂-t₁)*(t₁*c₁+t₂*c₂)/2 - end - #Log trapezoidal auc - function logauc(t₁, t₂, c₁, c₂) - return (t₂-t₁)*(c₂-c₁)/log(c₂/c₁) - end - #Log trapezoidal aumc - function logaumc(t₁, t₂, c₁, c₂) - return (t₂-t₁) * (t₂*c₂-t₁*c₁) / log(c₂/c₁) - (t₂-t₁)^2 * (c₂-c₁) / log(c₂/c₁)^2 - end - #Intrapolation - #linear prediction bx from ax, a1 < ax < a2 - function linpredict(a₁, a₂, ax, b₁, b₂) - return (ax - a₁) / (a₂ - a₁)*(b₂ - b₁) + b₁ - end - - function logtpredict(c₁, c₂, cx, t₁, t₂) - return log(cx/c₁)/log(c₂/c₁)*(t₂-t₁)+t₁ - end - - function logcpredict(t₁, t₂, tx, c₁, c₂) - return exp(log(c₁) + (tx-t₁)/(t₂-t₁)*(log(c₂) - log(c₁))) - end - function cpredict(t₁, t₂, tx, c₁, c₂, calcm) - if calcm == :lint || c₂ >= c₁ - return linpredict(t₁, t₂, tx, c₁, c₂) - else - return logcpredict(t₁, t₂, tx, c₁, c₂) - end - end -#--------------------------------------------------------------------------- -#= - function aucpart(t₁, t₂, c₁, c₂, calcm, aftertmax) - if calcm == :lint - auc = linauc(t₁, t₂, c₁, c₂) # (data[i - 1, conc] + data[i, conc]) * (data[i, time] - data[i - 1, time])/2 - aumc = linaumc(t₁, t₂, c₁, c₂) # (data[i - 1, conc]*data[i - 1, time] + data[i, conc]*data[i, time]) * (data[i, time] - data[i - 1, time])/2 - elseif calcm == :logt - if aftertmax - if c₁ > 0 && c₂ > 0 - auc = logauc(t₁, t₂, c₁, c₂) - aumc = logaumc(t₁, t₂, c₁, c₂) - else - auc = linauc(t₁, t₂, c₁, c₂) - aumc = linaumc(t₁, t₂, c₁, c₂) - end - else - auc = linauc(t₁, t₂, c₁, c₂) - aumc = linaumc(t₁, t₂, c₁, c₂) - end - elseif calcm == :luldt - if aftertmax - if c₁ > c₂ > 0 - auc = logauc(t₁, t₂, c₁, c₂) - aumc = logaumc(t₁, t₂, c₁, c₂) - else - auc = linauc(t₁, t₂, c₁, c₂) - aumc = linaumc(t₁, t₂, c₁, c₂) - end - else - auc = linauc(t₁, t₂, c₁, c₂) - aumc = linaumc(t₁, t₂, c₁, c₂) - end - elseif calcm == :luld - if c₁ > c₂ > 0 - auc = logauc(t₁, t₂, c₁, c₂) - aumc = logaumc(t₁, t₂, c₁, c₂) - else - auc = linauc(t₁, t₂, c₁, c₂) - aumc = linaumc(t₁, t₂, c₁, c₂) - end - end - return auc, aumc - end -=# -function aucpart(t₁, t₂, c₁::T1, c₂::T2, calcm, aftertmax, dosetime) where T1 where T2 - - if calcm == :lint - auc = linauc(t₁, t₂, c₁, c₂) - aumc = linaumc(t₁ - dosetime, t₂ - dosetime, c₁, c₂) - elseif calcm == :logt && aftertmax && c₁ > zero(T1) && c₂ > zero(T2) - auc = logauc(t₁, t₂, c₁, c₂) - aumc = logaumc(t₁ - dosetime, t₂ - dosetime, c₁, c₂) - elseif calcm == :luldt && aftertmax && c₁ > c₂ > zero(T2) - auc = logauc(t₁, t₂, c₁, c₂) - aumc = logaumc(t₁ - dosetime, t₂ - dosetime, c₁, c₂) - elseif calcm == :luld && c₁ > c₂ > zero(T2) - auc = logauc(t₁, t₂, c₁, c₂) - aumc = logaumc(t₁ - dosetime, t₂ - dosetime, c₁, c₂) - else - auc = linauc(t₁, t₂, c₁, c₂) - aumc = linaumc(t₁ - dosetime, t₂ - dosetime, c₁, c₂) - end - return auc, aumc -end -#------------------------------------------------------------------------------- -""" - nca!(data::PKSubject; adm = :ev, calcm = :lint, intp = :lint, - verbose = false, warn = true, io::IO = stdout) - -Pharmacokinetics non-compartment analysis for one PK subject. - -`adm` - administration: - -* `:ev` - extravascular; -* `:iv` - intra vascular bolus; - -`calcm` - calculation method; - -- :lint - Linear trapezoidal everywhere; -- :logt - Log-trapezoidat rule after Tmax if c₁ > 0 and c₂ > 0, else Linear trapezoidal used; -- :luld - Linear Up - Log Down everywhere if c₁ > c₂ > 0, else Linear trapezoidal used; -- :luldt - Linear Up - Log Down after Tmax if c₁ > c₂ > 0, else Linear trapezoidal used; - -`intp` - interpolation rule; -- :lint - linear interpolation; -- :logt - log interpolation; - - -`verbose` - print to out stream if "true"; - -- true; -- false. - -`warn` - warnings enabled; - -`io` - out stream (stdout by default). -""" -function nca!(data::PKSubject; adm = :ev, calcm = :lint, intp = :lint, verbose = false, warn = true, io::IO = stdout) - result = Dict() - #= - result = Dict(:Obsnum => 0, :Tmax => 0, :Tmaxn => 0, :Tlast => 0, - :Cmax => 0, :Cdose => NaN, :Clast => 0, :Ctau => NaN, :Ctaumin => NaN, :Cavg => NaN, - :Swing => NaN, :Swingtau => NaN, :Fluc => NaN, :Fluctau => NaN, - :AUClast => 0, :AUCall => 0, :AUCtau => 0, :AUMClast => 0, :AUMCall => 0, :AUMCtau => 0, - :MRTlast => 0, :Kel => NaN, :HL => NaN, - :Rsq => NaN, :ARsq => NaN, :Rsqn => NaN, - :AUCinf => NaN, :AUCpct => NaN) - =# - if warn - if !allunique(data.time) - @warn "Not all time values is unique!" - end - end - #if calcm == :logt / :luld / :luldt calculation method used - log interpolation using - if calcm == :logt || calcm == :luld || calcm == :luldt intp = :logt end - - #time = data.time - #obs = data.obs - time, obs = dropbeforedosetime(data.time, data.obs, data.dosetime) - #dropnan - #nanlist = findall(isnan, obs) - #if length(nanlist) > 0 - # deleteat!(time, nanlist) - # deleteat!(obs, nanlist) - #end - result[:Obsnum] = length(obs) - rmn = length(data.obs) - length(obs) - if result[:Obsnum] < 2 - return PKPDProfile(data, result; method = calcm) - end - - result[:Cmax] = maximum(obs) - doseaucpart = 0.0 - doseaumcpart = 0.0 - daucpartl = 0.0 - daumcpartl = 0.0 - - for i = result[:Obsnum]:-1:1 - if obs[i] > 0 - result[:Tlast] = time[i] - result[:Clast] = obs[i] - break - end - end - result[:Cmax], result[:Tmax], result[:Tmaxn] = ctmax(time, obs) - #result[:Tmaxn] += rmn - #----------------------------------------------------------------------- - if verbose - println(io, "Non-compartmental Pharmacokinetic Analysis") - printsortval(io, data.sort) - println(io, "Settings:") - println(io, "Method: $(calcm)") - end - #----------------------------------------------------------------------- - #Areas - aucpartl = Array{Number, 1}(undef, result[:Obsnum] - 1) - aumcpartl = Array{Number, 1}(undef, result[:Obsnum] - 1) - #Calculate all UAC part based on data - for i = 1:(result[:Obsnum] - 1) - aucpartl[i], aumcpartl[i] = aucpart(time[i], time[i + 1], obs[i], obs[i + 1], calcm, i + 1 > result[:Tmaxn], data.dosetime.time) - end - #pmask = Array{Bool, 1}(undef, result[:Obsnum] - 1) - #pmask .= true - pmask = trues(result[:Obsnum] - 1) - - ncas = nothing - ncae = nothing - - # If TAU set, calculates start and edt timepoints for AUCtau - if data.dosetime.tau > 0 - ncas, ncae = ncarange(time, data.dosetime.time, data.dosetime.tau) - result[:Ctaumin] = minimum(view(obs, ncas:ncae)) - end - - #Elimination - keldata = KelData() - if data.kelauto - if result[:Obsnum] - result[:Tmaxn] >= 3 - logconc = log.(obs) - cmask = trues(result[:Obsnum]) - kelexclrmn = data.kelrange.kelexcl .- rmn - filter!(x-> 0 < x <= result[:Obsnum], kelexclrmn) - view(cmask, kelexclrmn) .= false - kellast = findlast(x -> x, cmask) - if kellast - result[:Tmaxn] >= 3 - if adm == :iv tmnc = 0 else tmnc = 1 end - for i = result[:Tmaxn] + tmnc:kellast - 2 - cmask .= false - #cmask[i:result[:Obsnum]] .= true - view(cmask, i:kellast) .= true - if length(data.kelrange.kelexcl) > 0 - view(cmask, kelexclrmn) .= false - end - if sum(cmask) >= 3 - sl = slope(time[cmask], logconc[cmask]) - if sl[1] < 0 - push!(keldata, i + rmn, kellast + rmn, sl[1], sl[2], sl[3], sl[4]) - end - end - end - end - end - else - logconc = log.(obs) - #cmask = Array{Bool, 1}(undef, result[:Obsnum]) - #cmask .= false - cmask = falses(result[:Obsnum]) - cmask[data.kelrange.kelstart-rmn:data.kelrange.kelend-rmn] .= true - cmask[data.kelrange.kelexcl .- rmn] .= false - sl = slope(view(time, cmask), view(logconc, cmask)) - push!(keldata, data.kelrange.kelstart, data.kelrange.kelend, sl[1], sl[2], sl[3], sl[4]) - end - - #Calcalation partial areas (doseaucpart, doseaumcpart) - #Dosetime < first time - if data.dosetime.time == time[1] - result[:Cdose] = obs[1] - elseif data.dosetime.time < time[1] - if adm == :iv - if obs[1] > obs[2] > 0 - result[:Cdose] = logcpredict(time[1], time[2], data.dosetime.time, obs[1], obs[2]) - else - result[:Cdose] = obs[findfirst(x->x>0, obs)] - end - else - if data.dosetime.tau > 0 - result[:Cdose] = result[:Ctaumin] - else - result[:Cdose] = 0 - end - end - doseaucpart, doseaumcpart = aucpart(data.dosetime.time, time[1], result[:Cdose], obs[1], :lint, false, data.dosetime.time) #always :lint? - else - error("Some concentration before dose time after filtering!!!") - end - - #sum all full AUC parts and part before dose - result[:AUCall] = sum(view(aucpartl, pmask)) + doseaucpart - result[:AUMCall] = sum(view(aumcpartl, pmask)) + doseaumcpart - #----------------------------------------------------------------------- - #----------------------------------------------------------------------- - #Find last mesurable concentration (>0) from start, all after excluded - for i = result[:Tmaxn]:result[:Obsnum] - 1 - if obs[i+1] <= 0 pmask[i:end] .= false; break end - end - #= - #Find last mesurable concentration (>0) from end, all after excluded - for i = result[:Obsnum]:-1:1 - if obs[i] <= 0 pmask[i] = false else break end - end - =# - result[:AUClast] = sum(view(aucpartl, pmask)) + doseaucpart - result[:AUMClast] = sum(view(aumcpartl, pmask)) + doseaumcpart - result[:MRTlast] = result[:AUMClast] / result[:AUClast] - - - result[:Cllast] = data.dosetime.dose / result[:AUClast] - #----------------------------------------------------------------------- - #----------------------------------------------------------------------- - - if length(keldata) > 0 - result[:ARsq], result[:Rsqn] = findmax(keldata.ar) - data.kelrange.kelstart = keldata.s[result[:Rsqn]] - result[:Kelstart] = data.kelrange.kelstart - data.kelrange.kelend = keldata.e[result[:Rsqn]] - result[:Kelend] = data.kelrange.kelend - data.keldata = keldata - result[:Rsq] = keldata.r[result[:Rsqn]] - result[:Kel] = abs(keldata.a[result[:Rsqn]]) - result[:LZ] = keldata.a[result[:Rsqn]] - result[:LZint] = keldata.b[result[:Rsqn]] - result[:Clast_pred] = exp(result[:LZint] + result[:LZ]*result[:Tlast]) - result[:HL] = LOG2 / result[:Kel] - result[:AUCinf] = result[:AUClast] + result[:Clast] / result[:Kel] - result[:AUCinf_pred] = result[:AUClast] + result[:Clast_pred] / result[:Kel] - result[:AUCpct] = (result[:AUCinf] - result[:AUClast]) / result[:AUCinf] * 100.0 - result[:AUMCinf] = result[:AUMClast] + result[:Tlast] * result[:Clast] / result[:Kel] + result[:Clast] / result[:Kel] ^ 2 - result[:MRTinf] = result[:AUMCinf] / result[:AUCinf] - result[:Vzlast] = data.dosetime.dose / result[:AUClast] / result[:Kel] - result[:Vzinf] = data.dosetime.dose / result[:AUCinf] / result[:Kel] - result[:Clinf] = data.dosetime.dose / result[:AUCinf] - result[:Vssinf] = result[:Clinf] * result[:MRTinf] - - else - result[:Kel] = NaN - result[:HL] = NaN - result[:AUCinf] = NaN - result[:LZint] = NaN - end - #----------------------------------------------------------------------- - #Steady-state - tautime = data.dosetime.time + data.dosetime.tau - if data.dosetime.tau > 0 - eaucpartl = eaumcpartl = 0.0 - if tautime < time[end] - if tautime > result[:Tmax] aftertmax = true else aftertmax = false end - result[:Ctau] = cpredict(time[ncae], time[ncae + 1], tautime, obs[ncae], obs[ncae + 1], intp) - eaucpartl, eaumcpartl = aucpart(time[ncae], tautime, obs[ncae], result[:Ctau], calcm, aftertmax, data.dosetime.time) - #remoove part after tau - if ncae < result[:Obsnum] - 1 pmask[ncae:end] .= false end - elseif tautime > time[end] && result[:LZint] !== NaN - #extrapolation - result[:Ctau] = exp(result[:LZint] + result[:LZ] * tautime) - eaucpartl, eaumcpartl = aucpart(time[ncae], tautime, obs[ncae], result[:Ctau], calcm, true, data.dosetime.time) - else - result[:Ctau] = obs[end] - end - result[:AUCtau] = sum(view(aucpartl, pmask)) + eaucpartl + doseaucpart - result[:AUMCtau] = sum(view(aumcpartl, pmask)) + eaumcpartl + doseaumcpart - result[:Cavg] = result[:AUCtau]/data.dosetime.tau - if result[:Ctaumin] != 0 - result[:Swing] = (result[:Cmax] - result[:Ctaumin])/result[:Ctaumin] - else - result[:Swing] = NaN - end - if result[:Ctau] != 0 - result[:Swingtau] = (result[:Cmax] - result[:Ctau])/result[:Ctau] - else - result[:Swingtau] = NaN - end - result[:Fluc] = (result[:Cmax] - result[:Ctaumin])/result[:Cavg] - result[:Fluctau] = (result[:Cmax] - result[:Ctau])/result[:Cavg] - #If Kel calculated - if result[:Kel] !== NaN - result[:Accind] = 1 / (1 - (exp(-result[:Kel] * data.dosetime.tau))) - else - result[:Accind] = NaN - end - result[:MRTtauinf] = (result[:AUMCtau] + data.dosetime.tau * (result[:AUCinf] - result[:AUCtau])) / result[:AUCtau] - result[:Vztau] = data.dosetime.dose / result[:AUCtau] / result[:Kel] - result[:Cltau] = data.dosetime.dose / result[:AUCtau] - result[:Vsstau] = result[:Cltau] * result[:MRTtauinf] - end - if verbose - aucpartlsum = similar(aucpartl) - aumcpartlsum = similar(aumcpartl) - for i = 1:length(aucpartl) - aucpartlsum[i] = sum(view(aucpartl, 1:i)) - aumcpartlsum[i] = sum(view(aumcpartl, 1:i)) - end - astx = Vector{String}(undef, length(time)) - astx[1] = "" - for i = 1:length(pmask) - if pmask[i] astx[i+1] = "Yes" else astx[i+1] = "No" end - end - - aucpartlsum .+= doseaucpart - - aumcpartlsum .+= doseaumcpart - if doseaucpart > 0 || doseaumcpart > 0 astx[1] = "D" end - mx = hcat(time, obs, round.(vcat([doseaucpart], aucpartl), digits = 3), round.(vcat([doseaucpart], aucpartlsum), digits = 3), round.(vcat([doseaumcpart], aumcpartl), digits = 3), round.(vcat([doseaumcpart], aumcpartlsum), digits = 3), astx) - mx = vcat(permutedims(["Time", "Concentrtion", "AUC", "AUC (cumulate)", "AUMC", "AUMC (cumulate)", "Include"]), mx) - printmatrix(io, mx) - println(io, "") - println(io, "Cdose: $(result[:Cdose]), Dose time: $(data.dosetime.time)") - if data.dosetime.time > time[1] - println("Dose AUC part: $(doseaucpart)") - println("Dose AUMC part: $(doseaumcpart)") - end - println(io, "") - if tautime < time[end] && tautime > 0 - println(io, "Tau + dosetime is less then end time. Interpolation used.") - println(io, "Interpolation between: $(time[ncae]) - $( time[ncae + 1]), method: $(intp)") - println(io, "Ctau: $(result[:Ctau])") - println(io, "AUC final part: $(eaucpartl)") - println(io, "AUMC final part: $(eaumcpartl)") - end - end - result[:Tmaxn] += rmn #! - #----------------------------------------------------------------------- - return PKPDProfile(data, result; method = calcm) - end -""" - nca!(data::DataSet{PKSubject}; adm = :ev, calcm = :lint, intp = :lint, - verbose = false, warn = true, io::IO = stdout, sort = nothing) - -* `sort` - only for this type of subjects. - -Pharmacokinetics non-compartment analysis for PK subjects DataSet. -""" -function nca!(data::DataSet{PKSubject}; adm = :ev, calcm = :lint, intp = :lint, verbose = false, warn = true, io::IO = stdout, sort = nothing) - results = Array{PKPDProfile, 1}(undef, 0) - if isnothing(sort) - for i = 1:length(data.data) - push!(results, nca!(data.data[i]; adm = adm, calcm = calcm, intp = intp, verbose = verbose, warn = warn, io = io)) - end - else - for i = 1:length(data.data) - if sort ⊆ data.data[i].sort - push!(results, nca!(data.data[i]; adm = adm, calcm = calcm, intp = intp, verbose = verbose, warn = warn, io = io)) - end - end - end - return DataSet(results) -end - -#--------------------------------------------------------------------------- -""" - nca!(data::PDSubject; verbose = false, io::IO = stdout)::PKPDProfile{PDSubject} - -Pharmacodynamics non-compartment analysis for one PD subject. -""" -function nca!(data::PDSubject; verbose = false, io::IO = stdout)::PKPDProfile{PDSubject} - result = Dict(:TH => NaN, :AUCABL => 0, :AUCBBL => 0, :AUCBLNET => 0, :TABL => 0, :TBBL => 0) - result[:Obsnum] = length(data) - result[:TH] = data.th - result[:BL] = data.bl - result[:RMAX] = maximum(data.obs) - result[:RMIN] = minimum(data.obs) - for i = 2:length(data) #BASELINE - if data.obs[i - 1] <= result[:BL] && data.obs[i] <= result[:BL] - result[:TBBL] += data.time[i,] - data.time[i - 1] - result[:AUCBBL] += linauc(data.time[i - 1], data.time[i], result[:BL] - data.obs[i - 1], result[:BL] - data.obs[i]) - elseif data.obs[i - 1] <= result[:BL] && data.obs[i] > result[:BL] - tx = linpredict(data.obs[i - 1], data.obs[i], result[:BL], data.time[i - 1], data.time[i]) - result[:TBBL] += tx - data.time[i - 1] - result[:TABL] += data.time[i] - tx - result[:AUCBBL] += (tx - data.time[i - 1]) * (result[:BL] - data.obs[i - 1]) / 2 - result[:AUCABL] += (data.time[i] - tx) * (data.obs[i] - result[:BL]) / 2 #Ok - elseif data.obs[i - 1] > result[:BL] && data.obs[i] <= result[:BL] - tx = linpredict(data.obs[i - 1], data.obs[i], result[:BL], data.time[i - 1], data.time[i]) - result[:TBBL] += data.time[i] - tx - result[:TABL] += tx - data.time[i - 1] - result[:AUCBBL] += (data.time[i] - tx) * (result[:BL] - data.obs[i]) / 2 - result[:AUCABL] += (tx - data.time[i - 1]) * (data.obs[i - 1] - result[:BL]) / 2 - elseif data.obs[i - 1] > result[:BL] && data.obs[i] > result[:BL] - result[:TABL] += data.time[i] - data.time[i - 1] - result[:AUCABL] += linauc(data.time[i - 1], data.time[i], data.obs[i - 1] - result[:BL], data.obs[i] - result[:BL]) - end - end #BASELINE - - #THRESHOLD - if result[:TH] !== NaN - result[:AUCATH] = 0 - result[:AUCBTH] = 0 - result[:TATH] = 0 - result[:TBTH] = 0 - result[:AUCTHNET] = 0 - result[:AUCDBLTH] = 0 - for i = 2:length(data) - if data.obs[i - 1] <= result[:TH] && data.obs[i] <= result[:TH] - result[:TBTH] += data.time[i] - data.time[i - 1] - result[:AUCBTH] += linauc(data.time[i - 1], data.time[i], result[:TH] - data.obs[i - 1], result[:TH] - data.obs[i]) - elseif data.obs[i - 1] <= result[:TH] && data.obs[i] > result[:TH] - tx = linpredict(data.obs[i - 1], data.obs[i], result[:TH], data.time[i - 1], data.time[i]) - result[:TBTH] += tx - data.time[i - 1] - result[:TATH] += data.time[i] - tx - result[:AUCBTH] += (tx - data.time[i - 1]) * (result[:TH] - data.obs[i - 1]) / 2 - result[:AUCATH] += (data.time[i] - tx) * (data.obs[i] - result[:TH]) / 2 #Ok - elseif data.obs[i - 1] > result[:TH] && data.obs[i] <= result[:TH] - tx = linpredict(data.obs[i - 1], data.obs[i], result[:TH], data.time[i - 1], data.time[i]) - result[:TBTH] += data.time[i] - tx - result[:TATH] += tx - data.time[i - 1] - result[:AUCBTH] += (data.time[i] - tx) * (result[:TH] - data.obs[i]) / 2 - result[:AUCATH] += (tx - data.time[i - 1]) * (data.obs[i - 1] - result[:TH]) / 2 - elseif data.obs[i - 1] > result[:TH] && data.obs[i] > result[:TH] - result[:TATH] += data.time[i] - data.time[i - 1] - result[:AUCATH] += linauc(data.time[i - 1], data.time[i], data.obs[i - 1] - result[:TH], data.obs[i] - result[:TH]) - end - end - if result[:BL] > result[:TH] - result[:AUCDBLTH] = result[:AUCATH] - result[:AUCABL] - else - result[:AUCDBLTH] = result[:AUCABL] -result[:AUCATH] - end - result[:AUCTHNET] = result[:AUCATH] - result[:AUCBTH] - end - result[:AUCBLNET] = result[:AUCABL] - result[:AUCBBL] - return PKPDProfile(data, result; method = :lint) - end -""" - nca!(data::PDSubject; verbose = false, io::IO = stdout)::PKPDProfile{PDSubject} - -Pharmacodynamics non-compartment analysis for PD subjects set. -""" -function nca!(data::DataSet{PDSubject}; verbose = false, io::IO = stdout) - results = Array{PKPDProfile, 1}(undef, 0) - for i = 1:length(data.data) - push!(results, nca!(data.data[i]; verbose = verbose, io = io)) - end - return DataSet(results) -end - - -""" - nca!(data::UPKSubject; verbose = false, io::IO = stdout)::PKPDProfile{UPKSubject} - -Pharmacodynamics non-compartment analysis for urine data. -""" -function nca!(data::UPKSubject; verbose = false, io::IO = stdout)::PKPDProfile{UPKSubject} - - result = Dict() - #time = (data.stime .+ data.etime) / 2 - #exrate = (data.conc .* data.vol) ./ (data.etime .- data.stime) - result[:maxrate], result[:mTmax], result[:mTmaxn] = ctmax(data.time, data.obs) - result[:ar] = sum(data.conc .* data.vol) - result[:volume] = sum(data.vol) - - #result[:AUCrate] - #result[:tmaxrate] - #result[:arp] - #result[:Kel] - #result[:HL] - return PKPDProfile(data, result; method = :upk) -end -""" - nca!(data::DataSet{UPKSubject}; verbose = false, io::IO = stdout) - -Pharmacodynamics non-compartment analysis for PD subjects set. -""" -function nca!(data::DataSet{UPKSubject}; verbose = false, io::IO = stdout) - results = Vector{PKPDProfile}(undef, 0) - for i = 1:length(data) - push!(results, nca!(data[i]; verbose = verbose, io = io)) - end - return DataSet(results) -end -#= -stime::Vector -etime::Vector -obs::Vector -vol::Vector -=# -#------------------------------------------------------------------------------- - - -function getdatai(data, sort, cols, func; sortby = nothing) - sortlist = unique(data[!, sort]) - for i = 1:size(sortlist, 1) - inds = Vector{Int}(undef, 0) - for c = 1:size(data, 1) #For each line in data - if data[c, sort] == sortlist[i,:] - push!(inds, c) - end - end - datai = data[inds, cols] - if sortby !== nothing - sort!(datai, sortby) - end - func(datai, sortlist[i,:]) - end -end - -#------------------------------------------------------------------------------- - -function tryfloatparse!(x) - for i = 1:length(x) - if typeof(x[i]) <: AbstractString - try - x[i] = parse(Float64, x[i]) - catch - x[i] = NaN - end - else - try - if x[i] === missing - x[i] = NaN - else - x[i] = float(x[i]) - end - catch - x[i] = NaN - end - end - end - return convert(Vector{Float64}, x) -end - -""" - pkimport(data, sort::Array, rule::LimitRule; conc::Symbol, time::Symbol) - -Pharmacokinetics data import. - -- data - sourece data; -- sort - sorting columns; -- rule - applied LimitRule. - -- conc - concentration column; -- time - time column. -""" -function pkimport(data, sort::Array, rule::LimitRule; conc::Symbol, time::Symbol)::DataSet - results = Array{PKSubject, 1}(undef, 0) - if !(eltype(data[!, conc]) <: Real) throw(ArgumentError("Type of concentration data not <: Real!"))end - getdatai(data, sort, [conc, time], (x, y) -> push!(results, PKSubject(x[!, time], x[!, conc], sort = Dict(sort .=> collect(y)))); sortby = time) - results = DataSet(results) - applyncarule!(results, rule) - return results -end -""" - pkimport(data, sort::Array; time::Symbol, conc::Symbol) - -Pharmacokinetics data import. - -- data - sourece data; -- sort - sorting columns. - -- conc - concentration column; -- time - time column. -""" -function pkimport(data, sort::Array; time::Symbol, conc::Symbol)::DataSet - rule = LimitRule() - return pkimport(data, sort, rule; conc = conc, time = time) -end -""" - pkimport(data, sort::Array; time::Symbol, conc::Symbol) - -Pharmacokinetics data import for one subject. - -- data - sourece data; -- rule - applied LimitRule. - -- conc - concentration column; -- time - time column. -""" -function pkimport(data, rule::LimitRule; time::Symbol, conc::Symbol)::DataSet - #rule = LimitRule() - datai = ncarule!(copy(data[!,[time, conc]]), conc, time, rule) - return DataSet([PKSubject(datai[!, time], datai[!, conc])]) -end -""" - pkimport(data, sort::Array; time::Symbol, conc::Symbol) - -Pharmacokinetics data import for one subject. - -- data - sourece data; - -- conc - concentration column; -- time - time column. -""" -function pkimport(data; time::Symbol, conc::Symbol)::DataSet - datai = sort(data[!,[time, conc]], time) - return DataSet([PKSubject(datai[!, time], datai[!, conc])]) -end - -""" - pkimport(time::Vector, conc::Vector; sort = Dict())::PKSubject - -Pharmacokinetics data import. - -- time - time vector; -- conc - concentration vector. - -""" -function pkimport(time::Vector, conc::Vector; sort = Dict())::PKSubject - PKSubject(time, conc; sort = sort) -end - #--------------------------------------------------------------------------- -""" - pdimport(data, sort::Array; resp::Symbol, time::Symbol, - bl::Real = 0, th::Real = NaN) - -Pharmacodynamics data import. - -- data - sourece data; -- sort - sorting columns; - -- resp - responce column; -- time - time column -- bl - baseline; -- th - treashold. -""" -function pdimport(data, sort::Array; resp::Symbol, time::Symbol, bl::Real = 0, th::Real = NaN)::DataSet - sortlist = unique(data[:, sort]) - results = Array{PDSubject, 1}(undef, 0) - getdatai(data, sort, [resp, time], (x, y) -> push!(results, PDSubject(x[!, time], x[!, resp], bl, th, sort = Dict(sort .=> collect(y)))); sortby = time) - return DataSet(results) -end - -function pdimport(data; resp::Symbol, time::Symbol, bl = 0, th = NaN) - datai = sort(data[!,[time, resp]], time) - return DataSet([PDSubject(datai[!, time], datai[!, resp], bl, th)]) -end -#------------------------------------------------------------------------------- -function upkimport(data, sort::Array; stime::Symbol, etime::Symbol, conc::Symbol, vol::Symbol)::DataSet - results = Array{UPKSubject, 1}(undef, 0) - getdatai(data, sort, [stime, etime, conc, vol], (x, y) -> push!(results, UPKSubject(x[!, stime], x[!, etime], x[!, conc], x[!, vol]; sort =Dict(sort .=> collect(y)))); sortby = stime) - results = DataSet(results) - return results -end - -#------------------------------------------------------------------------------- -""" - applyncarule!(data::PKSubject, rule::LimitRule) - -Apply rule to PK subject . - -STEP 1 (NaN step): replace all NaN values with nan reyword value (if nan !== NaN); -STEP 2 (LLOQ step): replace values below lloq with btmax value if this value befor Tmax or with atmax if this value after Tmax (if lloq !== NaN); -STEP 3 (remove NaN): rm == true, then remove all NaN values; -""" -function applyncarule!(data::PKSubject, rule::LimitRule) - cmax, tmax, tmaxn = ctmax(data) - #NaN Rule - obsn = length(data) - if rule.nan !== NaN - for i = 1:length(data) - if data.obs[i] === NaN || data.obs[i] === missing - data.obs[i] = rule.nan - end - end - end - #LLOQ rule - if rule.lloq !== NaN - for i = 1:obsn - if data.obs[i] <= rule.lloq - if i <= tmaxn - data.obs[i] = rule.btmax - else - data.obs[i] = rule.atmax - end - end - end - end - #NaN Remove rule - if rule.rm - inds = findall(isnanormissing, data.obs) - deleteat!(data.time, inds) - deleteat!(data.obs, inds) - #filterv = data.obs .!== NaN - #data.time = data.time[filterv] - #data.obs = data.obs[filterv] - end -end - -function applyncarule!(data::DataSet{PKSubject}, rule::LimitRule; sort = nothing) - if isnothing(sort) - for i in data - applyncarule!(i, rule) - end - else - for i in data - if sort ⊆ i.sort - applyncarule!(i, rule) - end - end - end - data -end - -#--------------------------------------------------------------------------- -""" - setelimrange!(data::PKSubject, range::ElimRange; kelauto = false) - -Set range for elimination parameters calculation for subject. - -data - PK subject; -range - ElimRange object; -kelauto - set kelauto value. - -Set kelauto `false` if kelend > kelstart > 0. -""" -function setelimrange!(data::PKSubject, range::ElimRange; kelauto = false) - if range.kelend > length(data) throw(ArgumentError("Kel endpoint out of range")) end - if range.kelend > range.kelstart > 0 setkelauto!(data, kelauto) end - data.kelrange = range -end -""" - setelimrange!(data::DataSet{PKSubject}, range::ElimRange) - -Set range for elimination parameters calculation for DataSet. - -data - DataSet of PK subject; -range - ElimRange object. -""" -function setelimrange!(data::DataSet{PKSubject}, range::ElimRange) - for i = 1:length(data) - setelimrange!(data[i], range) - end - data -end -""" - setelimrange!(data::DataSet{PKSubject}, range::ElimRange, subj::Array{Int,1}) - -Set range for elimination parameters calculation for DataSet. - -data - DataSet of PK subject; -range - ElimRange object; -subj - subjects. -""" -function setelimrange!(data::DataSet{PKSubject}, range::ElimRange, subj::Array{Int,1}) - for i in subj - setelimrange!(data[i], range) - end - data -end -""" - setelimrange!(data::DataSet{PKSubject}, range::ElimRange, subj::Int) - -Set range for elimination parameters calculation for DataSet. - -data - DataSet of PK subject; -range - ElimRange object; -subj - subject. -""" -function setelimrange!(data::DataSet{PKSubject}, range::ElimRange, subj::Int) - setelimrange!(data[subj], range) - data -end -""" - setelimrange!(data::DataSet{PKSubject}, range::ElimRange, sort::Dict; kelauto = false) - -Set range for elimination parameters calculation for DataSet. - -data - DataSet of PK subject; -range - ElimRange object; -sort - subject sort. -""" -function setelimrange!(data::DataSet{PKSubject}, range::ElimRange, sort::Dict; kelauto = false) - for i = 1:length(data) - if sort ⊆ data[i].sort setelimrange!(data[i], range; kelauto = kelauto) end - end - data -end -#------------------------------------------------------------------------------- -""" - setkelauto!(data::PKSubject, kelauto::Bool) - -Set range for elimination parameters calculation for subject. - -data - PK subject; -kelauto - value. -""" -function setkelauto!(data::PKSubject, kelauto::Bool) - data.kelauto = kelauto -end - -function setkelauto!(data::DataSet{PKSubject}, kelauto::Bool) - for i = 1:length(data) - setkelauto!(data[i], kelauto) - end - data -end -function setkelauto!(data::DataSet{PKSubject}, kelauto::Bool, subj::Array{Int,1}) - for i in subj - setkelauto!(data[i], kelauto) - end - data -end -function setkelauto!(data::DataSet{PKSubject}, kelauto::Bool, subj::Int) - setkelauto!(data[subj], kelauto) - data -end - -#------------------------------------------------------------------------------- -function getkelauto(data::PKSubject) - return data.kelauto -end - -#------------------------------------------------------------------------------- -""" - applyelimrange!(data::PKPDProfile{PKSubject}, range::ElimRange) - -Set range for elimination parameters calculation for profile's subject and recalculate NCA parameters. - -data - PK subject; -range - ElimRange object. -""" -function applyelimrange!(data::PKPDProfile{PKSubject}, range::ElimRange) - setelimrange!(data.subject, range) - data.result = nca!(data.subject, calcm = data.method).result - data -end -function applyelimrange!(data::DataSet{PKPDProfile}, range::ElimRange) - for i = 1:length(data) - applyelimrange!(data[i], range) - end - data -end -function applyelimrange!(data::DataSet{PKPDProfile}, range::ElimRange, subj::Array{Int,1}) - for i = 1:length(data) - if i ∈ subj applyelimrange!(data[i], range) end - end - data -end -function applyelimrange!(data::DataSet{PKPDProfile}, range::ElimRange, subj::Int) - applyelimrange!(data[subj], range) - data -end -function applyelimrange!(data::DataSet{PKPDProfile}, range::ElimRange, sort::Dict) - for i = 1:length(data) - if sort ⊆ data[i].subject.sort - applyelimrange!(data[i], range) - end - end - data -end -#--------------------------------------------------------------------------- -""" - setdosetime!(data::PKSubject, dosetime::DoseTime) - -Set dosing time and tau for subject. - -data - PK subject; -dosetime - DoseTime object. -""" -function setdosetime!(data::PKSubject, dosetime::DoseTime) - data.dosetime = dosetime - data -end -""" - setdosetime!(data::DataSet{PKSubject}, dosetime::DoseTime) - -Set dosing time and tau for all subjects. - -data - dataset of PK subjects; -dosetime - DoseTime object. -""" -function setdosetime!(data::DataSet{PKSubject}, dosetime::DoseTime) - for i = 1:length(data) - setdosetime!(data[i], dosetime) - end - data -end -""" - setdosetime!(data::DataSet{PKSubject}, dosetime::DoseTime, sort::Dict) - -Set dosing time and tau for subjects in sort dict. - -data - dataset of PK subjects; -dosetime - DoseTime object. -""" -function setdosetime!(data::DataSet{PKSubject}, dosetime::DoseTime, sort::Dict) - for i = 1:length(data) - if sort ⊆ data[i].sort setdosetime!(data[i], dosetime) end - end - data -end -""" - findfirst(sort::Dict, data::DataSet{PKSubject}) - -The first item in the dataset that satisfies the condition - sort dictionary. -""" -function findfirst(sort::Dict, data::DataSet{PKSubject}) - for i = 1:length(data) - if sort ⊆ data[i].sort - return data[i] - end - end -end - -""" - dosetime(data::PKSubject) - -Return dosing time and tau for subject. - -data - PK subject. -""" -function dosetime(data::PKSubject) - data.dosetime -end - -""" - setth!(data::PDSubject, th::Real) - -Set treashold for subject. - -data - PK subject; -th - treashold. -""" -function setth!(data::PDSubject, th::Real) - data.th = th - data -end -function setth!(data::DataSet{PDSubject}, th::Real) - for i = 1:length(data) - setth!(data[i], th) - end - data -end -function setth!(data::DataSet{PDSubject}, th::Real, sort::Dict) - for i = 1:length(data) - if sort ⊆ data[i].sort setth!(data[i], th) end - end - data -end - -""" - setbl!(data::PDSubject, bl::Real) - -Set baseline for subject. - -data - PK subject; -bl - baseline. -""" -function setbl!(data::PDSubject, bl::Real) - data.bl = bl - data -end -function setbl!(data::DataSet{PDSubject}, bl::Real) - for i = 1:length(data) - setbl!(data[i], bl) - end - data -end -function setbl!(data::DataSet{PDSubject}, bl::Real, sort::Dict) - for i = 1:length(data) - if sort ⊆ data[i].sort setbl!(data[i], bl) end - end - data -end -#------------------------------------------------------------------------------- -""" - keldata(data::PKPDProfile{PKSubject}) - -Return elimination data for PK subject in profile. -""" -function keldata(data::PKPDProfile{PKSubject}) - return data.subject.keldata -end -function keldata(data::DataSet{PKPDProfile{PKSubject}}, i::Int) - return data[i].subject.keldata -end - - -#------------------------------------------------------------------------------- -# Param functions - -function cmax(data::PKSubject) - return maximum(data.obs) -end - -function auc() -end -#---- -function datatable(data; stack = true) - if stack return datatable_st(data) else return datatable_unst(data) end -end -function datatable_st(data) where T - - sortcols = Vector{Any}(undef, 0) - sorttype = Dict() - for i = 1:length(data) - if length(data[i].sort) > 0 - for k in keys(data[i].sort) - if k ∉ sortcols push!(sortcols, k) end - if k ∉ keys(sorttype) sorttype[k] = typeof(data[i].sort[k]) else sorttype[k] = promote_type(sorttype[k], typeof(data[i].sort[k])) end - end - end - end - pname = :Parameter - vname = :Value - - for i = 1:length(sortcols) - if sortcols[i] == pname pname = Symbol(string(pname)*"_") end - if sortcols[i] == vname pname = Symbol(string(vname)*"_") end - end - t = Vector{Vector}(undef, length(sortcols)) - for i = 1:length(t) - t[i] = Vector{sorttype[sortcols[i]]}(undef, 0) - end - nt = NamedTuple{Tuple(append!(copy(sortcols), [pname, vname]))}(Tuple(append!(copy(t), [Vector{Symbol}(undef, 0), Vector{Float64}(undef, 0)]))) - for i = 1:length(data) - sort = Dict() - for c = 1:length(sortcols) - if sortcols[c] ∈ keys(data[i].sort) sort[sortcols[c]] = data[i].sort[sortcols[c]] else sort[sortcols[c]] = missing end - end - for k in keys(data[i].result) - for s in keys(sort) - push!(nt[s], sort[s]) - end - push!(nt[pname], k) - push!(nt[vname], data[i].result[k]) - end - end - return nt -end -function datatable_unst(data) -end - -""" - auc_sparse(data::PKSubject) - -AUC for sparse data. - -```math -w_1 = (t_2 - t_1) / 2 -w_j = (t_{j+1} - t_{j-1}) / 2 (2 \\leq j \\leq J - 1) -w_J = (t_J - t_{J-1}) / 2 - -AUC = \\sum_{j=1}^J \\mu_j w_j -``` - -where `math \\mu_j` is the mean drug concentration at time t. -""" -function auc_sparse(data::PKSubject) - wts = Vector{Float64}(undef, length(data)) - wts[1] = (data.time[2] - data.time[1])/2.0 - wts[end] = (data.time[end] - data.time[end-1])/2.0 - if length(data) > 2 - for i = 2:length(data)-1 - wts[i] = (data.time[i+1] - data.time[i-1])/2.0 - end - end - return data.obs'*wts -end diff --git a/src/plots.jl b/src/plots.jl deleted file mode 100644 index 65c83f0..0000000 --- a/src/plots.jl +++ /dev/null @@ -1,426 +0,0 @@ -# Clinical Trial Utilities -# Author: Vladimir Arnautov aka PharmCat -# Copyright © 2019 Vladimir Arnautov aka PharmCat (mail@pharmcat.net) - -struct PKPlotStyle - linestyle::Symbol - linecolor::Symbol - markershape::Symbol - markercolor::Symbol -end - -const PKPLOTSTYLE = [ -PKPlotStyle(:solid, :blue, :circle, :blue), -PKPlotStyle(:solid, :red, :utriangle, :red), -PKPlotStyle(:solid, :green, :diamond, :green), -PKPlotStyle(:solid, :magenta, :pentagon, :magenta), -PKPlotStyle(:solid, :purple, :heptagon, :purple), -PKPlotStyle(:solid, :indigo, :octagon, :indigo), -PKPlotStyle(:solid, :gold, :star, :gold), -PKPlotStyle(:solid, :yellow, :rect, :yellow), -PKPlotStyle(:solid, :gray, :xcross, :gray), -PKPlotStyle(:solid, :cyan, :cross, :cyan), -PKPlotStyle(:dot, :blue, :utriangle, :blue), -PKPlotStyle(:dot, :red, :circle, :red), -PKPlotStyle(:dot, :green, :rect, :green), -PKPlotStyle(:dot, :yellow, :diamond, :yellow), -PKPlotStyle(:dot, :gray, :cross, :gray), -PKPlotStyle(:dot, :cyan, :xcross, :cyan), -PKPlotStyle(:dot, :gold, :pentagon, :gold), -PKPlotStyle(:dot, :magenta, :star, :magenta), -PKPlotStyle(:dot, :purple, :octagon, :purple), -PKPlotStyle(:dot, :indigo, :heptagon, :indigo), -] - -function randomstyle(rng) - linestyle = ([:solid, :dash, :dot, :dashdot, :dashdotdot])[sample(rng, 1:5)] - linecolor = ([:blue, :red, :green, :yellow, :gray, :cyan, :gold, :magenta, :purple, :indigo])[sample(rng, 1:10)] - markershape = ([:circle, :rect, :star5, :diamond, :hexagon, :cross, :xcross, :utriangle, :dtriangle, :rtriangle, :ltriangle, :pentagon, :heptagon, :octagon, :star4, :star6, :star7, :star8, :vline, :hline])[sample(rng, 1:20)] - markercolor = ([:blue, :red, :green, :yellow, :gray, :cyan, :gold, :magenta, :purple, :indigo])[sample(rng, 1:10)] - return PKPlotStyle(linestyle, linecolor, markershape, markercolor) -end - -@userplot PKPlot -@userplot PKElimpPlot - -function luceil(x) - fl = Int(floor(log10(x))) - if fl < 0 fl = 0 end - ceil(x/10^fl)*10^fl -end - -@recipe function f(subj::PKPlot; lcd = 5, tcd = 6) - x, y = subj.args - - if isa(lcd, Real) - lc = luceil(maximum(x->isnan(x) ? -Inf : x, y)/lcd) - yt = 0:lc:lc*lcd - elseif isa(lcd, StepRange) - yt = lcd - elseif lcd == :all - yt = y - end - - if isa(tcd, Real) - tc = luceil(maximum(x)/tcd) - xt = 0:tc:tc*tcd - elseif isa(tcd, StepRange) - xt = tcd - elseif tcd == :all - xt = x - end - - widen --> true - seriestype --> :line - xguide --> "Time" - link --> :both - legend --> true - grid --> true - gridstyle --> :auto - #ticks := [nothing :auto nothing] - #xlims --> (minimum(x), maximum(x)*1.1) - #ylims --> (0, maximum(y)*1.1) - - if !isa(lcd, Symbol) || lcd != :auto - yticks --> yt - end - if !isa(tcd, Symbol) || tcd != :auto - xticks --> xt - end - - seriescolor --> :blue - markershape --> :circle - markersize --> 3 - markercolor --> :match - markerstrokealpha --> 0 - (x, y) -end - -@recipe function f(subj::PKElimpPlot) - x, y = subj.args - seriestype --> :line - legend --> false - markersize --> 0 - markerstrokealpha --> 0 - (x, y) -end - -function plotlabel(d) - title = "" - if length(d) > 0 - for (k, v) in d - title *= "$(k) = $(v); " - end - title = title[1:length(title)] - end - return title -end - -""" - pkplot(subj::T; plotstyle::PKPlotStyle = PKPLOTSTYLE[1], kwargs...) where T <: AbstractSubject - -Plot for subject - -* subj - subject; -* plotstyle - styles for plots. - -""" -function pkplot(subj::T; plotstyle::PKPlotStyle = PKPLOTSTYLE[1], ls = false, elim = false, xticksn = 6, yticksn = 5, kwargs...) where T <: AbstractSubject - time = subj.time - obs = subj.obs - kwargs = Dict{Symbol, Any}(kwargs) - k = keys(kwargs) - - if !(:title in k) - kwargs[:title] = plotlabel(subj.sort) - end - if !(:xlims in k) - kwargs[:xlims] = (minimum(subj.time), maximum(subj.time)) - end - if !(:legend in k) - kwargs[:legend] = true - end - if !(:ylabel in k) - kwargs[:ylabel] = "Concentration" - end - if :yscale in k - if kwargs[:yscale] in [:ln, :log, :log2, :log10] ls = true end - end - if ls == true - inds = findall(x->x>0, subj.obs) - time = subj.time[inds] - obs = log.(subj.obs[inds]) - if (:ylims in k) - kwargs[:ylims] = (0,log(kwargs[:ylims][2])) - end - end - - p = pkplot(time, obs; lcd = yticksn, tcd = xticksn, - linestyle = plotstyle.linestyle, - linecolor = plotstyle.linecolor, - markershape = plotstyle.markershape, - markercolor = plotstyle.markercolor, - kwargs... - ) - if elim - if length(subj.keldata) > 0 - rsq, rsqn = findmax(subj.keldata.ar) - lz = subj.keldata.a[rsqn] - lzint = subj.keldata.b[rsqn] - ts = subj.time[subj.kelrange.kelstart] - te = subj.time[subj.kelrange.kelend] - - if ls true - x = [ts,te] - y = lzint .+ lz .* x - else - x = collect(ts:(te-ts)/100:te) - y = exp.(lzint .+ lz .* x) - end - pkelimpplot!(p, x, y) - end - end - return p - -end - -function pkplot!(subj::T; plotstyle::PKPlotStyle = PKPLOTSTYLE[1], ls = false, elim = false, xticksn = :auto, yticksn = :auto, kwargs...) where T <: AbstractSubject - time = subj.time - obs = subj.obs - kwargs = Dict{Symbol, Any}(kwargs) - k = keys(kwargs) - - if !(:xlims in k) - kwargs[:xlims] = (minimum(subj.time), maximum(subj.time)) - end - if !(:legend in k) - kwargs[:legend] = true - end - if :yscale in k - if kwargs[:yscale] in [:ln, :log, :log2, :log10] ls = true end - end - if ls == true - inds = findall(x->x>0, subj.obs) - time = subj.time[inds] - obs = log.(subj.obs[inds]) - if (:ylims in k) - kwargs[:ylims] = (0,log(kwargs[:ylims][2])) - end - end - - - p = pkplot!(time, obs; lcd = yticksn, tcd = xticksn, - linestyle = plotstyle.linestyle, - linecolor = plotstyle.linecolor, - markershape = plotstyle.markershape, - markercolor = plotstyle.markercolor, - kwargs... - ) - if :blline in k - if kwargs[:blline] && isa(subj, PDSubject) - hline!(p, [subj.bl], color = :red, label = "bl") - end - end - if :thline in k - if kwargs[:thline] && isa(subj, PDSubject) - hline!(p, [subj.th], color = :blue, label = "th") - end - end - return p -end - -function usort(data::DataSet{T}, list) where T <: AbstractData - dl = Vector{Dict}(undef, 0) - for i in data - subd = Dict(k => i.sort[k] for k in list) - if subd ∉ dl push!(dl, subd) end - end - dl -end - -function pageplot(pagedatatmp, styledict, typesort, utypes; xticksn = 6, yticksn = 5, kwargs...) - kwargs = Dict{Symbol, Any}(kwargs) - k = keys(kwargs) - - if !(:title in k) - kwargs[:title] = "" - end - if !(:xlims in k) - kwargs[:xlims] = (minimum(subj.time), maximum(subj.time)) - end - if !(:legend in k) - kwargs[:legend] = true - end - if !(:ylabel in k) - kwargs[:ylabel] = "Concentration" - end - - #utypes = keys(styledict) - - labels = Vector{String}(undef, 0) - fst = true - p = nothing - ylm = ylm2 = 0.0 - for u in utypes - for d in pagedatatmp - if u ⊆ d.sort - label = "AUTO" - if typesort !== nothing - tempsd = Dict(k => d.sort[k] for k in typesort) - style = styledict[tempsd] - label = plotlabel(tempsd) - if label ∉ labels - push!(labels, label) - else - label = "" - end - end - - kwargs[:label] = label - if fst - ylm = maximum(x->isnan(x) ? -Inf : x, d.obs) - #MNIMUM - p = pkplot(d; plotstyle = style, xticksn = xticksn, yticksn = yticksn, ylims = (0, ylm*1.2), kwargs...) - fst = false - else - ylm2 = maximum(x->isnan(x) ? -Inf : x, d.obs) - #MNIMUM - if ylm2 > ylm ylm = ylm2 end - pkplot!( d; plotstyle = style, xticksn = xticksn, yticksn = yticksn, ylims = (0, ylm*1.2), kwargs...) - end - end - end - end - if :vline in k - vline!(p, kwargs[:vline]) - end - if :hline in k - hline!(p, kwargs[:hline]) - end - p -end -function pageplot(pagedatatmp; elim = false, xticksn = 6, yticksn = 5, kwargs...) - kwargs = Dict{Symbol, Any}(kwargs) - k = keys(kwargs) - - if !(:title in k) - kwargs[:title] = "" - end - if !(:xlims in k) - kwargs[:xlims] = (minimum(subj.time), maximum(subj.time)) - end - if !(:legend in k) - kwargs[:legend] = true - end - if !(:ylabel in k) - kwargs[:ylabel] = "Concentration" - end - - fst = true - p = nothing - if length(pagedatatmp) > 1 elim = false end - ylm = ylm2 = 0.0 - for d in pagedatatmp - kwargs[:label] = "AUTO" - if fst - ylm = maximum(x->isnan(x) ? -Inf : x, d.obs) - #MNIMUM - p = pkplot(d; plotstyle = PKPLOTSTYLE[1], elim = elim, xticksn = xticksn, yticksn = yticksn, ylims = (0, ylm*1.2), kwargs...) - fst = false - else - ylm2 = maximum(x->isnan(x) ? -Inf : x, d.obs) - #MNIMUM - if ylm2 > ylm ylm = ylm2 end - pkplot!( d; plotstyle = PKPLOTSTYLE[1], xticksn = xticksn, yticksn = yticksn, ylims = (0, ylm*1.2), kwargs...) - end - end - if :vline in k - vline!(p, kwargs[:vline]) - end - if :hline in k - hline!(p, kwargs[:hline]) - end - p -end - -""" - pkplot(data::DataSet{T}; pagesort::Union{Nothing, Vector{Symbol}} = nothing, typesort::Union{Nothing, Vector{Symbol}} = nothing, kwargs...) where T <: AbstractSubject - -Plot for subjects in dataset. - -* data - subjects dataset; -* pagesort - subject page groupping; -* typesort - subject sorting within page; -""" -function pkplot(data::DataSet{T}; pagesort::Union{Nothing, Vector{Symbol}} = nothing, typesort::Union{Nothing, Vector{Symbol}} = nothing, elim = false, xticksn = :auto, yticksn = :auto, kwargs...) where T <: AbstractSubject - # Style types - kwargs = Dict{Symbol, Any}(kwargs) - k = keys(kwargs) - - if !(:title in k) - kwargs[:title] = :AUTO - end - if !(:xlims in k) - kwargs[:xlims] = nothing - end - if !(:legend in k) - kwargs[:legend] = true - end - if !(:ylabel in k) - kwargs[:ylabel] = "Concentration" - end - - styledict = nothing - local utypes - #style = PKPLOTSTYLE[1] - if pagesort !== nothing - if !isa(pagesort, Array) pagesort = [pagesort] end - end - if typesort !== nothing - styledict = Dict() - if !isa(typesort, Array) typesort = [typesort] end - utypes = usort(data, typesort) - for i = 1:length(utypes) - if i <= 20 - styledict[utypes[i]] = PKPLOTSTYLE[i] - else - styledict[utypes[i]] = randomstyle(MersenneTwister(34534)) - end - end - end - autotitle = false - if kwargs[:title] == :AUTO autotitle = true end - #--------------------------------------------------------------------------- - # PAGES - plots = Vector{Any}(undef, 0) - p = nothing - if pagesort !== nothing - pagedict = usort(data, pagesort) - for i in pagedict - if autotitle - kwargs[:title] = plotlabel(i) - end - pagedatatmp = Vector{eltype(data)}(undef, 0) - for d in data - if i ⊆ d.sort - push!(pagedatatmp, d) - end - end - if typesort !== nothing - p = pageplot(pagedatatmp, styledict, typesort, utypes; xticksn = xticksn, yticksn = yticksn, kwargs...) - else - p = pageplot(pagedatatmp; elim = elim, xticksn = xticksn, yticksn = yticksn, kwargs...) - end - push!(plots, p) - end - else - if kwargs[:title] == :AUTO kwargs[:title] = "" end - if typesort !== nothing - p = pageplot(data, styledict, typesort, utypes; xticksn = xticksn, yticksn = yticksn, kwargs...) - else - p = pageplot(data; elim = elim, xticksn = xticksn, yticksn = yticksn, kwargs...) - end - return p - end - plots -end diff --git a/test/csv/nanpk.csv b/test/csv/nanpk.csv deleted file mode 100644 index a60d43e..0000000 --- a/test/csv/nanpk.csv +++ /dev/null @@ -1,22 +0,0 @@ -Concentration,Time,Subject,Formulation -,0,1,T -0.2,1,1,T -0.3,2,1,T -0.4,3,1,T -0.3,4,1,T -,5,1,T -,6,1,T -0,0,2,T -0.3,1,2,T -0.4,2,2,T -0.5,3,2,T -,4,2,T -0.1,5,2,T -0.05,6,2,T -0,0,1,R -,1,1,R -0.9,2,1,R -0.8,3,1,R -0.3,4,1,R -0,5,1,R -0.1,6,1,R diff --git a/test/csv/pddata.csv b/test/csv/pddata.csv deleted file mode 100644 index e549dac..0000000 --- a/test/csv/pddata.csv +++ /dev/null @@ -1,14 +0,0 @@ -time;effect -0;0 -1;1 -2;4 -2.5;7 -3;5 -3.3;4 -3.6;3 -4;2 -5;8 -6;1 -7;2 -8;1 -9;1 diff --git a/test/csv/pkdata.csv b/test/csv/pkdata.csv deleted file mode 100644 index e163d80..0000000 --- a/test/csv/pkdata.csv +++ /dev/null @@ -1,22 +0,0 @@ -Concentration,Time,Subject,Formulation -0,0,1,T -0.2,1,1,T -0.3,2,1,T -0.4,3,1,T -0.3,4,1,T -0.2,5,1,T -0.1,6,1,T -0,0,2,T -0.3,1,2,T -0.4,2,2,T -0.5,3,2,T -0.3,4,2,T -0.1,5,2,T -0.05,6,2,T -0.1,0,1,R -0.2,1,1,R -0.9,2,1,R -0.8,3,1,R -0.3,4,1,R -0.2,5,1,R -0.1,6,1,R diff --git a/test/csv/pkdata2.csv b/test/csv/pkdata2.csv deleted file mode 100644 index 6a40996..0000000 --- a/test/csv/pkdata2.csv +++ /dev/null @@ -1,161 +0,0 @@ -Subject Formulation Time Concentration -1 T 0 0 -1 T 0.5 178.949 -1 T 1 190.869 -1 T 1.5 164.927 -1 T 2 139.962 -1 T 2.5 129.59 -1 T 3 131.369 -1 T 4 150.854 -1 T 5 121.239 -1 T 6 139.229 -1 T 8 128.52 -1 T 10 143.243 -1 T 12 144.964 -1 T 24 133.16 -1 T 48 137.271 -1 T 72 112.846 -2 R 0 0 -2 R 0.5 62.222 -2 R 1 261.177 -2 R 1.5 234.063 -2 R 2 234.091 -2 R 2.5 222.881 -2 R 3 213.896 -2 R 4 196.026 -2 R 5 199.634 -2 R 6 196.037 -2 R 8 213.352 -2 R 10 200.088 -2 R 12 196.035 -2 R 24 160.338 -2 R 48 110.28 -2 R 72 85.241 -3 R 0 0 -3 R 0.5 49.849 -3 R 1 77.367 -3 R 1.5 105.345 -3 R 2 100.943 -3 R 2.5 72.746 -3 R 3 69.985 -3 R 4 93.565 -3 R 5 91.981 -3 R 6 82.71 -3 R 8 84.205 -3 R 10 85.342 -3 R 12 76.027 -3 R 24 81.259 -3 R 48 70.107 -3 R 72 67.901 -4 R 0 0 -4 R 0.5 52.421 -4 R 1 208.542 -4 R 1.5 188.923 -4 R 2 165.177 -4 R 2.5 146.996 -4 R 3 152.701 -4 R 4 154.345 -4 R 5 128.398 -4 R 6 149.807 -4 R 8 151.066 -4 R 10 136.819 -4 R 12 132.257 -4 R 24 141.247 -4 R 48 129.138 -4 R 72 97.625 -5 T 0 0 -5 T 0.5 0 -5 T 1 9.545 -5 T 1.5 153.964 -5 T 2 152.34 -5 T 2.5 151.452 -5 T 3 161.312 -5 T 4 169.334 -5 T 5 162.907 -5 T 6 166.651 -5 T 8 168.668 -5 T 10 155.103 -5 T 12 154.066 -5 T 24 162.974 -5 T 48 109.814 -5 T 72 110.778 -6 T 0 0 -6 T 0.5 57.882 -6 T 1 100.498 -6 T 1.5 138.651 -6 T 2 147.287 -6 T 2.5 154.648 -6 T 3 122.316 -6 T 4 132.857 -6 T 5 126.067 -6 T 6 140.466 -6 T 8 115.542 -6 T 10 102.16 -6 T 12 113.751 -6 T 24 101.049 -6 T 48 92.55 -6 T 72 69.501 -7 R 0 0 -7 R 0.5 19.95 -7 R 1 128.405 -7 R 1.5 136.807 -7 R 2 113.109 -7 R 2.5 153.254 -7 R 3 123.606 -7 R 4 142.655 -7 R 5 112.347 -7 R 6 139.919 -7 R 8 105.513 -7 R 10 134.408 -7 R 12 123.37 -7 R 24 110.511 -7 R 48 90.291 -7 R 72 58.051 -8 R 0 0 -8 R 0.5 136.91 -8 R 1 126.646 -8 R 1.5 118.5 -8 R 2 134.926 -8 R 2.5 113.213 -8 R 3 130.896 -8 R 4 138.327 -8 R 5 22.724 -8 R 6 53.774 -8 R 8 55.107 -8 R 10 102.871 -8 R 12 134.133 -8 R 24 108.021 -8 R 48 98.466 -8 R 72 74.437 -9 T 0 0 -9 T 0.5 113.362 -9 T 1 128.273 -9 T 1.5 125.395 -9 T 2 146.933 -9 T 2.5 140.559 -9 T 3 167.347 -9 T 4 157.504 -9 T 5 141.35 -9 T 6 140.282 -9 T 8 105.438 -9 T 10 164.843 -9 T 12 135.58 -9 T 24 117.125 -9 T 48 109.745 -9 T 72 93.44 -10 R 0 0 -10 R 0.5 13.634 -10 R 1 62.561 -10 R 1.5 112.655 -10 R 2 125.482 -10 R 2.5 116.255 -10 R 3 112.674 -10 R 4 116.986 -10 R 5 119.81 -10 R 6 107.557 -10 R 8 120.479 -10 R 10 124.171 -10 R 12 106.476 -10 R 24 116.508 -10 R 48 45.204 -10 R 72 42.191 diff --git a/test/csv/upkdata.csv b/test/csv/upkdata.csv deleted file mode 100644 index 636a2e9..0000000 --- a/test/csv/upkdata.csv +++ /dev/null @@ -1,6 +0,0 @@ -subj conc st et vol -1 1 0 1 1 -1 2 1 2 2 -1 2 2 6 3 -1 1 6 12 3 -1 1 12 18 2 diff --git a/test/internal.jl b/test/internal.jl index 532c7fe..6d53445 100644 --- a/test/internal.jl +++ b/test/internal.jl @@ -83,26 +83,4 @@ end @test ci[3] == 2 @test ci[4] == 0.05 end - -@testset " DataSet " begin - pkds = ClinicalTrialUtilities.pkimport(pkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - @test eltype(pkds) == eltype(pkds.data) - a = pkds[:Subject => 1] - b = pkds[Dict(:Subject => 1)] - @test a === b - - b = findall(pkds, Dict(:Formulation => "T")) - @test length(b) == 2 - - deleteat!(pkds, Dict(:Formulation => "T")) - @test length(pkds) == 1 - - pkds = ClinicalTrialUtilities.pkimport(pkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - deleteat!(pkds, 1) - @test length(pkds) == 2 - - pkds = ClinicalTrialUtilities.pkimport(pkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - deleteat!(pkds, [1,2]) - @test length(pkds) == 1 -end end diff --git a/test/pktest.jl b/test/pktest.jl deleted file mode 100644 index 7a85546..0000000 --- a/test/pktest.jl +++ /dev/null @@ -1,1264 +0,0 @@ -println(" ---------------------------------- ") -@testset "#10 Pharmacokinetics " begin - - # Basic tests - pkds = ClinicalTrialUtilities.pkimport(pkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - pk = ClinicalTrialUtilities.nca!(pkds) - @test pk[1, :AUCinf] ≈ 1.63205 atol=1E-5 - @test pk[1, :Cmax] ≈ 0.4 atol=1E-5 - @test pk[1, :MRTlast] ≈ 3.10345 atol=1E-5 - @test pk[1, :Tmax] ≈ 3.0 atol=1E-5 - - pk = ClinicalTrialUtilities.nca!(pkds; calcm = :logt) - @test pk[1, :AUClast] ≈ 1.43851 atol=1E-5 - @test pk[1, :AUMClast] ≈ 4.49504 atol=1E-5 - - pk = ClinicalTrialUtilities.nca!(pkds; calcm = :luld) - @test pk[1, :AUClast] ≈ 1.43851 atol=1E-5 - @test pk[1, :AUMClast] ≈ 4.49504 atol=1E-5 - - #setkelauto! getkelauto - ClinicalTrialUtilities.setkelauto!(pkds, false, [1,2]) - ClinicalTrialUtilities.setkelauto!(pkds, true, 1) - @test ClinicalTrialUtilities.getkelauto(pkds[1]) == true - - ClinicalTrialUtilities.applyelimrange!(pk, ClinicalTrialUtilities.ElimRange(4,7)) - ClinicalTrialUtilities.applyelimrange!(pk, ClinicalTrialUtilities.ElimRange(5,7), [1,2,3]) - ClinicalTrialUtilities.applyelimrange!(pk, ClinicalTrialUtilities.ElimRange(4,7), 1) - ClinicalTrialUtilities.applyelimrange!(pk, ClinicalTrialUtilities.ElimRange(5,7), Dict(:Subject => 2)) - @test pk[1].subject.kelrange.kelstart == 4 - @test pk[2].subject.kelrange.kelstart == 5 - - pkds = ClinicalTrialUtilities.pkimport(pkdata[1:7,:], [:Subject, :Formulation]; conc = :Concentration, time = :Time) - pk = ClinicalTrialUtilities.nca!(pkds) - @test pk[1, :AUCinf] ≈ 1.63205 atol=1E-5 - @test pk[1, :Cmax] ≈ 0.4 atol=1E-5 - @test pk[1, :MRTlast] ≈ 3.10345 atol=1E-5 - @test pk[1, :Tmax] ≈ 3.0 atol=1E-5 - - pkds = ClinicalTrialUtilities.pkimport(pkdata[1:7,2], pkdata[1:7,1]) - pk = ClinicalTrialUtilities.nca!(pkds) - @test pk[:AUClast] ≈ 1.4499999999999997 - -#Program: - #AUClast - #AUCinf - #AUCpct - #Tmax - #Cmax - #Clast - #Adjusted R sq - #Kel - #HL - #AUMClast - #MRTlast - #Clast_pred - #Clinf / Cl_F_Obs - #Vzinf / Vz_F_Pred - - # -------------------------------------------------------------------------- - # Linear Trapezoidal Linear Interpolation - # calcm = :lint, intp = :lint - # DS1 - pkds = ClinicalTrialUtilities.pkimport(pkdata2, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 100, time = 0)) - pk = ClinicalTrialUtilities.nca!(pkds) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUC last - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9585.4218 - 10112.176 - 5396.5498 - 9317.8358 - 9561.26 - 6966.598 - 7029.5735 - 7110.6745 - 8315.0803 - 5620.8945], sigdigits = 6) - #AUCinf - @test round.(df[!, :AUCinf], sigdigits = 6) == round.([42925.019 - 16154.93 - 26026.183 - 22004.078 - 25820.275 - 16001.76 - 11688.953 - 15446.21 - 24865.246 - 8075.3242], sigdigits = 6) - #AUCpct - @test round.(df[!, :AUCpct], sigdigits = 5) == round.([77.669383 - 37.405019 - 79.26492 - 57.65405 - 62.969953 - 56.463551 - 39.861391 - 53.964925 - 66.559429 - 30.394194], sigdigits = 5) - #Tmax - @test df[!, :Tmax] == [1 - 1 - 1.5 - 1 - 4 - 2.5 - 2.5 - 4 - 3 - 2] - #Cmax - @test df[!, :Cmax] == [190.869 - 261.177 - 105.345 - 208.542 - 169.334 - 154.648 - 153.254 - 138.327 - 167.347 - 125.482] - #Clast - @test df[!, :Clast] == [112.846 - 85.241 - 67.901 - 97.625 - 110.778 - 69.501 - 58.051 - 74.437 - 93.44 - 42.191] - # Adjusted R sq - @test round.(df[!, :ARsq], digits = 6) == round.([0.71476928 - 0.99035145 - 0.77630678 - 0.83771737 - 0.82891994 - 0.92517856 - 0.96041642 - 0.92195356 - 0.92130684 - 0.86391165], digits = 6) - - # R sq - @test round.(df[!, :Rsq], digits = 6) == round.([0.78607696 - 0.99276359 - 0.81358898 - 0.91885869 - 0.85335995 - 0.95011904 - 0.97031231 - 0.94796904 - 0.94753789 - 0.88092269], digits = 6) - #Kel - @test round.(df[!, :Kel], sigdigits = 6) == round.([0.0033847439 - 0.014106315 - 0.0032914304 - 0.0076953442 - 0.0068133279 - 0.0076922807 - 0.012458956 - 0.0089300798 - 0.0056458649 - 0.017189737], sigdigits = 6) - #HL - @test round.(df[!, :HL], sigdigits = 5) == round.([204.78571 - 49.137367 - 210.59148 - 90.073577 - 101.73401 - 90.10945 - 55.634451 - 77.619371 - 122.77077 - 40.323315], sigdigits = 5) - #AUMC last - @test round.(df[!, :AUMClast], sigdigits = 6) == round.([333582.48 - 298701.39 - 186032.06 - 313955.9 - 315181.56 - 226977.06 - 219797.71 - 240526.05 - 277613.98 - 154893.06], sigdigits = 6) - #MRTlast - @test round.(df[!, :MRTlast], digits = 6) == round.([34.801023 - 29.538786 - 34.472406 - 33.69408 - 32.964438 - 32.58076 - 31.267574 - 33.826053 - 33.386807 - 27.556657], digits = 6) - #Clast_pred - @test round.(df[!, :Clast_pred], sigdigits = 6) == round.([117.30578 - 82.53669 - 66.931057 - 100.76793 - 105.29832 - 71.939942 - 61.172702 - 75.604277 - 93.761762 - 38.810857], sigdigits = 6) - - @test round.(df[!, :Clinf], sigdigits = 6) == round.([0.0023296437 - 0.0061900608 - 0.0038422846 - 0.0045446122 - 0.0038729255 - 0.0062493127 - 0.0085550864 - 0.0064740799 - 0.0040216775 - 0.012383404], sigdigits = 6) - - @test round.(df[!, :Vzinf], sigdigits = 6) == round.([0.68827768 - 0.43881487 - 1.1673601 - 0.59056646 - 0.56843374 - 0.8124135 - 0.68666158 - 0.72497447 - 0.71232266 - 0.72039519], sigdigits = 6) - -ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 120, time = 0, tau = 12)) -pk = ClinicalTrialUtilities.nca!(pkds; adm = :iv) -df = DataFrame(pk; unst = true) -sort!(df, :Subject) - -#Cltau / CLss -@test round.(df[!, :Cltau], sigdigits = 5) == round.([0.07185191 -0.050414459 -0.12240579 -0.070132959 -0.06902661 -0.085106504 -0.083532913 -0.10859036 -0.073251565 -0.092756742], sigdigits = 5) - -#AUMCtau -@test round.(df[!, :AUMCtau], sigdigits = 6) == round.([9984.8168 -14630.069 -6024.4953 -10299.721 -11466.123 -8467.3568 -9003.0193 -6457.0058 -10095.818 -8367.3005], sigdigits = 6) - -#AUCinf -@test round.(df[!, :AUCinf], sigdigits = 6) == round.([42925.019 -16154.93 -26026.183 -22004.078 -25714.393 -16001.76 -11688.953 -15446.21 -24865.246 -8171.1624], sigdigits = 6) - -#MRTtauinf / MRTinf -@test round.(df[!, :MRTtauinf], sigdigits = 6) == round.([302.40303 -75.590599 -312.72083 -148.34069 -172.0933 -130.19061 -91.908297 -161.57402 -176.30461 -70.260736], sigdigits = 6) - - -@test round.(df[!, :Vsstau], sigdigits = 4) == round.([21.728235 -3.8108592 -38.278842 -10.403572 -11.879017 -11.080068 -7.6773678 -17.545381 -12.914588 -6.517157], sigdigits = 4) - - - # -------------------------------------------------------------------------- - # Linear Log Trapezoidal - # Log-trapezoidat rule after Tmax if c₁ > 0 and c₂ > 0, else Linear trapezoidal used; - - - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :logt, intp = :logt, io = io, verbose = true) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9572.8582 - 10054.037 - 5391.5322 - 9296.2179 - 9518.6531 - 6948.5757 - 6987.0645 - 7064.7816 - 8298.9634 - 5485.6538], sigdigits = 6) - #AUCinf - @test round.(df[!, :AUCinf], sigdigits = 6) == round.([42912.456 - 16096.791 - 26021.165 - 21982.46 - 25777.668 - 15983.737 - 11646.444 - 15400.317 - 24849.129 - 7940.0834], sigdigits = 6) - #AUCpct - @test round.(df[!, :AUCpct], sigdigits = 5) == round.([77.692122 - 37.540119 - 79.280204 - 57.710748 - 63.074033 - 56.527216 - 40.006884 - 54.12574 - 66.602599 - 30.911888], sigdigits = 5) - #Tmax - @test df[!, :Tmax] == [1 - 1 - 1.5 - 1 - 4 - 2.5 - 2.5 - 4 - 3 - 2] - #Cmax - @test df[!, :Cmax] == [190.869 - 261.177 - 105.345 - 208.542 - 169.334 - 154.648 - 153.254 - 138.327 - 167.347 - 125.482] - #Clast - @test df[!, :Clast] == [112.846 - 85.241 - 67.901 - 97.625 - 110.778 - 69.501 - 58.051 - 74.437 - 93.44 - 42.191] - #Adjusted R sq - @test round.(df[!, :ARsq], digits = 6) == round.([0.71476928 - 0.99035145 - 0.77630678 - 0.83771737 - 0.82891994 - 0.92517856 - 0.96041642 - 0.92195356 - 0.92130684 - 0.86391165], digits = 6) - #Kel - @test round.(df[!, :Kel], sigdigits = 6) == round.([0.003384744 - 0.014106315 - 0.00329143 - 0.007695344 - 0.006813328 - 0.007692281 - 0.012458956 - 0.00893008 - 0.005645865 - 0.017189737], sigdigits = 6) - #HL - @test round.(df[!, :HL], sigdigits = 5) == round.([204.78571 - 49.137367 - 210.59148 - 90.073577 - 101.73401 - 90.10945 - 55.634451 - 77.619371 - 122.77077 - 40.323315], sigdigits = 5) - - #Clast_pred - @test round.(df[!, :Clast_pred], sigdigits = 6) == round.([117.30578 - 82.53669 - 66.931057 - 100.76793 - 105.29832 - 71.939942 - 61.172702 - 75.604277 - 93.761762 - 38.810857], sigdigits = 6) - - - # -------------------------------------------------------------------------- - # Linear Up Log Down - # Linear Up Log Down everywhere if c₁ > c₂ > 0, else Linear trapezoidal used; - - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :luld, intp = :logt, io = io, verbose = true) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9573.8106 - 10054.286 - 5392.4572 - 9297.0963 - 9519.1809 - 6948.9856 - 6988.7726 - 7073.0922 - 8303.3586 - 5486.8389], sigdigits = 6) - #AUCinf - @test round.(df[!, :AUCinf], sigdigits = 6) == round.([42913.408 - 16097.041 - 26022.09 - 21983.338 - 25778.196 - 15984.147 - 11648.152 - 15408.628 - 24853.524 - 7941.2686], sigdigits = 6) - #AUCpct - #Tmax - @test round.(df[!, :Tmax], sigdigits = 6) == round.([1 - 1 - 1.5 - 1 - 4 - 2.5 - 2.5 - 4 - 3 - 2], sigdigits = 6) - #Cmax - @test round.(df[!, :Cmax], sigdigits = 6) == round.([190.869 - 261.177 - 105.345 - 208.542 - 169.334 - 154.648 - 153.254 - 138.327 - 167.347 - 125.482], sigdigits = 6) - #Clast - @test round.(df[!, :Clast], sigdigits = 6) == round.([112.846 - 85.241 - 67.901 - 97.625 - 110.778 - 69.501 - 58.051 - 74.437 - 93.44 - 42.191], sigdigits = 6) - #Adjusted R sq - @test round.(df[!, :ARsq], sigdigits = 6) == round.([0.71476928 - 0.99035145 - 0.77630678 - 0.83771737 - 0.82891994 - 0.92517856 - 0.96041642 - 0.92195356 - 0.92130684 - 0.86391165], sigdigits = 6) - #Kel - @test round.(df[!, :Kel], sigdigits = 6) == round.([0.003384744 - 0.014106315 - 0.00329143 - 0.007695344 - 0.006813328 - 0.007692281 - 0.012458956 - 0.00893008 - 0.005645865 - 0.017189737], sigdigits = 6) - #HL - @test round.(df[!, :HL], sigdigits = 5) == round.([204.78571 - 49.137367 - 210.59148 - 90.073577 - 101.73401 - 90.10945 - 55.634451 - 77.619371 - 122.77077 - 40.323315], sigdigits = 5) - #AUMClast - #MRTlast - #Clast_pred - @test round.(df[!, :Clast_pred], sigdigits = 6) == round.([117.30578 - 82.53669 - 66.931057 - 100.76793 - 105.29832 - 71.939942 - 61.172702 - 75.604277 - 93.761762 - 38.810857], sigdigits = 6) - - # -------------------------------------------------------------------------- - # Linear Trapezoidal Linear/Log Interpolation. - # calcm = :lint, intp = :logt - # Linear trapezoidal everywhere + log interpolation - - pkds = ClinicalTrialUtilities.pkimport(pkdata2, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - pk = ClinicalTrialUtilities.nca!(pkds; intp = :logt) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9585.4218 - 10112.176 - 5396.5498 - 9317.8358 - 9561.26 - 6966.598 - 7029.5735 - 7110.6745 - 8315.0803 - 5620.8945], sigdigits = 6) - - - # -------------------------------------------------------------------------- - # Not present in Phoenix - # Linear Up Log Down after Tmax if c₁ > c₂ > 0, else Linear trapezoidal used; - - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :luldt, intp = :logt, io = io, verbose = true) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9573.810558691312 - 10054.286478059563 - 5392.457219413793 - 9297.096334450325 - 9519.181808199797 - 6948.985621117448 - 6988.960344867885 - 7073.306755718137 - 8303.373085532965 - 5486.838889441992], sigdigits = 6) - #AUCinf - @test round.(df[!, :AUCinf], sigdigits = 6) == round.([42913.407881300096 - 16097.04111262767 - 26022.090028134975 - 21983.33845321829 - 25778.19670343543 - 15984.14736468632 - 11648.33951823218 - 15408.84256127436 - 24853.53879891218 - 7941.268553852982], sigdigits = 6) - - - # -------------------------------------------------------------------------- - # NOT USED - # Log-trapezoidal Linear Interpolation - # if calcm == :logt / :luld / :luldt calculation method used - log interpolation used when possible - #= - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 100, time = 0, tau = 0)) - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :logt, intp = :lint) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9572.8582 - 10054.037 - 5391.5322 - 9296.2179 - 9518.6531 - 6948.5757 - 6987.0645 - 7064.7816 - 8298.9634 - 5485.6538], sigdigits = 6) - #AUCinf - @test round.(df[!, :AUCinf], sigdigits = 6) == round.([42912.456 - 16096.791 - 26021.165 - 21982.46 - 25777.668 - 15983.737 - 11646.444 - 15400.317 - 24849.129 - 7940.0834], sigdigits = 6) - #AUCpct - @test round.(df[!, :AUCpct], sigdigits = 6) == round.([77.692122 - 37.540119 - 79.280204 - 57.710748 - 63.074033 - 56.527216 - 40.006884 - 54.12574 - 66.602599 - 30.911888], sigdigits = 6) - #Tmax - @test round.(df[!, :Tmax], sigdigits = 6) == round.([1 - 1 - 1.5 - 1 - 4 - 2.5 - 2.5 - 4 - 3 - 2], sigdigits = 6) - #Cmax - @test round.(df[!, :Cmax], sigdigits = 6) == round.([190.869 - 261.177 - 105.345 - 208.542 - 169.334 - 154.648 - 153.254 - 138.327 - 167.347 - 125.482], sigdigits = 6) - #Clast - @test round.(df[!, :Clast], sigdigits = 6) == round.([112.846 - 85.241 - 67.901 - 97.625 - 110.778 - 69.501 - 58.051 - 74.437 - 93.44 - 42.191], sigdigits = 6) - #Adjusted R sq :ARsq - @test round.(df[!, :ARsq], sigdigits = 6) == round.([0.71476928 - 0.99035145 - 0.77630678 - 0.83771737 - 0.82891994 - 0.92517856 - 0.96041642 - 0.92195356 - 0.92130684 - 0.86391165], sigdigits = 6) - #Kel - @test round.(df[!, :Kel], sigdigits = 6) == round.([0.003384744 - 0.014106315 - 0.00329143 - 0.007695344 - 0.006813328 - 0.007692281 - 0.012458956 - 0.00893008 - 0.005645865 - 0.017189737], sigdigits = 6) - #HL - @test round.(df[!, :HL], sigdigits = 5) == round.([204.78571 - 49.137367 - 210.59148 - 90.073577 - 101.73401 - 90.10945 - 55.634451 - 77.619371 - 122.77077 - 40.323315], sigdigits = 5) - =# - - - # -------------------------------------------------------------------------- - # -------------------------------------------------------------------------- - # Steady state or Tau defined - # -------------------------------------------------------------------------- - # Linear-trapezoidal, Log Interpolation - - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 120, time = 0, tau = 12)) - pk = ClinicalTrialUtilities.nca!(pkds, intp = :logt) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUCtau - #Ctau - #Cavg - #Swingtau - #Fluctau - #Accind - - #AUCtau - @test round.(df[!, :AUCtau], sigdigits = 6) == round.([1670.1018 - 2380.2695 - 980.34575 - 1711.0358 - 1738.46 - 1409.998 - 1436.5595 - 1105.0705 - 1638.1903 - 1293.7065], sigdigits = 6) - - #Ctau - @test round.(df[!, :Ctau], sigdigits = 6) == round.([144.964 -196.035 -76.027 -132.257 -154.066 -113.751 -123.37 -134.133 -135.58 -106.476], sigdigits = 6) - - #Cavg - @test round.(df[!, :Cavg], sigdigits = 6) == round.([139.17515 -198.35579 -81.695479 -142.58631 -144.87167 -117.49983 -119.71329 -92.089208 -136.51585 -107.80888], sigdigits = 6) - - #Cavg - @test round.(df[!, :Fluctau], sigdigits = 2) .* 100 == round.([32.983619 -32.840987 -35.886931 -53.500928 -10.538983 -34.806007 -24.962976 -4.5542796 -23.269825 -17.629346], sigdigits = 2) - - #Swingtau - @test round.(df[!, :Swingtau], sigdigits = 6) == round.([0.31666483 -0.3322978 -0.38562616 -0.57679367 -0.099100386 -0.35953091 -0.24223069 -0.031267473 -0.23430447 -0.17850032 -], sigdigits = 6) - - #Accind - @test round.(df[!, :Accind], sigdigits = 6) == round.([25.123662 -6.4216193 -25.821565 -11.336753 -12.737742 -11.341063 -7.2010832 -9.8406852 -15.26571 -5.3650315], sigdigits = 6) - - #Cltau / CLss_F - @test round.(df[!, :Cltau], sigdigits = 6) == round.([0.07185191 -0.050414459 -0.12240579 -0.070132959 -0.06902661 -0.085106504 -0.083532913 -0.10859036 -0.073251565 -0.092756742], sigdigits = 6) - -#Vztau / Vz_F -@test round.(df[!, :Vztau], sigdigits = 6) == round.([21.228167 -3.5738929 -37.18924 -9.113687 -10.131115 -11.063884 -6.7046479 -12.160066 -12.974374 -5.3960536], sigdigits = 6) - - # -------------------------------------------------------------------------- - # Linear-trapezoidal Linear Interpolation - # TAU 2 - 10 - - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 120, time = 2, tau = 10)) - pk = ClinicalTrialUtilities.nca!(pkds) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - # Cmax - @test round.(df[!, :Cmax], sigdigits = 6) == round.([150.854 - 234.091 - 100.943 - 165.177 - 169.334 - 154.648 - 153.254 - 138.327 - 167.347 - 125.482], sigdigits = 6) - - #AUCtau - @test round.(df[!, :AUCtau], sigdigits = 6) == round.([ - 1367.7388 - 2043.0158 - 838.8295 - 1444.7985 - 1618.6205 - 1224.6608 - 1265.7013 - 880.311 - 1417.942 - 1167.911], sigdigits = 6) - - # -------------------------------------------------------------------------- - # Linear Up Log Down, Log Interpolation - # TAU 0 - 36 - - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 100, time = 0, tau = 36)) - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :luld, intp = :logt) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9573.810559 - 10054.28648 - 5392.457219 - 9297.096334 - 9519.180874 - 6948.985621 - 6988.772632 - 7073.092214 - 8303.358586 - 5486.838889], sigdigits = 6) - #AUCtau - @test round.(df[!, :AUCtau], sigdigits = 6) == round.([4947.052768 - 6265.637822 - 2863.392745 - 5008.644865 - 5415.246398 - 3882.203934 - 4096.937951 - 3802.30601 - 4531.779637 - 3744.768624], sigdigits = 6) - #Ctau - @test round.(df[!, :Ctau], sigdigits = 6) == round.([135.2155 - 132.9739623 - 75.47731257 - 135.0568587 - 133.7790224 - 96.70617845 - 99.89068375 - 103.1329035 - 113.3749669 - 72.57153458], sigdigits = 6) - #Cavg - @test round.(df[!, :Cavg], sigdigits = 6) == round.([137.4181324 - 174.0454951 - 79.53868737 - 139.129024 - 150.423511 - 107.8389982 - 113.803832 - 105.6196114 - 125.8827677 - 104.0213507], sigdigits = 6) - #Swingtau - @test round.(df[!, :Swingtau], sigdigits = 6) == round.([0.411591127 - 0.964121363 - 0.39571742 - 0.544105216 - 0.265773938 - 0.599153255 - 0.534217149 - 0.341249934 - 0.476048942 - 0.729080151], sigdigits = 6) - #Fluctau - @test round.(df[!, :Fluctau] * 100, sigdigits = 6) == round.([40.49938608 - 73.660647 - 37.5511445 - 52.81798086 - 23.6365827 - 53.72993308 - 46.8906146 - 33.32155462 - 42.87483828 - 50.8650052], sigdigits = 6) - #Accind - @test round.(df[!, :Accind], sigdigits = 6) == round.([8.716910716 - 2.511311393 - 8.949296376 - 4.132742748 - 4.597396017 - 4.1341712 - 2.766795121 - 3.637329819 - 5.43694762 - 2.167194353], sigdigits = 6) - - # Steady state, partial areas - # -------------------------------------------------------------------------- - # Linear Up Log Downl, Log Interpolation - # TAU 0.25 - 9 - - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 100, time = 0.25, tau = 9)) - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :luld, intp = :logt, io = io, verbose = true) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9566.596809 - 10054.28648 - 5392.457219 - 9297.096334 - 9519.180874 - 6948.985621 - 6988.772632 - 7058.818964 - 8302.368086 - 5486.838889], sigdigits = 6) - #AUCtau - @test round.(df[!, :AUCtau], sigdigits = 6) == round.([1268.275631 - 1831.820528 - 754.6493604 - 1336.480932 - 1310.904516 - 1114.240351 - 1079.366854 - 766.620245 - 1219.631864 - 970.3062692], sigdigits = 6) - # AUMCtau - @test round.(df[!, :AUMCtau], sigdigits = 6) == round.([5477.20423544297 - 8367.57088170951 - 3455.34643479800 - 6014.64604481587 - 6609.78830163090 - 5064.72384740413 - 4976.96365993911 - 2863.00517022791 - 5386.88322025614 - 4713.47970846693], sigdigits = 6) - #Ctau - @test round.(df[!, :Ctau], sigdigits = 6) == round.([137.721875 - 204.9625055 - 84.915625 - 141.9969549 - 160.0570781 - 106.9862605 - 123.572375 - 84.9595 - 142.566125 - 122.7865], sigdigits = 6) - #Cavg - @test round.(df[!, :Cavg], sigdigits = 6) == round.([140.9195145 - 203.5356142 - 83.84992893 - 148.4978814 - 145.6560573 - 123.8044835 - 119.9296505 - 85.18002722 - 135.5146515 - 107.8118077], sigdigits = 6) - #Swing tau - @test round.(df[!, :Swingtau], sigdigits = 6) == round.([0.38590184 - 0.27426721 - 0.240584404 - 0.468637128 - 0.057960085 - 0.445494022 - 0.240196282 - 0.628152237 - 0.173820219 - 0.021952739], sigdigits = 6) - #Fluctuation%_Tau - @test round.(df[!, :Fluctau] * 100, sigdigits = 6) == round.([37.71452462 - 27.61899665 - 24.36421266 - 44.8121175 - 6.369060133 - 38.4975876 - 24.74919661 - 62.65259796 - 18.28649133 - 2.500189968], sigdigits = 6) - #Accind - @test round.(df[!, :Accind], sigdigits = 6) == round.([33.3295745 - 8.38726982 - 34.26016611 - 14.94451581 - 16.81301567 - 14.95026394 - 9.427514161 - 12.94903929 - 20.18432092 - 6.976692409], sigdigits = 6) - - #= - #Kel - @test round.(df[!, :Kel], sigdigits = 6) == round.([0.003384744 - 0.014106315 - 0.00329143 - 0.007695344 - 0.006813328 - 0.007692281 - 0.012458956 - 0.00893008 - 0.005645865 - 0.017189737], sigdigits = 6) - =# - - # TAU 0.0 - 100 - - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 100, time = 0.0, tau = 100)) - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :luld, intp = :logt) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - #AUClast - @test round.(df[!, :AUClast], sigdigits = 6) == round.([9573.810559 - 10054.28648 - 5392.457219 - 9297.096334 - 9519.180874 - 6948.985621 - 6988.772632 - 7073.092214 - 8303.358586 - 5486.838889], sigdigits = 6) - #AUCtau - @test round.(df[!, :AUCtau], sigdigits = 6) == round.([12646.63632 - 11996.6718 - 7195.902904 - 11794.11692 - 12274.83395 - 8729.151856 - 8395.400098 - 8930.999936 - 10727.4135 - 6389.420453], sigdigits = 6) - #Ctau - @test round.(df[!, :Ctau], sigdigits = 6) == round.([106.6989373 - 55.60460942 - 61.0383917 - 81.23535402 - 87.01010737 - 58.00027625 - 43.15724396 - 58.8781831 - 80.05171762 - 23.98401112], sigdigits = 6) - @test round.(df[!, :Cavg], sigdigits = 6) == round.([126.4663632 - 119.966718 - 71.95902904 - 117.9411692 - 122.7483395 - 87.29151856 - 83.95400098 - 89.30999936 - 107.274135 - 63.89420453], sigdigits = 6) - #Swing tau - @test round.(df[!, :Swingtau], sigdigits = 6) == round.([0.78885568 - 3.697038658 - 0.725880992 - 1.567133516 - 0.94614172 - 1.666332128 - 2.551060864 - 1.349376165 - 1.090486063 - 4.231902178], sigdigits = 6) - - - # If TAU = Tlast then AUCtau = AUClast - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(dose = 100, time = 0, tau = 72)) - pk = ClinicalTrialUtilities.nca!(pkds, calcm = :luld, intp = :logt) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - @test df[!, :AUCtau] == df[!, :AUClast] - - - - # Utilities - #--------------------------------------------------------------------------- - # Elimination range check - ClinicalTrialUtilities.setkelauto!(pkds, false) - ClinicalTrialUtilities.setelimrange!(pkds, ClinicalTrialUtilities.ElimRange(11, 15)) - pk = ClinicalTrialUtilities.nca!(pkds) - df = DataFrame(pk; unst = true) - sort!(df, :Subject) - - @test round.(df[!, :AUCinf], sigdigits = 6) == round.([655061.1994219155 - 15395.482644599213 - 22227.779049490928 - 53919.75612131339 - 21428.55291240885 - 21778.845190342177 - 15100.132319764292 - 25511.20156654014 - 27923.624963124363 - 7386.834076155846], sigdigits = 6) - - # glucose2 - - pkds = ClinicalTrialUtilities.pkimport(glucose2, [:Subject, :Date]; conc = :glucose, time = :Time) - pk = ClinicalTrialUtilities.nca!(pkds) - df = DataFrame(pk; unst = true) - p = ClinicalTrialUtilities.pkplot(pkds; pagesort = [:Date], typesort = [:Subject]) - @test length(p) == 2 - print(io, pk[1].subject.keldata) - p = ClinicalTrialUtilities.pkplot(pkds; pagesort = [:Date], typesort = [:Subject], elim = true) - p = ClinicalTrialUtilities.pkplot(pkds; pagesort = nothing, typesort = [:Subject]) - @test length(p) == 1 - p = ClinicalTrialUtilities.pkplot(pkds; pagesort = nothing, typesort = nothing, legend = false, xlabel = "L1", ylabel = "L2", xlims = (0,12)) - @test length(p) == 1 - p = ClinicalTrialUtilities.pkplot(pkds; pagesort = [:Date], typesort = nothing) - @test length(p) == 2 - p = ClinicalTrialUtilities.pkplot(pkds[1]) - - # The Theoph dataset: 132 observations from 12 subjects - # see also https://github.com/asancpt/NonCompart-tests/blob/master/docs/validation_0.4.4.pdf - # https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6989226/ - - pkds = ClinicalTrialUtilities.pkimport(theodata, [:Subject]; conc = :conc, time = :Time) - pk = ClinicalTrialUtilities.nca!(pkds; calcm = :luld) - df = DataFrame(pk; unst = true) - @test round.(df[!, :AUClast], sigdigits = 6) == round.([ 147.2347485 - 88.7312755 - 95.8781978 - 102.6336232 - 118.1793538 - 71.6970150 - 87.9692274 - 86.8065635 - 83.9374360 - 135.5760701 - 77.8934723 - 115.2202082], sigdigits = 6) - #= - @test round.(df[!, :AUCinf], sigdigits = 6) == round.([214.9236316 - 97.3779346 - 106.1276685 - 114.2162046 - 136.3047316 - 81.74333453119011 #82.1758833 #6 Rsq& - 100.9876292 - 102.1533003 - 97.52000174742838 #97.5200039 - 167.8600307 - 86.9026173 - 125.8315397], sigdigits = 6) - =# - # NaN PK LimitRule test - nanpkdata.Concentration = ClinicalTrialUtilities.tryfloatparse!(nanpkdata.Concentration) - pkds = ClinicalTrialUtilities.pkimport(nanpkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - rule = ClinicalTrialUtilities.LimitRule(0, 0, NaN, 0, true) - ClinicalTrialUtilities.applyncarule!(pkds, rule) - @test length(pkds[1]) == 5 - @test pkds[3].obs[1] == 0.0 - @test pkds[2].obs[5] == 0.1 - pkds = ClinicalTrialUtilities.pkimport(nanpkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - rule = ClinicalTrialUtilities.LimitRule(0, 0, NaN, -1, false) - ClinicalTrialUtilities.applyncarule!(pkds, rule) - @test pkds[3].obs[6] === NaN - - ClinicalTrialUtilities.setdosetime!(pkds, ClinicalTrialUtilities.DoseTime(100,0), Dict(:Subject => 1, :Formulation => "R")) - @test pkds[3].dosetime.dose == 100 - ClinicalTrialUtilities.setelimrange!(pkds, ClinicalTrialUtilities.ElimRange(3,7), 1) - ClinicalTrialUtilities.setelimrange!(pkds, ClinicalTrialUtilities.ElimRange(3,7), [1,2]) - ClinicalTrialUtilities.setelimrange!(pkds, ClinicalTrialUtilities.ElimRange(3,7), Dict(:Subject => 1, :Formulation => "R")) - pks = ClinicalTrialUtilities.findfirst(Dict(:Subject => 1, :Formulation => "R"), pkds) - #experimental - pke = ClinicalTrialUtilities.nca!(pkds; sort = Dict(:Formulation => "R"), adm = :iv) - pkresulttable = ClinicalTrialUtilities.datatable_st(pke) -end - -println(" ---------------------------------- ") -@testset "#11 Urine PK " begin - upk = ClinicalTrialUtilities.upkimport(upkdata, [:subj]; stime = :st, etime = :et, conc = :conc, vol = :vol) - unca = ClinicalTrialUtilities.nca!(upk[1]) - unca = ClinicalTrialUtilities.nca!(upk) - @test unca[1][:maxrate] ≈ 4.0 - @test unca[1][:mTmax] ≈ 1.5 - @test unca[1][:ar] ≈ 16.0 - @test unca[1][:volume] ≈ 11.0 - p = ClinicalTrialUtilities.pkplot(upk, ylabel = "Excretion") -end -println(" ---------------------------------- ") -@testset "#12 Pharmacodynamics " begin - - - #PD - pdds = ClinicalTrialUtilities.pdimport(pddata; time=:time, resp=:effect, bl = 3.0) - pd = ClinicalTrialUtilities.nca!(pdds) - @test pd[1,:AUCABL] ≈ 7.38571428571429 atol=1E-5 - @test pd[1,:AUCBBL] ≈ 8.73571428571429 atol=1E-5 - ClinicalTrialUtilities.setth!(pdds, 1.5) - pd = ClinicalTrialUtilities.nca!(pdds) - @test pd[1,:AUCATH] ≈ 13.9595238095238 atol=1E-5 - @test pd[1,:AUCBTH] ≈ 1.80952380952381 atol=1E-5 - @test pd[1,:TABL] ≈ 3.48095238095238 atol=1E-5 - @test pd[1,:TBBL] ≈ 5.51904761904762 atol=1E-5 - @test pd[1,:TATH] ≈ 5.76190476190476 atol=1E-5 - @test pd[1,:TBTH] ≈ 3.23809523809524 atol=1E-5 - @test pd[1,:AUCBLNET] ≈ -1.35 atol=1E-5 - @test pd[1,:AUCTHNET] ≈ 12.15 atol=1E-5 - - pdds = ClinicalTrialUtilities.pdimport(pkdata, [:Formulation, :Subject]; time=:Time, resp=:Concentration, bl=0.2, th=0.3) - pd = ClinicalTrialUtilities.nca!(pdds) - @test pd[2, :AUCDBLTH] ≈ 0.3416666666666665 atol=1E-5 - ClinicalTrialUtilities.setbl!(pdds, 0.3) - ClinicalTrialUtilities.setth!(pdds, 0.2) - pd = ClinicalTrialUtilities.nca!(pdds) - @test pd[3, :AUCDBLTH] ≈ 0.3428571428571429 atol=1E-5 - - print(io, pdds) - print(io, pdds[1]) - print(io, pd) - print(io, pd[1]) -end -println(" ---------------------------------- ") -@testset "#13 Sparse PK " begin - pks = ClinicalTrialUtilities.PKSubject(sparse_pk; conc = :Concentration, time = :Time) - auc = ClinicalTrialUtilities.auc_sparse(pks) - @test auc ≈ 1.35 atol=1E-5 -end diff --git a/test/test.jl b/test/test.jl index fba98ec..623917a 100644 --- a/test/test.jl +++ b/test/test.jl @@ -229,52 +229,9 @@ println(" ---------------------------------- ") @test rdf[:Group] == [1, 2, 1, 2, 2, 1, 1, 2, 2, 1] end - -#PK -include("pktest.jl") - - include("sim.jl") -println(" ---------------------------------- ") -@testset " Scenario " begin - #Pharmacokinetics statistics - pkds = ClinicalTrialUtilities.pkimport(pkdata2, [:Subject, :Formulation]; time = :Time, conc = :Concentration) - pk = ClinicalTrialUtilities.nca!(pkds) - df = ClinicalTrialUtilities.DataFrame(pk; unst = true) - ds = ClinicalTrialUtilities.descriptive(df, stats = [:n, :mean, :sd], sort = [:Formulation]) - df = ClinicalTrialUtilities.DataFrame(ds; unst = true) - Base.show(io, pkds) - Base.show(io, pkds[1]) - Base.show(io, pk) - Base.show(io, pk[1]) - Base.show(io, ds) - Base.show(io, ds[1]) - - - @test ds[Dict(:Variable=>:AUClast,:Formulation=>"R")][:mean] ≈ 7431.283916666667 - - pdds = ClinicalTrialUtilities.pdimport(pkdata2, [:Subject, :Formulation]; time = :Time, resp = :Concentration) - pd = ClinicalTrialUtilities.nca!(pdds) - df = ClinicalTrialUtilities.DataFrame(pd; unst = true) - ds = ClinicalTrialUtilities.descriptive(df, stats = :default, sort = [:Formulation]) - - @test ds[Dict(:Variable=>:AUCBLNET,:Formulation=>"T")][:mean] ≈ 8607.09 - @test ds[(:Variable=>:AUCBLNET,:Formulation=>"R")][:sd] ≈ 1919.7954123986283 - -end - -println(" ---------------------------------- ") -@testset " HTML export " begin - - pkds = ClinicalTrialUtilities.pkimport(pkdata2, [:Subject, :Formulation]; time = :Time, conc = :Concentration) - pk = ClinicalTrialUtilities.nca!(pkds) - df = ClinicalTrialUtilities.DataFrame(pk; unst = true) - ClinicalTrialUtilities.htmlexport(df; io = io, sort = :Subject, rspan=:all, title = "PK Parameters", dict = :pk) - @test true -end - println(" ---------------------------------- ") @testset " Errors " begin @@ -306,13 +263,6 @@ println(" ---------------------------------- ") @test_throws ArgumentError ci = ClinicalTrialUtilities.rrpropci(30, 100, 40, 90; alpha=0.05, method=:mnnn) @test_throws ArgumentError ci = ClinicalTrialUtilities.orpropci(30, 100, 40, 90; alpha=0.05, method=:awoolfff) - #DS - pkds = ClinicalTrialUtilities.pkimport(pkdata, [:Subject, :Formulation]; conc = :Concentration, time = :Time) - @test_throws BoundsError pkds[20] - @test_throws ArgumentError pkds[:Subject => 20] - @test_throws ArgumentError pkds[Dict(:Subject => 20)] - - #Type check ClinicalTrialUtilities.besamplen(;theta0=0, theta1=-1, theta2=1, sd=2, logscale = false) end diff --git a/test/testdata.jl b/test/testdata.jl index f25f5e7..5bb2884 100644 --- a/test/testdata.jl +++ b/test/testdata.jl @@ -1,12 +1,4 @@ -#Simple PK DataFrame -pkdata = CSV.File(path*"/csv/pkdata.csv") |> DataFrame -sparse_pk = pkdata[2:7, :] -#Simple PK DataFrame -nanpkdata = CSV.File(path*"/csv/nanpk.csv") |> DataFrame -#Simple PD DataFrame -pddata = CSV.File(path*"/csv/pddata.csv") |> DataFrame -# Multiple subjects PK DataFrame -pkdata2 = CSV.File(path*"/csv/pkdata2.csv") |> DataFrame + #Glucose2 #Pinheiro, J. C. and Bates, D. M. (2000), Mixed-Effects Models in S and S-PLUS, Springer, New York. (Appendix A.10) #Hand, D. and Crowder, M. (1996), Practical Longitudinal Data Analysis, Chapman and Hall, London. @@ -14,8 +6,7 @@ glucose2 = CSV.File(path*"/csv/glucose2.csv") |> DataFrame #theo #https://cran.r-project.org/web/packages/PKNCA/vignettes/Example-theophylline.html theodata = CSV.File(path*"/csv/theo.csv") |> DataFrame -#Simple urine PK DataFrame -upkdata = CSV.File(path*"/csv/upkdata.csv") |> DataFrame + #Simple frequency dataset freqdat = CSV.File(path*"/csv/freqdat.csv") |> DataFrame #Small dataset with negative and zero value diff --git a/test/validation_report.jmd b/test/validation_report.jmd deleted file mode 100644 index a84c75c..0000000 --- a/test/validation_report.jmd +++ /dev/null @@ -1,142 +0,0 @@ ---- -title: ClinicalTrialUtilities:NCA validation report -author: -date: `j import Dates; Dates.Date(Dates.now())` -mainfont: romanuni.ttf -sansfont: NotoSans-Regular.ttf -monofont: NotoSansMono-Regular.ttf -mathfont: texgyredejavu-math.otf ---- - -```julia; echo = false -using Dates -``` - -\pagebreak - -# Introduction and package description - - -Stable documentation: []() - -## Validation purpose and scope - - -## Requirements - - * Julia 1.5.* (or higher) installed - * Julia packages from dependence list installed (see [Project.toml]()) - -## Developer software life cycle - - * Development stage - * Testing procedures development - * Performing testing procedures on local machine - * Push to master branch - * Performing testing procedures with GitHub Actions - * Make pull request to the official registry of general Julia packages (if nessesary) - * Make release (if previous completed) - -### Versions - - * X.Y.Z - patch release (no breaking changes) - * X.Y.# - minor release (may include breaking changes for versions 0.Y.#) - * X.#.# - major release (breaking changes, changes in public API) - * 0.#.# - no stable public API - * ≥1.#.# - stable public API - - -## Build support - -### Tier 1 - - * julia-version: 1.5, 1.6 - * julia-arch: x64 - * os: ubuntu-18.04, macos-10.15, windows-2019 - -\pagebreak - -# Installation - -## System information - - * Julia version: `j Sys.VERSION` - * Current machine: `j Sys.MACHINE` - -## Installation method - - -```julia; eval = false -import Pkg; Pkg.add("") -``` - -## Version check - -The installation process is checking within each testing job via GitHub Actions. -Also GitHub Action [chek](https://github.com/JuliaRegistries/General/blob/master/.github/workflows/automerge.yml) -performed before merging into JuliaRegistries/General repository -(see [Automatic merging of pull requests](https://github.com/JuliaRegistries/General#automatic-merging-of-pull-requests)). - -```julia; echo = false; results = "hidden" -using Pkg, -pkgversion(m::Module) = Pkg.TOML.parsefile(joinpath(dirname(string(first(methods(m.eval)).file)), "..", "Project.toml"))["version"] -ver = pkgversion(Metida) -``` - -Current package version: -```julia; echo = false; results = "tex" -ver -``` - -\pagebreak - -# Operation qualification - -This part of validation based on testing procedures entails running software products under known conditions with defined inputs and -documented outcomes that can be compared to their predefined expectations. All documented public API included in testing procedures and part of -critical internal methods. - -## Coverage - -Code coverage report available on [Codecov.io](). -Test procedures include all public API methods check. - -* Coverage goal: ≥ 90.0% - -## Data - -For operation checks generated data used. For any purpose, this data available in the repository and included in the package. - -Datasets: - - * - -## Testing results - -```julia -Pkg.test("") -``` - -\pagebreak - -# Performance qualification - -Purpose of this testing procedures to demonstrate performance for critical tasks. - - -\pagebreak - -# Glossary - - * Installation qualification (IQ) - Establishing confidence that process equipment and ancillary systems are compliant with appropriate codes and approved design intentions, and that manufacturer's recommendations are suitably considered. - * Operational qualification (OQ) Establishing confidence that process equipment and sub-systems are capable of consistently operating within established limits and tolerances. - * Product performance qualification (PQ) - Establishing confidence through appropriate testing that the finished product produced by a specified process meets all release requirements for functionality and safety. - * Repository - GitHub repository: https://github.com/PharmCat/Metida.jl - * Master branch - main branch on GitHub ([link](https://github.com/PharmCat/Metida.jl/tree/master)). - * Current machine - pc that used for validation report generating. - -# Reference - -* [General Principles of Software Validation; Final Guidance for Industry and FDA Staff](https://www.fda.gov/media/73141/download) -* [Guidance for Industry Process Validation: General Principles and Practices](https://www.fda.gov/files/drugs/published/Process-Validation--General-Principles-and-Practices.pdf) -* [Glossary of Computer System Software Development Terminology](https://www.fda.gov/inspections-compliance-enforcement-and-criminal-investigations/inspection-guides/glossary-computer-system-software-development-terminology-895) diff --git a/test/weave.jl b/test/weave.jl deleted file mode 100644 index b2d3f1c..0000000 --- a/test/weave.jl +++ /dev/null @@ -1,5 +0,0 @@ -using Weave, ClinicalTrialUtilities -weave(joinpath(dirname(pathof(ClinicalTrialUtilities)), "..", "test", "validation_report.jmd"); -doctype = "pandoc2pdf", -out_path = :pwd, -pandoc_options=["--toc", "-V colorlinks=true" , "-V linkcolor=blue", "-V urlcolor=red", "-V toccolor=gray"]) From 378e40cbb180ddfb0ea0bbf0a57d99031a4c198e Mon Sep 17 00:00:00 2001 From: PharmCat <13901158+PharmCat@users.noreply.github.com> Date: Mon, 28 Nov 2022 16:38:03 +0300 Subject: [PATCH 2/2] fix --- src/ClinicalTrialUtilities.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/ClinicalTrialUtilities.jl b/src/ClinicalTrialUtilities.jl index 42bd94e..f39aec2 100644 --- a/src/ClinicalTrialUtilities.jl +++ b/src/ClinicalTrialUtilities.jl @@ -14,7 +14,7 @@ __precompile__(true) module ClinicalTrialUtilities -using Distributions, Random, Roots, QuadGK, RecipesBase +using Distributions, Random, Roots, QuadGK using StatsBase