BACKGROUND
glucocorticoids  are a class of steroid hormones present in almost every animal cell, playing a central role in a wide range of physiological responses  <cit> .
because of their potent anti-inflammatory and immunosuppressive effects, synthetic glucocorticoids referred as corticosteroids   have been used widely in pharmacology as a therapeutic option for a wide range of autoimmune and inflammatory diseases  <cit> .
however, beneficial effects are derived from magnifying the physiological actions of endogenous glucocorticoids, causing a variety of side effects following long-term treatment with this class of drugs e.g.
hyperglycemia, dyslipidemia, arteriosclerosis, muscle wasting, and osteoporosis  <cit> .
the physiological and pharmacological effects of corticosteroids are complex and manifest themselves with expression changes of many genes across multiple tissues  <cit> .
it has been observed that even in a single tissue different dosing regimens of cs administration can induce different patterns of expression  <cit> .
as such genes with similar expression profiles under acute cs administration may not exhibit similar expression patterns during continuous infusion, pointing to the possibility of alternative regulatory mechanisms.
therefore, a better understanding of corticosteroid pharmacogenomic effects from multiple dosing regimens are very valuable not only to reveal the transcriptional dynamics under different patterns of input perturbations but also to provide an insight into the underlying molecular mechanisms of action, for both beneficial and detrimental effects, and thus for the optimization of clinical therapies.
it has been noted that genes affected by cs include both immunosuppressive genes, mostly associated with therapeutic effects, and metabolic genes often associated with adverse effects whose regulation is mainly controlled by glucocorticoid receptor gene mediated pathways  <cit> .
unbound cs binds with cytosolic free glucocorticoid receptors  releasing it from the heat shock complex allowing dimerization and translocation into the nucleus where it binds to glucocorticoid response element  of the target genes, leading to enhancement or inhibition of the target gene expression.
as a result, long-term treatment with corticosteroids results in sustained up- or down-regulation of numerous genes, leading to a new steady state which might be the basis for occurrence of adverse effects.
however, it has also been noted that chronic infusion of cs causes a sustained down-regulation of the receptor   <cit> .
while several alternative mechanisms have been proposed  <cit>  it is still not understood why drug effects remain strong although gr mrna is down-regulated to the point of almost being eliminated.
a plausible explanation is that besides direct regulation through gre binding sites in the 5' regulatory regions of genes, there are changes in expression that are also the indirect results of effects due to changes in expression of transcription factors  that act as secondary biosignals directly or indirectly modulating the transcription of genes  <cit> .
thus, along with identification of expression patterns, predicted regulatory control structures are also an essential source of information towards understanding corticosteroid effects.
in this study we address the question as to whether  significant transcriptional modules coupled with complex patterns of mrna changes across multiple dosing regimens of corticosteroids and  their common regulatory controls can be computationally identified.
hypothetically, transcriptional modules that are significantly coexpressed under different dosing regimens will be important gene clusters for further analysis towards better understanding of both beneficial and adverse effects of corticosteroids, especially the metabolic side-effects since these patterns are survived under a long-term treatment of corticosteroids.
the hypothesis explored here is that if two or more genes have the same temporal expression profiles in response to different dosing regimens, they are more likely to share some common regulatory mechanisms.
the liver was selected because of its major role in both the physiological efficacious and adverse effects of corticosteroids e.g.
altering the expression of serum proteins that regulate immune/inflammatory responses  <cit> , enhancing the expression of liver enzymes involved in metabolic effects   <cit> .
however, rich in vivo datasets of pharmacological time-series across multiple dosing regimens are often obtained from different microarray platforms and time-sets  <cit> , leading to a problematical issue for computational analysis  <cit> .
as an example, in a study comparing normal and chronic lymphocytic leukemia b-cells, wang et al.
<cit>  identified only  <dig> differentially expressed genes across all three studies, when combining results from three different platforms, while there are at least  <dig>  differentially expressed genes in each individual platform.
in general, there are two important issues relevant to the analysis of data derived from different platforms:  genes may be present in one platform but not in the other, and  genes present on both platforms may not be represented by the same probes.
since different microarray platforms do not contain the same probesets, and even do not have a similar hardware design and sample processing protocols, standard analyses may not yield comparable expression level quantifications across platforms, leading to many challenges for computational models aiming at the analysis of microarray data from heterogeneous sources  <cit> .
a number of approaches have been proposed and are generally classified into two main categories:  integrate raw expression profiles from different studies into one dataset so that available computational models can be directly applied, and  develop and/or utilize a unitless statistic as a primary analysis and then combine the result through a meta-level analysis.
the former category can be further divided into two sub-classes, namely combining raw data through a normalization and/or transformation procedure  <cit>  and pooling raw information from common probes that can be mapped to the same unigene clusters or full-length mrna transcripts  <cit> .
however, these approaches are not general enough to make data from different platforms fully compatible  <cit> .
since combining data across different platforms remains a serious challenge, meta-analysis - the second category - has been identified as a more popular technique in order to combine results, and thus data, from a number of independent studies  <cit> .
the assumption here is that while the raw expression levels from different platforms may not be comparable, the results of the primary analysis should be.
however, almost all prior studies has focused on the discovery of genes that are differentially expressed in conjunction with standard models such as effect size models  <cit> , bayesian models  <cit> .
consequently, in order to identify significant clusters of genes that share common expression patterns across multiple dosing regimens, we extend our prior study  <cit>  in the aspect of  producing an agreement matrix  that describes the agreement levels of co-expression of genes across multiple conditions and  successively searching clusterable subsets to infer all such gene clusters.
the approach follows the concept of meta-analysis to avoid the limitation of incompatible data across multiple datasets from different platforms .
the unitless statistic, expressing the confidence level of co-expression is the agreement level of cluster assignments drawn from multiple clustering runs.
there remain a number of open critical issues associated with a single clustering run , and thus consensus clustering coupled with the examination of am distribution has been designed with the aims of reducing aforementioned limitations  <cit> .
once the am is obtained for each condition independently , an average agreement matrix is calculated to estimate the confidence levels of coexpression between genes across multiple conditions, thus combining data from different datasets into a single input for the next analysis.
for the analysis at the meta-level, we extend the selection and clustering processes  to identify all possible clusters of genes that are highly coexpressed with the average am above as the input.
as such these clusters of genes will share common patterns of expression across multiple dosing regimens.
additionally, due to the selection of all possible patterns of expression several clusters may have similar expression patterns and thus we also provide a heuristic as an optional procedure to merge such similar clusters based on a criterion of maximizing the total homogeneity and separation of selected clusters.
subsequently, we analyze promoter regions of genes in every cluster in order to predict putative transcriptional regulators, aiming at providing another source of information towards better understanding those complex patterns of expression and the underlying regulatory mechanisms of corticosteroid effects.
our results demonstrate that the proposed computational approach is highly effective on both synthetic and real data.
when applying the approach to real time-series datasets , selected patterns of transcriptional responses are enriched in a biological sense with relevant putative-regulatory controls and significant metabolic pathways in each pattern.
computational results are further validated predicated upon literature evidence.
methods
<dig>  datasets
synthetic data
a number of synthetic datasets from the open literature are utilized to assess our approach for finding common sets of genes that are highly coexpressed across multiple conditions.
specifically, we used a series of four high-noise 20-timepoint sine-format synthetic datasets with different number of replicates at each time-point  from  <cit> .
each dataset contain five separate sets with  <dig> genes allocated equally in  <dig> classes, each of which contains the same list of genes but has different patterns across five conditions.
for each set, in the first step the data are generated according to an artificial pattern Φ which shows the values for gene i at time-point t in cluster l; four of six clusters follow the sine function i.e.
Φ  = sin , and the other two follow the non-periodic linear function  = t/ <dig> and Φ = -t/20), i =  <dig> ..., <dig>  t =  <dig> ..., <dig>  l =  <dig> ..., <dig>  in the second step, let x be the error-added value for gene i, time-point t and repeat; x is randomly drawn from a normal distribution n where μ is the value of the synthetic pattern Φ and σ is equal to λσit .
high-noise synthetic data are generated with λ =  <dig>  <cit> .
datasets are downloaded from the links provided in synthetic data in additional file  <dig>
acute corticosteroid administration
forty-seven male adx wistar rats weighting from  <dig> to  <dig> g underwent right jugular vein cannulation under light ether anesthesia  <dig> day before the study  <cit> .
forty-three rats were injected with a single intravenous bolus dose of methylprednisolone  of  <dig> mg/kg.
animals were sacrificed by exsanguinations under anesthesia and liver samples were harvested at  <dig> ,  <dig> ,  <dig> ,  <dig>   <dig>   <dig>   <dig>   <dig> ,  <dig>   <dig>   <dig>   <dig>   <dig>   <dig>   <dig>  and  <dig> after dosing.
the sampling time points were selected based on preliminary studies describing gr dynamics and enzyme induction in liver.
four untreated rats were sacrificed at random times and nominally considered as  <dig> h controls.
the gene expression was obtained via the affymetrix rg-u34a array which consists of  <dig>  probesets.
the data are publicly available through the geo omnibus database under the accession number gds <dig>
chronic corticosteroid administration
in a similar experiment model, forty rats were given  <dig>  mg/kg/hr infusions of mpl over  <dig> h via an azlet pump  <cit> .
the pump drug solutions were prepared for each rat based on its predose body weight.
animals were sacrificed at various times up to  <dig> days; specifically the time-points included are  <dig>   <dig>   <dig>   <dig>   <dig>   <dig>   <dig>   <dig>   <dig>  and  <dig> h. a control group of four animals was implanted with a saline-filled pump and killed at various times throughout the 7-day study period.
unlike the previous experiment, the microarray platform for this dataset is the rae230a which consists of  <dig>  probesets.
the data are publicly available through the geo omnibus database under the accession number gds <dig>
all protocols followed the principles of laboratory animal care  and were approved by the university at buffalo institutional animal care and use committee.
<dig>  identifying critical transcriptional modules
the general computational problem can be briefly defined as follows.
we are given a set of n genes g={gi}i=1n and k conditions.
for each condition k, every gene is characterized by one or more time-series expression profiles with rki corresponding probesets over tk time-points gk={gki}i=1n,gki={gkir,r∈rki},gkir={gkitr}t=1tk,k= <dig> ...,k. the question then becomes to search for clusters of genes that are highly coexpressed across all k conditions with a confidence level δ. the term 'highly coexpressed' is used in the sense that ∀gi,gj∈c,1k∑k=1kpk≥δ where c denotes a, yet to be determined, cluster and pk is the confidence level that two gene profiles i and j are clustered together in condition k; a gene profile includes sets of corresponding probesets rki of gene i in condition k, k =  <dig> ...,k. the subscripts {i, j}, t, k, r indicate the {gene id}, time, condition, and probesets respectively.
it should be also noted that in this work, we used three different terms to refer to the same object : 'cluster' when designing the algorithm, 'pattern' when exhibiting the expression changes, and 'module' when charactering the biological function.
the framework contains several step displayed in figure  <dig>
the pre-processing step
each dataset is pre-filtered to identify differentially expressed probesets.
since we would like to identify gene clusters with common expression patterns across multiple conditions, input datasets must contain the same set of genes.
thus using the respective platform information, probesets in each dataset are mapped to a list of genes and then the intersection across those gene lists is evaluated to extract a common set of genes which are differentially expressed across multiple conditions .
however, genes are sometimes characterized by multiple probesets whose expression profiles may be similar or sometimes different, but not identical.
these probesets can be considered as replicates of expression profiles for a gene and thus taking an average expression profiles across all these probesets to characterize for the expression profile of the gene may lose useful potential information.
therefore, from the common set of genes we re-map genes to corresponding probesets in each dataset with the respective platform before starting the analysis.
construction of the agreement matrix
the agreement matrix  quantifies the likelihood that two objects  are assigned to the same cluster .
if m clustering runs are performed on the data, each entry  will show the frequency with which two objects are assigned to the same cluster over 'm' clustering runs.
the am entries are defined as follows:
mxy=1m∑h=1mmwhere m={1 if x and y are clustered together when running method m0 otherwise
in addition to the various clustering methods that were utilized, different distance metrics  are also explored in order to attenuate the biases associated with individual distance metrics.
in our implementation, we are using hierarchical clustering , divisive analysis clustering , fuzzy analysis clustering , partitioning around medoid  with pearson correlation and manhattan metric, k-means , fuzzy c-means , self-organizing map , and model-based clustering  with euclidean metric as the core clustering methods   <cit> .
since clustering results are highly dependent on the input number of clusters , the sensitivity of the am as a function of nc was examined to find a 'suggestive' number of clusters  for each particular dataset.
after identifying nc* based on the procedure in our prior work  <cit> , all clustering runs are repeated with nc* to produce the final am for further analysis .
if two probesets  are clustered together, it is implied that their expression profiles are similar under a specific condition k. therefore, the fraction of times  they are clustered together over multiple clustering runs can be considered as the confidence level that they are coexpressed since mxy, by construction, aims at eliminating method-specific biases and assumptions.
subsequently, we calculate the average agreement levels between sets of corresponding probesets of any two genes to estimate the confidence level that those two genes are coexpressed in a specific condition.
the am entries in condition k is re-estimated as follows
amij=1|rki||rkj|∑x∈rki∑y∈rkjmxy,i,j= <dig> ...,n
with the assumption that the unitless statistics, i.e.
the confidence level of co-expression, is comparable across multiple conditions and different platforms  <cit> , we estimate the confidence level of co-expression between two genes across multiple conditions by taking the average.
while combining raw data remain challenges, the estimation of a unitless statistics provides a simple but efficient combination of heterogeneous data for further analyses.
amij=1k∑k=1kamij,i,j= <dig> ...,n
as a result, we obtain an agreement matrix whose entries exhibit a quantity that shows how confident genes are coexpressed.
this will be the input for the selection and clustering process.
selection and clustering
with the hypothesis that the more clusterable the data is the more biologically relevant it is, we applied our previously proposed procedure to select a more 'hypothetically clusterable' subset from the entire set of genes  <cit> .
the main hypothesis underlying the selection is that am entries associated with genes at the 'hypothetical' core of an expression pattern  will be consistently grouped together over multiple clustering runs.
this should be manifested by high corresponding values in the am, whereas genes belonging to the 'hypothetical' core of two clearly distinct clusters are associated with consistently low am entries.
on the contrary, cluster assignments associated with genes at cluster boundaries or between clusters would be very sensitive to the method used and thus they would have relatively moderate agreement levels with other genes.
as a result, with a user-defined confidence level δ genes associated with moderate am entries <δ) are eliminated to produce a more 'clusterable' subset of genes .
the process starts removing genes associated with the highest number of moderate am entries and then updates the am for the next loop until no moderate am entry exists.
the corresponding subset of genes is now considered as a 'hypothetically clusterable' subset since any two genes are highly coexpressed or non-coexpressed with the confidence level at least δ. subsequently, we used the concept of consensus clustering  <cit>  to divide the subset of genes into a number of clusters by applying the hierarchical clustering with the selected am as input data.
the algorithm starts with every gene filling a cluster and then grouping two clusters into a new one for each loop so that any pair of genes belonging to a new cluster always has an agreement level greater than or equal to δ. the iteration is stopped when no more new cluster is formed .
however, since there should be existed clusters of genes located closely to other clusters in the data and the input number of clusters for the core analysis is only a suggestive one, those clusters may not be completely separated.
as a result, although genes that belong to those clusters are identified as highly coexpressed, their relationship to genes in other clusters cannot be uniquely determined.
therefore, some significant clusters may be not included in the selected subset due to the constraint of 'clusterable' selection.
since we would like to obtain all significant patterns of expression, the procedure of selection and clustering is repeated on the removed domain.
the removed domain consists of a set of unselected genes whose coexpression levels are still high as quantified agreement levels in the original am.
after extracting the sub-agreement matrix corresponding to the set of unselected genes, the entire process of selection and clustering is applied again with the same confidence level δ as before.
the procedure is reiterated until no more clusters of genes are recognized.
figure  <dig> presents the pseudo-code of the procedure and an example to illustrate the process.
furthermore, due to the nature of clustering, trivial clusters may be identified in the final results.
in order to exclude such trivial clusters, each cluster c is assigned with a simple hypothetical quantity called 'cluster significance' which represents how large the cluster is in this study.
we then create the distribution of cluster significance on random data to estimate the cluster significance threshold corresponding to a user-defined threshold p-value for cluster selection.
without loss of generality we select k =  <dig> for the random data and assume that each probeset in the mapped dataset d correspond to a gene.
the suggestive number of clusters nc* for d is searched with the process in  <cit> .
subsequently, d is randomly resampled  a number of times , for each of which the entire process starting from building the am with the same nc* to extracting clusters of highly coexpressed genes is repeated and the resulting random clusters are returned.
after that, the procedure estimates the cluster significance, which is simply the cluster-size in this study, for these random clusters and constructs a distribution of cluster significance.
the corresponding p-value of cluster significance cs is defined as the number of random clusters whose significance is greater than cs over the total number of clusters identified in nr resamples: pvalue =∑clusters with cluster significance≥cs∑random clusters.
as a result, given a threshold p-value , the corresponding cluster significance threshold is inferred  and only clusters with significance greater than this threshold are selected.
merging similar patterns
because of the nature of the approach, it is quite reasonable to expect that the clustering process can break out patterns of expression into several sub-patterns.
thus, we repeat the process on the eliminated domains to extract all possible significant clusters, resulting in that several clusters may have a similar expression pattern but are separated into two or more clusters.
because cluster homogeneity reflects the similarity of expression profiles within a given cluster and cluster separation quantifies how effectively expression profiles are discriminated, we provide an optional heuristic in order to merge similar patterns together according to the criterion of maximizing sum of homogeneity and separation of all final output clusters.
starting with all significantly selected clusters, the procedure searches for a combination of two similar patterns so that their combination will generate a maximal increase of the sum of homogeneity and separation of all current clusters after merging those two patterns.
the process is repeated until no more combinations are found i.e.
any new combination always reduces the sum of homogeneity and separation.
eventually, a list of significant expression patterns that characterize the underlying transcriptional responses is generated.
the metric used during the evaluation of the heuristic is quantified as follows:
maxc{homogeneity+separation}=maxc{1k n∑k∑phk+2k n∑k∑p<qsk};hk=∑gi,gj∈cp;gi≠gjsim||cp||;sk=∑gi∈ck;gj∈cldis||cp||  ||cq||;sim=1‖rki‖‖rkj‖similarity;dis=1‖rki‖ ‖rkj‖dissimilarity;p= <dig> ...,n,q= <dig> ...,n,k= <dig> ...,k,r∈rki,r'∈rkj
where c is the current set of selected clusters c={cp}p=1n and n is the current number of clusters; hk  is the homogeneity of cluster cp in condition k and sk  is the separation between cluster cp and cq in condition k; sim and dis are the average similarity and dissimilarity  respectively between all probesets of gene i and gene j in condition k. similarity is measured by the pearson correlation coefficient and dissimilarity is estimated by the pearson correlation distance.
<dig>  predicting putative transcriptional regulators
promoter identification
promoters of genes are extracted from a rich database of promoter information with a default length  if there is no experimentally defined length as suggested by genomatix  <cit> .
in order to identify putative transcriptional regulators, we explore the basic underlying assumption of comparative genomics which states that functional regions evolve in a constrained fashion and thus at a lower rate than non-functional regions  <cit> .
it implies that conserved regions in a set of orthologous sequences survive due to their special functional implications i.e.
tfbss located on these conserved regions will be more promising as functional binding sites and thus associated tf families are more relevant to our context.
therefore, each promoter is characterized by a set of promoters from orthologous genes of other vertebrate species, if available .
to be consistent in the search for conserved regions on promoter sequences in order to identify putative transcription factor binding sites  we eliminate those that do not consist of more than two orthologous promoters.
putative functional binding sites
in order to identify conserved regions for each promoter dialign  <cit>  was used to perform multiple sequence alignments with the input sequences including each sequence as well as its orthologous promoters.
dialign was selected because it has many applications in comparative genomics  <cit> .
also, a recent benchmark study for the alignment of non-coding dna sequences has concluded that it can produce alignments with high sensitivity as well as specificity to detect constrained sites  <cit> .
following the alignment among orthologous promoter sequences, we relied on the conserved scores returning from dialign  to locate conserved regions which are defined as sub-sequences that are longer than  <dig> bp and continuously scored greater than the average score of all the alignment's conserved scores .
we next apply matinspector  <cit>  to scan for all physical tfbss and only those that overlap with the conserved regions selected above are kept for further analysis.
we used a common core similarity  <dig>  and utilized the optimal matrix similarity threshold for each position weight matrix  suggested from matbase, genomatix  <cit>  which ensure that a minimum number of matches are found in non-regulatory sequences i.e.
the false positive matches is minimized.
however, a gene may have multiple alternative promoters  <cit>  and virtually in all cases, it is not known which promoter of the gene is activated.
therefore, all putative tfbss detected from all alternative promoters of a gene are considered as candidates to infer putative transcriptional regulators for the gene.
subsequently, we estimate the common level of each candidate above in each corresponding module and select those tfbss present more than a common threshold  .
associated tf families with those selected tfbss are inferred and considered as transcriptional regulators for corresponding transcriptional modules.
RESULTS
method evaluation on synthetic data
in order to evaluate the effectiveness of the proposed approach, we use synthetic data with known class structure as described earlier.
the process of evaluation is repeated four times with four different datasets that are created with different number of replicates for each time-point .
in each time, we use five high-noise sets as the data for five conditions , each of which has  <dig> genes distributed across  <dig> clusters; each cluster has different patterns across five conditions but has the same set of gene ids.
we set the same parameters for all evaluation in this study and also for the analysis on real time-series datasets, specifically the confidence level of coexpression δ = 70% and p-value =  <dig>  for the selection of significant clusters.
furthermore, the testing process on synthetic data is done without the merging option.
without loss of generality, we assume that each gene has only one probeset in this evaluation.
the performance of the approach is assessed through its ability to recover the number of cluster structures and the list of gene ids identified in each cluster.
we use the adjusted rand index  <cit>  which is a statistic that measures the extent of concurrence between the clustering results and the underlying known class structure to evaluate the approach's performance in identifying gene clusters that are coexpressed across multiple conditions.
the larger the rand index is, the higher the agreement between the results and prior knowledge of cluster structure.
the number of selected genes, recovered cluster structures and the accuracy on the selected domain are listed in table  <dig>
*: the accuracy is only estimated on the selected domain
due to the fact that these datasets are high-noise synthetic data, some cluster structures may be missed when there is only one measurement at each time-point.
however, when the number of replicates is increased, the number of cluster structures is recovered.
as discussed in our previous study  <cit> , this is a reasonable observation due to the effect of gene expression replicates on clustering performance.
additionally, we also examine an alternative approach which is more intuitive in identifying gene clusters that are coexpressed across multiple conditions.
instead of performing a meta-analysis to avoid the limitation of incompatible data across different platforms, we can separately identify significant clusters of genes that are coexpressed in each condition  and then obtain the intersection among these gene clusters across all conditions.
in this experiment, we used pam  <cit> , mclust  <cit> , and consensus clustering   <cit>  as standard single clustering methods to identify clusters in each set of data, for which nc* =  <dig> is the input number of clusters.
we then simply took the intersection between clusters from set to set and only keep those clusters that contain more than  <dig> genes as significant clusters for the final estimation of accuracy.
the number of selected genes, number of clusters, and accuracy on the selected domain are listed in table  <dig>
+: only clusters with more than  <dig> genes;*: the accuracy is estimated on the selected domain; the number of clusters and the accuracy are formatted as 'before merging' | 'after merging'
in general, this approach selects a smaller number of genes with an equal or greater number of cluster structures, resulting in lower accuracy.
as an example, in each set of dataset  <dig> there are two cluster structures that are not clearly identifiable.
as a result, a single clustering methods  may fail to properly separate them in each set, leading to the situation where the intersection between clusters from set to set divides those cluster structures into many sub-clusters with a small number of genes.
on the contrary, by taking the average of the co-expression levels across multiple sets, the relationship of whether two genes are coexpressed across multiple conditions can be recovered.
consequently, our proposed approach is more advantage, resulting in a final highly correct classification as illustrated in table  <dig>  furthermore, since this simpler alternative approach produces many resulting clusters, we also attempted to apply the proposed merging process to reduce the number of clusters as well as improve the accuracy if applicable.
however, its testing performance is still not as high as that of our proposed approach although we do not apply the merging process for the proposed approach in this test.
additionally, the alternative approach is highly sensitive with the initial number of clusters.
for instance, when we constantly set nc* =  <dig> and test on dataset  <dig>  without the merging option our approach still recovers the correct number of cluster structures with high accuracy:  =  whereas 'pam' approach yields , 'mclust' provides , and 'consclust' does .
since this information is not available for real datasets, the more sensitive with it the less robust the approach is.
therefore, by taking the average of the co-expression levels between two genes across multiple datasets, our proposed approach provides more robust results.
acute vs. chronic cs administration
for the analysis of corticosteroid administration, the pre-processing step  is performed to provide corresponding mapped datasets.
the datasets are first filtered for differentially expressed probesets using anova technique  implemented in r  <cit>  and also customized by our previous work for easy uses  <cit> .
<dig>  probesets in the acute and  <dig>  probesets in the chronic are selected for further analysis.
to obtain the common set of genes across two conditions, these probesets are mapped into sets of genes based on the corresponding platform information.
<dig>  differentially expressed probesets in the acute are mapped into a set of  <dig>  genes and  <dig>  probesets in the chronic are mapped into another set of  <dig>  genes.
the intersection of these two gene-sets yields  <dig> genes in common for both dosing regimens.
from this common gene set, the re-mapping process subsequently returns a corresponding set of  <dig>  probesets for the acute and a set of  <dig>  probesets for the chronic data.
all datasets  are pre-processed with the model in our previous study to estimate the 'true' expression profiles that are integrated with potential information in replicates instead of simply taking the average expression profiles  <cit> .
the suggestive number of clusters nc* for both datasets is  <dig>
subsequently, we apply the proposed approach with the merging option to the intersection set of  <dig> genes that are affected by corticosteroid administration across the two dosing regimens.
we obtain  <dig> significant clusters with  <dig> genes in total.
these clusters are hypothesized to be transcriptional modules which share common regulatory mechanisms since they consist of genes that exhibit similar expression patterns in both acute and chronic dosing regimen.
table  <dig> shows the distribution of these  <dig> genes over six modules and also briefly describes how the pattern of expression changes.
although genes may exhibit simple or complex patterns of expression during corticosteroid administration, we crudely classify those patterns into up- or down- with one or two phases of regulation.
*: patterns consist of one-phase regulation , two-phase regulation  or simply up- or down-regulation.
details are listed in expression patterns in additional file  <dig>
a detailed description for patterns of these transcriptional modules is shown in figure  <dig> with the average expression patterns of all probesets clustered in each module following acute and chronic dosing.
in brief, transcriptional module  <dig>  is characterized by one-phase regulation in acute but two phases in chronic dosing.
genes in this module exhibited a fast and robust decline in mrna, which reached its peak between  <dig> h and  <dig> h, and returned to the baseline after about  <dig> h. however, when mpl is infused  this set of genes shows a more complex pattern involving both enhanced and suppressed regulation.
although a strong down regulation is observed at the beginning, it is subsequently followed by a sharp induction with the maximum around  <dig> h and then gradually returned to the baseline indicating some kind of possible tolerance.
the second transcriptional module  shows a similar pattern of expression in both acute and chronic regimen with two phases of regulation.
genes in this module exhibit an early up-regulation and reached their corresponding peaks at around  <dig> h in the acute and  <dig> h in the chronic.
subsequently, both profiles denote a clear down-regulation  and possible slight fluctuation before returning to base line.
an interesting dynamics is observed in the  <dig> genes of transcriptional module  <dig>  in the acute dosing, the genes in this module clearly exhibit an expression pattern with two phases of regulation .
yet, in chronic administration they exhibited an early transient decline in mrna followed by robust, sustained, up-regulation.
similar to module  <dig> is the transcriptional dynamics exhibited by transcriptional module  <dig>  characterized by an early induction with a maximum at  <dig>  h in the acute and  <dig> h in the chronic.
a typical pattern with down regulation for both acute and chronic administration is illustrated by transcriptional module  <dig> .
however, genes in the acute regimen exhibited a fluctuated repression with a maximum at around  <dig> h and then followed by an induction to return to the baseline as late as  <dig> h. meanwhile, genes in the chronic regimen characterized a pattern with a slightly transient up-regulation followed by a sustained down-regulation and eventual convergence to a new steady state in the presence of the drug.
the last transcriptional module  has a similar acute pattern of expression with two phases of regulation as that of module  <dig>  however, in the chronic regimen after falling to a value below the baseline  this set of genes was further sustained a slight suppression.
while comparing these expression patterns, we observe that modules  <dig>   <dig>  and  <dig> have similar expression patterns in acute  or chronic  with a slight difference in the other dosing regimen .
although the difference is not large enough to be intuitively recognized, the merging process could not merge them together, implying that the difference is significant.
furthermore, the separation of these expression patterns is also reinforced with different functional characteristics which will be illustrated below.
in summary, selected transcriptional modules exhibit a number of typical expression patterns under corticosteroid administration.
the pattern can be simply expressed as an up- or down- regulation or as a more complex one with two phases of regulation plus some fluctuation .
putative transcriptional regulators of critical transcriptional modules
it has been widely accepted that after corticosteroids bind to cytosolic glucocorticoid receptors , the activated steroid-receptor complex is rapidly translocated into the nucleus where it can alter the expression of target genes.
however, the drug seems to be cleared within about  <dig> h following a bolus injection, suggesting that the mrna levels of cs-target genes will return to the base line after that  <cit> .
in the contrast, the drug will reach and remain to a stable steady state after  <dig> h in the chronic administration.
yet, the gr is greatly diminished in response to corticosteroids  <cit> , suggesting that the mrna levels of cs-target genes in the chronic regimen should also return to the base line.
this mechanism is corresponding to the first-phase regulation of target genes.
however, almost all chronic patterns involve two phases of regulation and some  are only half-phase patterns i.e.
persistent up or down without returning to the baseline.
these complexities in expression patterns of cs-target genes can be explained by a number of possibilities previous studies have shown  <cit> , including multiple gr isoforms, multiple gres with different affinities to the drug receptor complex, or some other receptors that can mediate the effect of corticosteroids and thus affected genes in this case can reach a new steady state in the presence of the drug .
however, another possibility is a mechanism that results in the regulation of secondary biosignals which transcription factors are the most potential factors.
after affected by corticosteroids, they in turn further modulate the expression of glucocorticoid-regulated genes as a continuing cascade of events that were initiated by the drug.
as a result, this possibility suggests a possible interpretation of the complexities in expression changes of multiple cs-target genes with the second phase of regulation .
in order to reveal some underlying regulatory mechanism of these selected transcriptional modules, we start analyzing the promoter regions of genes to search for significant putative transcriptional regulators as well as possible relationships of regulation.
the hypothesis we explore here is that if two or more genes have similar temporal profiles in response to multiple dosing regimens, they are more likely to share some common regulatory mechanisms.
for the  <dig> genes in six transcriptional modules, we extract  <dig> rattus norvegicus's promoter sequences, of which we only keep  <dig> genes with  <dig> promoters that include sufficient information of orthologous promoters for further analysis.
figure  <dig> shows the identified putative regulation between tf families and transcriptional modules.
this finding highlights the possibility that secondary biosignals are involved in the regulatory complexities of expression changes for cs-affected genes.
almost all suggested tf families do consist of transcription factor members that are recognized as differentially expressed genes in one or both dosing regimens .
since transcription factors are characterized by pleiotropic effects, it is reasonable to observe a significant overlap across various transcriptional modules  <cit> .
while comparing these regulatory combinations, we observe that some tf families seem to be common regulators for all modules  whereas some are very specific to particular modules .
this could possibly explain the difference among the expression patterns of these modules.
it is likely that the more similar the expression pattern of clusters the more likely they share a larger fraction of common regulators, e.g.
tf families in this case.
for example, there are a large number of transcriptional regulators that are common between modules  <dig>   <dig> and  <dig> but it seems little overlap exists between the transcriptional regulators of modules  <dig> and  <dig>   <dig> and  <dig>   <dig> and  <dig>  except common regulators on the top of the figure.
furthermore, one tf family consists of several tf members whose expression patterns may be similar or different across multiple dosing regimens.
figure  <dig> shows two examples , each of which is represented by two tf members whose expression patterns are similar  or different  between the acute and the chronic corticosteroid administration.
provided that theses tfs are directly involved in the regulation of module  <dig>   <dig>  and  <dig>  their combinations can produce different ways of regulation, leading to different patterns of expression on target genes among different modules and even in a module but different dosing regimens.
as a result, such observations provide a complex, but insightful, perspective into the regulation of these transcriptional modules or the complexities of their expression patterns across multiple dosing regimens.
however, it has been widely recognized that genes affected by cs include both immunosuppressive genes and metabolic genes.
upon the identification of putative transcriptional regulators, their relevance to immune response is demonstrated based on current literature evidence.
specifically, nine among the  <dig> recognized ets transcription factors are known to regulate genes involved in immunity  <cit> ; forkhead transcription factors  play a major role in the control of apoptosis  <cit> ; and especially creb has been showed as an essential factor for interactions of glucocorticoid receptors to mediate gene expression  <cit> .
a number of others are overlapped with earlier in silico studies e.g.
e2ff, egrf, hoxf, nkxh, sp1f  <cit> .
however, given that the experiment of corticosteroid administration has been studied on normal rats, the relevance to adverse effects may be more important than the relevance to immune response.
in fact, almost all enriched functions  in these transcriptional modules are relevant to metabolic side-effects .
also, due to this reason nfkb and ap- <dig> families widely considered as factors involved in inflammation are not present as direct transcriptional regulators for these sets of genes.
furthermore, we identify a number of transcriptional regulators known to be critical factors in metabolic syndrome including obesity, dyslipidemia, hypertension, insulin resistance, etc.
e.g.
rxrf  <cit> , fkhd  <cit> , sp1f  <cit> .
for instance, the deletion of rxr in mouse liver results in abnormalities of all metabolic pathways regulated by retinoid x receptors heterodimers  <cit> ; foxos, members of fkhd family, are able to increase hepatic glucose production, decrease insulin secretion, and affect glucose or lipid metabolism  <cit> .
functional characterization of critical transcriptional modules
since selected transcriptional modules consist of sets of genes that are coexpressed across all dosing regimens, we hypothesize that these genes are more likely involved in critical functions following the drug treatment.
consequently, we search for enriched functions in these modules to explore the functional effects of corticosteroids on target genes as well as evaluate the importance of the selected modules.
using arraytrack  <cit> , we first identify the gene ontology terms  that are significant in each transcriptional module .
we then classify them into super-categories  based on the branch of molecular function and biological process in the go tree.
table  <dig> lists the distribution of main functions across selected transcriptional modules.
in general, all modules are involved in metabolic processes and binding category .
some modules seem to share almost all main functions e.g.
module  <dig>   <dig> and  <dig> whereas others seem to share less e.g.
module  <dig> and  <dig>   <dig> and  <dig>  or  <dig> and  <dig>  however, they are shown to have different roles with specific functions in those main categories.
for example although module  <dig> and  <dig> are involved in metabolic processes and binding, module  <dig> is associated with rnas and nucleotides whereas module  <dig> is specialized in proteins and macro-molecules.
these functional differences  can be linked to the similarities/differences in their corresponding expression patterns, strengthening the phenomenon that they are classified as distinct transcriptional modules although their expression patterns are not intuitively separated.
however, the most important conclusion drawn from this analysis is that all these transcriptional modules consist of components that participate in metabolic processes, implying that they include genes that experience metabolic effects under corticosteroid administration.
*: details are listed in functional characterization in additional file  <dig>
+: glucocorticoid receptor element - gre binding sites; √: gre binding sites are present on the promoters of almost all genes in the corresponding function group; x: possibly because of not enough promoter information to be considered.
details are listed in functional characterization in additional file  <dig>
using arraytrack, we also searched for enriched pathways in these transcriptional modules .
a large proportion of significant pathways selected in each module are metabolic pathways of amino acid metabolism or biosynthesis, providing another support that selected transcriptional modules are critical and able to capture metabolic side effects for further analysis.
table  <dig> shows significant pathways in each transcription module.
it is generally accepted that expression levels of many cs-affected genes are mediated through the binding motifs, called gres - glucocorticoid response elements, on their control regions.
we thus examine the presence of this binding site on the promoter of genes in each of the enriched pathways in order to assess the possible effect of gre of metabolic functions.
however, such gres are short  and fairly degenerate, leading to matches occurring by chance alone thus not implying any kind of functionality.
in order to address this issue, after extracting gene promoters from the genomatix database we identified conserved regions across sets of orthologous promoters.
as a result, those matches located on these conserved regions would be more reliable estimates of functional binding sites.
although it is currently believed that gres are composed of two hexamers with a three-nucleotide random-hinge region in between, the general consensus is that towards one hexamer, namely tgttct  <cit> .
we therefore search for this motif on conserved promoter regions across orthologous promoters of the selected genes.
the results are shown in table  <dig> and detailed information is provided in additional files in functional characterization in additional file  <dig>  in general, almost all metabolic pathways contain genes with the gre binding sites, implying that these genes are more likely to be directly regulated by the complex between corticosteroids and glucocorticoid receptors.
additionally, we also examine how frequently the gre binding sites are present on the control regions of all selected genes .
furthermore, we determined that given a background set of  <dig>  randomly selected genes, the frequency of gres in a set of genes is similar to that in the random set , implying that not all genes in those modules are directly regulated by the drug and that the presence of gre binding sites on the control regions of genes in enriched pathways is very significant and not random.
CONCLUSIONS
in summary, we have proposed a systematic computational approach that can identify critical transcriptional modules coupled with their common regulatory controls under the cs administration.
the approach provides a framework to handle challenging issues related to different platforms, time-grids, genes with multiple probesets, and also different tissues if applicable.
even if the datasets across multiple conditions are present on the same platform, time-grid and tissue, the approach is still useful since genes contain multiple probesets and estimation of a single gene profile by taking the average across these probesets may lose some useful information.
however, the analysis may be limited due to the small common set of genes across different platforms.
the computational effectiveness of the approach has been demonstrated on synthetic data.
when applying to real time-series datasets, the approach not only yields critical transcriptional modules but also provides an insight into the complexities of regulation of expression patterns.
these complexities are further analyzed by techniques in promoter analysis and functional analysis to deduce useful information of transcriptional regulators and enriched metabolic pathways, providing a better understanding towards regulatory mechanisms and adverse pharmacogenomic effects of corticosteroids.
authors' contributions
ttn designed the algorithms and experiments, devised and implemented the algorithms.
rra, dcd and wjj reviewed the material and contributed to the discussion.
ipa organized the activities and structured the approach.
all authors read and approved the final manuscript.
supplementary material
additional file 1
provide links to sources of synthetic datasets used in this study.
click here for file
additional file 2
provide detailed results of six significant transcriptional modules, including gene ids, cluster ids, probeset ids and corresponding 'true' expression profiles of identified probesets in response to acute and chronic corticosteroid administration respectively.
click here for file
additional file 3
provide detailed results of functional analysis, including gene ontology, pathway enrichment, binding information, and corresponding expression profiles of transcription factors found in the study.
click here for file
acknowledgements
ttn and ipa acknowledge financial support from the nih under grant gm <dig> and the epa under grant gad r 832721- <dig>  rra, dcd and wjj acknowledge financial support from the nih under grant gm  <dig>

 glucocorticoids are a class of steroid hormones present in almost every animal cell, playing a central role in a wide range of physiological responses. because of their potent anti-inflammatory and immunosuppressive effects, synthetic glucocorticoids referred as corticosteroids have been used widely in pharmacology as a therapeutic option for a wide range of autoimmune and inflammatory diseases. however, beneficial effects are derived from magnifying the physiological actions of endogenous glucocorticoids, causing a variety of side effects following long-term treatment with this class of drugs e.g. hyperglycemia, dyslipidemia, arteriosclerosis, muscle wasting, and osteoporosis. the physiological and pharmacological effects of corticosteroids are complex and manifest themselves with expression changes of many genes across multiple tissues. it has been observed that even in a single tissue different dosing regimens of cs administration can induce different patterns of expression. as such genes with similar expression profiles under acute cs administration may not exhibit similar expression patterns during continuous infusion, pointing to the possibility of alternative regulatory mechanisms. therefore, a better understanding of corticosteroid pharmacogenomic effects from multiple dosing regimens are very valuable not only to reveal the transcriptional dynamics under different patterns of input perturbations but also to provide an insight into the underlying molecular mechanisms of action, for both beneficial and detrimental effects, and thus for the optimization of clinical therapies. it has been noted that genes affected by cs include both immunosuppressive genes, mostly associated with therapeutic effects, and metabolic genes often associated with adverse effects whose regulation is mainly controlled by glucocorticoid receptor gene mediated pathways. unbound cs binds with cytosolic free glucocorticoid receptors releasing it from the heat shock complex allowing dimerization and translocation into the nucleus where it binds to glucocorticoid response element of the target genes, leading to enhancement or inhibition of the target gene expression. as a result, long-term treatment with corticosteroids results in sustained up- or down-regulation of numerous genes, leading to a new steady state which might be the basis for occurrence of adverse effects. however, it has also been noted that chronic infusion of cs causes a sustained down-regulation of the receptor. while several alternative mechanisms have been proposed it is still not understood why drug effects remain strong although gr mrna is down-regulated to the point of almost being eliminated. a plausible explanation is that besides direct regulation through gre binding sites in the 5' regulatory regions of genes, there are changes in expression that are also the indirect results of effects due to changes in expression of transcription factors that act as secondary biosignals directly or indirectly modulating the transcription of genes. thus, along with identification of expression patterns, predicted regulatory control structures are also an essential source of information towards understanding corticosteroid effects. in this study we address the question as to whether significant transcriptional modules coupled with complex patterns of mrna changes across multiple dosing regimens of corticosteroids and their common regulatory controls can be computationally identified. hypothetically, transcriptional modules that are significantly coexpressed under different dosing regimens will be important gene clusters for further analysis towards better understanding of both beneficial and adverse effects of corticosteroids, especially the metabolic side-effects since these patterns are survived under a long-term treatment of corticosteroids. the hypothesis explored here is that if two or more genes have the same temporal expression profiles in response to different dosing regimens, they are more likely to share some common regulatory mechanisms. the liver was selected because of its major role in both the physiological efficacious and adverse effects of corticosteroids e.g. altering the expression of serum proteins that regulate immune/inflammatory responses, enhancing the expression of liver enzymes involved in metabolic effects. however, rich in vivo datasets of pharmacological time-series across multiple dosing regimens are often obtained from different microarray platforms and time-sets, leading to a problematical issue for computational analysis. as an example, in a study comparing normal and chronic lymphocytic leukemia b-cells, wang et al. identified only differentially expressed genes across all three studies, when combining results from three different platforms, while there are at least differentially expressed genes in each individual platform. in general, there are two important issues relevant to the analysis of data derived from different platforms: genes may be present in one platform but not in the other, and genes present on both platforms may not be represented by the same probes. since different microarray platforms do not contain the same probesets, and even do not have a similar hardware design and sample processing protocols, standard analyses may not yield comparable expression level quantifications across platforms, leading to many challenges for computational models aiming at the analysis of microarray data from heterogeneous sources. a number of approaches have been proposed and are generally classified into two main categories: integrate raw expression profiles from different studies into one dataset so that available computational models can be directly applied, and develop and/or utilize a unitless statistic as a primary analysis and then combine the result through a meta-level analysis. the former category can be further divided into two sub-classes, namely combining raw data through a normalization and/or transformation procedure and pooling raw information from common probes that can be mapped to the same unigene clusters or full-length mrna transcripts. however, these approaches are not general enough to make data from different platforms fully compatible. since combining data across different platforms remains a serious challenge, meta-analysis - the second category - has been identified as a more popular technique in order to combine results, and thus data, from a number of independent studies. the assumption here is that while the raw expression levels from different platforms may not be comparable, the results of the primary analysis should be. however, almost all prior studies has focused on the discovery of genes that are differentially expressed in conjunction with standard models such as effect size models, bayesian models. consequently, in order to identify significant clusters of genes that share common expression patterns across multiple dosing regimens, we extend our prior study in the aspect of producing an agreement matrix that describes the agreement levels of co-expression of genes across multiple conditions and successively searching clusterable subsets to infer all such gene clusters. the approach follows the concept of meta-analysis to avoid the limitation of incompatible data across multiple datasets from different platforms. the unitless statistic, expressing the confidence level of co-expression is the agreement level of cluster assignments drawn from multiple clustering runs. there remain a number of open critical issues associated with a single clustering run, and thus consensus clustering coupled with the examination of am distribution has been designed with the aims of reducing aforementioned limitations. once the am is obtained for each condition independently, an average agreement matrix is calculated to estimate the confidence levels of coexpression between genes across multiple conditions, thus combining data from different datasets into a single input for the next analysis. for the analysis at the meta-level, we extend the selection and clustering processes to identify all possible clusters of genes that are highly coexpressed with the average am above as the input. as such these clusters of genes will share common patterns of expression across multiple dosing regimens. additionally, due to the selection of all possible patterns of expression several clusters may have similar expression patterns and thus we also provide a heuristic as an optional procedure to merge such similar clusters based on a criterion of maximizing the total homogeneity and separation of selected clusters. subsequently, we analyze promoter regions of genes in every cluster in order to predict putative transcriptional regulators, aiming at providing another source of information towards better understanding those complex patterns of expression and the underlying regulatory mechanisms of corticosteroid effects. our results demonstrate that the proposed computational approach is highly effective on both synthetic and real data. when applying the approach to real time-series datasets, selected patterns of transcriptional responses are enriched in a biological sense with relevant putative-regulatory controls and significant metabolic pathways in each pattern. computational results are further validated predicated upon literature evidence. methods datasets synthetic data a number of synthetic datasets from the open literature are utilized to assess our approach for finding common sets of genes that are highly coexpressed across multiple conditions. specifically, we used a series of four high-noise 20-timepoint sine-format synthetic datasets with different number of replicates at each time-point from. each dataset contain five separate sets with genes allocated equally in classes, each of which contains the same list of genes but has different patterns across five conditions. for each set, in the first step the data are generated according to an artificial pattern Φ which shows the values for gene i at time-point t in cluster l; four of six clusters follow the sine function i.e. Φ = sin, and the other two follow the non-periodic linear function = t/ and Φ = -t/20), i =..., t =..., l =..., in the second step, let x be the error-added value for gene i, time-point t and repeat; x is randomly drawn from a normal distribution n where μ is the value of the synthetic pattern Φ and σ is equal to λσit. high-noise synthetic data are generated with λ =. datasets are downloaded from the links provided in synthetic data in additional file acute corticosteroid administration forty-seven male adx wistar rats weighting from to g underwent right jugular vein cannulation under light ether anesthesia day before the study. forty-three rats were injected with a single intravenous bolus dose of methylprednisolone of mg/kg. animals were sacrificed by exsanguinations under anesthesia and liver samples were harvested at,,,, and after dosing. the sampling time points were selected based on preliminary studies describing gr dynamics and enzyme induction in liver. four untreated rats were sacrificed at random times and nominally considered as h controls. the gene expression was obtained via the affymetrix rg-u34a array which consists of probesets. the data are publicly available through the geo omnibus database under the accession number gds chronic corticosteroid administration in a similar experiment model, forty rats were given mg/kg/hr infusions of mpl over h via an azlet pump. the pump drug solutions were prepared for each rat based on its predose body weight. animals were sacrificed at various times up to days; specifically the time-points included are and h. a control group of four animals was implanted with a saline-filled pump and killed at various times throughout the 7-day study period. unlike the previous experiment, the microarray platform for this dataset is the rae230a which consists of probesets. the data are publicly available through the geo omnibus database under the accession number gds all protocols followed the principles of laboratory animal care and were approved by the university at buffalo institutional animal care and use committee. identifying critical transcriptional modules the general computational problem can be briefly defined as follows. we are given a set of n genes g={gi}i=1n and k conditions. for each condition k, every gene is characterized by one or more time-series expression profiles with rki corresponding probesets over tk time-points gk={gki}i=1n,gki={gkir,r∈rki},gkir={gkitr}t=1tk,k=...,k. the question then becomes to search for clusters of genes that are highly coexpressed across all k conditions with a confidence level δ. the term'highly coexpressed' is used in the sense that ∀gi,gj∈c,1k∑k=1kpk≥δ where c denotes a, yet to be determined, cluster and pk is the confidence level that two gene profiles i and j are clustered together in condition k; a gene profile includes sets of corresponding probesets rki of gene i in condition k, k =...,k. the subscripts {i, j}, t, k, r indicate the {gene id}, time, condition, and probesets respectively. it should be also noted that in this work, we used three different terms to refer to the same object:'cluster' when designing the algorithm,'pattern' when exhibiting the expression changes, and'module' when charactering the biological function. the framework contains several step displayed in figure the pre-processing step each dataset is pre-filtered to identify differentially expressed probesets. since we would like to identify gene clusters with common expression patterns across multiple conditions, input datasets must contain the same set of genes. thus using the respective platform information, probesets in each dataset are mapped to a list of genes and then the intersection across those gene lists is evaluated to extract a common set of genes which are differentially expressed across multiple conditions. however, genes are sometimes characterized by multiple probesets whose expression profiles may be similar or sometimes different, but not identical. these probesets can be considered as replicates of expression profiles for a gene and thus taking an average expression profiles across all these probesets to characterize for the expression profile of the gene may lose useful potential information. therefore, from the common set of genes we re-map genes to corresponding probesets in each dataset with the respective platform before starting the analysis. construction of the agreement matrix the agreement matrix quantifies the likelihood that two objects are assigned to the same cluster. if m clustering runs are performed on the data, each entry will show the frequency with which two objects are assigned to the same cluster over'm' clustering runs. the am entries are defined as follows: mxy=1m∑h=1mmwhere m={1 if x and y are clustered together when running method m0 otherwise in addition to the various clustering methods that were utilized, different distance metrics are also explored in order to attenuate the biases associated with individual distance metrics. in our implementation, we are using hierarchical clustering, divisive analysis clustering, fuzzy analysis clustering, partitioning around medoid with pearson correlation and manhattan metric, k-means, fuzzy c-means, self-organizing map, and model-based clustering with euclidean metric as the core clustering methods. since clustering results are highly dependent on the input number of clusters, the sensitivity of the am as a function of nc was examined to find a'suggestive' number of clusters for each particular dataset. after identifying nc* based on the procedure in our prior work, all clustering runs are repeated with nc* to produce the final am for further analysis. if two probesets are clustered together, it is implied that their expression profiles are similar under a specific condition k. therefore, the fraction of times they are clustered together over multiple clustering runs can be considered as the confidence level that they are coexpressed since mxy, by construction, aims at eliminating method-specific biases and assumptions. subsequently, we calculate the average agreement levels between sets of corresponding probesets of any two genes to estimate the confidence level that those two genes are coexpressed in a specific condition. the am entries in condition k is re-estimated as follows amij=1|rki||rkj|∑x∈rki∑y∈rkjmxy,i,j=...,n with the assumption that the unitless statistics, i.e. the confidence level of co-expression, is comparable across multiple conditions and different platforms, we estimate the confidence level of co-expression between two genes across multiple conditions by taking the average. while combining raw data remain challenges, the estimation of a unitless statistics provides a simple but efficient combination of heterogeneous data for further analyses. amij=1k∑k=1kamij,i,j=...,n as a result, we obtain an agreement matrix whose entries exhibit a quantity that shows how confident genes are coexpressed. this will be the input for the selection and clustering process. selection and clustering with the hypothesis that the more clusterable the data is the more biologically relevant it is, we applied our previously proposed procedure to select a more'hypothetically clusterable' subset from the entire set of genes. the main hypothesis underlying the selection is that am entries associated with genes at the'hypothetical' core of an expression pattern will be consistently grouped together over multiple clustering runs. this should be manifested by high corresponding values in the am, whereas genes belonging to the'hypothetical' core of two clearly distinct clusters are associated with consistently low am entries. on the contrary, cluster assignments associated with genes at cluster boundaries or between clusters would be very sensitive to the method used and thus they would have relatively moderate agreement levels with other genes. as a result, with a user-defined confidence level δ genes associated with moderate am entries <δ) are eliminated to produce a more'clusterable' subset of genes. the process starts removing genes associated with the highest number of moderate am entries and then updates the am for the next loop until no moderate am entry exists. the corresponding subset of genes is now considered as a'hypothetically clusterable' subset since any two genes are highly coexpressed or non-coexpressed with the confidence level at least δ. subsequently, we used the concept of consensus clustering to divide the subset of genes into a number of clusters by applying the hierarchical clustering with the selected am as input data. the algorithm starts with every gene filling a cluster and then grouping two clusters into a new one for each loop so that any pair of genes belonging to a new cluster always has an agreement level greater than or equal to δ. the iteration is stopped when no more new cluster is formed. however, since there should be existed clusters of genes located closely to other clusters in the data and the input number of clusters for the core analysis is only a suggestive one, those clusters may not be completely separated. as a result, although genes that belong to those clusters are identified as highly coexpressed, their relationship to genes in other clusters cannot be uniquely determined. therefore, some significant clusters may be not included in the selected subset due to the constraint of'clusterable' selection. since we would like to obtain all significant patterns of expression, the procedure of selection and clustering is repeated on the removed domain. the removed domain consists of a set of unselected genes whose coexpression levels are still high as quantified agreement levels in the original am. after extracting the sub-agreement matrix corresponding to the set of unselected genes, the entire process of selection and clustering is applied again with the same confidence level δ as before. the procedure is reiterated until no more clusters of genes are recognized. figure presents the pseudo-code of the procedure and an example to illustrate the process. furthermore, due to the nature of clustering, trivial clusters may be identified in the final results. in order to exclude such trivial clusters, each cluster c is assigned with a simple hypothetical quantity called'cluster significance' which represents how large the cluster is in this study. we then create the distribution of cluster significance on random data to estimate the cluster significance threshold corresponding to a user-defined threshold p-value for cluster selection. without loss of generality we select k = for the random data and assume that each probeset in the mapped dataset d correspond to a gene. the suggestive number of clusters nc* for d is searched with the process in. subsequently, d is randomly resampled a number of times, for each of which the entire process starting from building the am with the same nc* to extracting clusters of highly coexpressed genes is repeated and the resulting random clusters are returned. after that, the procedure estimates the cluster significance, which is simply the cluster-size in this study, for these random clusters and constructs a distribution of cluster significance. the corresponding p-value of cluster significance cs is defined as the number of random clusters whose significance is greater than cs over the total number of clusters identified in nr resamples: pvalue =∑clusters with cluster significance≥cs∑random clusters. as a result, given a threshold p-value, the corresponding cluster significance threshold is inferred and only clusters with significance greater than this threshold are selected. merging similar patterns because of the nature of the approach, it is quite reasonable to expect that the clustering process can break out patterns of expression into several sub-patterns. thus, we repeat the process on the eliminated domains to extract all possible significant clusters, resulting in that several clusters may have a similar expression pattern but are separated into two or more clusters. because cluster homogeneity reflects the similarity of expression profiles within a given cluster and cluster separation quantifies how effectively expression profiles are discriminated, we provide an optional heuristic in order to merge similar patterns together according to the criterion of maximizing sum of homogeneity and separation of all final output clusters. starting with all significantly selected clusters, the procedure searches for a combination of two similar patterns so that their combination will generate a maximal increase of the sum of homogeneity and separation of all current clusters after merging those two patterns. the process is repeated until no more combinations are found i.e. any new combination always reduces the sum of homogeneity and separation. eventually, a list of significant expression patterns that characterize the underlying transcriptional responses is generated. the metric used during the evaluation of the heuristic is quantified as follows: maxc{homogeneity+separation}=maxc{1k n∑k∑phk+2k n∑k∑p<qsk};hk=∑gi,gj∈cp;gi≠gjsim||cp||;sk=∑gi∈ck;gj∈cldis||cp||  ||cq||;sim=1‖rki‖‖rkj‖similarity;dis=1‖rki‖ ‖rkj‖dissimilarity;p=...,n,q=...,n,k=...,k,r∈rki,r'∈rkj where c is the current set of selected clusters c={cp}p=1n and n is the current number of clusters; hk is the homogeneity of cluster cp in condition k and sk is the separation between cluster cp and cq in condition k; sim and dis are the average similarity and dissimilarity respectively between all probesets of gene i and gene j in condition k. similarity is measured by the pearson correlation coefficient and dissimilarity is estimated by the pearson correlation distance. predicting putative transcriptional regulators promoter identification promoters of genes are extracted from a rich database of promoter information with a default length if there is no experimentally defined length as suggested by genomatix. in order to identify putative transcriptional regulators, we explore the basic underlying assumption of comparative genomics which states that functional regions evolve in a constrained fashion and thus at a lower rate than non-functional regions. it implies that conserved regions in a set of orthologous sequences survive due to their special functional implications i.e. tfbss located on these conserved regions will be more promising as functional binding sites and thus associated tf families are more relevant to our context. therefore, each promoter is characterized by a set of promoters from orthologous genes of other vertebrate species, if available. to be consistent in the search for conserved regions on promoter sequences in order to identify putative transcription factor binding sites we eliminate those that do not consist of more than two orthologous promoters. putative functional binding sites in order to identify conserved regions for each promoter dialign was used to perform multiple sequence alignments with the input sequences including each sequence as well as its orthologous promoters. dialign was selected because it has many applications in comparative genomics. also, a recent benchmark study for the alignment of non-coding dna sequences has concluded that it can produce alignments with high sensitivity as well as specificity to detect constrained sites. following the alignment among orthologous promoter sequences, we relied on the conserved scores returning from dialign to locate conserved regions which are defined as sub-sequences that are longer than bp and continuously scored greater than the average score of all the alignment's conserved scores. we next apply matinspector to scan for all physical tfbss and only those that overlap with the conserved regions selected above are kept for further analysis. we used a common core similarity and utilized the optimal matrix similarity threshold for each position weight matrix suggested from matbase, genomatix which ensure that a minimum number of matches are found in non-regulatory sequences i.e. the false positive matches is minimized. however, a gene may have multiple alternative promoters and virtually in all cases, it is not known which promoter of the gene is activated. therefore, all putative tfbss detected from all alternative promoters of a gene are considered as candidates to infer putative transcriptional regulators for the gene. subsequently, we estimate the common level of each candidate above in each corresponding module and select those tfbss present more than a common threshold. associated tf families with those selected tfbss are inferred and considered as transcriptional regulators for corresponding transcriptional modules.  method evaluation on synthetic data in order to evaluate the effectiveness of the proposed approach, we use synthetic data with known class structure as described earlier. the process of evaluation is repeated four times with four different datasets that are created with different number of replicates for each time-point. in each time, we use five high-noise sets as the data for five conditions, each of which has genes distributed across clusters; each cluster has different patterns across five conditions but has the same set of gene ids. we set the same parameters for all evaluation in this study and also for the analysis on real time-series datasets, specifically the confidence level of coexpression δ = 70% and p-value = for the selection of significant clusters. furthermore, the testing process on synthetic data is done without the merging option. without loss of generality, we assume that each gene has only one probeset in this evaluation. the performance of the approach is assessed through its ability to recover the number of cluster structures and the list of gene ids identified in each cluster. we use the adjusted rand index which is a statistic that measures the extent of concurrence between the clustering results and the underlying known class structure to evaluate the approach's performance in identifying gene clusters that are coexpressed across multiple conditions. the larger the rand index is, the higher the agreement between the results and prior knowledge of cluster structure. the number of selected genes, recovered cluster structures and the accuracy on the selected domain are listed in table *: the accuracy is only estimated on the selected domain due to the fact that these datasets are high-noise synthetic data, some cluster structures may be missed when there is only one measurement at each time-point. however, when the number of replicates is increased, the number of cluster structures is recovered. as discussed in our previous study, this is a reasonable observation due to the effect of gene expression replicates on clustering performance. additionally, we also examine an alternative approach which is more intuitive in identifying gene clusters that are coexpressed across multiple conditions. instead of performing a meta-analysis to avoid the limitation of incompatible data across different platforms, we can separately identify significant clusters of genes that are coexpressed in each condition and then obtain the intersection among these gene clusters across all conditions. in this experiment, we used pam, mclust, and consensus clustering as standard single clustering methods to identify clusters in each set of data, for which nc* = is the input number of clusters. we then simply took the intersection between clusters from set to set and only keep those clusters that contain more than genes as significant clusters for the final estimation of accuracy. the number of selected genes, number of clusters, and accuracy on the selected domain are listed in table +: only clusters with more than genes;*: the accuracy is estimated on the selected domain; the number of clusters and the accuracy are formatted as'before merging' |'after merging' in general, this approach selects a smaller number of genes with an equal or greater number of cluster structures, resulting in lower accuracy. as an example, in each set of dataset there are two cluster structures that are not clearly identifiable. as a result, a single clustering methods may fail to properly separate them in each set, leading to the situation where the intersection between clusters from set to set divides those cluster structures into many sub-clusters with a small number of genes. on the contrary, by taking the average of the co-expression levels across multiple sets, the relationship of whether two genes are coexpressed across multiple conditions can be recovered. consequently, our proposed approach is more advantage, resulting in a final highly correct classification as illustrated in table furthermore, since this simpler alternative approach produces many resulting clusters, we also attempted to apply the proposed merging process to reduce the number of clusters as well as improve the accuracy if applicable. however, its testing performance is still not as high as that of our proposed approach although we do not apply the merging process for the proposed approach in this test. additionally, the alternative approach is highly sensitive with the initial number of clusters. for instance, when we constantly set nc* = and test on dataset without the merging option our approach still recovers the correct number of cluster structures with high accuracy: = whereas'pam' approach yields,'mclust' provides, and'consclust' does. since this information is not available for real datasets, the more sensitive with it the less robust the approach is. therefore, by taking the average of the co-expression levels between two genes across multiple datasets, our proposed approach provides more robust results. acute vs. chronic cs administration for the analysis of corticosteroid administration, the pre-processing step is performed to provide corresponding mapped datasets. the datasets are first filtered for differentially expressed probesets using anova technique implemented in r and also customized by our previous work for easy uses. probesets in the acute and probesets in the chronic are selected for further analysis. to obtain the common set of genes across two conditions, these probesets are mapped into sets of genes based on the corresponding platform information. differentially expressed probesets in the acute are mapped into a set of genes and probesets in the chronic are mapped into another set of genes. the intersection of these two gene-sets yields genes in common for both dosing regimens. from this common gene set, the re-mapping process subsequently returns a corresponding set of probesets for the acute and a set of probesets for the chronic data. all datasets are pre-processed with the model in our previous study to estimate the'true' expression profiles that are integrated with potential information in replicates instead of simply taking the average expression profiles. the suggestive number of clusters nc* for both datasets is subsequently, we apply the proposed approach with the merging option to the intersection set of genes that are affected by corticosteroid administration across the two dosing regimens. we obtain significant clusters with genes in total. these clusters are hypothesized to be transcriptional modules which share common regulatory mechanisms since they consist of genes that exhibit similar expression patterns in both acute and chronic dosing regimen. table shows the distribution of these genes over six modules and also briefly describes how the pattern of expression changes. although genes may exhibit simple or complex patterns of expression during corticosteroid administration, we crudely classify those patterns into up- or down- with one or two phases of regulation. *: patterns consist of one-phase regulation, two-phase regulation or simply up- or down-regulation. details are listed in expression patterns in additional file a detailed description for patterns of these transcriptional modules is shown in figure with the average expression patterns of all probesets clustered in each module following acute and chronic dosing. in brief, transcriptional module is characterized by one-phase regulation in acute but two phases in chronic dosing. genes in this module exhibited a fast and robust decline in mrna, which reached its peak between h and h, and returned to the baseline after about h. however, when mpl is infused this set of genes shows a more complex pattern involving both enhanced and suppressed regulation. although a strong down regulation is observed at the beginning, it is subsequently followed by a sharp induction with the maximum around h and then gradually returned to the baseline indicating some kind of possible tolerance. the second transcriptional module shows a similar pattern of expression in both acute and chronic regimen with two phases of regulation. genes in this module exhibit an early up-regulation and reached their corresponding peaks at around h in the acute and h in the chronic. subsequently, both profiles denote a clear down-regulation and possible slight fluctuation before returning to base line. an interesting dynamics is observed in the genes of transcriptional module in the acute dosing, the genes in this module clearly exhibit an expression pattern with two phases of regulation. yet, in chronic administration they exhibited an early transient decline in mrna followed by robust, sustained, up-regulation. similar to module is the transcriptional dynamics exhibited by transcriptional module characterized by an early induction with a maximum at h in the acute and h in the chronic. a typical pattern with down regulation for both acute and chronic administration is illustrated by transcriptional module. however, genes in the acute regimen exhibited a fluctuated repression with a maximum at around h and then followed by an induction to return to the baseline as late as h. meanwhile, genes in the chronic regimen characterized a pattern with a slightly transient up-regulation followed by a sustained down-regulation and eventual convergence to a new steady state in the presence of the drug. the last transcriptional module has a similar acute pattern of expression with two phases of regulation as that of module however, in the chronic regimen after falling to a value below the baseline this set of genes was further sustained a slight suppression. while comparing these expression patterns, we observe that modules and have similar expression patterns in acute or chronic with a slight difference in the other dosing regimen. although the difference is not large enough to be intuitively recognized, the merging process could not merge them together, implying that the difference is significant. furthermore, the separation of these expression patterns is also reinforced with different functional characteristics which will be illustrated below. in summary, selected transcriptional modules exhibit a number of typical expression patterns under corticosteroid administration. the pattern can be simply expressed as an up- or down- regulation or as a more complex one with two phases of regulation plus some fluctuation. putative transcriptional regulators of critical transcriptional modules it has been widely accepted that after corticosteroids bind to cytosolic glucocorticoid receptors, the activated steroid-receptor complex is rapidly translocated into the nucleus where it can alter the expression of target genes. however, the drug seems to be cleared within about h following a bolus injection, suggesting that the mrna levels of cs-target genes will return to the base line after that. in the contrast, the drug will reach and remain to a stable steady state after h in the chronic administration. yet, the gr is greatly diminished in response to corticosteroids, suggesting that the mrna levels of cs-target genes in the chronic regimen should also return to the base line. this mechanism is corresponding to the first-phase regulation of target genes. however, almost all chronic patterns involve two phases of regulation and some are only half-phase patterns i.e. persistent up or down without returning to the baseline. these complexities in expression patterns of cs-target genes can be explained by a number of possibilities previous studies have shown, including multiple gr isoforms, multiple gres with different affinities to the drug receptor complex, or some other receptors that can mediate the effect of corticosteroids and thus affected genes in this case can reach a new steady state in the presence of the drug. however, another possibility is a mechanism that results in the regulation of secondary biosignals which transcription factors are the most potential factors. after affected by corticosteroids, they in turn further modulate the expression of glucocorticoid-regulated genes as a continuing cascade of events that were initiated by the drug. as a result, this possibility suggests a possible interpretation of the complexities in expression changes of multiple cs-target genes with the second phase of regulation. in order to reveal some underlying regulatory mechanism of these selected transcriptional modules, we start analyzing the promoter regions of genes to search for significant putative transcriptional regulators as well as possible relationships of regulation. the hypothesis we explore here is that if two or more genes have similar temporal profiles in response to multiple dosing regimens, they are more likely to share some common regulatory mechanisms. for the genes in six transcriptional modules, we extract rattus norvegicus's promoter sequences, of which we only keep genes with promoters that include sufficient information of orthologous promoters for further analysis. figure shows the identified putative regulation between tf families and transcriptional modules. this finding highlights the possibility that secondary biosignals are involved in the regulatory complexities of expression changes for cs-affected genes. almost all suggested tf families do consist of transcription factor members that are recognized as differentially expressed genes in one or both dosing regimens. since transcription factors are characterized by pleiotropic effects, it is reasonable to observe a significant overlap across various transcriptional modules. while comparing these regulatory combinations, we observe that some tf families seem to be common regulators for all modules whereas some are very specific to particular modules. this could possibly explain the difference among the expression patterns of these modules. it is likely that the more similar the expression pattern of clusters the more likely they share a larger fraction of common regulators, e.g. tf families in this case. for example, there are a large number of transcriptional regulators that are common between modules and but it seems little overlap exists between the transcriptional regulators of modules and and and except common regulators on the top of the figure. furthermore, one tf family consists of several tf members whose expression patterns may be similar or different across multiple dosing regimens. figure shows two examples, each of which is represented by two tf members whose expression patterns are similar or different between the acute and the chronic corticosteroid administration. provided that theses tfs are directly involved in the regulation of module and their combinations can produce different ways of regulation, leading to different patterns of expression on target genes among different modules and even in a module but different dosing regimens. as a result, such observations provide a complex, but insightful, perspective into the regulation of these transcriptional modules or the complexities of their expression patterns across multiple dosing regimens. however, it has been widely recognized that genes affected by cs include both immunosuppressive genes and metabolic genes. upon the identification of putative transcriptional regulators, their relevance to immune response is demonstrated based on current literature evidence. specifically, nine among the recognized ets transcription factors are known to regulate genes involved in immunity; forkhead transcription factors play a major role in the control of apoptosis; and especially creb has been showed as an essential factor for interactions of glucocorticoid receptors to mediate gene expression. a number of others are overlapped with earlier in silico studies e.g. e2ff, egrf, hoxf, nkxh, sp1f. however, given that the experiment of corticosteroid administration has been studied on normal rats, the relevance to adverse effects may be more important than the relevance to immune response. in fact, almost all enriched functions in these transcriptional modules are relevant to metabolic side-effects. also, due to this reason nfkb and ap- families widely considered as factors involved in inflammation are not present as direct transcriptional regulators for these sets of genes. furthermore, we identify a number of transcriptional regulators known to be critical factors in metabolic syndrome including obesity, dyslipidemia, hypertension, insulin resistance, etc. e.g. rxrf, fkhd, sp1f. for instance, the deletion of rxr in mouse liver results in abnormalities of all metabolic pathways regulated by retinoid x receptors heterodimers; foxos, members of fkhd family, are able to increase hepatic glucose production, decrease insulin secretion, and affect glucose or lipid metabolism. functional characterization of critical transcriptional modules since selected transcriptional modules consist of sets of genes that are coexpressed across all dosing regimens, we hypothesize that these genes are more likely involved in critical functions following the drug treatment. consequently, we search for enriched functions in these modules to explore the functional effects of corticosteroids on target genes as well as evaluate the importance of the selected modules. using arraytrack, we first identify the gene ontology terms that are significant in each transcriptional module. we then classify them into super-categories based on the branch of molecular function and biological process in the go tree. table lists the distribution of main functions across selected transcriptional modules. in general, all modules are involved in metabolic processes and binding category. some modules seem to share almost all main functions e.g. module and whereas others seem to share less e.g. module and and or and however, they are shown to have different roles with specific functions in those main categories. for example although module and are involved in metabolic processes and binding, module is associated with rnas and nucleotides whereas module is specialized in proteins and macro-molecules. these functional differences can be linked to the similarities/differences in their corresponding expression patterns, strengthening the phenomenon that they are classified as distinct transcriptional modules although their expression patterns are not intuitively separated. however, the most important conclusion drawn from this analysis is that all these transcriptional modules consist of components that participate in metabolic processes, implying that they include genes that experience metabolic effects under corticosteroid administration. *: details are listed in functional characterization in additional file +: glucocorticoid receptor element - gre binding sites; √: gre binding sites are present on the promoters of almost all genes in the corresponding function group; x: possibly because of not enough promoter information to be considered. details are listed in functional characterization in additional file using arraytrack, we also searched for enriched pathways in these transcriptional modules. a large proportion of significant pathways selected in each module are metabolic pathways of amino acid metabolism or biosynthesis, providing another support that selected transcriptional modules are critical and able to capture metabolic side effects for further analysis. table shows significant pathways in each transcription module. it is generally accepted that expression levels of many cs-affected genes are mediated through the binding motifs, called gres - glucocorticoid response elements, on their control regions. we thus examine the presence of this binding site on the promoter of genes in each of the enriched pathways in order to assess the possible effect of gre of metabolic functions. however, such gres are short and fairly degenerate, leading to matches occurring by chance alone thus not implying any kind of functionality. in order to address this issue, after extracting gene promoters from the genomatix database we identified conserved regions across sets of orthologous promoters. as a result, those matches located on these conserved regions would be more reliable estimates of functional binding sites. although it is currently believed that gres are composed of two hexamers with a three-nucleotide random-hinge region in between, the general consensus is that towards one hexamer, namely tgttct. we therefore search for this motif on conserved promoter regions across orthologous promoters of the selected genes. the results are shown in table and detailed information is provided in additional files in functional characterization in additional file in general, almost all metabolic pathways contain genes with the gre binding sites, implying that these genes are more likely to be directly regulated by the complex between corticosteroids and glucocorticoid receptors. additionally, we also examine how frequently the gre binding sites are present on the control regions of all selected genes. furthermore, we determined that given a background set of randomly selected genes, the frequency of gres in a set of genes is similar to that in the random set, implying that not all genes in those modules are directly regulated by the drug and that the presence of gre binding sites on the control regions of genes in enriched pathways is very significant and not random.  in summary, we have proposed a systematic computational approach that can identify critical transcriptional modules coupled with their common regulatory controls under the cs administration. the approach provides a framework to handle challenging issues related to different platforms, time-grids, genes with multiple probesets, and also different tissues if applicable. even if the datasets across multiple conditions are present on the same platform, time-grid and tissue, the approach is still useful since genes contain multiple probesets and estimation of a single gene profile by taking the average across these probesets may lose some useful information. however, the analysis may be limited due to the small common set of genes across different platforms. the computational effectiveness of the approach has been demonstrated on synthetic data. when applying to real time-series datasets, the approach not only yields critical transcriptional modules but also provides an insight into the complexities of regulation of expression patterns. these complexities are further analyzed by techniques in promoter analysis and functional analysis to deduce useful information of transcriptional regulators and enriched metabolic pathways, providing a better understanding towards regulatory mechanisms and adverse pharmacogenomic effects of corticosteroids. authors' contributions ttn designed the algorithms and experiments, devised and implemented the algorithms. rra, dcd and wjj reviewed the material and contributed to the discussion. ipa organized the activities and structured the approach. all authors read and approved the final manuscript. supplementary material additional file 1 provide links to sources of synthetic datasets used in this study. click here for file additional file 2 provide detailed results of six significant transcriptional modules, including gene ids, cluster ids, probeset ids and corresponding'true' expression profiles of identified probesets in response to acute and chronic corticosteroid administration respectively. click here for file additional file 3 provide detailed results of functional analysis, including gene ontology, pathway enrichment, binding information, and corresponding expression profiles of transcription factors found in the study. click here for file acknowledgements ttn and ipa acknowledge financial support from the nih under grant gm and the epa under grant gad r 832721- rra, dcd and wjj acknowledge financial support from the nih under grant gm 
the study addresses two issues, including identifying significant transcriptional modules coupled with dynamic expression patterns and predicting relevant common transcriptional controls to better understand the underlying mechanisms of corticosteroid adverse effects. following the orientation of meta-analysis, an extended computational approach that explores the concept of agreement matrix from consensus clustering has been proposed with the aims of identifying gene clusters that share common expression patterns across multiple dosing regimens as well as handling challenges in the analysis of microarray data from heterogeneous sources, e.g. relevant putative transcriptional regulators are also predicted to provide another source of information towards better understanding the complexities of expression patterns and the underlying regulatory mechanisms of those modules. therefore, rich in vivo datasets of pharmacological time-series with two dosing regimens sampled from rat liver are examined for temporal patterns of changes in gene expression and their regulatory commonalities. nevertheless, even in a single tissue different methods of corticosteroid administration can induce different patterns of expression and regulatory control structures. functional analysis reveals that virtually all enriched functions in these modules are shown to be related to metabolic processes, implying the importance of these modules in adverse effects under the administration of corticosteroids. comprehensively understanding corticosteroid pharmacogenomic effects is an essential step towards an insight into the underlying molecular mechanisms for both beneficial and detrimental clinical effects. we have proposed a framework to identify significant coexpressed clusters of genes across multiple conditions experimented from different microarray platforms, time-grids, and also tissues if applicable. six significant transcriptional modules coupled with typical patterns of expression have been identified. analysis on rich in vivo datasets of corticosteroid time-series yielded significant insights into the pharmacogenomic effects of corticosteroids, especially the relevance to metabolic side-effects. this has been illustrated through enriched metabolic functions in those transcriptional modules and the presence of gre binding motifs in those enriched pathways, providing significant modules for further analysis on pharmacogenomic corticosteroid effects. different platforms and time-grids in this study.
 comprehensively understanding corticosteroid pharmacogenomic effects is an essential step towards an insight into the underlying molecular mechanisms for both beneficial and detrimental clinical effects. nevertheless, even in a single tissue different methods of corticosteroid administration can induce different patterns of expression and regulatory control structures. therefore, rich in vivo datasets of pharmacological time-series with two dosing regimens sampled from rat liver are examined for temporal patterns of changes in gene expression and their regulatory commonalities.  the study addresses two issues, including identifying significant transcriptional modules coupled with dynamic expression patterns and predicting relevant common transcriptional controls to better understand the underlying mechanisms of corticosteroid adverse effects. following the orientation of meta-analysis, an extended computational approach that explores the concept of agreement matrix from consensus clustering has been proposed with the aims of identifying gene clusters that share common expression patterns across multiple dosing regimens as well as handling challenges in the analysis of microarray data from heterogeneous sources, e.g. different platforms and time-grids in this study. six significant transcriptional modules coupled with typical patterns of expression have been identified. functional analysis reveals that virtually all enriched functions in these modules are shown to be related to metabolic processes, implying the importance of these modules in adverse effects under the administration of corticosteroids. relevant putative transcriptional regulators are also predicted to provide another source of information towards better understanding the complexities of expression patterns and the underlying regulatory mechanisms of those modules.  we have proposed a framework to identify significant coexpressed clusters of genes across multiple conditions experimented from different microarray platforms, time-grids, and also tissues if applicable. analysis on rich in vivo datasets of corticosteroid time-series yielded significant insights into the pharmacogenomic effects of corticosteroids, especially the relevance to metabolic side-effects. this has been illustrated through enriched metabolic functions in those transcriptional modules and the presence of gre binding motifs in those enriched pathways, providing significant modules for further analysis on pharmacogenomic corticosteroid effects. 