mvboxcox implements an interpretable nonlinear extension
of logistic regression for binary outcomes and positive continuous
predictors, described in Bivariate logistic Box-Cox regression for
interpretable nonlinear exposure-response modeling by Shiyu Xu and
Xuekui Zhang. It is particularly useful when predictor-outcome
relationships may be nonlinear, including applications involving
right-skewed exposures. The package extends the univariate logistic
Box-Cox model of Xing et al. (2021) [2] to two positive
predictors, each with its own data-driven Box-Cox transformation
parameter, estimated by K-fold cross-validation with adaptive grid
refinement and thin-plate spline (TPS) refinement.
Install the package with:
install.packages("mvboxcox")For a positive predictor (x > 0), the Box-Cox transformation is (x^{()} = (x^- 1)/) for (> 0), and ((x)) for (= 0). The bivariate model is
[ {(Y_i = 1)} = _0 + 1 X{i1}^{(_1)} + 2 X{i2}^{(_2)} + ^{T}Z_i, ]
where (X_{i1}, X_{i2}) are positive exposures (formulaA)
and (Z_i) holds covariates left untransformed (formulaB).
(_j) describes the shape of the exposure-response relationship; (_j)
describes the association on the transformed scale. Because the
magnitude of (_j) depends on that scale, coefficients from different
transformation shapes should not be compared directly.
Let (m_j) denote the marginal median of a positive exposure (X_j). The median effect relative to a reference transformation (q) is
[ _j(q) = _j m_j^{_j-q}. ]
At (q=1), this is the effective linear slope at the median exposure and provides a common scale for comparing effects across transformation shapes. No parametric distribution is required. If (X_j) is log-normal with location parameter (_j), then (m_j=(_j)), and the expression reduces to (_j(q)=_j{(_j-q)_j}). Numerical comparisons between different predictors remain dependent on their measurement units.
The shape parameter for an exposure is identified through its association with the outcome. When the corresponding coefficient is close to zero, the shape parameter can be weakly identified. A flat cross-validation surface, a boundary estimate, or variability in the selected lambda may therefore reflect limited exposure-outcome information rather than numerical failure.
library(mvboxcox)
sim <- mvbc.simulator(
vLambda = c(0.5, 1.5), vBeta = c(-2.2, -0.4, -0.2, -0.005),
vMean = c(-0.08, -0.01, 50), vSd = c(0.93, 0.8, 18.12),
vNames = c("mercury", "lead", "age"), n = 2000, seed = 1
)
sim_test <- mvbc.simulator(simModel = sim)
fit <- suppressWarnings(
mvbc.train(Ybin ~ mercury + lead, ~ age, sim$data, depth = 2)
)
p_hat <- mvbc.predict(fit, sim_test$data)
mvbc.trainer.ssr(sim_test$data$Ybin, p_hat)
mvbc.median.effect(fit, sim$data, q = 1)mvbc.train(), mvbc.trainer(), and
mvbc.optimizer() also accept an observation
weights vector; with survey = TRUE they fit
each candidate model with a weighted glm() instead of an
unweighted glm(), and weight the cross-validation residual
criteria accordingly. See
vignette("introduction", package = "mvboxcox") for a full
walkthrough, including a weighted example.
The survey-weighted implementation incorporates observation weights but does not accept survey strata or primary sampling-unit identifiers. Results should therefore be interpreted as sampling-weighted rather than as complete design-based survey estimates.
The bundled depress data frame (8,893 rows, 6 columns)
is the analytic sample from the paper’s NHANES 2005-2006 and 2007-2008
application. It contains depression status, blood
mercury and blood_lead, age,
gender (1 = male, 0 = female), and the combined-cycle Day 1
dietary sampling weight.
data(depress, package = "mvboxcox")
summary(depress)The paper’s sampling-weighted BLBC model can be fitted with:
fit_nhanes <- mvbc.train(
depression ~ mercury + blood_lead,
~ age + factor(gender),
data = depress,
weights = depress$weight,
survey = TRUE
)
mvbc.median.effect(
fit_nhanes,
depress,
q = 1,
weights = depress$weight
)| Function | Purpose |
|---|---|
mvbc.simulator() |
Generate data from a logistic Box-Cox model |
ld() |
Apply the Box-Cox transformation |
mvbc.trainer() |
Evaluate lambda tuples by K-fold cross-validation |
mvbc.optimizer() |
Refine a trained lambda surface using TPS and L-BFGS-B |
mvbc.train() |
Run the trainer and optimizer as an end-to-end pipeline |
mvbc.predict() |
Predict probabilities from an mvbc.train.model |
mvbc.median.effect() |
Compute empirical or sampling-weighted median effects |
mvbc.tester() |
Evaluate repeated simulated train/test experiments |
mvbc.trainer.ssr() |
Calculate Pearson and deviance residual criteria |
Use help(package = "mvboxcox") for the complete
reference.
[1] Box, G. E. P., & Cox, D. R. (1964). An analysis of transformations. Journal of the Royal Statistical Society: Series B (Methodological), 26(2), 211-243.
[2] Xing, L., Zhang, X., Burstyn, I., & Gustafson, P. (2021). On logistic Box-Cox regression for flexibly estimating the shape and strength of exposure-disease relationships. Canadian Journal of Statistics, 49(3), 808-825.
[3] Xu, S., & Zhang, X. Bivariate logistic Box-Cox regression for interpretable nonlinear exposure-response modeling. Manuscript.