| reversmean {upndown} | R Documentation |
Reversal-anchored averaging estimators for Up-and-Down
Description
Dose-averaging target estimation for Up-and-Down experiments, historically the most popular approach, but not recommended as primary nowadays. Provided for completeness.
Usage
reversmean(
x,
y,
rstart = 3,
all = TRUE,
before = FALSE,
maxExclude = 1/3,
full = FALSE
)
reversals(y)
Arguments
x |
numeric vector: sequence of administered doses, treatments, stimuli, etc. |
y |
numeric vector: sequence of observed responses. Must be same length as |
rstart |
the reversal point from which the averaging begins. Default 3, considered a good compromise between performance and robustness. See Details. |
all |
logical: from the starting point onwards, should all values of |
before |
logical: whether to start the averaging one step earlier than the starting reversal point. Default |
maxExclude |
a fraction in |
full |
logical: should more detailed information be returned, or only the estimate? (default |
Details
Up-and-Down designs (UDDs) allocate doses in a random walk centered nearly symmetrically around a balance point. Therefore, a modified average of allocated doses could be a plausible estimate of the balance point's location.
During UDDs' first generation, a variety of dose-averaging estimators was developed, with the one proposed by Wetherill et al. (1966) eventually becoming the most popular. This estimator uses only doses observed at reversal points: points with a negative response following a positive one, or vice versa. More recent research (Kershaw 1985, 1987; Oron et al. 2022, supplement) strongly indicates that in fact it is better to use all doses starting from some cut-point, rather than skip most of them and choose only reversals.
The reversals() utility identifies reversal points, whereas reversmean() produces a dose-averaging estimate whose starting cut-point is determined by a reversal. User can choose whether to use all doses from that cut-point onwards, or only the reversals as in the older approaches. A few additional options make the estimate even more flexible.
More broadly, dose-averaging despite some advantages is not very robust, and also lacks an interval estimate with reliable coverage. Therefore, reversmean() provides neither a confidence interval nor a standard error.
For UDD target estimation we recommend using centered isotonic regression, available via udest, an up-and-down adapted wrapper to cir::quickInverse(). See Oron et al. 2022 (both article and supplement) for further information, as well as the cir package vignette.
Value
For reversals(), the indices of reversal points. For reversmean(), if full=FALSE returns the point estimate and otherwise returns a data frame with the estimate as well, as the index of the cutoff point used to start the averaging.
Author(s)
Assaf P. Oron <assaf.oron.at.gmail.com>
References
Kershaw CD: A comparison of the estimators of the ED50 in up-and-down experiments. J Stat Comput Simul 1987; 27:175–84.
Oron AP, Souter MJ, Flournoy N. Understanding Research Methods: Up-and-down Designs for Dose-finding. Anesthesiology 2022; 137:137–50. See in particular the open-access Supplement.
Wetherill GB, Chen H, Vasudeva RB: Sequential estimation of quantal response curves: A new method of estimation. Biometrika 1966; 53:439–54
See Also
-
udest, the recommended estimation method for up-and-down targets. -
adaptmean, an unpublished but arguably better approach to dose-averaging (this is not the recommended method though; that would beudestreferenced above).
Examples
#' **An up-and-down experiment that has generated some controversy**
#'
#' Van Elstraete, AC et al. The Median Effective Dose of Preemptive Gabapentin
#' on Postoperative Morphine Consumption After Posterior Lumbar Spinal Fusion.
#' *Anesthesia & Analgesia* 2008, 106: 305-308.
# It was a classical median-finding up-and-down study.
doses = c(4:7, 6:13, 12:19, 18:21, 20, 19:23, 22, 21:23, 22:19, 20:23,
22:24, 23, 22, 23, 22:25, 24:22, rep(23:24,2), 23, 22)
# With U&D, responses (except the last one) can be read off the doses:
responses = c( (1 - sign(diff(doses)))/2, 0 )
### Let us plot the dose-allocation time series.
# Saving current settings as now required by the CRAN powers-that-be :0
op <- par(no.readonly = TRUE)
par(mar=c(4,4,4,1), mgp=c(2.5,0.8,0), cex.axis = 0.7, las = 1)
udplot(doses, responses, main='Van Elstraete et al. 2008 Study',
xtitle = "Patient Number", ytitle = 'Gabapentin (mg/kg)')
#' Overlay the ED50 reported in the article (21.7 mg/kg):
abline(h = 21.7)
#' The authors cite a little-known 1991 article by Dixon as the method source.
#' However, in their author rejoinder they claim to have used the Dixon-Mood (1948) estimate.
# Our package does include the Dixon-Mood point estimate.
# (w/o the CIs, because we do not endorse this estimation approach)
# Does it reproduce the article estimate?
dixonmood(doses, responses)
# Not at all! Let us overlay this one in red
abline(h = dixonmood(doses, responses), col=2)
# We have found that many articles claiming to use Dixon-Mood (or Dixon-Massey) actually
# Do something else. For example, in this article they report that
# "it is necessary to reject sequences with three to six identical results".
# Nothing like this appears in the original Dixon-Mood article, where the estimation method
# involves identifying the less-common response (either 0 or 1), and using only x values
# associated with these responses; obviating the need to exclude specific sequences.
#
# More generally, these historical estimates have long passed their expiry dates.
# Their foundation is not nearly as solid as, e.g., linear regression,
# and it's time to stop using them.
# That said, our package does offer two more types of dose-averaging estimates.
# Both are able to take advantage of the "n+1" dose-allocation, which is determined by
# the last dose and response:
n = length(doses)
dosePlus1 = doses[n] + ifelse(responses[n]==0, 1, -1)
reversmean(c(doses, dosePlus1), responses)
# Interestingly, in this particular case the answer is very similar to the Dixon-Mood estimate.
# The `reversmean()` default averages all doses from the 3rd reversal point onwards.
# By the way, at what point did the third reversal happen?
# It'll be the 3rd number in this vector:
reversals(responses)
# Far more commonly in literature, particularly in sensory studies,
# one encounters the 1960s-era approach (led by Wetherill) of taking *only doses
# at reversal points, usually starting from the first one. `reversmean()` can do that too:
wetherill = reversmean(c(doses, dosePlus1), responses, all = FALSE, rstart = 1)
wetherill
# This one gives an even lower result than the previous ones.
abline(h = wetherill, col = 3)
# There's another approach to dose-averaging, although it is not in use anywhere that we know of.
# It does not require the y values at all. The underlying assumption is that the dose
# sequence has done enough meandering around the true balance point, to provide information
# about where (approximately) the starting-dose effect is neutralized.
adaptmean(c(doses, dosePlus1))
# Again a bit curiously, this relatively recent approach gives a result similar to what
# the authors reported (but not similar to the original Dixon-Mood).
# This is not too surprising, since here `adaptmean()` excludes the first one-third of doses,
# which is approximately what happened if indeed the authors excluded all those long dose-increase
# sequences at the start.
# All this shows how dicey dose-averaging, at face value a simple and effective method, can become.
# The sample size here is rather large for up-and-down studies, and yet because of the unlucky
# choice of starting point (which in many studies, due to safety concerns cannot be evaded)
# there is really no good option of which observations to exclude.
# This is one reason why we strongly recommend using Centered Isotonic Regression as default:
defest = udest(doses, responses, target = 0.5)
abline(h = defest$point, col = 'purple')
# For this dataset, it is the highest of all the estimates.
legend('bottomright', col = c(1:3, 'purple'),
legend = c("Article's estimate", 'Dixon-Mood', 'Reversals (Wetherill)', 'Standard (CIR)'),
lty = 1, bty='n', cex = 0.8)
par(op) # Back to business as usual ;)