The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.
library(graphicalMCP)
library(lrstat)
library(gMCP)
There are two R packages that cover graphical multiple comparison
procedures (MCPs): gMCP
(Rohmeyer
and Klinglmueller 2024) and lrstat
(Lu 2023). The development of
graphicalMCP
benefited from these two packages. Here we
provide some comparisons between graphicalMCP
and other
packages with respect to key functions.
graphicalMCP::graph_generate_weights()
gMCP::generateWeights()
graphicalMCP::graph_test_shortcut()
gMCP::gMCP()
graphicalMCP::graph_calculate_power()
gMCP::calcPower()
graphicalMCP::graph_test_closure()
gMCP::gMCP()
graphicalMCP::graph_calculate_power()
gMCP::calcPower()
graphicalMCP::graph_test_closure()
gMCP::gMCP()
graphicalMCP::graph_calculate_power()
gMCP::calcPower()
A random graph for five hypotheses will be generated and used for the
comparison. Weighting strategies from the following two functions will
be compared: graphicalMCP::graph_generate_weights()
and
gMCP::generateWeights()
. This process is repeated 1000
times. Weighting strategies are matched for all 1000 cases.
set.seed(1234)
<- NULL
identical for (i in 1:1000) {
<- random_graph(5)
graph <- graphicalMCP::graph_generate_weights(graph)
graphicalmcp_weights dimnames(graphicalmcp_weights) <- list(NULL, NULL)
<-
gmcp_weights ::generateWeights(graph$transitions, graph$hypotheses)
gMCP<- gmcp_weights[nrow(gmcp_weights):1, ] # Reorder rows
gmcp_weights <- c(
identical
identical,all.equal(gmcp_weights, graphicalmcp_weights, tolerance = 1e-7)
)
}all(identical)
#> [1] TRUE
A random graph for five hypotheses will be generated and used for the
comparison. A set of p-values is randomly generated to be used for the
graphical MCP. Adjusted p-values are calculated and compared using the
following functions: graphicalMCP::graph_test_shortcut()
and gMCP::gMCP()
. This process is repeated 10000 times.
Adjusted p-values are matched for all 10000 cases.
set.seed(1234)
<- 0.025
alpha <- NULL
identical for (i in 1:10000) {
<- random_graph(5)
graph <- runif(5, 0, alpha)
p <-
graphicalmcp_test_shortcut graph_test_shortcut(graph, p, alpha = alpha)$outputs$adjusted_p
<-
gmcp_test_shortcut gMCP(as_graphMCP(graph), p, alpha = alpha)@adjPValues
<- c(
identical
identical,all.equal(graphicalmcp_test_shortcut, gmcp_test_shortcut, tolerance = 1e-7)
)
}all(identical)
#> [1] TRUE
A random graph for five hypotheses will be generated and used for the
comparison. A set of marginal power (without multiplicity adjustment) is
randomly generated. Local power (with multiplicity adjustment) is
calculated and compared using the following functions:
graphicalMCP::graph_calculate_power()
and
gMCP::calcPower()
. Since different simulation methods are
used, results are slightly different. The maximum absolute difference in
local power is 0.0051 (0.51%) among 1000 cases, which is relatively
small.
set.seed(1234)
<- 0.025
alpha <- matrix(.5, 5, 5)
sim_corr diag(sim_corr) <- 1
<- NULL
graphicalmcp_power <- NULL
gmcp_power for (i in 1:1000) {
<- random_graph(5)
graph <- runif(5, 0.5, 0.9)
marginal_power <-
noncentrality_parameter qnorm(1 - 0.025, lower.tail = TRUE) -
qnorm(1 - marginal_power, lower.tail = TRUE)
set.seed(1234 + i - 1)
<- rbind(
graphicalmcp_power
graphicalmcp_power,graph_calculate_power(
graph,alpha = alpha,
power_marginal = marginal_power,
sim_corr = sim_corr,
sim_n = 2^17
$power$power_local
)
)
set.seed(1234 + i - 1)
<- rbind(
gmcp_power
gmcp_power,calcPower(
$hypotheses,
graphalpha = alpha,
$transitions,
graphmean = noncentrality_parameter,
corr.sim = sim_corr,
n.sim = 2^17
$LocalPower
)
)
}
<- data.frame(
diff rbind(graphicalmcp_power, gmcp_power),
procedure = rep(c("graphicalMCP", "gMCP"), each = nrow(graphicalmcp_power))
)
write.csv(
diff,::here("vignettes/cache/comparisons_power_shortcut.csv"),
hererow.names = FALSE
)
<- read.csv(here::here("vignettes/cache/comparisons_power_shortcut.csv"))
diff <- subset(diff, procedure == "graphicalMCP")
graphicalmcp_power <- subset(diff, procedure == "gMCP")
gmcp_power round(
max(
abs(
-ncol(graphicalmcp_power)] -
graphicalmcp_power[, -ncol(gmcp_power)]
gmcp_power[,
)
),4
# Maximum difference in local power among 1000 cases )
A successive graph with two primary and two secondary hypotheses will
be generated and used for the comparison. A set of p-values is randomly
generated to be used for the graphical MCP. Adjusted p-values are
calculated and compared using the following functions:
graphicalMCP::graph_test_closure()
and
gMCP::gMCP()
. Parametric tests are used for two primary
hypotheses. This process is repeated 10000 times. Adjusted p-values are
matched for all 10000 cases.
<- c(0.5, 0.5, 0, 0)
hypotheses <- rbind(
transitions c(0, 0.5, 0.5, 0),
c(0.5, 0, 0, 0.5),
c(0, 1, 0, 0),
c(1, 0, 0, 0)
)<- graph_create(hypotheses, transitions)
graph
set.seed(1234)
<- 0.025
alpha <- NULL
identical <- rbind(
test_corr c(1, 0.5, NA, NA),
c(0.5, 1, NA, NA),
c(NA, NA, 1, NA),
c(NA, NA, NA, 1)
)for (i in 1:10000) {
<- runif(4, 0, alpha)
p <- graph_test_closure(
graphicalmcp_test_parametric
graph,
p,alpha = alpha,
test_groups = list(1:2, 3:4),
test_types = c("parametric", "bonferroni"),
test_corr = list(test_corr[1:2, 1:2], NA)
$outputs$adjusted_p
)<- gMCP(
gmcp_test_parametric as_graphMCP(graph),
p,alpha = 0.025,
correlation = test_corr
@adjPValues
)<- c(
identical
identical,all.equal(graphicalmcp_test_parametric, gmcp_test_parametric, tolerance = 1e-7)
)
}all(identical)
#> [1] TRUE
A successive graph with two primary and two secondary hypotheses will
be generated and used for the comparison. A set of marginal power
(without multiplicity adjustment) is randomly generated. Local power
(with multiplicity adjustment) is calculated and compared using the
following functions: graphicalMCP::graph_calculate_power()
and gMCP::calcPower()
. Parametric tests are used for two
primary hypotheses. This process is repeated 100 times. Since different
simulation methods are used, results are slightly different. The maximum
absolute difference in local power is 0.0142 (1.42%) among 100 cases,
which is small.
<- c(0.5, 0.5, 0, 0)
hypotheses <- rbind(
transitions c(0, 0.5, 0.5, 0),
c(0.5, 0, 0, 0.5),
c(0, 1, 0, 0),
c(1, 0, 0, 0)
)<- graph_create(hypotheses, transitions)
graph <- rbind(
test_corr c(1, 0.5, NA, NA),
c(0.5, 1, NA, NA),
c(NA, NA, 1, NA),
c(NA, NA, NA, 1)
)<- matrix(0.5, 4, 4)
sim_corr diag(sim_corr) <- 1
set.seed(1234)
<- 0.025
alpha <- NULL
graphicalmcp_power_parametric <- NULL
gmcp_power_parametric for (i in 1:100) {
<- runif(4, 0.5, 0.9)
marginal_power <-
noncentrality_parameter qnorm(1 - 0.025, lower.tail = TRUE) -
qnorm(1 - marginal_power, lower.tail = TRUE)
set.seed(1234 + i - 1)
<- rbind(
graphicalmcp_power_parametric
graphicalmcp_power_parametric,graph_calculate_power(
graph,alpha = alpha,
test_groups = list(1:2, 3:4),
test_types = c("parametric", "bonferroni"),
test_corr = list(test_corr[1:2, 1:2], NA),
power_marginal = marginal_power,
sim_corr = sim_corr,
sim_n = 2^14
$power$power_local
)
)
set.seed(1234 + i - 1)
<- rbind(
gmcp_power_parametric
gmcp_power_parametric,calcPower(
$hypotheses,
graphalpha = alpha,
$transitions,
graphcorr.test = test_corr,
mean = noncentrality_parameter,
corr.sim = sim_corr,
n.sim = 2^14
$LocalPower
)
)
}
<- data.frame(
diff rbind(graphicalmcp_power_parametric, gmcp_power_parametric),
procedure = rep(c("graphicalMCP", "gMCP"), each = nrow(gmcp_power_parametric))
)
write.csv(
diff,::here("vignettes/cache/comparisons_power_parametric.csv"),
hererow.names = FALSE
)
<- read.csv(here::here("vignettes/cache/comparisons_power_parametric.csv"))
diff <- subset(diff, procedure == "graphicalMCP")
graphicalmcp_power <- subset(diff, procedure == "gMCP")
gmcp_power round(
max(
abs(
-ncol(graphicalmcp_power_parametric)] -
graphicalmcp_power_parametric[, -ncol(gmcp_power)]
gmcp_power_parametric[,
)
),4
# Maximum difference in local power among 100 cases )
A successive graph with two primary and two secondary hypotheses will
be generated and used for the comparison. A set of p-values is randomly
generated to be used for the graphical MCP. Adjusted p-values are
calculated and compared using the following functions:
graphicalMCP::graph_test_closure()
and
lrstat::fadjpsim()
. Simes tests are used for two primary
hypotheses. This process is repeated 10000 times. Adjusted p-values are
matched for all 10000 cases.
<- c(0.5, 0.5, 0, 0)
hypotheses <- 0.0001
eps <- rbind(
transitions c(0, 1 - eps, eps, 0),
c(1 - eps, 0, 0, eps),
c(0, 1, 0, 0),
c(1, 0, 0, 0)
)<- graph_create(hypotheses, transitions)
graph
set.seed(1234)
<- 0.025
alpha <- NULL
identical <- rbind(
family c(1, 1, 0, 0),
c(0, 0, 1, 0),
c(0, 0, 0, 1)
)for (i in 1:10000) {
<- runif(4, 0, alpha)
p <- graph_test_closure(
graphicalmcp_test_simes
graph,
p,alpha = alpha,
test_groups = list(1:2, 3:4),
test_types = c("simes", "bonferroni")
$outputs$adjusted_p
)names(graphicalmcp_test_simes) <- NULL
<-
lrstat_test_simes fadjpsim(
fwgtmat(graph$hypotheses, graph$transitions),
p,
family
)
<- c(
identical
identical,all.equal(graphicalmcp_test_simes, lrstat_test_simes, tolerance = 1e-7)
)
}all(identical)
#> [1] TRUE
Power simulations are not available in lrstat
. Thus a
comparison could be done to compare
graphicalMCP::graph_calculate_power()
and a manual
repetition of lrstat::fadjpsim()
for many sets of marginal
power. This process is the same as the above comparison of adjusted
p-values, and thus omitted.
These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.