dBH_lm.RddBH_lm computes the rejection set of \(dBH_\gamma(\alpha)\) and \(dBH^2_\gamma(\alpha)\), as well as
\(dSU_{\gamma, \Delta}(\alpha)\) and \(dSU^2_{\gamma, \Delta}(\alpha)\) a broad class of step-up procedures for
linear models \(y = X\beta + \epsilon\) where \(X\) is a fixed-design matrix and \(\epsilon\) has i.i.d.
components drawn from \(N(0, \sigma^2)\) for some unknonw variance \(\sigma^2\). dBH_lm can handle
both one-sided tests \(H_i: \beta_i \le 0\) or \(H_i: \beta_i \ge 0\), and two-sided tests \(H_i: \beta_i = 0\).
dBH_lm( y, X, subset = 1:ncol(X), intercept = TRUE, side = c("right", "left", "two"), alpha = 0.05, gamma = NULL, tautype = "QC", niter = 1, avals = NULL, avals_type = c("BH", "geom", "bonf", "manual"), geom_fac = 2, eps = 0.05, qcap = 2, gridsize = 20, exptcap = 0.9, is_safe = NULL, verbose = FALSE )
| y | the outcome vector. |
|---|---|
| X | the design matrix. The intercept term should be added manually into |
| subset | a subset of |
| intercept | a logical indicating whether an intercept is included in the linear regression. By default,
the intercept term is not tested. If it needs to be tested together with other variables, manually add it into
|
| side | a string that takes values in {"right", "left", "two"}, with "right" for one-sided tests \(H_i: \mu_i \le 0\), "left" for one-sided tests, \(H_i: \mu_i \ge 0\), and "two" for two-sided tests \(H_i: \mu_i = 0\). |
| alpha | the target FDR level. |
| gamma | the parameter for the dBH and dSU procedures. The default is |
| tautype | the type of tau function. In the current version, only "QC" (q-value-calibration) is supported with \(\tau(c; X) = cR_{BH}(c) / m\). |
| niter | the number of iterations. In the current version it can only be 1, for dBH/dSU, or 2, for \(dBH^2\)/\(dSU^2\). |
| avals | the a-values in the step-up procedures defined in Appendix C.1. The default is NULL, in which case
|
| avals_type | a string that takes values in {"BH", "geom", "bonf", "manual"}, which determines the type of thresholds in the step-up procedures. See Details. |
| geom_fac | a real number that is larger than 1. This is the growing rate of thresholds when
|
| eps | a real number in [0, 1], which is used to determine \(t_{hi}\) discussed in Section 4.2. |
| qcap | a real number that is larger than 1. It is used to filter out hypotheses with q-values above
|
| gridsize | an integer for the size of the grid used when |
| exptcap | a real number in [0, 1]. It is used by \(dBH^2\)/\(dSU^2\) to filter out hypotheses with
\(g_i*(q_i | S_i)\) below |
| is_safe | a logical or NULL indicating whether the procedure is taken as safe. DON'T set |
| verbose | a logical indicating whether a progress bar is shown. |
a list with the following attributes
rejs: the indices of rejected hypotheses (after the randomized pruning step if any);
initrejs: the indices of rejected hypotheses (before the randomized pruning step if any);
cand: the set of candidate hypotheses for which \(g_i^{*}(q_i | S_i)\) is evaluated;
expt: \(g_i^{*}(q_i | S_i)\) for each hypothesis in cand;
safe: TRUE iff the procedure is safe;
secBH: TRUE iff the randomzied pruning step (a.k.a. the secondary BH procedure) is invoked;
secBH_fac: a vector that gives \(\hat{R}_i / R_{+}\) that is defined in Section 2.2. It only shows up
in the output if secBH = TRUE.
dBH_lm can handle all dSU procedures that are defined in Appendix C.1. with
$$\Delta_{\alpha}(r) = \frac{\alpha a_{\ell}}{m}, r\in [a_{\ell}, a_{\ell + 1}), \ell = 0, 1, \ldots, L$$
for any set of integer a-values \(1\le a_1 < \ldots < a_L\le m\) (with \(a_0 = 0, a_{L+1} = m+1\)). There are two-ways to input the a-values.
(Recommended) use avals_type while leave avals as its default. Three types of avals_type
are supported.
When avals_type = "BH", the BH procedure is used (i.e. \(a_\ell = \ell, \ell = 1, \ldots, m\));
When avals_type = "geom", the geometrically increasing a-values that are defined in Appendix C.1.
are used with the growing rate specified by geom_fac, whose default is 2;
When avals_type = "bonf", the Bonferroni procedure is used (i.e. \(a_1 = 1, L = 1\)).
(Not recommended) use avals while set avals_type = "manual". This option allows any set of
increasing a-values. However, it can be much slower than the recommended approach.
# Generate beta n <- 150 p <- 100 beta1 <- 0.5 nalt <- 10 beta <- c(rep(beta1, nalt), rep(0, p - nalt)) # Generate X and y set.seed(1) X <- matrix(rnorm(n * p), nrow = n) y <- X %*% beta + rnorm(n) # Run dBH_1(\alpha) for one-sided tests alpha <- 0.05 res <- dBH_lm(y = y, X = X, side = "right", alpha = alpha, gamma = 1, niter = 1, avals_type = "BH") # Run dBH_1(\alpha) for two-sided tests res <- dBH_lm(y = y, X = X, side = "right", alpha = alpha, gamma = 1, niter = 1, avals_type = "BH") # Run dBH^2_1(\alpha) for one-sided tests alpha <- 0.05 res <- dBH_lm(y = y, X = X, side = "right", alpha = alpha, gamma = 1, niter = 2, avals_type = "BH") # Run dBH^2_1(\alpha) for one-sided tests res <- dBH_lm(y = y, X = X, side = "right", alpha = alpha, gamma = 1, niter = 2, avals_type = "BH") # Run dSU_1(\alpha) with the geometrically increasing a-values for one-sided tests res <- dBH_lm(y = y, X = X, side = "right", alpha = alpha, gamma = 1, niter = 1, avals_type = "geom", geom_fac = 2) # Run dBY(\alpha) for one-sided tests res <- dBH_lm(y = y, X = X, side = "right", alpha = alpha, gamma = NULL, niter = 1, avals_type = "BH") # Run dBY^2(\alpha) for one-sided tests res <- dBH_lm(y = y, X = X, side = "right", alpha = alpha, gamma = NULL, niter = 2, avals_type = "BH")