BACKGROUND
affymetrix microarrays are high throughput assays for measuring the expression levels of thousands of gene transcripts simultaneously.
this type of microarray measures the expression of each transcript multiple times through a set of "probe pairs".
since the advent of the affymetrix microarray, numerous methods have been proposed for producing numerical expression summaries for each transcript based on the probe pair data.
several systematic studies have appeared comparing a number of methods on a common basis .
these studies rely heavily on calibration data sets derived from spike-in, dilution series, and mixture experiments for comparing methods.
our goal here was to carry out a comparative study of affymetrix array processing methods using data sets from typical biological experiments seeking differentially expressed genes in human tissue samples.
the following seven methods are considered here: dchip  <cit> , gcrma-eb and gcrma-mle  <cit> , mas <dig>  <cit> , pdnn  <cit> , rma  <cit> , and tm .
while not every popular method is included in our study, several highly distinctive and original approaches are studied.
for example, dchip was one of the first approaches to attempt to learn probe weights directly from the probe data, and rma pioneered the approach of disregarding the control mismatch probes.
pdnn uses physical modeling to determine probe weights, while the two gcrma methods use gc content of the probe sequences to reduce variance in the mismatch  probe levels.
the mas <dig> method is the current default method provided by affymetrix.
in addition to the six methods cited previously, we also include a method designated tm .
this is a simple method that has been used in a number of published investigations , but has not been considered in any previous systematic comparison of affymetrix processing methods.
to produce the probe-set summary score, the pm-mm differences are rank ordered, and the brightest 20% and dimmest 20% of values are deleted.
the mean of the remaining values is used as the summary score.
the scores for all probe-sets are then quantile normalized to a reference array using a piecewise linear spline with  <dig> knots.
an important feature of this study is the use of false discovery rate  to quantify the sensitivity of a processing method in terms of its ability to distinguish differentially expressed genes from genes having invariant expression.
this is a highly relevant property, as differential expression analysis is the most common application of microarray data.
a key advantage of using fdr to compare processing methods is that fdr values can be calculated accurately using real disease profiling data where the identities of differentially expressed genes are uncertain.
in contrast, most previous systematic comparisons of array processing methods have focused on calibration data sets in which concentrations of certain genes were experimentally manipulated.
when it is highly likely that at least one gene is differentially expressed, false discovery rate may be defined as the expected ratio of the number of false positive calls to the total number of positive calls in a differential expression analysis between two groups of samples  <cit> .
if the groups are biologically distinct, a sensitive processing method should result in many genes with low fdr.
thus to compare the performances of different array processing methods, we looked at two datasets in which a verified biological characteristic divided the samples into two classes, and compared the methods based on the number of genes having fdr smaller than various thresholds.
for this to be a valid basis for comparison, the fdr values must be estimated with reasonable accuracy.
following other recent work , we used a permutation approach for this estimation, arguing that there is no reason that this approach favors or disfavors any particular array processing method.
a small fdr is due either to a small numerator, a large denominator, or both.
the denominator of the fdr depends on the actual data distribution, so variation in this value may be due to factors such as accuracy in modeling the physical and chemical nature of probe binding.
variation in the fdr numerator, however, depends only on the distribution of values produced for randomized data, a purely statistical quantity reflecting the tendency of the method to incorrectly produce test statistic outliers.
our results suggest that both factors are important in determining sensitivity.
the best methods produce many large test statistic values in the actual data, and also produce consistently small test statistic values for randomized data.
poor performance of one method can be directly explained by the tendency of the method to produce outlier expression values, leading to greater numbers of incorrectly large test statistics.
for overall comparison, we evaluated every pair of methods on the basis of whether the first method is expected to call at least one truly differentially expressed gene that is not also called by the second method.
if this is not expected to occur, the second method is said to strongly outperform the first.
based on this comparison, two of the methods considered are clearly favored, two are inferior, and results for the other three methods are mixed.
RESULTS
sensitivity differences
our primary basis for comparison is sensitivity – the number of genes detected at a given fdr <dig> level, where fdr <dig> is a rescaled fdr .
figure  <dig> shows the key sensitivity results, using both the t-test statistic and the rank-sum statistic to assess differential expression.
setting aside at first differences between the seven processing methods, we note two findings.
first, in the colon data, analysis using the rank-sum statistic is substantially more sensitive than analysis using the t-test statistic.
for the ovary data, where the sample sizes would not naturally suggest a robust analysis, there is no harm to sensitivity in using the rank-sum statistic.
second, the ovary curves are substantially higher overall than the colon curves.
this may be due to a greater number of true positives in the ovary data, or it may be that the small sample size for the msi group makes it difficult to attain high evidence levels for differential expression in the colon data.
in any case, both data sets have many genes with small fdr values, supporting the biological relevance of the tumor groupings for both colon and ovary samples.
the more challenging colon set distinguishes the seven processing methods to a greater extent than the ovary set.
using fdr <dig> =  <dig>  as a reference point, there is roughly 7-fold variation across the seven methods in the number of detected genes in colon data using t-test statistics, while for rank-sum statistics the range is roughly 2-fold.
in the ovary data, the range is around  <dig> -fold for both statistics.
also notable is that variation in sensitivity due to the choice of test statistic  is smaller than variation in sensitivity due to the seven processing methods.
no single method stands out as having the best or worst performance in every case.
however some methods generally perform better than others.
the dchip and tm methods perform consistently well, while the gcrma-eb and mas <dig> methods consistently perform poorly.
pdnn performs well on the ovary data, but poorly on the colon data, and results for the other methods are mixed.
level of agreement between methods
identities of probe sets falling below a given fdr <dig> threshold vary across the methods.
figure  <dig> summarizes this variation.
the ratio of the number of probe sets falling below various fdr <dig> thresholds in k or more of the seven methods to the number of probe sets falling below the threshold for at least one method is plotted against the fdr <dig> threshold, for k =  <dig>   <dig>   <dig>   <dig>   <dig>  in the ovary data there is a very high level of agreement in this measure.
for the rank-sum analysis, almost 90% of called genes are called by at least four methods, and more than 70% of called genes are called by all seven methods.
for the t-test analysis, the agreement is slightly higher yet.
for the colon data, the methods are much more inconsistent.
for the rank-sum analysis, three of the methods agree on up to 90% of genes, but all seven methods only agree on around 30% of genes.
the t-test analysis is even worse, with only around 10% of genes common to all seven methods.
turning to pairwise agreement, table  <dig> shows the percentage of genes called by both members of a pair of methods out of the genes called by at least one of the two methods.
in the ovary data, mas <dig> shares the fewest calls with the other methods for both t-test and rank-sum analysis, while gcrma-eb has relatively weak agreement for the t-test analysis.
in the colon data, the gcrma-eb method is highly inconsistent, with less than a quarter of calls in common with four of the six other methods for t-tests.
a notable similarity is that the dchip and tm methods have at least 90% agreement in all analyses.
complementing comparison of the statistical tests, we also compared the expression levels produced by the seven processing methods.
for each pair of methods, and for each pair of samples within one of the two data sets, we calculated pearson correlation coefficients of expression levels over all genes.
these values were summarized by taking the median over all pairs of samples within a data set, shown in table  <dig>  interestingly, methods calling similar genes as differentially expressed do not exhibit particularly strong correlation in expression levels.
for example, tm and dchip perform very similarly in terms of which genes are identified as significant, but the pairwise correlation between expression levels for these two methods is less than the average.
on the other hand, the tm and mas <dig> methods are generally at the extreme high and low ends of the sensitivity scale respectively, but their expression levels are the most strongly correlated of any pair of methods.
calibration
variation in fdr across the seven methods is due to two factors – variation in the number of transcripts with large test statistics, and variation in the expected number of transcripts with large test statistics when there is no real differential expression.
here we investigate the second factor, which is driven by the tendency of each method to produce outlier expression values.
the numerator of the fdr aims to correct for variation in the number of false positives, so that a method claiming large numbers of differentially expressed genes is not considered superior unless it also produces relatively small numbers of false positives.
this can be viewed as a calibration, in which for each method, the test statistic must reach a certain threshold in order that the proportion of false positives is no greater than a specified value.
calibration results are summarized in figure  <dig>  for each method, the threshold test statistic value required to obtain fdr <dig> less than f was calculated, and plotted against f. for example, to achieve any fdr <dig> value between  <dig>  and  <dig>  in the colon rank-sum data, gcrma-mle requires the lowest test statistics, rma requires a rank-sum statistic  <dig>  units larger than that of gcrma-mle, and mas <dig> requires a rank-sum statistic  <dig>  units larger than that of gcrma-mle.
since calibration depends only on randomized data, it should be possible to trace variation in thresholds across the processing methods to statistical properties of the expression levels.
for example, if one method produces expression levels with heavier tails, it is easier to get a large t-test statistic value by chance, particularly for the colon data with small sample sizes.
this would necessitate a higher threshold.
to quantify this, let  denote the log <dig> expression level of transcript i in sample j for method k, where k =  <dig>  ...,  <dig> denotes the seven processing methods, and let
where  is the pth quantile of , and med is the median value.
this is an affine-invariant measure of the size of the right tail of the expression values.
values of bk for the seven methods and two data sets are shown in table  <dig>  for reference, a gaussian distribution has a b value of  <dig>  when the sample size is as in the ovary data, and  <dig>  when the sample size is as in the colon data.
the gcrma-eb method is seen to have a much greater propensity for producing extreme expression values, explaining its low sensitivity, poor agreement with other methods, and conservative calibration.
variation in observed test statistics
in addition to calibration differences, fdr variation is also influenced by the observed test statistic values.
this is summarized in figure  <dig>  for each method, and for a range of test statistic values t, the number of probe sets for which the observed test statistic value exceeds t was calculated and plotted against t. for example, in the colon rank-sum data, pdnn had the smallest test statistics, with mas <dig> having around  <dig> more probe sets meeting a log test statistic threshold of  <dig> compared to pdnn.
the dchip and tm methods have over a thousand more probe sets meeting this threshold.
variation in test statistic values across the methods is greater in colon than in ovary data, and generally tracks with sensitivity.
however note that in the colon rank-sum data, dchip has substantially larger test statistics than gcrma-mle, even while gcrma-mle has better sensitivity , due to its less stringent calibration .
identification of genes with large fold changes
an interesting possibility that can not always be excluded is that the intergroup differences are so vast that nearly every gene is affected to a small degree.
if this were the case, the fdr values for the t-test and wilcoxon statistics would converge to zero for every gene as the number of samples grows, making fdr values difficult to interpret.
to further investigate this issue, we repeated the analysis using t-statistics truncated to zero when the fold change is less than  <dig>  as test statistics for fdr analysis.
the corresponding fdr values remain bounded away from zero for genes having true fold change smaller than  <dig> , while genes with true fold change exceeding  <dig>  have fdr values converging to  <dig>  thus the statistic identifies a meaningful subset of genes even when all genes are differentially expressed to some degree.
results for this analysis are shown in figure  <dig>  in the ovary data, the gcrma-eb method performs best, with gcrma-mle, mas <dig>  and tm slightly inferior.
several of the methods, specifically pdnn, dchip, and rma exhibit flat curves indicating that only a limited number of genes meet the 50% change criterion.
in the colon data, gcrma-mle and tm are nearly tied as the best performers.
overall, variation in sensitivity across the methods exists at a similar level to that found in the t-test and wilcoxon analyses.
only the gcrma-mle and tm methods give consistently good performances in the two data sets for this analysis.
strong outperformance
thus far we have focused on sensitivity as a criterion for comparing methods.
however even if one method is less sensitive than another, if the overlap in the called gene sets is not too great then the less sensitive method may still contribute to our understanding of which genes are differentially expressed.
suppose two methods denoted  <dig> and  <dig> give n <dig> and n <dig> genes respectively at a given fdr level.
then nk = ·nk estimates the expected number of truly differentially expressed genes called by method k. now suppose that i is the number of genes called by both methods.
then nk - i is an estimated lower bound for the expected number of genes correctly called by method k but not by the other method.
we will say that method  <dig> strongly outperforms method  <dig> if n <dig> - i ≥  <dig> but n <dig> - i <  <dig>  this means that in terms of differential expression, method  <dig> is not expected to contribute any true positives that were not called by method  <dig>
discussion
impact of processing method choice
the choice of processing method for affymetrix array data evidently has a major impact on the ability to confidently report the results of differential expression analysis.
the effect is greater, for example, than the choice of using a robust or a non-robust analysis, even in the colon data where robust analysis results in substantial improvements.
differences among processing methods are much greater in the more challenging colon data set compared to the ovary data, yet it should be noted that the sample sizes in the colon data are not atypical in real investigations.
while results from two data sets can never conclusively determine the optimal method, it is notable that across both data sets, using both t-statistic and rank-sum analyses, there is a high degree of similarity in the rank ordering of the methods from the best to the worst performer.
the trimmed mean  and dchip methods consistently perform as well or better than any of the other methods.
a possible explanation for this is that the weights used by the dchip may tend to downweight the least and greatest pm-mm differences, just as the tm method excludes these differences.
interpretation of fdr comparisons
when comparing array processing methods using experimental data in which the identities of differentially expressed genes are unknown, great care must be taken to ensure that apparent differences in sensitivity are not due to other factors.
one critical point is that the null distribution providing the expected number of false positives at a given test statistic threshold  must fairly reflect the statistical behavior of null genes.
permutation approaches have been extensively used to produce empirical p-values  and were used by efron et al.
<cit>  to estimate fdr values.
although permutation approaches are known to be slightly biased for estimating the fdr, the size of the bias  can not explain the magnitude of differences found here.
in addition, for a comparative analysis, as carried out here, it is more crucial that the biases be relatively constant across the methods.
however, since permutation approaches may not be highly accurate when the sample size is small, it is important to check performance on multiple data sets before conclusions about performance are drawn.
while we have focused on fdr as the basis of comparison, the pursuit of small fdr values is not the only desirable operating characteristic of an array processing method, and other reports have also emphasized the accuracy of estimating the precise size of concentration differences.
however to the extent that most actual studies seek to find differential expression between groups, the use of small fdr values seems more instrumental as the basis for judging methods.
variation due to choice of test statistic
although our primary aim was to investigate variation in sensitivity due to the seven processing methods, all analysis was carried out independently for two test statistics.
the t-statistic is widely used in practice, but is well-known to be sensitive to outliers, particularly when the sample size is small.
we found that certain processing methods, particularly eb-gcrma, had a tendency to produce outlier expression values in the colon data set.
thus the combination of using the eb-gcrma method with t-statistics in the colon data led to particularly poor performance.
variation due to log transform and array normalization
in practice, the approach used for array normalization and for forming log-transformed expression values may be equally or more influential than the method used for producing probe set summaries  <cit> .
in this study, we used implementations of the seven processing methods as prepared by their developers, and thus array normalization and and log-transforms were applied in a method-specific fashion.
this provides a comparative analysis of the various methods as they are used in practice, which is most directly relevant since few investigators will override the default normalization and log-transform methods provided by the developers of each method.
nevertheless it remains of interest whether these routine processing steps are the determining factor of performance.
in a future study it will be important to investigate this question further by modifying the implementations of the processing methods so that uniform log transforms and array normalizations are applied.
comparison of methods using data from disease profiling data sets
a key point that we advocate in this work is that false discovery rates in actual disease profiling data constitute a valuable complement to benchmarking results obtained from spike-in, dilution series, and mixture experiments .
the primary obstacle that must be overcome is that proper null sampling distributions are essential to ensure that the methods are compared on a common basis.
since numerous data sets covering a wide range of affymetrix platforms are available, to the extent that multiple data sets are in agreement about relative performances it is unlikely that the randomization procedure used to calculate fdr values is systematically biased against a particular method.
in spite of the statistical challenges in using disease profiling data for benchmarking, we argue that these data sets also offer some unique advantages.
calibration data sets are relatively few in number and are not available for all platforms.
newer platforms in particular are under-represented.
therefore overtraining to the available calibration data through manipulation of the many tuning parameters in the more complicated processing methods is an unavoidable concern.
in addition, the calibration data sets likely do not represent the same degree of challenge as disease profiling data in that reproducibility of fold changes for affected and unaffected genes is quite high compared to data from, say, human tissues where a large number of uncontrolled sources of variability are present.
CONCLUSIONS
performances of multiple array processing methods on disease profiling data sets vary widely across the seven methods studied here, but results are generally consistent between the two data sets studied.
results of our analysis generally do not parallel results obtained using calibration data sets  <cit> , suggesting that such comparisons may not completely capture the most relevant aspects of performance.
a major determinant of sensitivity is test statistic variability for randomized data.
such variability will affect false discovery rates as well as empirical p-values, which are an often-used alternative approach for identifying differentially expressed genes .
therefore it will be important in future work to seek a better understanding of statistical sampling properties of array processing methods.
a particular focus should be the way that sampling variance in probe masking and probe weighting is controlled.
methods seeking to incorporate mechanistic information about the dynamics of probe binding, such as the two gcrma methods and pdnn, should in principal outperform more generic approaches such as the tm method.
our results, particularly in the colon data, suggest that in medium-sized data sets this potential is not yet reached.
in this comparative analysis we did not seek to draw definitive conclusions about the "best" or "worst" methods.
such conclusions may be made after investigating a greater number of data sets, including disease profiling data, data from controlled experiments, and calibration data.
moreover, it may be that the correct choice of method may depend on the scientific question being asked.
the key message of this work is that the wide range of data sets collected in actual scientific investigations may be used for comparison of processing methods, and that in at least the two data sets considered here, similar results were obtained in the rank ordering of the methods.

 affymetrix microarrays are high throughput assays for measuring the expression levels of thousands of gene transcripts simultaneously. this type of microarray measures the expression of each transcript multiple times through a set of "probe pairs". since the advent of the affymetrix microarray, numerous methods have been proposed for producing numerical expression summaries for each transcript based on the probe pair data. several systematic studies have appeared comparing a number of methods on a common basis. these studies rely heavily on calibration data sets derived from spike-in, dilution series, and mixture experiments for comparing methods. our goal here was to carry out a comparative study of affymetrix array processing methods using data sets from typical biological experiments seeking differentially expressed genes in human tissue samples. the following seven methods are considered here: dchip, gcrma-eb and gcrma-mle, mas, pdnn, rma, and tm. while not every popular method is included in our study, several highly distinctive and original approaches are studied. for example, dchip was one of the first approaches to attempt to learn probe weights directly from the probe data, and rma pioneered the approach of disregarding the control mismatch probes. pdnn uses physical modeling to determine probe weights, while the two gcrma methods use gc content of the probe sequences to reduce variance in the mismatch probe levels. the mas method is the current default method provided by affymetrix. in addition to the six methods cited previously, we also include a method designated tm. this is a simple method that has been used in a number of published investigations, but has not been considered in any previous systematic comparison of affymetrix processing methods. to produce the probe-set summary score, the pm-mm differences are rank ordered, and the brightest 20% and dimmest 20% of values are deleted. the mean of the remaining values is used as the summary score. the scores for all probe-sets are then quantile normalized to a reference array using a piecewise linear spline with knots. an important feature of this study is the use of false discovery rate to quantify the sensitivity of a processing method in terms of its ability to distinguish differentially expressed genes from genes having invariant expression. this is a highly relevant property, as differential expression analysis is the most common application of microarray data. a key advantage of using fdr to compare processing methods is that fdr values can be calculated accurately using real disease profiling data where the identities of differentially expressed genes are uncertain. in contrast, most previous systematic comparisons of array processing methods have focused on calibration data sets in which concentrations of certain genes were experimentally manipulated. when it is highly likely that at least one gene is differentially expressed, false discovery rate may be defined as the expected ratio of the number of false positive calls to the total number of positive calls in a differential expression analysis between two groups of samples. if the groups are biologically distinct, a sensitive processing method should result in many genes with low fdr. thus to compare the performances of different array processing methods, we looked at two datasets in which a verified biological characteristic divided the samples into two classes, and compared the methods based on the number of genes having fdr smaller than various thresholds. for this to be a valid basis for comparison, the fdr values must be estimated with reasonable accuracy. following other recent work, we used a permutation approach for this estimation, arguing that there is no reason that this approach favors or disfavors any particular array processing method. a small fdr is due either to a small numerator, a large denominator, or both. the denominator of the fdr depends on the actual data distribution, so variation in this value may be due to factors such as accuracy in modeling the physical and chemical nature of probe binding. variation in the fdr numerator, however, depends only on the distribution of values produced for randomized data, a purely statistical quantity reflecting the tendency of the method to incorrectly produce test statistic outliers. our results suggest that both factors are important in determining sensitivity. the best methods produce many large test statistic values in the actual data, and also produce consistently small test statistic values for randomized data. poor performance of one method can be directly explained by the tendency of the method to produce outlier expression values, leading to greater numbers of incorrectly large test statistics. for overall comparison, we evaluated every pair of methods on the basis of whether the first method is expected to call at least one truly differentially expressed gene that is not also called by the second method. if this is not expected to occur, the second method is said to strongly outperform the first. based on this comparison, two of the methods considered are clearly favored, two are inferior, and results for the other three methods are mixed.  sensitivity differences our primary basis for comparison is sensitivity – the number of genes detected at a given fdr level, where fdr is a rescaled fdr. figure shows the key sensitivity results, using both the t-test statistic and the rank-sum statistic to assess differential expression. setting aside at first differences between the seven processing methods, we note two findings. first, in the colon data, analysis using the rank-sum statistic is substantially more sensitive than analysis using the t-test statistic. for the ovary data, where the sample sizes would not naturally suggest a robust analysis, there is no harm to sensitivity in using the rank-sum statistic. second, the ovary curves are substantially higher overall than the colon curves. this may be due to a greater number of true positives in the ovary data, or it may be that the small sample size for the msi group makes it difficult to attain high evidence levels for differential expression in the colon data. in any case, both data sets have many genes with small fdr values, supporting the biological relevance of the tumor groupings for both colon and ovary samples. the more challenging colon set distinguishes the seven processing methods to a greater extent than the ovary set. using fdr = as a reference point, there is roughly 7-fold variation across the seven methods in the number of detected genes in colon data using t-test statistics, while for rank-sum statistics the range is roughly 2-fold. in the ovary data, the range is around -fold for both statistics. also notable is that variation in sensitivity due to the choice of test statistic is smaller than variation in sensitivity due to the seven processing methods. no single method stands out as having the best or worst performance in every case. however some methods generally perform better than others. the dchip and tm methods perform consistently well, while the gcrma-eb and mas methods consistently perform poorly. pdnn performs well on the ovary data, but poorly on the colon data, and results for the other methods are mixed. level of agreement between methods identities of probe sets falling below a given fdr threshold vary across the methods. figure summarizes this variation. the ratio of the number of probe sets falling below various fdr thresholds in k or more of the seven methods to the number of probe sets falling below the threshold for at least one method is plotted against the fdr threshold, for k = in the ovary data there is a very high level of agreement in this measure. for the rank-sum analysis, almost 90% of called genes are called by at least four methods, and more than 70% of called genes are called by all seven methods. for the t-test analysis, the agreement is slightly higher yet. for the colon data, the methods are much more inconsistent. for the rank-sum analysis, three of the methods agree on up to 90% of genes, but all seven methods only agree on around 30% of genes. the t-test analysis is even worse, with only around 10% of genes common to all seven methods. turning to pairwise agreement, table shows the percentage of genes called by both members of a pair of methods out of the genes called by at least one of the two methods. in the ovary data, mas shares the fewest calls with the other methods for both t-test and rank-sum analysis, while gcrma-eb has relatively weak agreement for the t-test analysis. in the colon data, the gcrma-eb method is highly inconsistent, with less than a quarter of calls in common with four of the six other methods for t-tests. a notable similarity is that the dchip and tm methods have at least 90% agreement in all analyses. complementing comparison of the statistical tests, we also compared the expression levels produced by the seven processing methods. for each pair of methods, and for each pair of samples within one of the two data sets, we calculated pearson correlation coefficients of expression levels over all genes. these values were summarized by taking the median over all pairs of samples within a data set, shown in table interestingly, methods calling similar genes as differentially expressed do not exhibit particularly strong correlation in expression levels. for example, tm and dchip perform very similarly in terms of which genes are identified as significant, but the pairwise correlation between expression levels for these two methods is less than the average. on the other hand, the tm and mas methods are generally at the extreme high and low ends of the sensitivity scale respectively, but their expression levels are the most strongly correlated of any pair of methods. calibration variation in fdr across the seven methods is due to two factors – variation in the number of transcripts with large test statistics, and variation in the expected number of transcripts with large test statistics when there is no real differential expression. here we investigate the second factor, which is driven by the tendency of each method to produce outlier expression values. the numerator of the fdr aims to correct for variation in the number of false positives, so that a method claiming large numbers of differentially expressed genes is not considered superior unless it also produces relatively small numbers of false positives. this can be viewed as a calibration, in which for each method, the test statistic must reach a certain threshold in order that the proportion of false positives is no greater than a specified value. calibration results are summarized in figure for each method, the threshold test statistic value required to obtain fdr less than f was calculated, and plotted against f. for example, to achieve any fdr value between and in the colon rank-sum data, gcrma-mle requires the lowest test statistics, rma requires a rank-sum statistic units larger than that of gcrma-mle, and mas requires a rank-sum statistic units larger than that of gcrma-mle. since calibration depends only on randomized data, it should be possible to trace variation in thresholds across the processing methods to statistical properties of the expression levels. for example, if one method produces expression levels with heavier tails, it is easier to get a large t-test statistic value by chance, particularly for the colon data with small sample sizes. this would necessitate a higher threshold. to quantify this, let denote the log expression level of transcript i in sample j for method k, where k =..., denotes the seven processing methods, and let where is the pth quantile of, and med is the median value. this is an affine-invariant measure of the size of the right tail of the expression values. values of bk for the seven methods and two data sets are shown in table for reference, a gaussian distribution has a b value of when the sample size is as in the ovary data, and when the sample size is as in the colon data. the gcrma-eb method is seen to have a much greater propensity for producing extreme expression values, explaining its low sensitivity, poor agreement with other methods, and conservative calibration. variation in observed test statistics in addition to calibration differences, fdr variation is also influenced by the observed test statistic values. this is summarized in figure for each method, and for a range of test statistic values t, the number of probe sets for which the observed test statistic value exceeds t was calculated and plotted against t. for example, in the colon rank-sum data, pdnn had the smallest test statistics, with mas having around more probe sets meeting a log test statistic threshold of compared to pdnn. the dchip and tm methods have over a thousand more probe sets meeting this threshold. variation in test statistic values across the methods is greater in colon than in ovary data, and generally tracks with sensitivity. however note that in the colon rank-sum data, dchip has substantially larger test statistics than gcrma-mle, even while gcrma-mle has better sensitivity, due to its less stringent calibration. identification of genes with large fold changes an interesting possibility that can not always be excluded is that the intergroup differences are so vast that nearly every gene is affected to a small degree. if this were the case, the fdr values for the t-test and wilcoxon statistics would converge to zero for every gene as the number of samples grows, making fdr values difficult to interpret. to further investigate this issue, we repeated the analysis using t-statistics truncated to zero when the fold change is less than as test statistics for fdr analysis. the corresponding fdr values remain bounded away from zero for genes having true fold change smaller than, while genes with true fold change exceeding have fdr values converging to thus the statistic identifies a meaningful subset of genes even when all genes are differentially expressed to some degree. results for this analysis are shown in figure in the ovary data, the gcrma-eb method performs best, with gcrma-mle, mas and tm slightly inferior. several of the methods, specifically pdnn, dchip, and rma exhibit flat curves indicating that only a limited number of genes meet the 50% change criterion. in the colon data, gcrma-mle and tm are nearly tied as the best performers. overall, variation in sensitivity across the methods exists at a similar level to that found in the t-test and wilcoxon analyses. only the gcrma-mle and tm methods give consistently good performances in the two data sets for this analysis. strong outperformance thus far we have focused on sensitivity as a criterion for comparing methods. however even if one method is less sensitive than another, if the overlap in the called gene sets is not too great then the less sensitive method may still contribute to our understanding of which genes are differentially expressed. suppose two methods denoted and give n and n genes respectively at a given fdr level. then nk = ·nk estimates the expected number of truly differentially expressed genes called by method k. now suppose that i is the number of genes called by both methods. then nk - i is an estimated lower bound for the expected number of genes correctly called by method k but not by the other method. we will say that method strongly outperforms method if n - i ≥ but n - i < this means that in terms of differential expression, method is not expected to contribute any true positives that were not called by method discussion impact of processing method choice the choice of processing method for affymetrix array data evidently has a major impact on the ability to confidently report the results of differential expression analysis. the effect is greater, for example, than the choice of using a robust or a non-robust analysis, even in the colon data where robust analysis results in substantial improvements. differences among processing methods are much greater in the more challenging colon data set compared to the ovary data, yet it should be noted that the sample sizes in the colon data are not atypical in real investigations. while results from two data sets can never conclusively determine the optimal method, it is notable that across both data sets, using both t-statistic and rank-sum analyses, there is a high degree of similarity in the rank ordering of the methods from the best to the worst performer. the trimmed mean and dchip methods consistently perform as well or better than any of the other methods. a possible explanation for this is that the weights used by the dchip may tend to downweight the least and greatest pm-mm differences, just as the tm method excludes these differences. interpretation of fdr comparisons when comparing array processing methods using experimental data in which the identities of differentially expressed genes are unknown, great care must be taken to ensure that apparent differences in sensitivity are not due to other factors. one critical point is that the null distribution providing the expected number of false positives at a given test statistic threshold must fairly reflect the statistical behavior of null genes. permutation approaches have been extensively used to produce empirical p-values and were used by efron et al. to estimate fdr values. although permutation approaches are known to be slightly biased for estimating the fdr, the size of the bias can not explain the magnitude of differences found here. in addition, for a comparative analysis, as carried out here, it is more crucial that the biases be relatively constant across the methods. however, since permutation approaches may not be highly accurate when the sample size is small, it is important to check performance on multiple data sets before conclusions about performance are drawn. while we have focused on fdr as the basis of comparison, the pursuit of small fdr values is not the only desirable operating characteristic of an array processing method, and other reports have also emphasized the accuracy of estimating the precise size of concentration differences. however to the extent that most actual studies seek to find differential expression between groups, the use of small fdr values seems more instrumental as the basis for judging methods. variation due to choice of test statistic although our primary aim was to investigate variation in sensitivity due to the seven processing methods, all analysis was carried out independently for two test statistics. the t-statistic is widely used in practice, but is well-known to be sensitive to outliers, particularly when the sample size is small. we found that certain processing methods, particularly eb-gcrma, had a tendency to produce outlier expression values in the colon data set. thus the combination of using the eb-gcrma method with t-statistics in the colon data led to particularly poor performance. variation due to log transform and array normalization in practice, the approach used for array normalization and for forming log-transformed expression values may be equally or more influential than the method used for producing probe set summaries. in this study, we used implementations of the seven processing methods as prepared by their developers, and thus array normalization and and log-transforms were applied in a method-specific fashion. this provides a comparative analysis of the various methods as they are used in practice, which is most directly relevant since few investigators will override the default normalization and log-transform methods provided by the developers of each method. nevertheless it remains of interest whether these routine processing steps are the determining factor of performance. in a future study it will be important to investigate this question further by modifying the implementations of the processing methods so that uniform log transforms and array normalizations are applied. comparison of methods using data from disease profiling data sets a key point that we advocate in this work is that false discovery rates in actual disease profiling data constitute a valuable complement to benchmarking results obtained from spike-in, dilution series, and mixture experiments. the primary obstacle that must be overcome is that proper null sampling distributions are essential to ensure that the methods are compared on a common basis. since numerous data sets covering a wide range of affymetrix platforms are available, to the extent that multiple data sets are in agreement about relative performances it is unlikely that the randomization procedure used to calculate fdr values is systematically biased against a particular method. in spite of the statistical challenges in using disease profiling data for benchmarking, we argue that these data sets also offer some unique advantages. calibration data sets are relatively few in number and are not available for all platforms. newer platforms in particular are under-represented. therefore overtraining to the available calibration data through manipulation of the many tuning parameters in the more complicated processing methods is an unavoidable concern. in addition, the calibration data sets likely do not represent the same degree of challenge as disease profiling data in that reproducibility of fold changes for affected and unaffected genes is quite high compared to data from, say, human tissues where a large number of uncontrolled sources of variability are present.  performances of multiple array processing methods on disease profiling data sets vary widely across the seven methods studied here, but results are generally consistent between the two data sets studied. results of our analysis generally do not parallel results obtained using calibration data sets, suggesting that such comparisons may not completely capture the most relevant aspects of performance. a major determinant of sensitivity is test statistic variability for randomized data. such variability will affect false discovery rates as well as empirical p-values, which are an often-used alternative approach for identifying differentially expressed genes. therefore it will be important in future work to seek a better understanding of statistical sampling properties of array processing methods. a particular focus should be the way that sampling variance in probe masking and probe weighting is controlled. methods seeking to incorporate mechanistic information about the dynamics of probe binding, such as the two gcrma methods and pdnn, should in principal outperform more generic approaches such as the tm method. our results, particularly in the colon data, suggest that in medium-sized data sets this potential is not yet reached. in this comparative analysis we did not seek to draw definitive conclusions about the "best" or "worst" methods. such conclusions may be made after investigating a greater number of data sets, including disease profiling data, data from controlled experiments, and calibration data. moreover, it may be that the correct choice of method may depend on the scientific question being asked. the key message of this work is that the wide range of data sets collected in actual scientific investigations may be used for comparison of processing methods, and that in at least the two data sets considered here, similar results were obtained in the rank ordering of the methods. 
performance of array processing methods in disease profiling and other realistic biological studies should be given greater consideration when comparing affymetrix processing methods. an emphasis is placed on understanding sensitivity for detecting differentially expressed genes in terms of two key statistical determinants: test statistic variability for non-differentially expressed genes, and test statistic size for truly differentially expressed genes. previously reported performance analyses using tightly controlled calibration data sets are not highly consistent with results reported here using data from human tissue samples. poor performance of one method was directly tied to a tendency to produce highly variable test statistic values under randomization. a critical step in processing oligonucleotide microarray data is combining the information in multiple probes to produce a single number that best captures the expression level of a rna transcript. choice of processing method has a major impact on differential expression analysis of microarray data. the best performing methods called up to 90% of the same genes differentially expressed, had less variable test statistics under randomization, and had a greater number of large test statistics in the experimental data. several systematic studies comparing multiple methods for array processing have used tightly controlled calibration data sets as the basis for comparison. here we compare performances for seven processing methods using two data sets originally collected for disease profiling studies. in the two data sets considered here, up to seven-fold variation across the processing methods was found in the number of genes detected at a given false discovery rate. two other methods are inferior, while results for the other three methods are mixed.
 a critical step in processing oligonucleotide microarray data is combining the information in multiple probes to produce a single number that best captures the expression level of a rna transcript. several systematic studies comparing multiple methods for array processing have used tightly controlled calibration data sets as the basis for comparison. here we compare performances for seven processing methods using two data sets originally collected for disease profiling studies. an emphasis is placed on understanding sensitivity for detecting differentially expressed genes in terms of two key statistical determinants: test statistic variability for non-differentially expressed genes, and test statistic size for truly differentially expressed genes.  in the two data sets considered here, up to seven-fold variation across the processing methods was found in the number of genes detected at a given false discovery rate. the best performing methods called up to 90% of the same genes differentially expressed, had less variable test statistics under randomization, and had a greater number of large test statistics in the experimental data. poor performance of one method was directly tied to a tendency to produce highly variable test statistic values under randomization. based on an overall measure of performance, two of the seven methods are superior in the two data sets considered here. two other methods are inferior, while results for the other three methods are mixed.  choice of processing method has a major impact on differential expression analysis of microarray data. previously reported performance analyses using tightly controlled calibration data sets are not highly consistent with results reported here using data from human tissue samples. performance of array processing methods in disease profiling and other realistic biological studies should be given greater consideration when comparing affymetrix processing methods. 