getting-started.RmdExposure-response and population PK analysis datasets have no CDISC
structure. Every study therefore invents one, and the pharmacometrician
receives something different each time. eradam produces
them to one structure and documents it.
Everything in this vignette runs on the fixtures shipped with the package: a fictional 24-subject study at three dose levels, dosed at 0 h and 24 h, with two analytes and a biomarker endpoint.
library(eradam)
pc <- eradam_example("pc") # concentrations
ex <- eradam_example("ex") # dosing
adsl <- eradam_example("adsl") # subject level
resp <- eradam_example("resp") # response endpoint
head(pc, 3)
#> # A tibble: 3 × 10
#> STUDYID USUBJID PARAMCD PARAM AVISIT NFRLT AFRLT AVAL AVALU LLOQ
#> <chr> <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <chr> <dbl>
#> 1 ERADAM01 ERADAM01-001 DRUGX DRUGX Plasma Concentration (ng… Day 1 0 0 NA ng/mL 0.5
#> 2 ERADAM01 ERADAM01-001 DRUGX DRUGX Plasma Concentration (ng… Day 1 0.5 0.58 230. ng/mL 0.5
#> 3 ERADAM01 ERADAM01-001 DRUGX DRUGX Plasma Concentration (ng… Day 1 1 1.01 348. ng/mL 0.5make_adpc() gives one record per subject, analyte and
timepoint, carrying nominal time (NFRLT), actual time
(AFRLT) and actual time since the most recent dose
(ARRLT) — the ADNCA v1.0 names.
adpc <- make_adpc(pc, adsl, ex, blq = "half")
adpc[
adpc$USUBJID == "ERADAM01-001" & adpc$PARAMCD == "DRUGX",
c("NFRLT", "AFRLT", "NRRLT", "ARRLT", "AVAL", "AVALC", "BLQFL", "DTYPE")
]
#>
#> ── eradam ADPC ─────────────────────────────────────────────────────────────────────────────────────
#> Warning: Unknown or uninitialised column: `USUBJID`.
#> Warning: Unknown or uninitialised column: `PARAMCD`.
#> ℹ 12 record(s) | 0 subject(s) | 0 analyte(s)
#> ℹ BLQ policy: half | 1 BLQ record(s) in output
#> # A tibble: 12 × 8
#> NFRLT AFRLT NRRLT ARRLT AVAL AVALC BLQFL DTYPE
#> <dbl> <dbl> <dbl> <dbl> <dbl> <chr> <chr> <chr>
#> 1 0 0 0 0 0.25 <0.5 Y LOQ
#> 2 0.5 0.58 0.5 0.58 230. 229.713 N NA
#> 3 1 1.01 1 1.01 348. 347.558 N NA
#> 4 2 2.03 2 2.03 378. 378.262 N NA
#> 5 4 4.03 4 4.03 355. 355.441 N NA
#> 6 8 7.99 8 7.99 139. 139.162 N NA
#> 7 12 12.0 12 12.0 83.1 83.057 N NA
#> 8 24 23.9 0 23.9 14.8 14.761 N NA
#> 9 25 25.0 1 1.03 343. 343.309 N NA
#> 10 28 28.0 4 3.98 298. 297.759 N NA
#> 11 36 36.0 12 12.0 94.8 94.821 N NA
#> 12 48 48.0 24 24.0 15.6 15.573 N NANote the eighth record. Its nominal time is 24 h, but the sample was
actually drawn at 23.9 h — before the second dose — so
ARRLT is 23.9, not 0. The ninth record, at 25 h, restarts
at ARRLT 1.0. That distinction is the whole reason both
times are carried.
The policy is an argument, never a hardcoded rule, because no
standard picks one for you. BLQFL records the judgement
whatever the policy, so the count is always recoverable from the
delivered dataset, and DTYPE = "LOQ" — the ADaM
controlled-terminology value — marks any record whose AVAL
was imputed.
policies <- c("half", "zero", "drop", "keep")
data.frame(
policy = policies,
records = vapply(policies, function(p) nrow(make_adpc(pc, adsl, ex, blq = p)), integer(1)),
n_blq = vapply(policies, function(p) sum(make_adpc(pc, adsl, ex, blq = p)$BLQFL == "Y"), integer(1)),
row.names = NULL
)
#> policy records n_blq
#> 1 half 576 48
#> 2 zero 576 48
#> 3 drop 528 0
#> 4 keep 576 48
expo <- derive_exposure_metrics(adpc[adpc$PARAMCD == "DRUGX", ], ex)
head(expo, 3)
#>
#> ── eradam exposure metrics ─────────────────────────────────────────────────────────────────────────
#> ℹ 3 group(s) | AUC method: linear
#> ℹ Analytes: DRUGX
#> ℹ CMAX: median 447.9 (range 378.3 to 491.3)
#> ℹ CTROUGH: median 53.56 (range 14.76 to 84.98)
#> ℹ AUCLST: median 11520 (range 6686 to 11980)
#> # A tibble: 3 × 13
#> USUBJID PARAMCD NOBS TFIRST TLAST CMAX TMAX CMIN CTROUGH AUCLST AUCMTH DOSE AUCLSTD
#> <chr> <chr> <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <chr> <dbl> <dbl>
#> 1 ERADAM01-001 DRUGX 12 0 48.0 378. 2.03 0.25 14.8 6686. linear 50 134.
#> 2 ERADAM01-002 DRUGX 12 0 48.0 491. 3.97 0.25 53.6 11525. linear 50 230.
#> 3 ERADAM01-003 DRUGX 12 0 48.0 448. 28 0.25 85.0 11978. linear 50 240.AUC is the linear trapezoidal rule by default.
method = "lin_up_log_down" uses the log trapezoid on
descending segments, which is the more accurate choice on a log-linear
elimination phase and is what non-compartmental analysis software does.
On a decaying profile it gives a slightly smaller number:
lin <- derive_exposure_metrics(adpc[adpc$PARAMCD == "DRUGX", ], ex)$AUCLST
lud <- derive_exposure_metrics(adpc[adpc$PARAMCD == "DRUGX", ], ex,
method = "lin_up_log_down")$AUCLST
summary(lin - lud)
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 333.6 475.8 963.4 975.7 1239.1 2475.4The arithmetic is worth being explicit about, so here it is on a
curve you can check by hand — concentrations 0, 10, 20, 10
at 0, 1, 2, 3 hours:
curve <- data.frame(
USUBJID = "X-001", PARAMCD = "DRUGX",
AFRLT = c(0, 1, 2, 3), AVAL = c(0, 10, 20, 10)
)
# linear: (0+10)/2 + (10+20)/2 + (20+10)/2 = 5 + 15 + 15 = 35
derive_exposure_metrics(curve)$AUCLST
#> [1] 35
# lin up / log down: the last segment becomes (20-10)/log(20/10) = 14.42695
derive_exposure_metrics(curve, method = "lin_up_log_down")$AUCLST
#> [1] 34.42695derive_exposure_metrics() reads numbers off the observed
points and stops there. Everything that needs a fitted terminal slope —
λz, half-life, AUC extrapolated to infinity — lives in
derive_terminal_phase(), and it returns the fit with the
estimate so you can judge it.
tp <- derive_terminal_phase(adpc[adpc$PARAMCD == "DRUGX", ], method = "lin_up_log_down")
tp[!is.na(tp$LAMZ), c("USUBJID", "LAMZ", "LAMZHL", "LAMZNPT",
"LAMZLL", "LAMZUL", "R2ADJ", "AUCIFO", "AUCPEO", "AUCPEFL")]
#> # A tibble: 14 × 10
#> USUBJID LAMZ LAMZHL LAMZNPT LAMZLL LAMZUL R2ADJ AUCIFO AUCPEO AUCPEFL
#> <chr> <dbl> <dbl> <int> <dbl> <dbl> <dbl> <dbl> <dbl> <chr>
#> 1 ERADAM01-001 0.148 4.69 3 28.0 48.0 0.999 6314. 1.67 N
#> 2 ERADAM01-002 0.101 6.85 3 28.1 48.0 0.996 11635. 5.48 N
#> 3 ERADAM01-005 0.114 6.09 3 28.0 48.0 0.997 6831. 3.96 N
#> 4 ERADAM01-008 0.0898 7.72 3 28.0 48.1 0.999 13492. 6.80 N
#> 5 ERADAM01-009 0.108 6.44 3 28.0 48.0 0.998 19682. 4.37 N
#> 6 ERADAM01-010 0.0727 9.53 3 28.0 48.0 0.987 33576. 10.3 N
#> 7 ERADAM01-011 0.105 6.61 3 28.0 48.0 0.979 21788. 4.61 N
#> 8 ERADAM01-012 0.120 5.78 3 28.0 48.1 0.993 22730. 3.12 N
#> 9 ERADAM01-014 0.132 5.25 3 28.0 48.0 0.997 18210. 2.48 N
#> 10 ERADAM01-015 0.163 4.25 3 28.0 48.0 0.997 16885. 1.12 N
#> 11 ERADAM01-017 0.112 6.18 3 28.0 48.1 0.999 51813. 4.10 N
#> 12 ERADAM01-018 0.105 6.60 3 28.0 48.1 0.997 35228. 4.61 N
#> 13 ERADAM01-019 0.101 6.89 3 28.0 48.0 0.999 42423. 5.23 N
#> 14 ERADAM01-021 0.0861 8.05 3 28.0 48.1 0.999 48696. 7.74 NLAMZNPT, LAMZLL, LAMZUL and
R2ADJ are the window the best-fit rule chose and how well
it fitted. AUCPEO is the percentage of AUCIFO
that was extrapolated rather than observed, and AUCPEFL
flags anything above the conventional 20% threshold — the value is still
returned, the flag says do not lean on it.
The interesting rows are the ones with no answer. This study doses
twice, at 0 h and 24 h, so for some subjects the peak of the
second dose is the profile’s TMAX and only two
points remain after it:
table(is.na(tp$LAMZ))
#>
#> FALSE TRUE
#> 14 10
unique(tp$LAMZNRS[!is.na(tp$LAMZNRS)])[1:2]
#> [1] "Only 2 quantifiable points after TMAX (28); 3 required."
#> [2] "Only 2 quantifiable points after TMAX (27.993); 3 required."Two points is enough to draw a line through and not enough to
believe, so the answer is NA with the reason attached
rather than a number. If you want λz per dosing interval, subset to the
interval before calling — this function treats everything in a
by group as one profile.
On a curve where the answer is known, it is exact.
C = 100 exp(-0.1 t) has λz = 0.1, half-life
log(2) / 0.1 and
AUC(0, ∞) = 100 / 0.1 = 1000:
t <- c(0, 0.5, 1, 2, 4, 8, 12, 18, 24, 36, 48)
mono <- data.frame(USUBJID = "X-001", PARAMCD = "DRUGX",
AFRLT = t, AVAL = 100 * exp(-0.1 * t))
m <- derive_terminal_phase(mono, method = "lin_up_log_down")
c(lambda_z = m$LAMZ, half_life = m$LAMZHL, aucinf = m$AUCIFO)
#> lambda_z half_life aucinf
#> 0.100000 6.931472 1000.000000derive_partial_auc() gives the area between two times
you state. The rule at the window edges is the part implementations
quietly disagree on, so it is written down:
method = "linear" interpolates linearly,
"lin_up_log_down" interpolates logarithmically on a
descending segment, and INTRPFL tells you whether an edge
had to be interpolated at all.
derive_partial_auc(adpc[adpc$PARAMCD == "DRUGX", ], 0, 24)[1:3, ]
#>
#> ── eradam partial AUC ──────────────────────────────────────────────────────────────────────────────
#> ℹ 3 profile(s) | window [0, 24] | AUC method: linear
#> ℹ AUCINT derived for 3/3 profile(s)
#> ℹ AUCINT: median 5514 (range 3304 to 5515)
#> ℹ 3 profile(s) needed an interpolated window edge
#> # A tibble: 3 × 8
#> USUBJID PARAMCD AUCINTLL AUCINTUL AUCINT NOBSINT INTRPFL AUCINTRS
#> <chr> <chr> <dbl> <dbl> <dbl> <int> <chr> <chr>
#> 1 ERADAM01-001 DRUGX 0 24 3304. 9 Y NA
#> 2 ERADAM01-002 DRUGX 0 24 5514. 9 Y NA
#> 3 ERADAM01-003 DRUGX 0 24 5515. 9 Y NAOver the whole observed window it reproduces AUCLST
exactly, which is the check that the edge handling is not quietly adding
or losing area:
one <- adpc[adpc$PARAMCD == "DRUGX" & adpc$USUBJID == "ERADAM01-001", ]
derive_partial_auc(one, min(one$AFRLT), max(one$AFRLT))$AUCINT
#> [1] 6686.057
derive_exposure_metrics(one)$AUCLST
#> [1] 6686.057A window that runs outside the observed range is refused, not guessed:
derive_partial_auc(one, 0, 96)$AUCINTRS
#> [1] "Window [0, 96] is not inside the observed range [0, 48.006]; no extrapolation is done."Back-extrapolating to C = 0 at time zero and
forward-extrapolating the tail with λz are both real conventions used
elsewhere. Picking one silently is how two AUC0-24 columns
end up disagreeing, so this function returns NA and the
reason. For the forward tail use derive_terminal_phase(),
which extrapolates on purpose and reports how much.
make_poppk() interleaves dosing and observation records
in time order per subject and names them the way an estimation tool
expects.
pk <- make_poppk(adpc, ex, adsl, covariates = c("AGE", "SEX", "WTBL", "CRCL"))
pk[pk$ID == 1, c("ID", "TIME", "TAD", "AMT", "DV", "MDV", "EVID", "CMT", "BLQ")][1:6, ]
#>
#> ── eradam population PK dataset ────────────────────────────────────────────────────────────────────
#> ℹ 6 record(s): 1 dosing, 5 observation
#> ℹ 1 subject(s) | TIME 0 to 1.013 h | 0 observation(s) with MDV=1
#> ℹ Compartments: 1 = dosing, 2 = DRUGX, 3 = DRUGXM1
#> # A tibble: 6 × 9
#> ID TIME TAD AMT DV MDV EVID CMT BLQ
#> <int> <dbl> <dbl> <dbl> <dbl> <int> <int> <int> <int>
#> 1 1 0 0 25 NA 1 1 1 0
#> 2 1 0 0 NA 0.25 0 0 2 1
#> 3 1 0 0 NA 0.1 0 0 3 1
#> 4 1 0.497 0.497 NA 79.6 0 0 3 0
#> 5 1 0.58 0.58 NA 230. 0 0 2 0
#> 6 1 1.01 1.01 NA 348. 0 0 2 0Records sort by ID, then TIME, then
EVID descending. That last key is what puts the dose at
time 0 before the two samples at time 0: a concentration
recorded at a given time reflects the dose already given. A sample
intended as a pre-dose trough must therefore carry a TIME
strictly less than the dose time — which is exactly what actual rather
than nominal sampling times give you.
Feeding make_poppk() the ADPC rather than the raw
PC is the recommended route: the BLQ policy has already
been applied and recorded, so imputed values arrive as DV
with MDV = 0, and BLQ = 1 still marks where
they came from.
table(BLQ = pk$BLQ, MDV = pk$MDV, EVID = pk$EVID)
#> , , EVID = 0
#>
#> MDV
#> BLQ 0 1
#> 0 528 0
#> 1 48 0
#>
#> , , EVID = 1
#>
#> MDV
#> BLQ 0 1
#> 0 0 48
#> 1 0 0
ader <- make_ader(expo, resp)
ader[1:3, c("USUBJID", "EXPANLT", "EXPPARM", "EXPVAL", "EXPQN", "EXPQ", "PARAMCD", "AVAL")]
#>
#> ── eradam exposure-response dataset ────────────────────────────────────────────────────────────────
#> ℹ 3 subject(s) | exposure: AUCLST | response: AVAL
#> ℹ Q1 (6680-11900): n = 2, mean response = -20.22
#> ℹ Q2 (11900-22400): n = 1, mean response = -26.38
#> # A tibble: 3 × 8
#> USUBJID EXPANLT EXPPARM EXPVAL EXPQN EXPQ PARAMCD AVAL
#> <chr> <chr> <chr> <dbl> <int> <ord> <chr> <dbl>
#> 1 ERADAM01-001 DRUGX AUCLST 6686. 1 Q1 (6680-11900) PCHGBIO -21.3
#> 2 ERADAM01-002 DRUGX AUCLST 11525. 1 Q1 (6680-11900) PCHGBIO -19.1
#> 3 ERADAM01-003 DRUGX AUCLST 11978. 2 Q2 (11900-22400) PCHGBIO -26.4
tapply(ader$AVAL, ader$EXPQ, mean)
#> Q1 (6680-11900) Q2 (11900-22400) Q3 (22400-32200) Q4 (32200-53000)
#> -17.89000 -27.94167 -32.53667 -41.66333EXPQ is an ordered factor labelled with the range of
each group, so a figure legend or a table row header is self-describing
without extra work.
Every variable, its meaning and its source, returned as data — ready for a define document or a data specification.
eradam_structure("ader")
#> # A tibble: 10 × 4
#> dataset variable label source
#> <chr> <chr> <chr> <chr>
#> 1 ader USUBJID Unique subject identifier (or the `by` variables supplied) `exposure` / `r…
#> 2 ader EXPANLT Analyte the exposure metric was derived for `exposure` colu…
#> 3 ader EXPPARM Name of the exposure metric used Argument `expos…
#> 4 ader EXPVAL Value of the exposure metric `exposure`
#> 5 ader EXPQN Exposure quantile group number (1 = lowest exposure) Derived
#> 6 ader EXPQ Exposure quantile group label, including the range of the group Derived
#> 7 ader PARAMCD Response parameter code `response` (kep…
#> 8 ader PARAM Response parameter description `response` (kep…
#> 9 ader AVAL Response value `response` (arg…
#> 10 ader <other> Any remaining columns of `exposure` and `response` `exposure` / `r…The same table is rendered into each function’s help page by roxygen, from this one object, so the code and the documentation cannot drift apart.
Where a real convention exists this package follows it:
NFRLT/AFRLT/
NRRLT/ARRLT and the negative pre-dose
ARRLT from the CDISC ADaM Implementation Guide for
Non-compartmental Analysis Input Data (ADNCA v1.0);
DTYPE = "LOQ", AVAL, AVALC,
PARAMCD, ANL01FL from ADaM; CMAX,
TMAX, CMIN, CTROUGH,
AUCLST, LAMZ, LAMZHL,
LAMZLL, LAMZUL, LAMZNPT,
R2, R2ADJ, AUCIFO,
AUCIFP, AUCPEO, AUCPEP,
CLST, TLST and AUCINT from the
CDISC PPTESTCD terminology; ID,
TIME, AMT, DV, MDV,
EVID, CMT from NONMEM. The best-fit λz rule,
the exclusion of TMAX, the three-point minimum,
t½ = ln(2)/λz and the 20% extrapolation threshold are the
accepted pharmacometric conventions.
Where none exists, the documentation says so rather than implying a
standard: the wide one-row-per-subject layouts, the CTROUGH
definition, the BLQFL and LLOQ names, the
compartment numbering, every EXP* variable, the
CLSTP / AUCPEFL / LAMZNRS /
AUCINTLL / AUCINTUL / NOBSINT /
INTRPFL / AUCINTRS names, returning
NA with a written reason instead of a doubtful fit, and
refusing a partial-AUC window outside the observed range are this
package’s own conventions. They are written down, which is the point —
the alternative is not a standard, it is a different undocumented
structure every study.
No steady-state parameters (AUCTAU, CAVG,
CLSS, accumulation, fluctuation), no clearance or volume
(CL, CLF, VZ, VZF,
MRT), no moment curves (AUMC,
VSS), no manual override of the λz window, and no
sparse-sampling variance estimation. The README’s Scope and
limits section carries the full list. This is an open
implementation of published rules, not validated NCA software.
Nothing above depends on the fixtures. Any PC-shaped,
EX-shaped and ADSL-shaped data frame works,
and the column names that are not fixed are arguments:
derive_exposure_metrics(my_pc, my_ex, conc = "PCSTRESN", time = "TIMEHR")
derive_terminal_phase(my_pc, conc = "PCSTRESN", time = "TIMEHR")
derive_partial_auc(my_pc, 0, 24, conc = "PCSTRESN", time = "TIMEHR")
make_adpc(my_pc, my_adsl, my_ex, lloq = 0.1, keep = c("AGE", "WTBL"))
make_poppk(my_adpc, my_ex, my_adsl, covariates = c("AGE", "WTBL", "CRCL"))Every return value is a plain tibble, so it hands straight on to
admiral, rtables or xportr.