BACKGROUND
dna microarray technology provides tools for studying the expression profiles of hundreds or thousands of distinct genes simultaneously.
a fundamental goal in microarray studies is to identify a subset of genes that are differentially expressed under experimental conditions of interest.
before conducting a microarray experiment, one important issue that needs to be determined is the number of arrays  required in order to have adequate power to identify differentially expressed genes.
many sample size estimation methods have been developed for various type i error specifications, such as family-wise error rate   <cit> , false discovery rate   <cit> , and the number of false positives  <cit> .
the sample size for a microarray study is commonly calculated as the number of arrays needed to achieve the specified power on average .
the power, the proportion of truly differentially expressed genes expected to be detected, is known as the sensitivity.
with the sample size estimate that is calculated to achieve a specified sensitivity on average, the proportion of truly differentially expressed genes detected would frequently be less than the average.
consequently, the sample size calculated tends to give an over-optimistic outcome.
alternatively, wang and chen  <cit> , tsai et al.
<cit>  and shao and tseng  <cit>  proposed an alternative formulation: the sample size is calculated to ensure detecting at least the specified sensitivity level with a specified probability.
this will be referred to as the  probability formulation.
when the sample size problem is formulated to achieve the specified sensitivity on average, we will show that the needed sample size can be simply calculated using the univariate sample size formula without considering dependency among genes.
on the other hand, if the problem is formulated to achieve a specified sensitivity with a specified probability, then it requires estimating a percentile of the distribution of sensitivity.
in this case, the dependency among genes needs to be taken into consideration.
tsai et al.
<cit>  presented an approach for controlling the comparison-wise error rate  under the model of independent or equi-correlated normal distribution with a constant power for all genes.
shao and tseng  <cit>  proposed a model-free procedure to estimate a general correlation matrix under the normal distribution.
they used a dataset of  <dig> samples to illustrate an estimation of the correlation matrix.
however, the size of pilot data is often small,  <dig> or fewer per group, and the estimated variances of the true positives are often negative  resulting in poor estimate of sample size in our simulation study.
tibshirani  <cit>  proposed a permutation method to estimate the fdr and average sensitivity for assessing a specific sample size.
tibshirani's method requires only a small number of pilot datasets and is completely model-free in the sense that no assumptions on the distribution, effect sizes, and correlations of the test statistics are required.
however, the standard deviation estimate  of a test statistic depends on the sample size.
a test statistic from a small sample size will have a larger variation than that from a larger sample size.
since the sample size of a pilot dataset is often small, the cutoff level based on a small pilot dataset often exceeds the true cutoff for needed samples and results in over-estimation of the needed sample size.
this paper presents an overview of the power and parameter specifications, and proposes a permutation procedure for sample size determination under the probability formulation .
the approach of tibshirani  <cit>  is improved to attain a more correct permutation distribution by incorporation of an adjustment factor.
the proposed method uses a small pilot dataset of  <dig> to  <dig> samples per group; the method requires fewer samples than the tibshirani  <cit>  method when the sample size for the pilot dataset is small relative to the needed sample size.
when the sample size for the pilot dataset is large, the proposed method and the tibshirani  <cit>  method are equivalent.
methods
let m denote the number of genes studied in an array of which m <dig> and m <dig> are the numbers of non-differentially and differentially expressed genes, respectively.
given the significance level α , the results of m tests can be summarized as a  <dig> ×  <dig> table .
v is the number of true null hypotheses that are falsely rejected;
u is the number of true alternative hypotheses that are correctly rejected;
s is the number of true null hypotheses that are correctly not rejected;
t is the number of alternative hypotheses that incorrectly not rejected;
r is the total number of null hypotheses rejected among the m tests.
v/m <dig> is the proportion of genes not differentially expressed that are declared significant, its expectation is the per comparison-wise error rate e/m <dig> = α. v/r is the proportion of declared significant genes among the total number of significances declared that are, in fact, not differentially expressed.
its expectation is the false discovery rate e = q, given r >  <dig>  u/m <dig> is the proportion of truly differentially expressed genes that are correctly declared.
in a diagnosis problem, this proportion is often referred as the true positive rate, or the sensitivity.
by taking expectation, we have the "average sensitivity", e/m <dig>  denoted by λ.
sample size estimation
in sample size estimation, m, m <dig>  and the  effect size δ =  for the differentially expressed genes are pre-specified by the investigator.
estimation of sample size needed to achieve the specified sensitivity λ <dig>  on average, is straightforward.
since m <dig> and λ <dig> are pre-specified, given a fdr level q* the corresponding significance level for per comparison-wise error rate α can easily be calculated.
setting α = /, the fdr will be controlled at q* for sufficiently large m <dig> and m <dig>
if δi = δ <dig> is constant for all i, then the comparison-wise power  of the univariate test is the same and exactly equal to λ <dig>  given α, δ <dig>  and  = λ <dig>  the sample size can be based on the univariate sample size calculation and is given as
where tα and tβ are the percentiles of a t-distribution.
if δi's are different, then βi = t- <dig> () from equation .
the sample size n* can be calculated from the following equation
the needed sample size is n = ⌈n*⌉, where ⌈n*⌉ is the smallest integer greater than or equal to n*.
given the sample size n as calculated, the outcome of a univariate test on a truly differentially expressed gene can be modeled by a bernoulli random variable with the success probability at least  since n ≥ n*.
the expected number of true detections is at least m <dig> λ <dig>  regardless of the correlation structure among genes and hence the desired sensitivity can be achieved on average.
most sample size estimation methods are either based on this approach or extensions  <cit> .
however, the sample size calculated under this formulation is inadequate; a simple demonstration under an independent model is shown below.
given m, π <dig> , a constant effect size δi = δ <dig>  q*, λ <dig>  and the calculated sample size n , under the independent model, the total number of truly differentially expressed genes detected u is a binomial random variable with success probability  .
the probability ϕλ <dig> of identifying at least λ <dig> fraction of m <dig> differentially expressed genes can be calculated as the sum of the binomial probabilities  <cit> :
the method of using equation  to estimate sample size is referred to as the univariate method.
column 3- <dig> of table  <dig> show the estimated sample size n, the average sensitivity λ and the probability ϕλ <dig> at λ <dig> =  <dig> ,  <dig> ,  <dig> ,  <dig> .
the parameters used in the calculation are: m =  <dig> , π <dig> = 5%, 10%, 20%, δ <dig> =  <dig> and q* =  <dig> .
it can be seen that the probability ϕλ <dig> can be less than 60%.
that is, using this formulation to calculate needed arrays may result in that an experiment will have the sensitivity less than the specified λ <dig> level with more than 40% probability.
a. estimated sample size n, average sensitivity λ and probability ϕλ <dig> for the specified sensitivity λ <dig> = 60%, 70%, 80%, 90%, under the independent model.
the parameters used in the calculation were: m =  <dig> , π <dig> = 5%, 10%, 20%, δ <dig> =  <dig> and q* =  <dig> .
b. sample size n is computed by the univariate method from equation  to achieve sensitivity λ <dig> on average.
c. sample size n is calculated using tsai et al.
<cit>  method to ensure the probability ϕλ <dig> of detecting at least λ <dig> fraction of differentially expressed genes is at least 95%.
alternatively, wang and chen  <cit>  formulated the problem as: the number of arrays needed to achieve the specified sensitivity λ  <dig> with a probability ϕ λ <dig>  in this formulation both λ  <dig> and ϕ λ <dig> need to be specified and not necessarily equal.
the ϕ λ <dig> is set at 95% since it is consistent with the common statistical practice of using the 95% confidence probability.
under this formulation, for specified λ  <dig> the needed number of arrays is calculated so that the average sensitivity is greater than λ  <dig> and the 5th percentile, λ  <dig>  of the distribution of the sensitivity u/m <dig> is greater than λ0:
in the independent and constant effect size model, tsai et al.
<cit>  used equations  and  to estimate the needed sample size which is referred to as the binomial method.
columns 6- <dig> of table  <dig> show the estimated sample size n, the average sensitivity λ, and the probability ϕλ <dig> for λ <dig> =  <dig> ,  <dig> ,  <dig> ,  <dig> .
the probabilities in column  <dig> are all higher than 95% due to n ≥ n*.
the procedure will ensure detecting the specified proportion of differentially expressed genes with at least 95% probability.
in table  <dig>  the theoretical results indicate that the two methods give quite close sample size estimates.
the difference of the estimates reflects the difference of the two formulations; when δ <dig> =  <dig>  the difference is up to  <dig>  for a given sensitivity, the needed sample size increases as the effect size δ <dig> decreasing, and the difference of the two formulations in the estimates is larger.
we calculated the sample sizes using the same parameters as table  <dig> for δ <dig> =  <dig>  the sample size differences increase at about four times those of table  <dig> .
permutation method for sample size estimation
tibshirani  <cit>  proposed a permutation method to account for both dependency and unequal effect sizes among genes using a pilot dataset for assessing sample size.
this method is applied here to estimate the required sample size.
because the sample size of the pilot data is typically smaller than the needed sample size, the null distributions generated from the pilot data have more variations; simply using the null distributions generated from a small pilot dataset can overestimate the needed sample size.
a procedure modified from the tibshirani  <cit>  method with adequate adjustment for sample size estimation is proposed below.
for simplicity, assume an equal sample size in each group, denoted as n = n <dig> = n <dig>  start with some pilot data with at least  <dig> samples per group, denoted as n0p and n1p for the control and treatment group, respectively.
for specified m, m <dig>  δ = , q*, and λ <dig>  the algorithm for a two sample t-test is described as follows.
algorithm: sample size estimation
<dig>  set α = /, use the method of tsai et al.
<cit>   to find the needed sample size as the initial sample size n.
<dig>  compute the adjustment factor f = f <dig> f <dig> where , , and tdf, p is the pth percentile of a t-distribution with df degrees of freedom.
<dig>  generate the b-th permutation samples.
<dig>  compute the t-statistics and sample standard deviations for the permutation samples for all genes.
<dig>  multiply each t-statistic by the factor f and add  to a set of randomly selected m <dig> t-statistic of differentially expressed genes to generate the permutation t-statistics sb = {s0b, s1b}, where s0b is the set for the non-differentially expressed genes, and s1b is the set for the differential expressed genes such that s0b = ft0b and s1b = ft1b + , where t0b and t1b are the vectors of the t-statistic, δ is a vector of the effect size and  is the vector of the sample standard deviation.
<dig>  store the permutation statistics sb.
<dig>  repeat 3- <dig> for all possible permutations, b =  <dig>   <dig>  ..., n, where n =  cn0p
<dig>  construct the null distribution by pooling all permutation statistics from the set of non-differentially expressed genes s <dig> = {s <dig>  s <dig>  ..., s0n}.
find the 100×th and 100×th percentiles as the critical values.
<dig>  compute the number of significances for the true positives ub for each statistic in s1b for each permutation sample b =  <dig>   <dig>  ..., n.
<dig>  order u <dig>  u <dig>  ..., un, and find the 5th percentile, denoted by u*.
<dig>  compare u* to m <dig> λ <dig>  if u* ≥ m1λ <dig>  stop and report n as the sample size estimate; otherwise, increase n by  <dig> and go to  <dig>
in the proposed algorithm, the permutation t-statistics of non-differentially expressed genes from all possible permutations were pooled to estimate the null distribution of the test statistics .
the number of true positives  was estimated for each permutation sample  since the set of differentially expressed genes in each permutation sample were known.
the distribution of the number of true positives u and its 5th percentile u* were estimated .
to reduce the excess variation of the permutation distribution, the proposed method includes the adjustment factor: f = f1f <dig>  the adjustment factor consists of two scale factors: f <dig> and f <dig>  the first factor, f <dig>  accounts for differential sample sizes between the pilot study and the planned study and the second scale factor, f <dig>  uses the maximum likelihood estimate of the t-statistic  <cit> .
when the sample size of pilot data is large, both factors f <dig> and f <dig> converge to  <dig> and the proposed and the tibshirani  <cit>  methods are equivalent.
since the permutation technique is used to estimate the critical value and the distribution of the sensitivity, no assumptions on the distribution of the t-statistic and the dependency among the statistics are made.
furthermore, the proposed method does not need to estimate the covariance matrix among all genes which can result in computation difficulty when the sample size of the pilot dataset is small.
RESULTS
two simulation analyses were conducted to evaluate the two formulations of sample size estimation described above.
the first analysis evaluated the two formulations under the independent and constant effect size model.
the theoretical results for the two formulations are shown in table  <dig>  the simulation analysis provides an empirical validation.
the second analysis evaluated the four methods under a correlated model: 1) the univariate method ; 2) the shao and tseng  <cit>  model-free method, 3) the tibshirani  <cit>  permutation method; and 4) the proposed permutation method.
the univariate method is designed for the average formulation, while the three other methods are considered with 95% probability with a use of a pilot dataset.
the same model parameters in table  <dig> were used in the evaluation.
the type i error rate was based on setting the fdr at q* =  <dig> .
note that there are many multiple testing fdr procedures with different strategies.
for example, the storey's fdr procedure  <cit>  involved an estimation of the number of non-differentially expressed genes m <dig>  however, to minimize the confounding effect brought by the variation in estimating m <dig>  we simply used the true m <dig> in our simulation analysis.
sample sizes were calculated for the given parameter values.
the empirical estimates of the fdr, average sensitivity λ and the probability ϕλ <dig> were then calculated and evaluated.
using the true m <dig> provides a direct validation of the proposed procedure with control of the fdr.
the purpose of the first simulation study was to validate the theoretical results of the sample size, sensitivity, and 95% probability for the two methods shown in table  <dig> under the independent model.
we generated  <dig>  simulation samples with sample sizes per group from the column  <dig> or column  <dig> of table  <dig>  for the null model, m <dig> = m ×  genes were generated from the independent standard normal n; for the alternative model, m <dig> = m × π <dig> genes were generated based on independent normal n. for each simulation sample set, the t-statistics and the correspondent p-values were computed, and the numbers of false positives and true positives at the fdr level q* =  <dig>  were recorded.
the empirical estimates of the fdr, average sensitivity λ and probability ϕλ <dig> were then calculated.
the estimate of ϕλ <dig> was the proportion of times out of the  <dig>  simulations that the number of true positives was not less than m <dig> × λ <dig>
a. empirical estimates of fdr q, average sensitivity λ, and probability ϕλ <dig> of the univariate method for the average formulation and of the binomial method for the 95% probability formulation.
the parameters used in the calculation were: m =  <dig> , δ <dig> =  <dig>  and q* =  <dig> .
b. sample size n is computed by the univariate method from equation  to achieve sensitivity λ <dig> on average.
c. sample size n is calculated using tsai et al.
<cit>  method to ensure sensitivity λ <dig> with 95% probability.
d. sample size n  is calculated using the proposed permutation method to ensure sensitivity λ <dig> with 95% probability with pilot study of group size  <dig> under the independent model.
for comparison purposes, the mean and standard deviation of the sample size estimates from the proposed permutation method using a pilot dataset of group size  <dig> are also provided in the last column of table  <dig>  the pilot data were randomly generated from the normal distribution in each simulation.
the proposed method tends to over-estimate the needed sample size by up to five arrays.
the second analysis was to evaluate the four methods, the univariate method , shao and tseng  <cit> , tibshirani  <cit> , and proposed permutation methods, under a correlated model using the well known colon cancer dataset  <cit> .
the colon cancer dataset  <cit>  consists of  <dig> normal and  <dig> colon tumor tissue samples with  <dig>  genes.
the analysis consisted of two steps.
the first step evaluated the sample size estimates obtained by the three 95% probability formulation methods based on a pilot dataset of sample size  <dig> and  <dig> per group.
the second step compared the sample sizes estimated by the proposed method from the first step with the estimates from the univariate method.
in the first step,  <dig> samples from the colon dataset were randomly selected without replacement from each group to form a pilot dataset.
the algorithm described above was used to estimate the sample size for the proposed method and the tibshirani  <cit>  method.
for example, for π <dig> = 5%, q* =  <dig>  and λ <dig> = 90%, the initial sample size was n =  <dig>  and α =  <dig> .
a constant effect size δi = δ <dig> =  <dig> was considered.
for the proposed permutation method, the initial adjustment factors for f were f <dig> =  <dig>  and f <dig> =  =  <dig> , while no adjustment was taken for the tibshirani  <cit>  method.
for the shao and tseng  <cit>  model-free method, a correlation matrix of t-statistics was estimated by using all possible permutation datasets from the pilot dataset.
however, the shao and tseng  <cit>  model-free method was found to have computational difficulty in most cases.
details are given later.
the procedure was repeated  <dig>  times to select different pilot datasets of size  <dig> from each group to account for the variation of pilot dataset.
the means and standard deviations of the sample size estimates from the tibshirani  <cit>  and proposed methods were calculated and are shown in columns  <dig> and  <dig> of table  <dig>  the univariate method is considered as the standard method, and the estimates are listed in column  <dig>  the needed sample size estimated from either the tibshirani  <cit>  or the proposed method is greater than that from the univariate method in each case.
the difference between the univariate method and the proposed method is less than  <dig> arrays per group in each case.
the mean and standard deviation estimates from the tibshirani  <cit>  method are much larger than the estimates from the proposed method.
the difference increases as λ <dig> increases or π <dig> decreases.
note that, under the independent model, the sample size and standard deviation estimates from the proposed method are smaller .
a. the sample size estimates are based on  <dig>  repetitions using the colon tumor data  <cit>  with  <dig> and  <dig> samples from each group as pilot dataset.
the parameters used in the calculation were: m =  <dig> , δ <dig> =  <dig> and q* =  <dig> .
b. the univariate method.
c. the proposed permutation method
d. the tibshirani  <cit>  permutation method.
e. the shao and tseng  <cit>  model-free method using the entire  <dig> samples.
the procedure was repeated with  <dig> samples for the initial pilot dataset.
the estimates are shown in columns  <dig> and  <dig>  the proposed procedure gives consistent results from the two pilot sample sizes; however, the results from the tibshirani  <cit>  method differ substantially.
the tibshirani approach does not adequately take the pilot sample size into consideration.
when the pilot sample size is much smaller than the needed sample size, the overestimation of the sample size by tibshirani  <cit>  method becomes severe.
as the pilot study size getting closer to the needed sample size, the tibshirani  <cit>  and the proposed methods will give similar results.
in our simulations, the algorithm b in shao and tseng  <cit>  couldn't successfully produce solutions for the pilot data of group size  <dig> in all  <dig>  replications.
when the group size increases to  <dig>  the algorithm works only when π <dig> = 20%, λ <dig> = 60% and 70%; the mean  of the sample size estimates are  <dig>  and  <dig> , respectively.
the estimated values appear too small to be correct.
this method does not appear to be applicable for small pilot sample sizes.
using the entire colon cancer dataset  <cit>  of  <dig> samples, the sample size estimates are shown in column  <dig>  the estimates generally need one or two more arrays than the univariate methods, but fewer than the proposed method.
since the tibshirani  <cit>  method gave larger estimates and the shao and tseng  <cit>  gave smaller estimates than the proposed method.
in the second step of analysis, the univariate method and the proposed method were evaluated.
comparison of the performance of the two methods is similar to that shown in table  <dig>  the data were sampled without replacement from the colon cancer dataset, instead of from the normal random variables under the independent model.
the sample sizes were based on column  <dig> or column  <dig> of table  <dig>  the data were then randomly permuted to remove the difference between two groups, and a common effect size δ <dig> =  <dig> was added to a set of randomly selected m <dig> genes in the tumor group.
for each re-sampled data set, the permutation test was used to generate a p-value and the numbers of false positives and true positives were computed using q* =  <dig> .
the number of repetitions to compute the permutation test was  <dig> .
the empirical estimates of fdr, λ and ϕλ <dig> were computed.
the entire procedure was repeated  <dig>  times.
the effect size of δ <dig> =  <dig>  was used to validate the proposed permutation method under a correlated model using the colon cancer dataset  <cit> .
in practice, the effect sizes can be much smaller.
we calculated the sample sizes using the same parameters as table  <dig> with an effect size δ <dig> =  <dig> for two pilot sample sizes  <dig> and  <dig>  the sample size estimates are shown in table  <dig>  the proposed procedure gives similar results for the two pilot sample sizes, which are consistent with the results for δ <dig> =  <dig> in table  <dig>  the difference between the univariate method and the proposed method is about  <dig> arrays per group.
the tibshirani  <cit>  method would require up to  <dig> and  <dig> extra arrays per group for  <dig> and  <dig> pilot samples, respectively.
the estimates for the shao and tseng  <cit>  method could be estimated only when the pilot study size is around or larger then the needed sample size.
a. the sample size estimates are based on  <dig>  repetitions using the colon tumor data  <cit>  with  <dig> and  <dig> samples from each group as pilot dataset.
the parameters used in the calculation were: m =  <dig> , δ <dig> =  <dig> and q* =  <dig> .
b. the univariate method.
c. the proposed permutation method
d. the tibshirani  <cit>  permutation method.
e. the shao and tseng  <cit>  model-free method using the entire  <dig> samples.
discussion and
CONCLUSIONS
determination of the needed sample size before conducting a microarray experiment is an important issue.
the sample size problem is commonly formulated as the number of arrays needed to achieve the specified sensitivity λ on average.
this paper demonstrates that the calculated sample size under this formulation may have the sensitivity λ at the specified level on average, but, the probability ϕλ that the specified sensitivity is achieved can be low  due to the variance in sensitivity distributions.
furthermore, under this formulation this paper shows that the sample size can be calculated by a univariate method, regardless of the correlation structure among the gene expression levels; the procedures to account for correlations, such as li et al.
<cit> , are not needed .
these findings agree with the results reported by jung  <cit>  and dobbin and simon  <cit> .
however, this paper provides a theoretical interpretation for this approach.
under the confidence probability formulation, consideration of the dependency among gene expressions is necessary in estimating the sample size since the percentile of the sensitivity distributions not only depends on the effect size of individual genes but also on their correlations.
we propose a permutation method based on the method proposed by tibshirani  <cit> , but with an inclusion of an adjustment factor and a requirement to achieve a specific sensitivity with 95% probability.
the adjustment factor provides more accurate estimates of the power and sample size.
shao and tseng  <cit>  also formulated the needed sample size in terms of confidence probability.
under the normality assumption, shao and tseng  <cit>  proposed algorithms for mild correlations among genes using a preliminary dataset.
they showed that their approach worked well for an example dataset of  <dig> samples.
however, using their algorithm b in our simulation for the colon dataset , the estimated variance of the true positives can be negative when the preliminary sample size is  <dig> or  <dig>  their procedure does not perform well for a small pilot dataset with small sample size.
in practice, sample sizes of pilot data are often small.
our simulation studies show that our procedure can work well with  <dig> to  <dig> samples per group.
however, our procedure seems to over-estimate the needed sample size when the correlations are very small, especially with small effect sizes.
in this situation, our simulation results indicate that the factor f <dig> may not be necessary .
the choice of a particular multiple testing procedure used for data analysis can affect the error rate and power in the sample size estimation.
using a conservative procedure in the data analysis may decrease the "power" of the study; sometimes, the calculated sample size may have sensitivity below the specified level.
for example, in this paper the calculation is based on the true number of non-differentially genes m <dig>  however, if the data analysis uses an overestimated m <dig> such as the benjamini and hochberg procedure  <cit> , then the power may be below the desired level.
an alternative is to use the total number of genes m instead of the number of non-differentially genes m <dig> to estimate the sample size.
this procedure is expected to generate an appropriate sample size to achieve the desired sensitivity with a specified probability, regardless of which multiple testing procedure is used for data analysis.
authors' contributions
jjc conceived the study and wrote the manuscript.
jjc and wjl developed the methodology and proved theoretical results.
wjl implemented the algorithms.
hmh improved the concepts of the average and 95% confidence probability formulations.
jjc, hmh and wjl performed the analysis.
all authors read and approved the final manuscript.
supplementary material
additional file 1
the software for the algorithm of the proposed method.
it provides software and an example for the algorithm of the proposed method.
click here for file
acknowledgements
huey-miin hsueh's research was done while visiting the nctr.
the authors are very grateful to reviewers for much helpful comments and suggestions for revising and improving this paper.
the views presented in this paper are those of the authors and do not necessarily represent those of the u.s. food and drug administration

 dna microarray technology provides tools for studying the expression profiles of hundreds or thousands of distinct genes simultaneously. a fundamental goal in microarray studies is to identify a subset of genes that are differentially expressed under experimental conditions of interest. before conducting a microarray experiment, one important issue that needs to be determined is the number of arrays required in order to have adequate power to identify differentially expressed genes. many sample size estimation methods have been developed for various type i error specifications, such as family-wise error rate, false discovery rate, and the number of false positives. the sample size for a microarray study is commonly calculated as the number of arrays needed to achieve the specified power on average. the power, the proportion of truly differentially expressed genes expected to be detected, is known as the sensitivity. with the sample size estimate that is calculated to achieve a specified sensitivity on average, the proportion of truly differentially expressed genes detected would frequently be less than the average. consequently, the sample size calculated tends to give an over-optimistic outcome. alternatively, wang and chen, tsai et al. and shao and tseng proposed an alternative formulation: the sample size is calculated to ensure detecting at least the specified sensitivity level with a specified probability. this will be referred to as the probability formulation. when the sample size problem is formulated to achieve the specified sensitivity on average, we will show that the needed sample size can be simply calculated using the univariate sample size formula without considering dependency among genes. on the other hand, if the problem is formulated to achieve a specified sensitivity with a specified probability, then it requires estimating a percentile of the distribution of sensitivity. in this case, the dependency among genes needs to be taken into consideration. tsai et al. presented an approach for controlling the comparison-wise error rate under the model of independent or equi-correlated normal distribution with a constant power for all genes. shao and tseng proposed a model-free procedure to estimate a general correlation matrix under the normal distribution. they used a dataset of samples to illustrate an estimation of the correlation matrix. however, the size of pilot data is often small, or fewer per group, and the estimated variances of the true positives are often negative resulting in poor estimate of sample size in our simulation study. tibshirani proposed a permutation method to estimate the fdr and average sensitivity for assessing a specific sample size. tibshirani's method requires only a small number of pilot datasets and is completely model-free in the sense that no assumptions on the distribution, effect sizes, and correlations of the test statistics are required. however, the standard deviation estimate of a test statistic depends on the sample size. a test statistic from a small sample size will have a larger variation than that from a larger sample size. since the sample size of a pilot dataset is often small, the cutoff level based on a small pilot dataset often exceeds the true cutoff for needed samples and results in over-estimation of the needed sample size. this paper presents an overview of the power and parameter specifications, and proposes a permutation procedure for sample size determination under the probability formulation. the approach of tibshirani is improved to attain a more correct permutation distribution by incorporation of an adjustment factor. the proposed method uses a small pilot dataset of to samples per group; the method requires fewer samples than the tibshirani method when the sample size for the pilot dataset is small relative to the needed sample size. when the sample size for the pilot dataset is large, the proposed method and the tibshirani method are equivalent. methods let m denote the number of genes studied in an array of which m and m are the numbers of non-differentially and differentially expressed genes, respectively. given the significance level α, the results of m tests can be summarized as a × table. v is the number of true null hypotheses that are falsely rejected; u is the number of true alternative hypotheses that are correctly rejected; s is the number of true null hypotheses that are correctly not rejected; t is the number of alternative hypotheses that incorrectly not rejected; r is the total number of null hypotheses rejected among the m tests. v/m is the proportion of genes not differentially expressed that are declared significant, its expectation is the per comparison-wise error rate e/m = α. v/r is the proportion of declared significant genes among the total number of significances declared that are, in fact, not differentially expressed. its expectation is the false discovery rate e = q, given r > u/m is the proportion of truly differentially expressed genes that are correctly declared. in a diagnosis problem, this proportion is often referred as the true positive rate, or the sensitivity. by taking expectation, we have the "average sensitivity", e/m denoted by λ. sample size estimation in sample size estimation, m, m and the effect size δ = for the differentially expressed genes are pre-specified by the investigator. estimation of sample size needed to achieve the specified sensitivity λ on average, is straightforward. since m and λ are pre-specified, given a fdr level q* the corresponding significance level for per comparison-wise error rate α can easily be calculated. setting α = /, the fdr will be controlled at q* for sufficiently large m and m if δi = δ is constant for all i, then the comparison-wise power of the univariate test is the same and exactly equal to λ given α, δ and = λ the sample size can be based on the univariate sample size calculation and is given as where tα and tβ are the percentiles of a t-distribution. if δi's are different, then βi = t- () from equation. the sample size n* can be calculated from the following equation the needed sample size is n = ⌈n*⌉, where ⌈n*⌉ is the smallest integer greater than or equal to n*. given the sample size n as calculated, the outcome of a univariate test on a truly differentially expressed gene can be modeled by a bernoulli random variable with the success probability at least since n ≥ n*. the expected number of true detections is at least m λ regardless of the correlation structure among genes and hence the desired sensitivity can be achieved on average. most sample size estimation methods are either based on this approach or extensions. however, the sample size calculated under this formulation is inadequate; a simple demonstration under an independent model is shown below. given m, π, a constant effect size δi = δ q*, λ and the calculated sample size n, under the independent model, the total number of truly differentially expressed genes detected u is a binomial random variable with success probability. the probability ϕλ of identifying at least λ fraction of m differentially expressed genes can be calculated as the sum of the binomial probabilities: the method of using equation to estimate sample size is referred to as the univariate method. column 3- of table show the estimated sample size n, the average sensitivity λ and the probability ϕλ at λ =,,,. the parameters used in the calculation are: m =, π = 5%, 10%, 20%, δ = and q* =. it can be seen that the probability ϕλ can be less than 60%. that is, using this formulation to calculate needed arrays may result in that an experiment will have the sensitivity less than the specified λ level with more than 40% probability. a. estimated sample size n, average sensitivity λ and probability ϕλ for the specified sensitivity λ = 60%, 70%, 80%, 90%, under the independent model. the parameters used in the calculation were: m =, π = 5%, 10%, 20%, δ = and q* =. b. sample size n is computed by the univariate method from equation to achieve sensitivity λ on average. c. sample size n is calculated using tsai et al. method to ensure the probability ϕλ of detecting at least λ fraction of differentially expressed genes is at least 95%. alternatively, wang and chen formulated the problem as: the number of arrays needed to achieve the specified sensitivity λ with a probability ϕ λ in this formulation both λ and ϕ λ need to be specified and not necessarily equal. the ϕ λ is set at 95% since it is consistent with the common statistical practice of using the 95% confidence probability. under this formulation, for specified λ the needed number of arrays is calculated so that the average sensitivity is greater than λ and the 5th percentile, λ of the distribution of the sensitivity u/m is greater than λ0: in the independent and constant effect size model, tsai et al. used equations and to estimate the needed sample size which is referred to as the binomial method. columns 6- of table show the estimated sample size n, the average sensitivity λ, and the probability ϕλ for λ =,,,. the probabilities in column are all higher than 95% due to n ≥ n*. the procedure will ensure detecting the specified proportion of differentially expressed genes with at least 95% probability. in table the theoretical results indicate that the two methods give quite close sample size estimates. the difference of the estimates reflects the difference of the two formulations; when δ = the difference is up to for a given sensitivity, the needed sample size increases as the effect size δ decreasing, and the difference of the two formulations in the estimates is larger. we calculated the sample sizes using the same parameters as table for δ = the sample size differences increase at about four times those of table. permutation method for sample size estimation tibshirani proposed a permutation method to account for both dependency and unequal effect sizes among genes using a pilot dataset for assessing sample size. this method is applied here to estimate the required sample size. because the sample size of the pilot data is typically smaller than the needed sample size, the null distributions generated from the pilot data have more variations; simply using the null distributions generated from a small pilot dataset can overestimate the needed sample size. a procedure modified from the tibshirani method with adequate adjustment for sample size estimation is proposed below. for simplicity, assume an equal sample size in each group, denoted as n = n = n start with some pilot data with at least samples per group, denoted as n0p and n1p for the control and treatment group, respectively. for specified m, m δ =, q*, and λ the algorithm for a two sample t-test is described as follows. algorithm: sample size estimation set α = /, use the method of tsai et al. to find the needed sample size as the initial sample size n. compute the adjustment factor f = f f where,, and tdf, p is the pth percentile of a t-distribution with df degrees of freedom. generate the b-th permutation samples. compute the t-statistics and sample standard deviations for the permutation samples for all genes. multiply each t-statistic by the factor f and add to a set of randomly selected m t-statistic of differentially expressed genes to generate the permutation t-statistics sb = {s0b, s1b}, where s0b is the set for the non-differentially expressed genes, and s1b is the set for the differential expressed genes such that s0b = ft0b and s1b = ft1b +, where t0b and t1b are the vectors of the t-statistic, δ is a vector of the effect size and is the vector of the sample standard deviation. store the permutation statistics sb. repeat 3- for all possible permutations, b =..., n, where n = cn0p construct the null distribution by pooling all permutation statistics from the set of non-differentially expressed genes s = {s s..., s0n}. find the 100×th and 100×th percentiles as the critical values. compute the number of significances for the true positives ub for each statistic in s1b for each permutation sample b =..., n. order u u..., un, and find the 5th percentile, denoted by u*. compare u* to m λ if u* ≥ m1λ stop and report n as the sample size estimate; otherwise, increase n by and go to in the proposed algorithm, the permutation t-statistics of non-differentially expressed genes from all possible permutations were pooled to estimate the null distribution of the test statistics. the number of true positives was estimated for each permutation sample since the set of differentially expressed genes in each permutation sample were known. the distribution of the number of true positives u and its 5th percentile u* were estimated. to reduce the excess variation of the permutation distribution, the proposed method includes the adjustment factor: f = f1f the adjustment factor consists of two scale factors: f and f the first factor, f accounts for differential sample sizes between the pilot study and the planned study and the second scale factor, f uses the maximum likelihood estimate of the t-statistic. when the sample size of pilot data is large, both factors f and f converge to and the proposed and the tibshirani methods are equivalent. since the permutation technique is used to estimate the critical value and the distribution of the sensitivity, no assumptions on the distribution of the t-statistic and the dependency among the statistics are made. furthermore, the proposed method does not need to estimate the covariance matrix among all genes which can result in computation difficulty when the sample size of the pilot dataset is small.  two simulation analyses were conducted to evaluate the two formulations of sample size estimation described above. the first analysis evaluated the two formulations under the independent and constant effect size model. the theoretical results for the two formulations are shown in table the simulation analysis provides an empirical validation. the second analysis evaluated the four methods under a correlated model: 1) the univariate method; 2) the shao and tseng model-free method, 3) the tibshirani permutation method; and 4) the proposed permutation method. the univariate method is designed for the average formulation, while the three other methods are considered with 95% probability with a use of a pilot dataset. the same model parameters in table were used in the evaluation. the type i error rate was based on setting the fdr at q* =. note that there are many multiple testing fdr procedures with different strategies. for example, the storey's fdr procedure involved an estimation of the number of non-differentially expressed genes m however, to minimize the confounding effect brought by the variation in estimating m we simply used the true m in our simulation analysis. sample sizes were calculated for the given parameter values. the empirical estimates of the fdr, average sensitivity λ and the probability ϕλ were then calculated and evaluated. using the true m provides a direct validation of the proposed procedure with control of the fdr. the purpose of the first simulation study was to validate the theoretical results of the sample size, sensitivity, and 95% probability for the two methods shown in table under the independent model. we generated simulation samples with sample sizes per group from the column or column of table for the null model, m = m × genes were generated from the independent standard normal n; for the alternative model, m = m × π genes were generated based on independent normal n. for each simulation sample set, the t-statistics and the correspondent p-values were computed, and the numbers of false positives and true positives at the fdr level q* = were recorded. the empirical estimates of the fdr, average sensitivity λ and probability ϕλ were then calculated. the estimate of ϕλ was the proportion of times out of the simulations that the number of true positives was not less than m × λ a. empirical estimates of fdr q, average sensitivity λ, and probability ϕλ of the univariate method for the average formulation and of the binomial method for the 95% probability formulation. the parameters used in the calculation were: m =, δ = and q* =. b. sample size n is computed by the univariate method from equation to achieve sensitivity λ on average. c. sample size n is calculated using tsai et al. method to ensure sensitivity λ with 95% probability. d. sample size n is calculated using the proposed permutation method to ensure sensitivity λ with 95% probability with pilot study of group size under the independent model. for comparison purposes, the mean and standard deviation of the sample size estimates from the proposed permutation method using a pilot dataset of group size are also provided in the last column of table the pilot data were randomly generated from the normal distribution in each simulation. the proposed method tends to over-estimate the needed sample size by up to five arrays. the second analysis was to evaluate the four methods, the univariate method, shao and tseng, tibshirani, and proposed permutation methods, under a correlated model using the well known colon cancer dataset. the colon cancer dataset consists of normal and colon tumor tissue samples with genes. the analysis consisted of two steps. the first step evaluated the sample size estimates obtained by the three 95% probability formulation methods based on a pilot dataset of sample size and per group. the second step compared the sample sizes estimated by the proposed method from the first step with the estimates from the univariate method. in the first step, samples from the colon dataset were randomly selected without replacement from each group to form a pilot dataset. the algorithm described above was used to estimate the sample size for the proposed method and the tibshirani method. for example, for π = 5%, q* = and λ = 90%, the initial sample size was n = and α =. a constant effect size δi = δ = was considered. for the proposed permutation method, the initial adjustment factors for f were f = and f = =, while no adjustment was taken for the tibshirani method. for the shao and tseng model-free method, a correlation matrix of t-statistics was estimated by using all possible permutation datasets from the pilot dataset. however, the shao and tseng model-free method was found to have computational difficulty in most cases. details are given later. the procedure was repeated times to select different pilot datasets of size from each group to account for the variation of pilot dataset. the means and standard deviations of the sample size estimates from the tibshirani and proposed methods were calculated and are shown in columns and of table the univariate method is considered as the standard method, and the estimates are listed in column the needed sample size estimated from either the tibshirani or the proposed method is greater than that from the univariate method in each case. the difference between the univariate method and the proposed method is less than arrays per group in each case. the mean and standard deviation estimates from the tibshirani method are much larger than the estimates from the proposed method. the difference increases as λ increases or π decreases. note that, under the independent model, the sample size and standard deviation estimates from the proposed method are smaller. a. the sample size estimates are based on repetitions using the colon tumor data with and samples from each group as pilot dataset. the parameters used in the calculation were: m =, δ = and q* =. b. the univariate method. c. the proposed permutation method d. the tibshirani permutation method. e. the shao and tseng model-free method using the entire samples. the procedure was repeated with samples for the initial pilot dataset. the estimates are shown in columns and the proposed procedure gives consistent results from the two pilot sample sizes; however, the results from the tibshirani method differ substantially. the tibshirani approach does not adequately take the pilot sample size into consideration. when the pilot sample size is much smaller than the needed sample size, the overestimation of the sample size by tibshirani method becomes severe. as the pilot study size getting closer to the needed sample size, the tibshirani and the proposed methods will give similar results. in our simulations, the algorithm b in shao and tseng couldn't successfully produce solutions for the pilot data of group size in all replications. when the group size increases to the algorithm works only when π = 20%, λ = 60% and 70%; the mean of the sample size estimates are and, respectively. the estimated values appear too small to be correct. this method does not appear to be applicable for small pilot sample sizes. using the entire colon cancer dataset of samples, the sample size estimates are shown in column the estimates generally need one or two more arrays than the univariate methods, but fewer than the proposed method. since the tibshirani method gave larger estimates and the shao and tseng gave smaller estimates than the proposed method. in the second step of analysis, the univariate method and the proposed method were evaluated. comparison of the performance of the two methods is similar to that shown in table the data were sampled without replacement from the colon cancer dataset, instead of from the normal random variables under the independent model. the sample sizes were based on column or column of table the data were then randomly permuted to remove the difference between two groups, and a common effect size δ = was added to a set of randomly selected m genes in the tumor group. for each re-sampled data set, the permutation test was used to generate a p-value and the numbers of false positives and true positives were computed using q* =. the number of repetitions to compute the permutation test was. the empirical estimates of fdr, λ and ϕλ were computed. the entire procedure was repeated times. the effect size of δ = was used to validate the proposed permutation method under a correlated model using the colon cancer dataset. in practice, the effect sizes can be much smaller. we calculated the sample sizes using the same parameters as table with an effect size δ = for two pilot sample sizes and the sample size estimates are shown in table the proposed procedure gives similar results for the two pilot sample sizes, which are consistent with the results for δ = in table the difference between the univariate method and the proposed method is about arrays per group. the tibshirani method would require up to and extra arrays per group for and pilot samples, respectively. the estimates for the shao and tseng method could be estimated only when the pilot study size is around or larger then the needed sample size. a. the sample size estimates are based on repetitions using the colon tumor data with and samples from each group as pilot dataset. the parameters used in the calculation were: m =, δ = and q* =. b. the univariate method. c. the proposed permutation method d. the tibshirani permutation method. e. the shao and tseng model-free method using the entire samples. discussion and  determination of the needed sample size before conducting a microarray experiment is an important issue. the sample size problem is commonly formulated as the number of arrays needed to achieve the specified sensitivity λ on average. this paper demonstrates that the calculated sample size under this formulation may have the sensitivity λ at the specified level on average, but, the probability ϕλ that the specified sensitivity is achieved can be low due to the variance in sensitivity distributions. furthermore, under this formulation this paper shows that the sample size can be calculated by a univariate method, regardless of the correlation structure among the gene expression levels; the procedures to account for correlations, such as li et al., are not needed. these findings agree with the results reported by jung and dobbin and simon. however, this paper provides a theoretical interpretation for this approach. under the confidence probability formulation, consideration of the dependency among gene expressions is necessary in estimating the sample size since the percentile of the sensitivity distributions not only depends on the effect size of individual genes but also on their correlations. we propose a permutation method based on the method proposed by tibshirani, but with an inclusion of an adjustment factor and a requirement to achieve a specific sensitivity with 95% probability. the adjustment factor provides more accurate estimates of the power and sample size. shao and tseng also formulated the needed sample size in terms of confidence probability. under the normality assumption, shao and tseng proposed algorithms for mild correlations among genes using a preliminary dataset. they showed that their approach worked well for an example dataset of samples. however, using their algorithm b in our simulation for the colon dataset, the estimated variance of the true positives can be negative when the preliminary sample size is or their procedure does not perform well for a small pilot dataset with small sample size. in practice, sample sizes of pilot data are often small. our simulation studies show that our procedure can work well with to samples per group. however, our procedure seems to over-estimate the needed sample size when the correlations are very small, especially with small effect sizes. in this situation, our simulation results indicate that the factor f may not be necessary. the choice of a particular multiple testing procedure used for data analysis can affect the error rate and power in the sample size estimation. using a conservative procedure in the data analysis may decrease the "power" of the study; sometimes, the calculated sample size may have sensitivity below the specified level. for example, in this paper the calculation is based on the true number of non-differentially genes m however, if the data analysis uses an overestimated m such as the benjamini and hochberg procedure, then the power may be below the desired level. an alternative is to use the total number of genes m instead of the number of non-differentially genes m to estimate the sample size. this procedure is expected to generate an appropriate sample size to achieve the desired sensitivity with a specified probability, regardless of which multiple testing procedure is used for data analysis. authors' contributions jjc conceived the study and wrote the manuscript. jjc and wjl developed the methodology and proved theoretical results. wjl implemented the algorithms. hmh improved the concepts of the average and 95% confidence probability formulations. jjc, hmh and wjl performed the analysis. all authors read and approved the final manuscript. supplementary material additional file 1 the software for the algorithm of the proposed method. it provides software and an example for the algorithm of the proposed method. click here for file acknowledgements huey-miin hsueh's research was done while visiting the nctr. the authors are very grateful to reviewers for much helpful comments and suggestions for revising and improving this paper. the views presented in this paper are those of the authors and do not necessarily represent those of the u.s. food and drug administration
this paper discusses some crucial issues in the problem formulation, parameter specifications, and approaches that are commonly proposed for sample size estimation in microarray experiments. before conducting a microarray experiment, one important issue that needs to be determined is the number of arrays required in order to have adequate power to identify differentially expressed genes. a permutation method using a small pilot dataset to estimate sample size is proposed. on the other hand, the needed sample size calculated by the proposed permutation method will ensure detecting at least the desired sensitivity with 95% probability. the proposed permutation method takes the correlation structure and effect size heterogeneity into consideration and works well using only a small pilot dataset. common methods for sample size estimation are formulated as the minimum sample size necessary to achieve a specified sensitivity on average at a specified false discovery rate level and specified expected proportion of the true differentially expression genes in the array. we formulate the sample size problem as the number of arrays needed to achieve a specified sensitivity with 95% probability at the specified significance level. a sample size estimate based on the common formulation, to achieve the desired sensitivity on average, can be calculated using a univariate method without taking the correlation among genes into consideration. this formulation of sample size problem is inadequate because the probability of detecting the specified sensitivity can be lower than 50%. we recommend that the sample size problem should be formulated to detect a specified proportion of differentially expressed genes with 95% probability. the method is shown to perform well for a real example dataset using a small pilot dataset with 4- samples per group. this formulation ensures finding the desired proportion of true positives with high probability.
 before conducting a microarray experiment, one important issue that needs to be determined is the number of arrays required in order to have adequate power to identify differentially expressed genes. this paper discusses some crucial issues in the problem formulation, parameter specifications, and approaches that are commonly proposed for sample size estimation in microarray experiments. common methods for sample size estimation are formulated as the minimum sample size necessary to achieve a specified sensitivity on average at a specified false discovery rate level and specified expected proportion of the true differentially expression genes in the array. unfortunately, the probability of detecting the specified sensitivity in such a formulation can be low. we formulate the sample size problem as the number of arrays needed to achieve a specified sensitivity with 95% probability at the specified significance level. a permutation method using a small pilot dataset to estimate sample size is proposed. this method accounts for correlation and effect size heterogeneity among genes.  a sample size estimate based on the common formulation, to achieve the desired sensitivity on average, can be calculated using a univariate method without taking the correlation among genes into consideration. this formulation of sample size problem is inadequate because the probability of detecting the specified sensitivity can be lower than 50%. on the other hand, the needed sample size calculated by the proposed permutation method will ensure detecting at least the desired sensitivity with 95% probability. the method is shown to perform well for a real example dataset using a small pilot dataset with 4- samples per group.  we recommend that the sample size problem should be formulated to detect a specified proportion of differentially expressed genes with 95% probability. this formulation ensures finding the desired proportion of true positives with high probability. the proposed permutation method takes the correlation structure and effect size heterogeneity into consideration and works well using only a small pilot dataset. 