mmcm.resamp {mmcm} | R Documentation |
The permuted modified maximum contrast method
Description
This function gives P
-value for the permuted modified maximum
contrast method.
Usage
mmcm.resamp(
x,
g,
contrast,
alternative = c("two.sided", "less", "greater"),
nsample = 20000,
abseps = 0.001,
seed = NULL,
nthread = 2
)
Arguments
x |
a numeric vector of data values |
g |
a integer vector giving the group for the corresponding elements of x |
contrast |
a numeric contrast coefficient matrix for permuted modified maximum contrast statistics |
alternative |
a character string specifying the alternative hypothesis, must be one of "two.sided" (default), "greater" or "less". You can specify just the initial letter. |
nsample |
specifies the number of resamples (default: 20000) |
abseps |
specifies the absolute error tolerance (default: 0.001) |
seed |
a single value, interpreted as an integer;
see |
nthread |
sthe number of threads used in parallel computing, or FALSE that means single threading (default: 2) |
Details
mmcm.resamp
performs the permuted modified maximum contrast
method that is detecting a true response pattern under the unequal sample size
situation.
Y_{ij} (i=1, 2, \ldots; j=1, 2, \ldots, n_i)
is an observed response for j
-th individual in
i
-th group.
\bm{C}
is coefficient matrix for permuted modified maximum contrast
statistics (i \times k
matrix, i
: No. of groups, k
:
No. of pattern).
\bm{C}=(\bm{c}_1, \bm{c}_2, \ldots, \bm{c}_k)^{\rm{T}}
\bm{c}_k
is coefficient vector of k
-th pattern.
\bm{c}_k=(c_{k1}, c_{k2}, \ldots, c_{ki})^{\rm{T}} \qquad (\textstyle \sum_i c_{ki}=0)
M_{\max}
is a permuted modified maximum contrast statistic.
\bar{Y}_i=\frac{\sum_{j=1}^{n_i} Y_{ij}}{n_{i}},
\bar{\bm{Y}}=(\bar{Y}_1, \bar{Y}_2, \ldots, \bar{Y}_i, \ldots, \bar{Y}_a)^{\rm{T}},
M_{k}=\frac{\bm{c}^t_k \bar{\bm{Y}}}{\sqrt{\bm{c}^t_k \bm{c}_k}},
M_{\max}=\max(M_1, M_2, \ldots, M_k).
Consider testing the overall null hypothesis H_0: \mu_1=\mu_2=\ldots=\mu_i
,
versus alternative hypotheses H_1
for response petterns
(H_1: \mu_1<\mu_2<\ldots<\mu_i, \mu_1=\mu_2<\ldots<\mu_i,
\mu_1<\mu_2<\ldots=\mu_i
).
The P
-value for the probability distribution of M_{\max}
under the overall null hypothesis is
P\mbox{-value}=\Pr(M_{\max}>m_{\max} \mid H_0)
m_{\max}
is observed value of statistics.
This function gives distribution of M_{\max}
by using the
permutation method, follow algorithm:
1. Initialize counting variable: COUNT = 0
.
Input parameters: NRESAMPMIN
(minimum resampling count,
we set 1000), NRESAMPMAX
(maximum resampling count), and
\epsilon
(absolute error tolerance).
2. Calculate m_{\max}
that is the observed value of the test
statistic.
3. Let y_{ij}^{(r)}
donate data, which are sampled without
replacement, and independently, form observed value y_{ij}
.
Where, (r)
is suffix of the resampling number
(r = 1,\, 2,\, \ldots)
.
4. Calculate m^{(r)}_{\max}
from y_{ij}^{(r)}
.
If m^{(r)}_{\max} > m_{\max}
, then increment the
counting variable: COUNT = COUNT + 1
. Calculate
approximate P-value \hat{p}^{(r)}=COUNT/r
, and
the simulation standard error SE(\hat{p}^{(r)})=\sqrt{\hat{p}^{(r)}(1-
\hat{p}^{(r)})/r}
.
5. Repeat 3–4, while r > 1000
and 3.5SE(\hat{p}^{(r)})
< \epsilon
(corresponding to 99% confidence
level), or NRESAMPMAX
times. Output the approximate P-value
\hat{p}^{(r)}
.
Value
statistic |
the value of the test statistic with a name describing it. |
p.value |
the p-value for the test. |
alternative |
a character string describing the alternative hypothesis. |
method |
the type of test applied. |
contrast |
a character string giving the names of the data. |
contrast.index |
a suffix of coefficient vector of the |
error |
estimated absolute error and, |
msg |
status messages. |
References
Nagashima, K., Sato, Y., Hamada, C. (2011). A modified maximum contrast method for unequal sample sizes in pharmacogenomic studies Stat Appl Genet Mol Biol. 10(1): Article 41. http://dx.doi.org/10.2202/1544-6115.1560
Sato, Y., Laird, N.M., Nagashima, K., et al. (2009). A new statistical screening approach for finding pharmacokinetics-related genes in genome-wide studies. Pharmacogenomics J. 9(2): 137–146. http://www.ncbi.nlm.nih.gov/pubmed/19104505
See Also
Examples
## Example 1 ##
# true response pattern: dominant model c=(1, 1, -2)
set.seed(136885)
x <- c(
rnorm(130, mean = 1 / 6, sd = 1),
rnorm( 90, mean = 1 / 6, sd = 1),
rnorm( 10, mean = -2 / 6, sd = 1)
)
g <- rep(1:3, c(130, 90, 10))
boxplot(
x ~ g,
width = c(length(g[g==1]), length(g[g==2]), length(g[g==3])),
main = "Dominant model (sample data)",
xlab = "Genotype", ylab="PK parameter"
)
# coefficient matrix
# c_1: additive, c_2: recessive, c_3: dominant
contrast <- rbind(
c(-1, 0, 1), c(-2, 1, 1), c(-1, -1, 2)
)
y <- mmcm.resamp(x, g, contrast, nsample = 20000,
abseps = 0.01, seed = 5784324)
y
## Example 2 ##
# for dataframe
# true response pattern:
# pos = 1 dominant model c=( 1, 1, -2)
# 2 additive model c=(-1, 0, 1)
# 3 recessive model c=( 2, -1, -1)
set.seed(3872435)
x <- c(
rnorm(130, mean = 1 / 6, sd = 1),
rnorm( 90, mean = 1 / 6, sd = 1),
rnorm( 10, mean = -2 / 6, sd = 1),
rnorm(130, mean = -1 / 4, sd = 1),
rnorm( 90, mean = 0 / 4, sd = 1),
rnorm( 10, mean = 1 / 4, sd = 1),
rnorm(130, mean = 2 / 6, sd = 1),
rnorm( 90, mean = -1 / 6, sd = 1),
rnorm( 10, mean = -1 / 6, sd = 1)
)
g <- rep(rep(1:3, c(130, 90, 10)), 3)
pos <- rep(c("rsXXXX", "rsYYYY", "rsZZZZ"), each = 230)
xx <- data.frame(pos = pos, x = x, g = g)
# coefficient matrix
# c_1: additive, c_2: recessive, c_3: dominant
contrast <- rbind(
c(-1, 0, 1), c(-2, 1, 1), c(-1, -1, 2)
)
y <- by(xx, xx$pos, function(x) mmcm.resamp(x$x, x$g,
contrast, abseps = 0.02, nsample = 10000))
y <- do.call(rbind, y)[,c(3,7,9)]
y