Dear all,
sorry, my last message contained a typo. Correct:
library(PowerTOST)
x??? <- 0.90
y??? <- 0.35
res? <- as.numeric(sampleN.TOST(theta0 = x, CV = y, design =
"2x2x4",
??????????????????????????????? method = "central", details = FALSE,
??????????????????????????????? print = FALSE)[7:8])
mesh <- 28
ys?? <- unique(sort(c(y, seq(y*.8, y*1.2, length.out = mesh))))
xs?? <- unique(sort(c(x, seq(x*0.95, 1, length.out = mesh))))
z??? <- matrix(nrow = length(ys), ncol = length(xs),
?????????????? dimnames = list(paste0("y.", signif(ys, 5)),
?????????????????????????????? paste0("x.", signif(xs, 5))))
for (i in seq_along(ys)) {
? for (j in seq_along(xs)) {
??? z[i, j] <- suppressMessages(
???????????????? power.TOST(CV = ys[i], theta0 = xs[j], design =
"2x2x4",
??????????????????????????? method = "central", n = res[1]))
? }
}
z??? <- z[nrow(z):1, ncol(z):1]????????? # reverse rows & columns
z[paste0("y.", y), paste0("x.", x)] == res[2] # should be
TRUE
contour(xs, ys, z, nlevels = 20, las = 1, labcex = 1,
??????? xlab = "x", ylab = "y", main = paste("n
=", res[1]))
abline(h = y, v = x, lty = 2)
points(x, y, col = "red", cex = 1.5)
text(x, y, labels = signif(z[paste0("y.", y), paste0("x.",
x)], 6),
???? col = "red", adj = c(-0.1, 1.5))
Helmut
--
Ing. Helmut Sch?tz
BEBAC?? Consultancy Services for
Bioequivalence and Bioavailability Studies
Neubaugasse 36/11
1070 Vienna, Austria
E helmut.schuetz at bebac.at
W https://bebac.at/
F https://forum.bebac.at/
You're doing a lot of manipulation of the z matrix; I haven't followed all of it, but that's where I'd look for problems. Generally if you keep your calculation of the z matrix very simple you are better off. For example, once you have xs and ys in the form you want, calculate z as z <- outer(x,y, function(x,y) ...) and then plot it using contour(x, y, z, ...) The only tricky issue here is that the function needs to be vectorized, so if your power.TOST function doesn't accept vectors for CV and theta0, you'll need to put it in a wrapper that does (perhaps using the Vectorize function). Duncan Murdoch On 28/09/2020 3:08 p.m., Helmut Sch?tz wrote:> Dear all, > > sorry, my last message contained a typo. Correct: > > library(PowerTOST) > x??? <- 0.90 > y??? <- 0.35 > res? <- as.numeric(sampleN.TOST(theta0 = x, CV = y, design = "2x2x4", > ??????????????????????????????? method = "central", details = FALSE, > ??????????????????????????????? print = FALSE)[7:8]) > mesh <- 28 > ys?? <- unique(sort(c(y, seq(y*.8, y*1.2, length.out = mesh)))) > xs?? <- unique(sort(c(x, seq(x*0.95, 1, length.out = mesh)))) > z??? <- matrix(nrow = length(ys), ncol = length(xs), > ?????????????? dimnames = list(paste0("y.", signif(ys, 5)), > ?????????????????????????????? paste0("x.", signif(xs, 5)))) > for (i in seq_along(ys)) { > ? for (j in seq_along(xs)) { > ??? z[i, j] <- suppressMessages( > ???????????????? power.TOST(CV = ys[i], theta0 = xs[j], design = "2x2x4", > ??????????????????????????? method = "central", n = res[1])) > ? } > } > z??? <- z[nrow(z):1, ncol(z):1]????????? # reverse rows & columns > z[paste0("y.", y), paste0("x.", x)] == res[2] # should be TRUE > contour(xs, ys, z, nlevels = 20, las = 1, labcex = 1, > ??????? xlab = "x", ylab = "y", main = paste("n =", res[1])) > abline(h = y, v = x, lty = 2) > points(x, y, col = "red", cex = 1.5) > text(x, y, labels = signif(z[paste0("y.", y), paste0("x.", x)], 6), > ???? col = "red", adj = c(-0.1, 1.5)) > > > Helmut >
Dear Duncan, Duncan Murdoch wrote on 2020-09-28 21:47:> You're doing a lot of manipulation of the z matrix; I haven't followed > all of it, but that's where I'd look for problems.? Generally if you > keep your calculation of the z matrix very simple you are better off. > For example, once you have xs and ys in the form you want, calculate z as > > z <- outer(x,y, function(x,y) ...) > > and then plot it using contour(x, y, z, ...) > > The only tricky issue here is that the function needs to be > vectorized, so if your power.TOST function doesn't accept vectors for > CV and theta0, > you'll need to put it in a wrapper that does (perhaps using the > Vectorize function).Here I'm lost. power.TOST(theta0, CV, ...) vectorizes properly for theta0 _or_ CV but no _both_. Hence library(PowerTOST) power.TOST(theta0 = c(0.9, 0.95, 1), CV = 0.25, n = 28) and power.TOST(theta0 = 0.95, CV = c(0.2, 0.25, 0.3), n = 28) work, whereas power.TOST(theta0 = c(0.9, 0.95, 1), 0.95, CV = c(0.2, 0.25, 0.3), n = 28) not. Of note, we will throw an error in the next release if both arguments are vectors. I tried f <- function(x, y) { ? power.TOST(theta0 = x, CV = y, n = 28) } x <- unique(sort(c(0.95, seq(0.95*0.95, 1, length.out = 28)))) y <- unique(sort(c(0.25, seq(0.25*0.8, 0.25*1.2, length.out = 28)))) Vectorize(f, c("x, y"), SIMPLIFY = "array") which is obviously not correct. Helmut -- Ing. Helmut Sch?tz BEBAC?? Consultancy Services for Bioequivalence and Bioavailability Studies Neubaugasse 36/11 1070 Vienna, Austria E helmut.schuetz at bebac.at W https://bebac.at/ F https://forum.bebac.at/