The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
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[@]dogado.de.
The regmedint function itself does not contain a
bootstrap standard error option, which may be perferred in some
settings. However, it is relatively easy to implmement in R using the
regmedint() function and the corresponding
coef() point estimate extraction method.
set.seed(138087069)
library(regmedint)
library(tidyverse)
## Prepare dataset
data(vv2015)
## Main fit
regmedint_obj <- regmedint(data = vv2015,
## Variables
yvar = "y",
avar = "x",
mvar = "m",
cvar = c("c"),
eventvar = "event",
## Values at which effects are evaluated
a0 = 0,
a1 = 1,
m_cde = 1,
c_cond = 0.5,
## Model types
mreg = "logistic",
yreg = "survAFT_weibull",
## Additional specification
interaction = TRUE,
casecontrol = FALSE)
coef(summary(regmedint_obj))## est se Z p lower upper
## cde 0.541070807 0.29422958 1.8389409 0.06592388 -0.03560858 1.11775019
## pnde 0.488930417 0.21049248 2.3227928 0.02019028 0.07637274 0.90148809
## tnie 0.018240025 0.03706111 0.4921608 0.62260566 -0.05439841 0.09087846
## tnde 0.498503455 0.21209540 2.3503737 0.01875457 0.08280410 0.91420281
## pnie 0.008666987 0.02730994 0.3173565 0.75097309 -0.04485951 0.06219348
## te 0.507170442 0.21090051 2.4047853 0.01618197 0.09381303 0.92052785
## pm 0.045436278 0.09119614 0.4982259 0.61832484 -0.13330488 0.22417743boot packageThe boot package is the classical way to perform
bootstrapping in R. It requires defining a wrapper function.
library(boot)
## Define a wrapper function
regmedint_boot <- function(data, ind) {
## Note the change in the data argument.
regmedint_obj <- regmedint(data = data[ind,],
## Variables
yvar = "y",
avar = "x",
mvar = "m",
cvar = c("c"),
eventvar = "event",
## Values at which effects are evaluated
a0 = 0,
a1 = 1,
m_cde = 1,
c_cond = 0.5,
## Model types
mreg = "logistic",
yreg = "survAFT_weibull",
## Additional specification
interaction = TRUE,
casecontrol = FALSE)
coef(regmedint_obj)
}
## Run bootstrapping
ncpus <- 1
## For parallization, use the following instead.
## ncpus <- parallel::detectCores()
boot_obj <- boot(data = vv2015, statistic = regmedint_boot, R = 1000,
## For palatalization
## See https://cran.r-project.org/package=boot
parallel = "multicore",
ncpus = ncpus)
## Confidence interval for the pm
boot.ci(boot_obj, type = "basic", index = 7)## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot_obj, type = "basic", index = 7)
##
## Intervals :
## Level Basic
## 95% (-0.4322, 0.4106 )
## Calculations and Intervals on Original Scalemodelr packageIn the tidyverse ecosystem, the modelr package can be
used to provide a potentially more flexible workflow in some
settings.
library(modelr)
library(future)
future::plan(sequential)
## For parallization, use the following instead.
## future::plan(multiprocess)
library(furrr)## Error in library(furrr): there is no package called 'furrr'## Bootstrapping
tib_obj <- vv2015 %>%
modelr::bootstrap(n = 1000) %>%
## Resampled data objects are in the list column named strap.
mutate(boot_fit = future_map(strap, function(strap) {
## Note the change in the data argument.
regmedint_obj <- regmedint(data = as_tibble(strap),
## Variables
yvar = "y",
avar = "x",
mvar = "m",
cvar = c("c"),
eventvar = "event",
## Values at which effects are evaluated
a0 = 0,
a1 = 1,
m_cde = 1,
c_cond = 0.5,
## Model types
mreg = "logistic",
yreg = "survAFT_weibull",
## Additional specification
interaction = TRUE,
casecontrol = FALSE)
## Trick to return a row tibble
mat <- t(matrix(coef(regmedint_obj)))
colnames(mat) <- names(coef(regmedint_obj))
as_tibble(mat)
})) %>%
select(-strap) %>%
unnest(cols = c(boot_fit))## Error in `mutate()`:
## ℹ In argument: `boot_fit = future_map(...)`.
## Caused by error in `future_map()`:
## ! could not find function "future_map"tib_obj2 <- tib_obj %>%
pivot_longer(-.id) %>%
mutate(name = factor(name, levels = names(coef(regmedint_obj)))) %>%
group_by(name) %>%
summarize(lower_boot = quantile(value, probs = 0.025),
upper_boot = quantile(value, probs = 0.975))## Error in eval(expr, envir, enclos): object 'tib_obj' not found## Error in eval(expr, envir, enclos): object 'tib_obj2' not foundtib_obj2 %>%
mutate(lower_delta = confint(regmedint_obj)[,"lower"],
upper_delta = confint(regmedint_obj)[,"upper"])## Error in eval(expr, envir, enclos): object 'tib_obj2' not foundThese 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.
Health stats visible at Monitor.