BACKGROUND
the rapid increase in the number of completely sequenced genomes in the past few years has generated much effort in functional genomics, particularly studies seeking to assign biological functions to dna sequences.
comparative gene expression profiling is widely used to study the functional role of genes.
the dna microarray assay provides an invaluable technique for large scale expression analysis.
in the two-channel dna microarray assay, rna from two samples is reverse transcribed to cdna and labeled with two distinct fluorescent dyes before being co-hybridized to immobilized dna strands on a microarray slide.
spotted arrays currently being used can be divided into two groups based on the nature of immobilized dna used: one in which the immobilized dna is comprised of both sense and antisense strands  and the other where the immobilized dna is single stranded consisting of only the sense strands .
during hybridization, the two fluorescently labeled cdna samples compete for hybridization to the immobilized strands.
hybridization reactions between complementary strands occur only between the labeled antisense strand and immobilized sense strand.
the ratio of the intensities of the two fluorescently labeled cdnas is used to quantify the relative levels of transcripts in the two samples  <cit> .
this method serves well for pair-wise comparison of transcript levels in two samples.
with over ten thousand different dna species immobilized on the microarray, the relative transcription level of all the corresponding genes in the two samples can be obtained in a single assay.
dna microarrays have found applications in gene discovery, disease diagnosis, pharmacogenomics and toxicology research.
they are increasingly used for a series of related samples, for which a comparison across all samples and all genes is desirable.
when a large number of samples are to be compared, a combinatorial approach pairing all possible pairs  is often taken.
this results in a large number of microarrays, requiring a large amount of each rna sample.
a 'loop design', where every sample is directly compared to two other samples to form a closed loop, has been proposed to overcome this problem  <cit> .
the ratios calculated using a loop design have variable levels of precision since some samples are more directly related than others  <cit> .
when a new sample is to be inserted into the earlier 'loop', rna for at least two of the previous samples is needed to pair with the new sample to form a new node in the closed loop.
another approach to tackle the issue of a large combinatorial pair-wise comparison is the 'reference design'  <cit>  in which a common reference sample is introduced with which all rna samples are hybridized.
two possible universal references are rna pooled from various samples and genomic dna  <cit> .
for a given set of samples, pooled rna provides an excellent reference.
however, if the experimental conditions change, the possibility arises that some new transcripts may not be represented in the initially-pooled rna.
genomic dna is an attractive candidate for use as a common reference.
it is isolated from cells or tissue and sheared to fragments in a narrow range of length.
it is easier to prepare, maintain and reproduce, as compared to rna.
it is especially useful for microorganisms, which lack repetitive sequences in their genome, and microarrays using genomic dna as a reference have been used to identify genes differentially expressed in various growth stages of mycobacterium tuberculosis  <cit> .
it has also been recently reported that the data obtained using genomic dna as a reference in microarray experiments with arabidopsis thaliana employing 70-mer oligonucleotide microarrays was in agreement with ratios obtained from direct hybridizations  <cit> .
genomic dna samples, isolated from stationary phase cultures where the chromosome is not being replicated, have the same representation of all genes as in the genome.
since the transcript level of each gene is being referenced to its own representation on the genome , the relative expression can be compared across different genes in the sample  as well as across different samples.
the use of microarrays using genomic dna in a range of applications including genomic diversity studies  <cit>  and aneuploidy detection using comparative genomic hybridization  <cit>  has been demonstrated.
in the conventional two-channel cdna hybridization  <cit> , both the cdna samples are antisense strands.
the probability of hybridization between strands in solution is very low.
on the other hand, using sheared genomic dna as a reference, hybridization between complementary sense and anti-sense strands can occur in solution between the complementary genomic dna strands, and between genomic dna sense strands and their cdna counterparts.
the number of strands lost to hybridization in the solution phase may differ for different rna species as well as for the two complementary strands of the same species.
this may result in decreased fidelity in the ratio of cdna to genomic dna as a representation of gene expression level.
with the complexity of hybridization in both the solution phase and the immobilized surface phase, and between double strands of genomic dna and single strands of cdna, it is difficult to assess the effect of using genomic dna as a common reference.
adapting a mathematical model we developed previously to assess diffusional constraints on dna microarray assay  <cit> , we have constructed a kinetic model for microarray hybridization and predicted the effectiveness as well as potential pitfalls of using genomic dna as a reference.
we examined the effect of the various parameters that may affect surface hybridization.
here we report the framework of the model and our findings.
mathematical model for microarray hybridization
a dna microarray consists of thousands of spots, each spot containing dna strands of known sequence immobilized on an impermeable surface.
we simulate hybridization reactions for one spot on the microarray, which is considered to be a sample chamber with two compartments; a spot phase and a solution phase .
the spot phase is the small volume in which the bound strands are assumed to be present at a uniform concentration.
the solution phase, which comprises a vast majority of the volume of the microarray chamber, contains only fluorescently labeled single strands at the beginning of the hybridization and double strands formed by hybridization of complementary mobile single strands during the course of the hybridization.
five kinds of single-stranded dna molecules are present in the system: labeled single stranded cdna  reverse transcribed from rna, denoted as s; the two identically labeled complementary  strands from genomic dna , denoted as r and r'; and the non-labeled anti-sense and sense strands immobilized on the array surface  denoted by b and b'.
the nomenclature used in this study is summarised in table  <dig>  the.
the double-stranded species are denoted by combining the constituent single strand symbols, for example rb' denotes the labeled complex formed by hybridization of genomic dna antisense strands to the sense bound strands.
we assume that each phase is well-mixed, and hence all mobile species are present at uniform concentrations within each phase.
we consider the case that pcr products are used for immobilization; thus, both sense and antisense strands of probe dna are immobilized.
we assume that the two complementary bound strands do not hybridize to each other  <cit>  due to the immobilization procedure used and neglect the formation of bb'.
the mobile species in the two compartments are considered to move across the phase boundary at a rate proportional to the difference in the concentrations of identical species present in the two phases.
the proportionality constant is the effective mass transfer coefficient for the transport of mobile dna strands and is estimated from the diffusivity as discussed below.
the model equations take the form of a mass balance on each component in each phase that accounts for the change in concentration due to reactions within that phase and transport between the two phases.
non-specific hybridization in both phases is neglected.
all these equations are of the form:
where  = r',  = b', a <dig> = r, a <dig> = b, a <dig> = s, subscript p denotes either solution or spot phase and p' is the other phase; kf is the forward rate of hybridization; kb is the backward rate of hybridization; kt is the rate of transport across the spot-phase solution-phase boundary; vt is the total volume of the sample chamber and vp is the volume of phase p.
model parameters
model geometry
the diameter of the spot phase is set to  <dig>  cm as seen in a typical microarray spotted on poly-lysine slides  <cit> , with a height of  <dig> × 10- <dig> cm .
this gives a spot volume of  <dig>  × 10- <dig> l. typically, labeled strands are resuspended in a volume of  <dig>  × 10- <dig> l, and this volume is applied under a cover slip of area  <dig>  cm ×  <dig>  cm.
for this geometry, the height of the microarray chamber is about  <dig>  × 10- <dig> cm.
the solution phase volume is the difference between the sample chamber volume  and the spot volume .
it has been shown before that only the mobile strands within a radius of  <dig>  cm from the spot are available for hybridization due to transport effects  <cit> .
hence, we have assumed the solution phase to have a diameter of  <dig>  cm and a height of  <dig>  × 10- <dig> cm corresponding to the height of the microarray chamber.
hybridization rate constants
the hybridization rate constants are assumed to be identical for dna strands from both genomic dna and cdna.
however, the rate constants for reaction between two mobile strands may differ from the rate constant for the reaction of a mobile strand and a bound strand.
a forward rate constant  for dna hybridizations in solution of  <dig> m- <dig> s- <dig>  <cit>  was used for simulations.
the rate constant of hybridization could be slower for hybridization of mobile strands to immobilized strands  on a solid surface and values in a range up to 100-fold lower than kf have been used in simulations.
we also discuss the effect of this reduced forward rate constant of hybridization between mobile and bound strands for a range between  <dig> m-1s- <dig> and  <dig> m-1s- <dig>  the backward rate constant for dissociation of dna double strands as calculated from equilibrium constants reported in literature ranges from 10- <dig> s- <dig> to 10- <dig> s- <dig>  <cit> .
simulations were carried out using a backward rate constant ranging from 10- <dig> s- <dig> to  <dig> s- <dig> .
in all simulations performed, the deviation from the true value was greatest when the hybridization reaction is set to be irreversible.
the results discussed below are for the case with the backward rate constant set to  <dig> s- <dig>  the results obtained when the reaction is set to be reversible are very close to the true value under all conditions tested.
rate constant for transport
the rate constant for transport  is estimated from the diffusion coefficient as , where d is the diffusion coefficient, h is the height of the spot phase and r is the radius of the spot phase.
the diffusivity of dna single strands in solution has been estimated to be 10- <dig> cm2/s  <cit>  which leads to an estimate of  <dig> s- <dig> for kt.
the simulations reported in the next section were carried out under both very fast transport  and very slow transport  conditions.
the transport rate has no effect on the performance of the assay in the range tested .
initial single-strand concentrations
from the yeast transcriptome data published by velculescu et al.
<cit> , the mass percentages of mrna that belong to the rare, intermediate and abundant classes are estimated to be  <dig> ,  <dig>  and  <dig> % respectively.
for an mrna sample of  <dig>  μg used in each microarray assay, the total amount of rare, intermediate and abundant species are  <dig>  μg,  <dig>  μg and  <dig>  × 10- <dig> μg respectively.
to convert those numbers to molar concentration we calculate the number of genes within each abundance class using the intensity data from microarray experiments of s. coelicolor.
the s. coelicolor transcripts were divided into three abundance classes using the following intensity cut-offs: intensity <  <dig> as rare,  <dig> < intensity <  <dig> as intermediate and > <dig> as abundant sequences.
normalized intensity values were used for this estimate.
<dig> %, 26% and  <dig> % of all genes were classified rare, intermediate and abundant respectively.
the corresponding number of rare, intermediate and abundant species are  <dig>   <dig> and  <dig> respectively.
assuming a sample volume of  <dig>  × 10- <dig> l, 100% reverse transcription efficiency, and an average strand length of  <dig> we calculate the rare species concentration as ~ <dig>  pm, intermediate species concentration ~ <dig> pm and abundant species concentration ~ <dig> pm.
previous experimental reports of the use of genomic dna as a reference for microarray hybridizations used genomic dna concentrations ranging from  <dig>  μg to  <dig> μg for m. tuberculosis  <cit> .
since m. tuberculosis has ~ <dig> genes, this translates to a concentration of each gene from ~ <dig> pm to  <dig> pm, which is the range of genomic dna concentrations used for simulations.
to calculate the concentration of the immobilized species, we assume the concentration of dna in the spotting solution to be  <dig>  g/l and that 10- <dig> l of the solution is spotted on the microarray.
we also assume that 75% of the dna thus spotted is washed away in the microarray post-processing steps.
<dig>  × 10- <dig> g dna remaining on the microarray is uniformly distributed in the  <dig>  × 10- <dig> l spot phase volume.
assuming an average dna strand length of  <dig> bp, this is approximately equal to an immobilized strand concentration in the spot phase of 10- <dig> m.
RESULTS
the mathematical model described above takes into account dna hybridization between single stranded cdna and double stranded genomic dna  in solution and immobilized double strands on a microarray surface.
this model considers hybridization only on one spot on the microarray.
the immobilized strands are distributed uniformly in the spot phase and the mobile strands are present both in the solution phase and spot phase.
hybridization between mobile and bound species occurs in the spot phase.
all concentrations described in the following sections are the concentrations in the spot phase.
the fluorescence intensity corresponding to hybridized cdna sample is expressed as is =  and the channel corresponding to the hybridized genomic dna reference as ir =  + .
the result obtained for hybridization with genomic dna used as a reference is a hybridization ratio γ = is/ir =  / .
since the concentration of all genes in a genomic dna sample is equal, the ratio  for different genes is an indication of the relative abundance of the transcript for those genes.
typically when genomic dna is used as the reference, the ratio  from one hybridization of cdna derived from sample one, is compared to another ratio  from sample two, to obtain the relative expression level of the transcripts in samples  <dig> and  <dig>  ideally this "ratio of ratios", i.e.
γ1/γ <dig>  should be equal to the ratio of the initial concentrations of the transcript in those samples.
we simulate this experimental process with two different initial cdna concentrations  and the same genomic dna concentration to obtain the ratios γ <dig> and γ <dig>  the accuracy of the microarray assay is quantified by the accuracy index α,
a value of one for α corresponds to a perfect assay, where the measured relative concentration of the transcript in the two samples  is exactly equal to the true relative concentration .
any deviation in α from unity is a measure of the error of the assay.
these results are also applicable to the comparison of expression levels of two genes in one cdna sample, as the model makes no distinction between hybridization on two spots on one microarray and hybridization to a spot corresponding to the same sequence in two different microarray experiments.
we systematically vary model parameters to investigate the effect of different hybridization conditions, transcript abundance levels, and degree of differential expression on the performance of the microarray assay.
effectiveness of using genomic dna as a reference
the effectiveness of using genomic dna as a reference in the microarray assay was predicted by simulations using the model and parameters described above.
figure 2a shows the variation in α with hybridization time for different rna abundance levels and differential expression ratios corresponding to biologically realistic scenarios as listed in table  <dig>  these cases represent a wide range of possible combinations of the three abundance classes and differential expression ratios .
a differential expression ratio of  <dig> is used as an example of a small change in expression level and  <dig> as an example of a large change in expression level.
the parameter values used are kf =  <dig> m-1s- <dig>  kf-bound =  <dig> m-1s- <dig>  kb =  <dig> s- <dig>  initial genomic dna concentration =  <dig> pm, initial bound strand concentration = 10- <dig> m, kt =  <dig> s- <dig>  for a wide range of initial cdna concentrations and differential expression ratios, the accuracy index α is within 5% from unity , indicating that the assay performance is robust to the initial concentration of single strands.
the only condition where α is significantly different from  <dig>  is when intermediate species are upregulated  <dig> fold, a situation not likely to happen frequently in cells.
furthermore, in most microarray assays, a 17% deviation from the true value is not considered large.
the model simulation predicts that for abundant species, the accuracy decreases with hybridization time, indicating short hybridization time will lead to better results.
however, the concentrations of the double strands formed as a result of hybridization of single strands in solution to the immobilized strands  increase with time for approximately the first  <dig> hours.
the fluorescent intensities detected when the microarray slide is scanned is proportional to the concentration of these double stranded species.
hence the hybridization time has to be long enough to obtain intensities sufficiently above background levels for accurate measurement.
the model prediction that using genomic dna as a reference can provide an accurate measurement in a microarray assay was verified experimentally.
the transcript levels of two samples were assayed using both direct cdna: cdna hybridization  as well as by hybridizing to genomic dna .
this ratio of the two samples  was obtained by dividing the ratios obtained from those two cdna samples /   microarray data can be found in the supplementary material.
figure  <dig> shows a scatter plot of the relative transcript level obtained by these two methods.
the ratio obtained from indirect comparison is within  <dig>  fold of that obtained from direct comparison for 91% of genes.
out of the remaining 9% genes, 81% have a ratio obtained from direct comparison within  <dig> -fold.
for  <dig> % of all the genes, the ratio obtained from indirect comparison is within  <dig> fold of that obtained from direct comparison.
70% of the remaining  <dig> % have a ratio obtained from direct comparison within  <dig> -fold.
thus, in general, the ratio obtained from hybridizations using genomic dna as a reference is consistent with those obtained from direct cdna: cdna hybridization.
as can be seen in figure  <dig>  this is true over a large range of differential expression .
also, since total rna samples containing a wide range of transcript abundance levels were used in this experiment, the dataset demonstrates that the accuracy of the assay is maintained over all mrna expression levels.
effect of genomic dna concentration
in previous investigations several concentrations of the genomic dna ranging from  <dig> to  <dig> pm have been used  <cit> .
to examine the effect of genomic dna concentrations on the microarray assay, we carried out simulations for a range of genomic dna  and cdna concentrations  with different degrees of downregulation .
α is within 6% of unity for all values of genomic dna concentrations examined .
thus, the assessment of the ratio between the transcript levels of a gene in two samples does not vary significantly with the genomic dna concentration used in the experiment.
however, in a microarray experiment, accurate assessment of the ratio of expression levels can be attained only if the fluorescent intensities can be accurately measured.
spots with lower intensities and closer to the background level are prone to increased noise interference and decreased accuracy.
at higher concentration of genomic dna, the intensity of the sample channel decreases.
this effect is more profound for rare species, compounding the problem of their low intensity.
this leads to the need to carefully select the genomic dna concentration for achieving intensities which are significantly above background levels.
effect of bound strand concentration
in the preparation of microarray slides the amount of dna immobilized on the slide varies.
many factors, including differences in the amount of dna deposited, spot morphology, dna retention in the spotting and post-processing procedure contribute to such variation.
in cdna:cdna hybridization, the ratio of transcript levels from two different samples is not affected by the amount of dna immobilized on the slides, since the same transcript species from both samples are affected to the same extent.
in contrast, when cdna:gdna is used, variation in the amount of immobilized dna may exert a different effect on different samples.
we simulated the effect on α of this changing bound strand concentration by using bound strand concentrations of 10- <dig> and 10- <dig> m for all cases listed in table  <dig>  this range spans values from 100-fold higher to 100-fold lower than the estimated bound strand concentration.
all other parameters are same as those used for the simulations in figure  <dig>  the findings are similar to the plot shown in figure 2a .
thus, the variation in the amount of dna immobilized on slides does not significantly affect the accuracy of the assay.
however, a lower bound strand concentration does result in reduced concentrations of the double strands formed, and hence lower intensities on the microarray.
as an illustration, when the bound strand concentration is lowered from 10- <dig> m to 10- <dig> m, the intensity of the sample channel  decreases 24-fold for rare species, 25-fold for intermediate species and 38-fold for abundant species.
accuracy affected by transcript abundance at low hybridization rate
an uncertainty in the simulation is the effect of reaction rate constants; especially, the effect of change in the relative magnitude of forward rate constants for hybridization between mobile and mobile, and mobile and bound strands, is further investigated.
to evaluate this, kf-bound was varied 100-fold from  <dig> m-1s- <dig> to  <dig> m-1s- <dig> and the results are shown in figure  <dig>  for all values of kf-bound, α is within 21% of unity for all the rna abundance classes and all differential expression ratios examined.
in general, the assay is more accurate for rare and intermediate species  as compared to abundant species and for lower differential expression ratios compared to higher ratios.
for the abundant species, as the forward rate constant for hybridization between mobile and bound species  is reduced, the accuracy of the assay decreases.
the model predicts that the highest error  will be observed for abundant species with high differential expression ratios.
the simulation results presented above illustrate that the accuracy index α is most sensitive to the rate constant of the forward reaction for hybridization between mobile and bound species  and the deviation of α from unity is highest for abundant species in the sample.
the ratio determined in the hybridization assay is .
ideally, this should reflect the true ratio .
the ratio /  is always close to s0/r <dig> since both s and r hybridize with the same species: r' in solution and b' on the surface.
this is because in our analysis, we assume that s and r have identical reaction kinetics, they are thus both stoichiometrically and kinetically indistinguishable from each other.
therefore, the deviation in /  from r0/r0'  governs the deviation in γ from the true ratio.
since the reaction is irreversible, for  to be equal to , the product   should be equal to  .
let εs be the amount of s reacted with r' and εr be the amount of r reacted with r'.
since the stoichiometric ratios for reactions of s with r' and r with r' are both  <dig>  r' =  -εr - εs and r = r <dig> - εr.
under conditions where s <dig> is low , εs is relatively small and r ≈ r'.
however, for abundant species, εs is significant and as a result, r/r' is greater than r0/.
also, since
s+b' → sb'
r+b' → rb'
r'+b → r'b
b' gets consumed more than b. again, for rare and intermediate abundance species, this difference in b and b' is insignificant since the bound strands are in excess and hence the assay works very well for rare and intermediate species.
the difference between b and b' is significant for abundant species  for which b/b' is greater than b0/b0'.
the exact deviation in b/b' and r/r' depends on the rate constants for hybridization and hence it is not surprising that the rate constants have a major effect on hybridization results.
this is illustrated by values for b, b', r and r' for s <dig> =  <dig> pm, r <dig> =  =  <dig> pm and kf-bound =  <dig> m- <dig> s- <dig> and  <dig> m- <dig> s- <dig>  at the end of  <dig> hours of hybridization time, presented in table  <dig>  the products  ,   and ratios /  and /  are also calculated.
thus, for the larger value of kf-bound,   ≈  , which is not the case for the lower values of kf-bound.
this means that factors, such as mixing conditions in the sample chamber, which affect reaction kinetics, will affect the accuracy of the assay.
one implication of the above discussion is that if single stranded species are used for immobilization , as are used in spotted oligo-arrays, the assay will be robust to a wide range of hybridization conditions.
this is also seen from our simulations where for all the conditions discussed in this paper, α is equal to one if only sense strands are immobilized on the microarray slide, a situation encountered during the use of oligonucleotide spotted arrays.
CONCLUSIONS
the use of genomic dna as a reference is useful to assess the expression levels of a large array of genes among different samples.
we have developed a kinetic model to predict the effect of using genomic dna in the microarray assay under a wide range of conditions.
the model predicts that the assay can accurately estimate the relative concentrations for a wide range of initial cdna concentrations and ratios from hybridizations using genomic dna as a reference will correspond to true ratios under all conditions if single stranded oligonucleotide microarrays are used.
the model also serves as a useful tool to predict the performance of such assays under varying conditions that are otherwise difficult to carry out experimentally.
despite a number of publications on its application, the use of genomic dna as a reference for microarray assay is still not wide spread.
we carried out this study to verify on a theoretical basis the validity of this approach and the results are indeed reassuring.
we expect that the use of genomic dna as reference will accelerate especially for comparative transcriptome analysis involving a wide range of samples from different sources.

 the rapid increase in the number of completely sequenced genomes in the past few years has generated much effort in functional genomics, particularly studies seeking to assign biological functions to dna sequences. comparative gene expression profiling is widely used to study the functional role of genes. the dna microarray assay provides an invaluable technique for large scale expression analysis. in the two-channel dna microarray assay, rna from two samples is reverse transcribed to cdna and labeled with two distinct fluorescent dyes before being co-hybridized to immobilized dna strands on a microarray slide. spotted arrays currently being used can be divided into two groups based on the nature of immobilized dna used: one in which the immobilized dna is comprised of both sense and antisense strands and the other where the immobilized dna is single stranded consisting of only the sense strands. during hybridization, the two fluorescently labeled cdna samples compete for hybridization to the immobilized strands. hybridization reactions between complementary strands occur only between the labeled antisense strand and immobilized sense strand. the ratio of the intensities of the two fluorescently labeled cdnas is used to quantify the relative levels of transcripts in the two samples. this method serves well for pair-wise comparison of transcript levels in two samples. with over ten thousand different dna species immobilized on the microarray, the relative transcription level of all the corresponding genes in the two samples can be obtained in a single assay. dna microarrays have found applications in gene discovery, disease diagnosis, pharmacogenomics and toxicology research. they are increasingly used for a series of related samples, for which a comparison across all samples and all genes is desirable. when a large number of samples are to be compared, a combinatorial approach pairing all possible pairs is often taken. this results in a large number of microarrays, requiring a large amount of each rna sample. a'loop design', where every sample is directly compared to two other samples to form a closed loop, has been proposed to overcome this problem. the ratios calculated using a loop design have variable levels of precision since some samples are more directly related than others. when a new sample is to be inserted into the earlier'loop', rna for at least two of the previous samples is needed to pair with the new sample to form a new node in the closed loop. another approach to tackle the issue of a large combinatorial pair-wise comparison is the'reference design' in which a common reference sample is introduced with which all rna samples are hybridized. two possible universal references are rna pooled from various samples and genomic dna. for a given set of samples, pooled rna provides an excellent reference. however, if the experimental conditions change, the possibility arises that some new transcripts may not be represented in the initially-pooled rna. genomic dna is an attractive candidate for use as a common reference. it is isolated from cells or tissue and sheared to fragments in a narrow range of length. it is easier to prepare, maintain and reproduce, as compared to rna. it is especially useful for microorganisms, which lack repetitive sequences in their genome, and microarrays using genomic dna as a reference have been used to identify genes differentially expressed in various growth stages of mycobacterium tuberculosis. it has also been recently reported that the data obtained using genomic dna as a reference in microarray experiments with arabidopsis thaliana employing 70-mer oligonucleotide microarrays was in agreement with ratios obtained from direct hybridizations. genomic dna samples, isolated from stationary phase cultures where the chromosome is not being replicated, have the same representation of all genes as in the genome. since the transcript level of each gene is being referenced to its own representation on the genome, the relative expression can be compared across different genes in the sample as well as across different samples. the use of microarrays using genomic dna in a range of applications including genomic diversity studies and aneuploidy detection using comparative genomic hybridization has been demonstrated. in the conventional two-channel cdna hybridization, both the cdna samples are antisense strands. the probability of hybridization between strands in solution is very low. on the other hand, using sheared genomic dna as a reference, hybridization between complementary sense and anti-sense strands can occur in solution between the complementary genomic dna strands, and between genomic dna sense strands and their cdna counterparts. the number of strands lost to hybridization in the solution phase may differ for different rna species as well as for the two complementary strands of the same species. this may result in decreased fidelity in the ratio of cdna to genomic dna as a representation of gene expression level. with the complexity of hybridization in both the solution phase and the immobilized surface phase, and between double strands of genomic dna and single strands of cdna, it is difficult to assess the effect of using genomic dna as a common reference. adapting a mathematical model we developed previously to assess diffusional constraints on dna microarray assay, we have constructed a kinetic model for microarray hybridization and predicted the effectiveness as well as potential pitfalls of using genomic dna as a reference. we examined the effect of the various parameters that may affect surface hybridization. here we report the framework of the model and our findings. mathematical model for microarray hybridization a dna microarray consists of thousands of spots, each spot containing dna strands of known sequence immobilized on an impermeable surface. we simulate hybridization reactions for one spot on the microarray, which is considered to be a sample chamber with two compartments; a spot phase and a solution phase. the spot phase is the small volume in which the bound strands are assumed to be present at a uniform concentration. the solution phase, which comprises a vast majority of the volume of the microarray chamber, contains only fluorescently labeled single strands at the beginning of the hybridization and double strands formed by hybridization of complementary mobile single strands during the course of the hybridization. five kinds of single-stranded dna molecules are present in the system: labeled single stranded cdna reverse transcribed from rna, denoted as s; the two identically labeled complementary strands from genomic dna, denoted as r and r'; and the non-labeled anti-sense and sense strands immobilized on the array surface denoted by b and b'. the nomenclature used in this study is summarised in table the. the double-stranded species are denoted by combining the constituent single strand symbols, for example rb' denotes the labeled complex formed by hybridization of genomic dna antisense strands to the sense bound strands. we assume that each phase is well-mixed, and hence all mobile species are present at uniform concentrations within each phase. we consider the case that pcr products are used for immobilization; thus, both sense and antisense strands of probe dna are immobilized. we assume that the two complementary bound strands do not hybridize to each other due to the immobilization procedure used and neglect the formation of bb'. the mobile species in the two compartments are considered to move across the phase boundary at a rate proportional to the difference in the concentrations of identical species present in the two phases. the proportionality constant is the effective mass transfer coefficient for the transport of mobile dna strands and is estimated from the diffusivity as discussed below. the model equations take the form of a mass balance on each component in each phase that accounts for the change in concentration due to reactions within that phase and transport between the two phases. non-specific hybridization in both phases is neglected. all these equations are of the form: where = r', = b', a = r, a = b, a = s, subscript p denotes either solution or spot phase and p' is the other phase; kf is the forward rate of hybridization; kb is the backward rate of hybridization; kt is the rate of transport across the spot-phase solution-phase boundary; vt is the total volume of the sample chamber and vp is the volume of phase p. model parameters model geometry the diameter of the spot phase is set to cm as seen in a typical microarray spotted on poly-lysine slides, with a height of × 10- cm. this gives a spot volume of × 10- l. typically, labeled strands are resuspended in a volume of × 10- l, and this volume is applied under a cover slip of area cm × cm. for this geometry, the height of the microarray chamber is about × 10- cm. the solution phase volume is the difference between the sample chamber volume and the spot volume. it has been shown before that only the mobile strands within a radius of cm from the spot are available for hybridization due to transport effects. hence, we have assumed the solution phase to have a diameter of cm and a height of × 10- cm corresponding to the height of the microarray chamber. hybridization rate constants the hybridization rate constants are assumed to be identical for dna strands from both genomic dna and cdna. however, the rate constants for reaction between two mobile strands may differ from the rate constant for the reaction of a mobile strand and a bound strand. a forward rate constant for dna hybridizations in solution of m- s- was used for simulations. the rate constant of hybridization could be slower for hybridization of mobile strands to immobilized strands on a solid surface and values in a range up to 100-fold lower than kf have been used in simulations. we also discuss the effect of this reduced forward rate constant of hybridization between mobile and bound strands for a range between m-1s- and m-1s- the backward rate constant for dissociation of dna double strands as calculated from equilibrium constants reported in literature ranges from 10- s- to 10- s-. simulations were carried out using a backward rate constant ranging from 10- s- to s-. in all simulations performed, the deviation from the true value was greatest when the hybridization reaction is set to be irreversible. the results discussed below are for the case with the backward rate constant set to s- the results obtained when the reaction is set to be reversible are very close to the true value under all conditions tested. rate constant for transport the rate constant for transport is estimated from the diffusion coefficient as, where d is the diffusion coefficient, h is the height of the spot phase and r is the radius of the spot phase. the diffusivity of dna single strands in solution has been estimated to be 10- cm2/s which leads to an estimate of s- for kt. the simulations reported in the next section were carried out under both very fast transport and very slow transport conditions. the transport rate has no effect on the performance of the assay in the range tested. initial single-strand concentrations from the yeast transcriptome data published by velculescu et al., the mass percentages of mrna that belong to the rare, intermediate and abundant classes are estimated to be, and % respectively. for an mrna sample of μg used in each microarray assay, the total amount of rare, intermediate and abundant species are μg, μg and × 10- μg respectively. to convert those numbers to molar concentration we calculate the number of genes within each abundance class using the intensity data from microarray experiments of s. coelicolor. the s. coelicolor transcripts were divided into three abundance classes using the following intensity cut-offs: intensity < as rare, < intensity < as intermediate and > as abundant sequences. normalized intensity values were used for this estimate. %, 26% and % of all genes were classified rare, intermediate and abundant respectively. the corresponding number of rare, intermediate and abundant species are and respectively. assuming a sample volume of × 10- l, 100% reverse transcription efficiency, and an average strand length of we calculate the rare species concentration as ~ pm, intermediate species concentration ~ pm and abundant species concentration ~ pm. previous experimental reports of the use of genomic dna as a reference for microarray hybridizations used genomic dna concentrations ranging from μg to μg for m. tuberculosis. since m. tuberculosis has ~ genes, this translates to a concentration of each gene from ~ pm to pm, which is the range of genomic dna concentrations used for simulations. to calculate the concentration of the immobilized species, we assume the concentration of dna in the spotting solution to be g/l and that 10- l of the solution is spotted on the microarray. we also assume that 75% of the dna thus spotted is washed away in the microarray post-processing steps. × 10- g dna remaining on the microarray is uniformly distributed in the × 10- l spot phase volume. assuming an average dna strand length of bp, this is approximately equal to an immobilized strand concentration in the spot phase of 10- m.  the mathematical model described above takes into account dna hybridization between single stranded cdna and double stranded genomic dna in solution and immobilized double strands on a microarray surface. this model considers hybridization only on one spot on the microarray. the immobilized strands are distributed uniformly in the spot phase and the mobile strands are present both in the solution phase and spot phase. hybridization between mobile and bound species occurs in the spot phase. all concentrations described in the following sections are the concentrations in the spot phase. the fluorescence intensity corresponding to hybridized cdna sample is expressed as is = and the channel corresponding to the hybridized genomic dna reference as ir = +. the result obtained for hybridization with genomic dna used as a reference is a hybridization ratio γ = is/ir = /. since the concentration of all genes in a genomic dna sample is equal, the ratio for different genes is an indication of the relative abundance of the transcript for those genes. typically when genomic dna is used as the reference, the ratio from one hybridization of cdna derived from sample one, is compared to another ratio from sample two, to obtain the relative expression level of the transcripts in samples and ideally this "ratio of ratios", i.e. γ1/γ should be equal to the ratio of the initial concentrations of the transcript in those samples. we simulate this experimental process with two different initial cdna concentrations and the same genomic dna concentration to obtain the ratios γ and γ the accuracy of the microarray assay is quantified by the accuracy index α, a value of one for α corresponds to a perfect assay, where the measured relative concentration of the transcript in the two samples is exactly equal to the true relative concentration. any deviation in α from unity is a measure of the error of the assay. these results are also applicable to the comparison of expression levels of two genes in one cdna sample, as the model makes no distinction between hybridization on two spots on one microarray and hybridization to a spot corresponding to the same sequence in two different microarray experiments. we systematically vary model parameters to investigate the effect of different hybridization conditions, transcript abundance levels, and degree of differential expression on the performance of the microarray assay. effectiveness of using genomic dna as a reference the effectiveness of using genomic dna as a reference in the microarray assay was predicted by simulations using the model and parameters described above. figure 2a shows the variation in α with hybridization time for different rna abundance levels and differential expression ratios corresponding to biologically realistic scenarios as listed in table these cases represent a wide range of possible combinations of the three abundance classes and differential expression ratios. a differential expression ratio of is used as an example of a small change in expression level and as an example of a large change in expression level. the parameter values used are kf = m-1s- kf-bound = m-1s- kb = s- initial genomic dna concentration = pm, initial bound strand concentration = 10- m, kt = s- for a wide range of initial cdna concentrations and differential expression ratios, the accuracy index α is within 5% from unity, indicating that the assay performance is robust to the initial concentration of single strands. the only condition where α is significantly different from is when intermediate species are upregulated fold, a situation not likely to happen frequently in cells. furthermore, in most microarray assays, a 17% deviation from the true value is not considered large. the model simulation predicts that for abundant species, the accuracy decreases with hybridization time, indicating short hybridization time will lead to better results. however, the concentrations of the double strands formed as a result of hybridization of single strands in solution to the immobilized strands increase with time for approximately the first hours. the fluorescent intensities detected when the microarray slide is scanned is proportional to the concentration of these double stranded species. hence the hybridization time has to be long enough to obtain intensities sufficiently above background levels for accurate measurement. the model prediction that using genomic dna as a reference can provide an accurate measurement in a microarray assay was verified experimentally. the transcript levels of two samples were assayed using both direct cdna: cdna hybridization as well as by hybridizing to genomic dna. this ratio of the two samples was obtained by dividing the ratios obtained from those two cdna samples / microarray data can be found in the supplementary material. figure shows a scatter plot of the relative transcript level obtained by these two methods. the ratio obtained from indirect comparison is within fold of that obtained from direct comparison for 91% of genes. out of the remaining 9% genes, 81% have a ratio obtained from direct comparison within -fold. for % of all the genes, the ratio obtained from indirect comparison is within fold of that obtained from direct comparison. 70% of the remaining % have a ratio obtained from direct comparison within -fold. thus, in general, the ratio obtained from hybridizations using genomic dna as a reference is consistent with those obtained from direct cdna: cdna hybridization. as can be seen in figure this is true over a large range of differential expression. also, since total rna samples containing a wide range of transcript abundance levels were used in this experiment, the dataset demonstrates that the accuracy of the assay is maintained over all mrna expression levels. effect of genomic dna concentration in previous investigations several concentrations of the genomic dna ranging from to pm have been used. to examine the effect of genomic dna concentrations on the microarray assay, we carried out simulations for a range of genomic dna and cdna concentrations with different degrees of downregulation. α is within 6% of unity for all values of genomic dna concentrations examined. thus, the assessment of the ratio between the transcript levels of a gene in two samples does not vary significantly with the genomic dna concentration used in the experiment. however, in a microarray experiment, accurate assessment of the ratio of expression levels can be attained only if the fluorescent intensities can be accurately measured. spots with lower intensities and closer to the background level are prone to increased noise interference and decreased accuracy. at higher concentration of genomic dna, the intensity of the sample channel decreases. this effect is more profound for rare species, compounding the problem of their low intensity. this leads to the need to carefully select the genomic dna concentration for achieving intensities which are significantly above background levels. effect of bound strand concentration in the preparation of microarray slides the amount of dna immobilized on the slide varies. many factors, including differences in the amount of dna deposited, spot morphology, dna retention in the spotting and post-processing procedure contribute to such variation. in cdna:cdna hybridization, the ratio of transcript levels from two different samples is not affected by the amount of dna immobilized on the slides, since the same transcript species from both samples are affected to the same extent. in contrast, when cdna:gdna is used, variation in the amount of immobilized dna may exert a different effect on different samples. we simulated the effect on α of this changing bound strand concentration by using bound strand concentrations of 10- and 10- m for all cases listed in table this range spans values from 100-fold higher to 100-fold lower than the estimated bound strand concentration. all other parameters are same as those used for the simulations in figure the findings are similar to the plot shown in figure 2a. thus, the variation in the amount of dna immobilized on slides does not significantly affect the accuracy of the assay. however, a lower bound strand concentration does result in reduced concentrations of the double strands formed, and hence lower intensities on the microarray. as an illustration, when the bound strand concentration is lowered from 10- m to 10- m, the intensity of the sample channel decreases 24-fold for rare species, 25-fold for intermediate species and 38-fold for abundant species. accuracy affected by transcript abundance at low hybridization rate an uncertainty in the simulation is the effect of reaction rate constants; especially, the effect of change in the relative magnitude of forward rate constants for hybridization between mobile and mobile, and mobile and bound strands, is further investigated. to evaluate this, kf-bound was varied 100-fold from m-1s- to m-1s- and the results are shown in figure for all values of kf-bound, α is within 21% of unity for all the rna abundance classes and all differential expression ratios examined. in general, the assay is more accurate for rare and intermediate species as compared to abundant species and for lower differential expression ratios compared to higher ratios. for the abundant species, as the forward rate constant for hybridization between mobile and bound species is reduced, the accuracy of the assay decreases. the model predicts that the highest error will be observed for abundant species with high differential expression ratios. the simulation results presented above illustrate that the accuracy index α is most sensitive to the rate constant of the forward reaction for hybridization between mobile and bound species and the deviation of α from unity is highest for abundant species in the sample. the ratio determined in the hybridization assay is. ideally, this should reflect the true ratio. the ratio / is always close to s0/r since both s and r hybridize with the same species: r' in solution and b' on the surface. this is because in our analysis, we assume that s and r have identical reaction kinetics, they are thus both stoichiometrically and kinetically indistinguishable from each other. therefore, the deviation in / from r0/r0' governs the deviation in γ from the true ratio. since the reaction is irreversible, for to be equal to, the product should be equal to. let εs be the amount of s reacted with r' and εr be the amount of r reacted with r'. since the stoichiometric ratios for reactions of s with r' and r with r' are both r' = -εr - εs and r = r - εr. under conditions where s is low, εs is relatively small and r ≈ r'. however, for abundant species, εs is significant and as a result, r/r' is greater than r0/. also, since s+b' → sb' r+b' → rb' r'+b → r'b b' gets consumed more than b. again, for rare and intermediate abundance species, this difference in b and b' is insignificant since the bound strands are in excess and hence the assay works very well for rare and intermediate species. the difference between b and b' is significant for abundant species for which b/b' is greater than b0/b0'. the exact deviation in b/b' and r/r' depends on the rate constants for hybridization and hence it is not surprising that the rate constants have a major effect on hybridization results. this is illustrated by values for b, b', r and r' for s = pm, r = = pm and kf-bound = m- s- and m- s- at the end of hours of hybridization time, presented in table the products, and ratios / and / are also calculated. thus, for the larger value of kf-bound, ≈, which is not the case for the lower values of kf-bound. this means that factors, such as mixing conditions in the sample chamber, which affect reaction kinetics, will affect the accuracy of the assay. one implication of the above discussion is that if single stranded species are used for immobilization, as are used in spotted oligo-arrays, the assay will be robust to a wide range of hybridization conditions. this is also seen from our simulations where for all the conditions discussed in this paper, α is equal to one if only sense strands are immobilized on the microarray slide, a situation encountered during the use of oligonucleotide spotted arrays.  the use of genomic dna as a reference is useful to assess the expression levels of a large array of genes among different samples. we have developed a kinetic model to predict the effect of using genomic dna in the microarray assay under a wide range of conditions. the model predicts that the assay can accurately estimate the relative concentrations for a wide range of initial cdna concentrations and ratios from hybridizations using genomic dna as a reference will correspond to true ratios under all conditions if single stranded oligonucleotide microarrays are used. the model also serves as a useful tool to predict the performance of such assays under varying conditions that are otherwise difficult to carry out experimentally. despite a number of publications on its application, the use of genomic dna as a reference for microarray assay is still not wide spread. we carried out this study to verify on a theoretical basis the validity of this approach and the results are indeed reassuring. we expect that the use of genomic dna as reference will accelerate especially for comparative transcriptome analysis involving a wide range of samples from different sources. 
to compare a large number of samples, the'reference design' is widely used, in which all rna samples are hybridized to a common reference. however, genomic dna, comprising of both the sense and anti-sense strands, is unlike the single stranded cdna usually used in microarray hybridizations. experimental results of dna microarray assay using genomic dna as a reference correlated well to those obtained by a direct hybridization between two cdna samples. for the case where a microarray with immobilized single strands is used, results from hybridizations using genomic dna as a reference will correspond to true ratios under all conditions. the model predicts that the assay can accurately estimate relative concentrations for a wide range of initial cdna concentrations. at low values of the rate constant for hybridization between immobilized and mobile strands, the assay performance varies with the hybridization time and initial cdna concentrations. genomic dna is an attractive candidate for use as a universal reference, especially for bacterial systems with a low percentage of non-coding sequences. simulation using the mathematical model, and the experimental study presented here show the potential utility of microarray assays using genomic dna as a reference. the presence of the antisense strand in the'reference' leads to reactions between complementary labeled strands in solution and may cause the assay result to deviate from true values. the model predicts that the initial concentrations of labeled genomic dna strands and immobilized strands, and the hybridization time do not significantly affect the assay performance. we have developed a mathematical model to predict the validity of using genomic dna as a reference in the microarray assay. we conclude that the use of genomic dna as reference dna should greatly facilitate comparative transcriptome analysis. dna microarray is an invaluable tool for gene expression explorations.
 dna microarray is an invaluable tool for gene expression explorations. in the two-dye microarray, fluorescence intensities of two samples, each labeled with a different dye, are compared after hybridization. to compare a large number of samples, the'reference design' is widely used, in which all rna samples are hybridized to a common reference. genomic dna is an attractive candidate for use as a universal reference, especially for bacterial systems with a low percentage of non-coding sequences. however, genomic dna, comprising of both the sense and anti-sense strands, is unlike the single stranded cdna usually used in microarray hybridizations. the presence of the antisense strand in the'reference' leads to reactions between complementary labeled strands in solution and may cause the assay result to deviate from true values.  we have developed a mathematical model to predict the validity of using genomic dna as a reference in the microarray assay. the model predicts that the assay can accurately estimate relative concentrations for a wide range of initial cdna concentrations. experimental results of dna microarray assay using genomic dna as a reference correlated well to those obtained by a direct hybridization between two cdna samples. the model predicts that the initial concentrations of labeled genomic dna strands and immobilized strands, and the hybridization time do not significantly affect the assay performance. at low values of the rate constant for hybridization between immobilized and mobile strands, the assay performance varies with the hybridization time and initial cdna concentrations. for the case where a microarray with immobilized single strands is used, results from hybridizations using genomic dna as a reference will correspond to true ratios under all conditions.  simulation using the mathematical model, and the experimental study presented here show the potential utility of microarray assays using genomic dna as a reference. we conclude that the use of genomic dna as reference dna should greatly facilitate comparative transcriptome analysis. 