Comparative analysis of Hox downstream genes in Drosophila

In animals, morphology along the various body axes is very diverse,requiring both a system that confers positional identity and a means to respond to these positional cues. It has long been known that the system specifying positional identity along the anteroposterior (AP) axis is based on an evolutionarily conserved set of regulators, the Hox genes(Carroll, 1995; Mann and Morata, 2000; McGinnis and Krumlauf, 1992). In Drosophila, Hox genes are expressed in defined domains along the AP axis, and their activity assigns distinct morphologies to the various body segments (McGinnis and Krumlauf,1992). In addition, Hox genes are very often expressed in overlapping domains and crossregulate each other(McGinnis and Krumlauf, 1992; Miller et al., 2001). Consequently, loss of function of one Hox gene frequently leads to the ectopic expression of neighboring Hox genes, which is one of the reasons for the drastic homeotic transformations of body segments initially identified by Ed Lewis (Lewis, 1978). Therefore, only a subset of Hox functions can be identified in loss-of-function mutants.

Hox genes code for transcription factors with a highly conserved DNA-binding domain, the homeodomain(McGinnis et al., 1984; Scott and Weiner, 1984), and it has been postulated that Hox proteins direct morphogenesis by regulating appropriate sets of downstream genes in a segment-specific manner(Graba et al., 1997; Hombria and Lovegrove, 2003). Although a wide range of strategies has been used to identify Hox downstream genes (Graba et al., 1997; Hombria and Lovegrove, 2003; Pradel and White, 1998), our knowledge of their nature is still limited. Initial attempts have focused on in vitro studies or on heterologous systems; however, Hox proteins acquire DNA-binding specificity mostly through interactions with various co-factors in vivo (Ebner et al., 2005; Gebelein et al., 2004; Mahaffey, 2005; Mann, 1995; Mann and Affolter, 1998). Therefore, most known Hox downstream genes have been identified by candidate gene approaches based on expression patterns or similar mutant phenotypes(Pearson et al., 2005),highlighting the power of in vivo strategies to identify Hox target genes. This notion is further supported by recent successful approaches combining loss- or gain-of-function alleles of Hox genes and microarray experiments to identify Hox downstream genes on a larger scale(Cobb and Duboule, 2005; Hedlund et al., 2004; Lei et al., 2005; Williams et al., 2005). Still,most previous efforts were biased toward the identification of direct Hox target genes, and, while knowledge of direct Hox targets is a prerequisite to understanding how Hox proteins acquire DNA-binding specificity in vivo, we need to know the entire Hox-dependent regulatory network with all its tiers of regulatory interactions to understand how Hox proteins control morphogenesis on a cellular level.

Most of the known Hox downstream genes code either for transcription factors or for signaling molecules(Hombria and Lovegrove, 2003; Pearson et al., 2005). These two classes represent the top tiers of regulatory cascades and are able to coordinate many downstream events. Hence, they are not informative for elucidating the role of Hox proteins in the specification of morphological properties on a cellular level per se. To this end, the functional analysis of the so-called realizators, which directly influence the morphology by regulating cytodifferentiation processes(Garcia-Bellido, 1975; Pradel and White, 1998), is required. Unfortunately, even though the concept of realizators was postulated more than 30 years ago, so far very few Hox realizator genes have been identified and studied mechanistically(Bello et al., 2003; Lohmann et al., 2002). One well-studied example of a realizator gene in Drosophila is the apoptosis-inducing gene reaper (rpr), which is expressed in the maxillary segment in Drosophila embryos and is directly controlled by the Hox protein Deformed (Dfd)(Lohmann et al., 2002). In addition, the Dfd-dependent expression of rpr and, consequently, the activation of apoptosis was shown to be necessary and sufficient for the maintenance of the boundary between the maxillary and mandibular segments of the embryonic head (Lohmann et al.,2002). This is one of the few examples demonstrating how a Hox protein can execute, via a single realizator gene, one specific aspect of segmental morphology on the cellular level. To understand and mechanistically link the many remaining Hox functions with morphogenetic outputs, we need to quantitatively identify Hox downstream genes. Functional analysis of this set will then allow us to elucidate all tiers of interactions within the Hox-regulatory network, and to establish links between Hox genes and realizator genes. This seems fundamental for a complete understanding of the role of Hox genes in development and evolution.

The wild-type strain used was Oregon-R. UAS-lacZ,UAS-Dfd, UAS-Antp and arm-GAL4 strains were from W. McGinnis, UAS-Scr, UAS-Ubx and UAS-Abd-B from F. Hirth (Hirth et al., 2001),UAS-Abd-A from A. M. Michelson(Michelson, 1994) and UAS-hep^act. from M. Mlodzik(Weber et al., 2000). For trans-heterozygous mutants the following alleles were used: Dfd^r11 and Dfd^w21 from W. McGinnis; Scr¹, Scr⁴, Abd-B^M2,wg^l-12 and wg^l-17 from the Bloomington Stock Center; and Abd-B^M5 from C. Nüsslein-Volhard(Tübingen Drosophila Stock Collection). Dfd mutant embryos for BrdU staining were Dfd^w21/TM3Sb[twi::GFP] crossed to Dfd^r11/TM3Sb[twi::GFP] and homozygous Dfd mutants (Dfd^w21/Dfd^r11)were identified by absence of GFP signal.

cDNAs were from the Drosophila Genomics Resource Center: CG5080 (LD34147), CG7447 (LD16414), disco(GH27656), Dll (LP01770), ImpL2 (SD07266), gt(RE29225), sage (RE59356), skl (RE14076), spz(SD07354), LysE (LP07339), CG8193 (GH07976), CG3097(RE43153), Mp20 (RE55741), CG17052 (LD43683), Ance(LD11258), Hsp23 (LD06759), sn (RH62992), mas(LP06006), pav (RE22456), wrapper (GH03113), wg(RE02607) and W (AT13267). prd cDNA was from W. McGinnis, Eip63E cDNA and predicted Dfd response elements tested by EMSA were PCR amplified, cloned and sequenced. Expression plasmids for Dfd and Ubx were obtained from W. McGinnis and S. Carroll, respectively.

In situ hybridization and immunochemistry were performed as described(Bergson and McGinnis, 1990; Tautz and Pfeifle, 1989), and BrdU labeling and scanning electron microscopy were done as described(Dolbeare and Selden, 1994; Lohmann et al., 2002). Hox protein expression was measured by the fluorescent intensity of a standardized area of individual nuclei using the Zeiss LSM 510 META confocal microscope. Twenty nuclei of four independent embryos were analyzed for each expression domain and genotype. Antibodies were: anti-Dfd, W. McGinnis; anti-Scr,anti-Antp, anti-Abd-B and anti-wg, Developmental Studies Hybridoma Bank (Iowa,University), anti-Ubx, R. White (Cambridge); anti-Abd-A, I. W. Duncan(Washington, University); anti-GFP, Torrey Pines Biolabs (Houston); anti-BrdU,Roche; anti-mouse AlexaFluor 488, anti-guinea pig AlexaFluor 488 and anti-rabbit AlexaFluor 488, Molecular Probes.

Microarray hybridizations were carried out as described(Schmid et al., 2003) in biological triplicates with RNA from pools of stage 11 or stage 12 embryos. Raw data were quantile normalized and expression estimates were calculated using gcRMA (Wu et al., 2004)implemented in R. Statistical testing for differential expression was carried out using LogitT (Lemon et al.,2003). Microarray data discussed here have been deposited with ArrayExpress database at the EBI(http://www.ebi.ac.uk/arrayexpress-old/;Accession number E-MEXP-879). For analysis of gene ontology categories, GO lists from FlyBase were used. Genes were sorted using a combination of molecular and biological GO terms. Genes containing the following description in their GO annotations were classified as realizators: apoptosis, cell death,cell adhesion, cell shape, cell cycle, mitosis, cell proliferation,cytoskeleton, proteolysis, peptidolysis, cytoskeleton, structural constituent of larval cuticle or peritrophic membrane.

Quantitative real-time PCR (qRT-PCR) was carried out in triplicates from RNA of pooled tissue as described (Schmid et al., 2003) using SYBR-green QPCR Master Mix (Invitrogen). Expression of β-Tubulin was used for crossexperiment normalization. Primer and probe sequences are available upon request.

For cluster identification the Cis-Analyst algorithm(Berman et al., 2004) was used with a Position-Frequency-Matrix (PFM) based on DNaseI footprint data(Bergman et al., 2005) and consensus sequences from the literature. The PFM shown in Fig. S1B in the supplementary material was generated by PATSER. To define enhancer parameters,such as length of enhancer, number of Dfd-binding sites per enhancer, distance between binding sites, known Dfd-dependent enhancers were analyzed. The parameters identified, as shown in Fig. S1B in the supplementary material,were used to predict clusters of Dfd-binding sites in the regulatory regions of selected genes in Drosophila melanogaster. To this end, intergenic and intronic sequences of D. melanogaster were aligned to a multiple sequence file, sorted and separated to segment-files, which included annotation information. The PATSER program used these segment-files as a template to generate P-values for each Dfd-binding site identified according to the PFM. Using this binding site information, clusters of Dfd-binding sites were predicted using the standalone version of cis-Analyst-helper. To validate this approach statistically and to optimize the parameters chosen, Dfd downstream genes identified in the microarray experiment were used. The logic of this approach is based on the assumption that direct Dfd target genes should be enriched among the Dfd downstream genes identified in the microarray experiment when compared with randomly selected genes. To identify Dfd clusters in other Drosophila species (D. simulans, D. yakuba, D. erecta, D. pseudoobscura), a NCBI BLAST search was performed. To consider clusters as being conserved, the following conservation criteria had to be fulfilled: (1) conservation of the enhancer in at least two other Drosophila species; (2) the length of the homologues enhancer had to be ≥50% of the enhancer length identified in D. melanogaster; (3) conservation of at least two Dfd-binding sites within the conserved enhancer elements; and (4) conserved enhancers with less then 50% of length conservation but more then two Dfd-binding sites conserved were treated as minor hits. Conserved enhancers were ranked according to the following parameters: (1) evolutionary distance of Drosophilaspecies; (2) overall sequence similarity of conserved enhancers; (3) numbers of binding sites present in conserved clusters; and (4) degree of variation in enhancer length. Binding site matches for other transcription factors located in the Hox response elements were identified by using rVISTA, Transfac and Jasper databases.

Electrophoretic mobility shift assays (EMSAs) were performed as described(Lohmann et al., 2002).

In order to systematically elucidate Hox-regulatory networks, we performed a comparative microarray screen using stage 11 and stage 12 embryos ubiquitously overexpressing six out of eight Hox genes - Dfd, Sex combs reduced (Scr), Antennapedia (Antp), Ultrabithorax (Ubx), abdominal A (abd-A)and Abdominal B (Abd-B) - by means of the UAS-GAL4 system(Brand and Perrimon, 1993)(Fig. 1A; for endogenous expression of Hox proteins see Fig. 1E). The Hox genes labial (lab) and proboscipedia (pb) were not included, as a similar study has been performed for lab (Leemans et al., 2001), and as pb mutant embryos do not exhibit any obvious defects (Wakimoto et al.,1984). We opted to use overexpression instead of a loss-of-function experiments for several reasons. First, the analysis of Hox mutants is complicated by the extensive crossregulatory interactions of Hox genes (Gould et al., 1997; Miller et al., 2001); loss of expression of one Hox gene often results in ectopic expression of other Hox genes, thereby obscuring the effects on downstream genes. Another important limitation of a loss-of-function approach in conjunction with microarray analysis is sensitivity; due to the small expression domains of many Hox genes(McGinnis and Krumlauf, 1992; Pearson et al., 2005), locally restricted differences in gene expression caused by Hox mutations will be diluted in RNA isolated from whole embryos, and therefore many downstream genes might not be detected. Isolation of cells expressing individual Hox genes by cell sorting (Wang et al.,2006) could provide a means to circumvent this problem; however,the required reporter genes that are expressed in specific Hox domains,although independent from Hox gene activity, currently do not exist. Conversely, although ectopic expression of one Hox gene affects the expression of other Hox genes, their ubiquitous overexpression should allow us to even detect genes whose expression is only weakly, or locally, affected, because the Hox expression domains are expanded manyfold and consequently their transcriptional output is amplified. To achieve ubiquitous Hox overexpression in the desired stages of development, we used an armadillo(arm)-GAL4 driver line (Sanson et al., 1996), which confers ubiquitous expression starting at stage 10, as judged from analyzing GFP activity in embryos carrying an additional UAS-2xEGFP transgene (Fig. 1A). Previous studies have shown that ubiquitous overexpression of Hox genes in UAS fly strains, which were also used in our study, is sufficient to induce ectopic differentiation of Hox-dependent structures without affecting the development of early embryonic stages in an unspecific manner(Li et al., 1999). Thus, a substantial part of Hox downstream genes seem to be responsive to Hox signaling even at ectopic locations and should be detectable by microarray analysis.

One concern with overexpression experiments, however, is that they might result in varying or even unphysiological concentrations of transcription factor proteins in the nucleus, which might unspecifically affect gene expression. Therefore, we confirmed similar RNA and protein expression levels in our overexpression lines by qRT-PCR and antibody stainings(Fig. 1B,E). In addition, we quantified Dfd and Abd-B protein levels in the nuclei of wild-type and the respective overexpression embryos and found that protein levels in transgenic embryos were on the lower end of the wild-type expression spectrum(Fig. 1C,D). Thus,non-physiological effects of Hox misexpression due to increased protein levels should not interfere with our experiment.

For the microarray experiment, total RNA from arm::lacZ(control), arm::Dfd, arm::Scr, arm::Antp,arm::Ubx, arm::abd-A and arm::Abd-Bembryos collected separately at the different stages was prepared and hybridized in biological triplicates to Affymetrix Drosophila Genome 1 arrays, which contains probe sets interrogating more than 13,500 genes. A combination of per-gene (Lemon et al.,2003) and common variance (>1.5-fold change) filtering was used to identify 1508 Hox responsive genes (Table S1, Table S2, Table S3, Table S4,Table S5 and Table S6 in the supplementary material). This list contained six of the 18 genes previously shown to be under Hox control during stages 11 and 12 (Fig. 2A), and in situ hybridizations for three of these transcripts confirmed the microarray results(Fig. 2A). The fact that we were only able to recover one-third of the known Hox targets can be explained by a number of differences between our experimental setup and the ones used before, such as detection method(Capovilla et al., 2001; Mahaffey et al., 1993), timing and level of overexpression (Feinstein et al., 1995), and use of mutants instead of overexpression(Mahaffey et al., 1993; Ryoo and Mann, 1999). Based on this observation, it is conceivable that the actual number of all Hox downstream genes is two- to threefold the number we have discovered in our study, which is still significantly lower than previously suggested(Liang and Biggin, 1998). The microarray data also showed that anterior Hox genes were repressed by those normally expressed more posteriorly (data not shown), a crossregulatory interaction known as posterior suppression(Miller et al., 2001). Again,we could confirm the microarray data by performing antibody stainings for all Hox proteins on embryos ubiquitously misexpressing either Dfd or Abd-B(Fig. 2B).

To verify differential expression of the newly identified genes at the cellular level, we carried out in situ hybridization on embryos misexpressing the various Hox genes (Fig. 3;see Fig. S1 in the supplementary material). Twenty-four of the 25 randomly selected genes that showed a specific in situ signal behaved as observed in the microarray experiment. In addition, for a selected subset of seven genes Hox-dependent regulation could also be shown in Hox mutants(Fig. 3), demonstrating the power of the initial microarray experiment. For example, three transcripts found to be induced by Dfd in the microarray experiment were sickle(skl), a known apoptosis activator(Wing et al., 2002), CG5080, a gene putatively involved in cytoskeletal regulation(Jasper et al., 2001), and CG7447, a gene of unknown function. In situ analysis confirmed strong and ectopic induction of all three genes in response to Dfd misexpression(Fig. 3B,F,J), and showed that their expression in the maxillary segment was lost in Dfd mutants(Fig. 3D,H,L). Similarly, mRNA levels of salivary gland-expressed bHLH (sage), a transcription factor gene exclusively expressed in the salivary gland primordium (Chandrasekaran and Beckendorf,2003), were increased in response to ectopic Scr activity(Fig. 3N). By contrast, sage expression was abolished in Scr mutants(Fig. 3P), consistent with Scr being a master regulator of salivary gland morphogenesis(Panzer et al., 1992). Among the genes that were induced by Abd-B were Ecdysone-inducible gene L2(ImpL2), putatively involved in cell adhesion(Garbe et al., 1993), and spätzle (spz), which encodes a Toll receptor ligand involved in embryonic axis specification(DeLotto et al., 2001). Again,we observed strong ectopic expression of ImpL2 and spz in arm::Abd-B embryos (Fig. 3R,V), whereas expression in the posterior end was lost in Abd-B mutants (Fig. 3T,X). In summary, in situ hybridization with probes for 24 randomly selected genes (Fig. 3; see Fig. S1 in the supplementary material) not only confirmed the microarray results, but also demonstrated that many of the identified Hox downstream genes responded in a converse manner in the respective Hox mutants.

To understand the logic of Hox-dependent morphogenesis, it is important to place the newly identified downstream genes within the underlying regulatory hierarchy. To this end, we developed a bioinformatics tool to detect direct Hox target genes, based on the identification of evolutionarily conserved clusters of Hox consensus binding sites in the genome (see Materials and methods and Fig. S2 in the supplementary material for detail). Using this approach, we were able to identify a large number of putative direct targets of Dfd. From the 240 genes found to be significantly regulated by Dfd, 75 had clusters of Dfd-binding sites (31% of all identified Dfd responsive genes),which was significantly more than expected by chance (P<0.001). In addition, 46 of these clusters were well conserved in at least two other Drosophila species (19% of all identified Dfd responsive genes) (see Tables S7 and S8 in the supplementary material). Most of the predicted Dfd response elements also contained binding sites for other transcription factors(data not shown), a known prerequisite for functional enhancer elements(Berman et al., 2004). We randomly selected six of the 75 predicted Dfd response elements and performed EMSA to test whether Dfd protein could bind to these elements. All enhancer elements tested were bound by Dfd in vitro(Fig. 4), whereas Ubx, a Hox protein specifying trunk identity, was not able to interact with these enhancers (see Fig. S2C in the supplementary material). In addition,competition experiments showed that Dfd specifically bound some, but not all,of the predicted Dfd-binding sites in these enhancers(Fig. 4A-D), demonstrating that the simple presence of a consensus binding site is not sufficient for Dfd binding in the context of these enhancers in vitro and/or that some of the predicted sites are not functional in vivo. Based on our results with Dfd, it seems likely that about 20 to 30% of the identified downstream genes are direct Hox targets. In sum, the combination of microarray analysis with bioinformatics approaches will allow us in the future to not only identify direct Hox target genes, but also to construct complete Hox-regulatory networks.

To assess the specificity of Hox gene regulation, the 1508 responsive genes were classified according to the number of Hox proteins that influenced their expression and the influence of the developmental stage. Remarkably, most downstream genes (1039, 68.9%) were affected by only a single Hox protein,with Abd-A having a very high proportion of unique response genes (two-thirds of its downstream genes were unique), whereas the fraction of unique response genes was smaller (18 to 36%) for the other Hox genes(Fig. 5A). The use of various statistical cut-offs showed that this result is not an artifact of arbitrary thresholding (data not shown). In addition, we were able to confirm the specificity of the Hox response by analyzing the expression of some of the unique downstream genes by in situ hybridizations in embryos misexpressing any of the Hox genes (see Fig. S3 in the supplementary material). About one-third of the identified downstream genes (449, 29.8%) were affected by several Hox proteins, and only 20 genes (1, 3%) responded to all Hox proteins,representing the classes of regional and common downstream genes, respectively(Fig. 5A). Even when we excluded the Abd-A experiment, which was performed slightly differently from the rest of the set and therefore could interfere with this type of analysis,the result did not change: 63% of the genes were uniquely regulated by only one Hox protein, 34.5% of the genes by some and 2.5% by all Hox proteins(Fig. 5A). Remarkably, among the predicted direct Dfd target genes the distribution of unique and regional Hox downstream genes was similar to their distribution among all identified Hox downstream genes (Fig. 5B). Taken together, these results indicate that the specific effects of Hox proteins on morphology are largely mediated by regulatory interactions with uniquely regulated downstream genes, and that despite the very similar DNA-binding sequences for all Hox proteins observed in vitro(Ekker et al., 1994), the overlap of commonly regulated genes in vivo is relatively small.

To evaluate the influence of co-factors on the regulation of Hox downstream genes, we analyzed the Hox response at the two developmental stages. The idea was that overexpression of Hox genes remained the same, while the regulatory environment changed during the progression of embryogenesis. We found that most of the downstream genes were Hox responsive at either embryonic stage 11 or stage 12 (Fig. 5C,D),consistent with the hypothesis that Hox proteins strongly change their transcriptional output through the interaction with differentially expressed co-factors (Gebelein et al.,2004; Mahaffey,2005; Merabet et al.,2005). To confirm stagespecific regulation of Hox downstream genes on a cellular level, we performed in situ hybridizations for some of the differentially expressed genes (Fig. 5D), and indeed found that most of these genes were Hox responsive primarily at one of the two stages (Fig. 5D).

To group the newly identified Hox downstream genes functionally, we used Gene Ontology (GO) annotations (Fig. 6A) and analyzed the distribution of GO categories within the Hox-responsive genes. We found that the two major groups of Hox downstream genes encoded metabolic and realizator functions, followed by the transcription or translation, signaling, transport, stress or defense response and DNA repair or replication classes (Fig. 6A). Supporting the significance of the realizator genes for the Hox response, we found that this group was the only class statistically over-represented among the Hox downstream genes by Fisher's exact test after Bonferoni correction (P<0.001). By contrast, the transcription or translation and transport classes were under-represented(P<0.001). Focusing on Dfd downstream genes, we also found that realizators were over-represented, whereas the transcription or translation functions were under-represented, albeit at weaker P-values. Intriguingly, even within the group of direct Dfd targets, the transcription or translation class was not over-represented. As sequence-specific transcription factors were never over-represented among any category tested,this suggests that the activation of regulatory proteins is not the preferred mode of Hox action. The fact that realizator processes at the same time are under direct as well as indirect Hox control demonstrates that Hox-dependent morphogenesis is indeed achieved by regulatory networks, rather than linear pathways. Furthermore, these results also argue against the idea that Hox genes mainly act on other transcription factors, a hypothesis that was based on the limited information on known direct Hox targets.

Since the realizator gene concept was postulated almost 30 years ago(Garcia-Bellido, 1975), but only a few such genes have previously been identified as Hox downstream genes in Drosophila, we decided to study this class of Hox response genes in more detail. Among the realizators, the largest subgroup comprises genes involved in proteolytic processes, followed by genes with cytoskeleton functions, a diverse group containing cuticle, chorion and peritrophic membrane genes, cell cycle or cell proliferation genes, apoptosis or cell death genes and cell adhesion genes (Fig. 6B). Genes within realizator subclasses are often coordinately regulated: most apoptotic (7/8) and cell cycle or proliferation genes (18/21)were activated, whereas almost all cell adhesion genes (12/14) and the majority of genes involved in proteolytic processes (56/75) were repressed by Hox proteins (Fig. 6C). Re-analyzing data from a more restricted microarray study, a similar trend can be identified for lab, another Hox gene(Leemans et al., 2001): one apoptotic gene and six cell cycle or cell proliferation genes were activated,whereas three cell adhesion genes were all repressed by lab. This suggests that a variety of cellular processes need to be regulated in a coordinated fashion in every segment in order to realize common aspects of segmental morphology. Support for this notion also comes from a previous analysis in Drosophila, showing that two Hox proteins, Dfd and Abd-B,locally activate the apoptosis gene rpr and thus the apoptotic machinery at segment boundaries for their maintenance(Lohmann et al., 2002).

Finally, we also wanted to analyze which role the identified Hox downstream genes play during the morphogenesis of segments. To this end, we first asked whether differences in morphology along the AP axis are reflected in the differential regulation of Hox responsive genes. One line of evidence supporting this idea is that many downstream genes responded only to a single Hox protein (Fig. 5A). In addition, we found that coordinated regulation of shared downstream genes was more frequent among Hox proteins specifying segments with similar morphologies than among those that specify diverse segmental morphologies(Fig. 6D). For example, both Dfd and Scr specify small lobe-shaped gnathal segments and only 4%(n=2) of their shared downstream genes are regulated in an opposite manner (Fig. 6D). One of the differentially regulated Dfd or Scr downstream genes is PAPS synthetase, which is repressed by Dfd and activated by Scr in our dataset(see Tables S1 and S2 in the supplementary material). Consistent with this observation, it had been shown that PAPS synthetase is specifically expressed in the salivary gland placodes in a Scr-dependent manner, while it is absent from Dfd-expressing cells in stage 11 and 12 embryos(Jullien et al., 1997). Another example for a strong correlation of differences in morphology and the differential regulation of shared downstream genes is the Scr-Abd-A pair. In this case, opposite regulation was found for 58% (n=55) of the targets shared by Scr and Abd-A, in line with the very different morphologies specified by these two Hox proteins (Fig. 6D). Again, two selected examples, pipe and PH4αSG2, both activated by Scr and repressed by Abd-A(see Table S2 and Table S5 in the supplementary material), are known to be expressed exclusively in the Scr-specified salivary glands at stages 11 and 12(Abrams and Andrew, 2002; Zhu et al., 2005). Interestingly, it has been shown only recently that pipe,differentially regulated by Dfd and Abd-A, and PAPS synthetase,differentially regulated by Scr and Dfd in our microarray analysis, are both necessary for the production of sulfated macromolecules in the salivary glands of Drosophila embryos (Zhu et al., 2005). Thus, it seems that the diversification of segments is achieved, on the one hand, through the regulation of unique downstream genes,and, on the other hand, through the differential regulation of shared downstream genes.

To analyze the morphogenetic function of Hox responsive genes in more detail, we focused on the potential role of several newly identified Dfd downstream genes during the development of the maxillary segment. It has long been known that Dfd is expressed in the maxillary and mandibular segments, and is necessary for the morphological specializations (mouth hooks, cirri,ventral organ) of these head segments(McGinnis et al., 1990). However, only a single cellular event necessary for the morphogenesis of the maxillary segment and under the control of Dfd has been explained mechanistically so far: the maintenance of the boundary between the maxillary and mandibular segments, which is dependent on Dfd-mediated activation of rpr expression in the anterior part of the maxillary segment(Lohmann et al., 2002). Consistently, rpr was found among the activated Dfd downstream genes in our microarray analysis (see Table S1 in the supplementary material). Another prominent feature of Dfd mutants is the displacement of maxillary and mandibular segments to a more dorsal position, caused by the accumulation of supernumerary cells at the ventral side of both segments,which had been observed almost 20 years ago(Fig. 7B)(Regulski et al., 1987). There are two alternative explanations for this defect: loss of cell death and/or overactivation of cell proliferation. Consistent with the former explanation,we observed reduced local expression of the apoptosis activator skl(Fig. 7K,P), one of the newly identified genes activated by Dfd, and a concurrent reduction in the number of apoptotic cells at the ventral side of the maxillary segment in Dfdmutants (data not shown) (Lohmann et al.,2002). Additionally, we were able to show that cell proliferation at the ventral region of the maxillary segment is increased in Dfdmutant embryos by performing BrdU labeling experiments(Fig. 7E,J). This might be attributed to the de-repression of two genes identified as repressed by Dfd in this study (see Table S1 in the supplementary material): Ecdysone-induced protein 63E (Eip63E), encoding a cyclin-dependent protein kinase(Stowers et al., 2000), and wingless (wg), encoding a signaling molecule known to play a role in cell proliferation (Giraldez and Cohen, 2003). In line with a role of these genes in shaping the maxillary segment, we found that both genes are ectopically expressed in the dorsal part of the maxillary segment in Dfd mutants(Fig. 7Q,R). Although the function of Eip63E during the morphogenesis of the maxillary segment could not be analyzed due to the lack of mutants, we could confirm that wg mutants have reduced gnathal lobes(Fig. 7D)(Rusch and Kaufman, 2000),suggesting an important role of wg in the regulation of cell proliferation in the maxillary segment. The third notable defect of Dfd mutants is the loss of the maxillary cirri primordium(Regulski et al., 1987). paired (prd), one of the transcription factor genes identified in our screen (see Table S1 in the supplementary material), is known to be important for development of cirri and the maxillary ventral organ(Vanario-Alonso et al., 1995). Because late prd expression is completely under the control of Dfd(Fig. 7N,S), we conclude that some aspects of ventral maxillary identity are specified by Dfd via prd regulation. Finally, we analyzed Dfd-dependent regulation of cell shape changes, because cells at ventral positions of wild-type maxillary segments are round (Fig. 7G),whereas in Dfd mutants many appeared elongated(Fig. 7H). The JNK pathway has been implicated in cell shape changes in Drosophila, for example during embryonic dorsal closure and adult thorax closure(Harden, 2002; Xia and Karin, 2004) and because we had identified several genes responsive to the JNK pathway(Jasper et al., 2001)[Ras-related protein (Rala), Angiotensin converting enzyme (Ance) and CG5080] (see Table S1 in the supplementary material) as Dfd downstream genes, we tested the contribution of the JNK pathway to the cell shape phenotype of Dfd mutants. After ubiquitous activation of the JNK pathway by overexpressing a constitutively active form of Hemipterous (Weber et al.,2000) using the arm-GAL4 driver, we observed elongated cells in the maxillary segment (Fig. 7I), as well as in other parts of the embryo (data not shown). As we could confirm for one of the JNK-responsive Dfd downstream genes, CG5080, implicated in the regulation of cytoskeletal dynamics(Jasper et al., 2001), strong upregulation by Dfd (Fig. 7T),we conclude that the JNK pathway plays a major role in organizing cell shapes in the maxillary segment.

More than 30 years ago Antonio Garcia-Bellido proposed that a hierarchy of three classes of genes, activators, selectors and realizators, accounts for cell differentiation during development, thereby providing a functional scheme for the control of morphogenetic processes. The key proposal was that, once activated in their appropriate territories by so-called activator genes,selector genes (he applied this name to homeotic genes) would not directly specify morphological differences between different segments, but would rather select a battery of subordinate downstream genes, the realizator genes,encoding cellular proteins directly required in cell differentiation processes(Garcia-Bellido, 1975). Until now much effort has gone into elucidating the nature and function of all three hierarchical levels, with a substantial amount of knowledge having accumulated at the activator and selector level. It is now well established that a genetic cascade comprising maternal and various classes of segmentation genes regulate the temporal and spatial expression of a unique combination of Hox genes in different segments, which subsequently specifies the identities of individual segments (McGinnis and Krumlauf,1992; St Johnston and Nusslein-Volhard, 1992). Additionally, the discovery that Hox proteins act as transcriptional regulators established the general view that each segment will enter a specific morphogenetic program and develop unique shape and function depending on Hox downstream genes, in particular the realizator genes. Although the question of Hox downstream gene identity and function is not a novel problem, and although a considerable amount of progress has been made in recent years, our knowledge of their nature and function is still far from complete, especially with regard to the realizator genes in the sense of Garcia-Bellido.

Many of the known Hox downstream genes coded either for transcription factors or signaling molecules, and only very few of them were realizators(Hombria and Lovegrove, 2003; Pearson et al., 2005; Pradel and White, 1998). This was puzzling, as the primary function of Hox proteins is to specify the morphology of different segments, thus one would have expected to find a bias toward realizators. Moreover, this finding established the view that most of the cellular responses mediated by Hox proteins, including realizator functions, are indirectly influenced through the action of intermediate regulatory molecules. Our analysis of Hox downstream genes in Drosophila, which was designed to allow for a quantitative identification of Hox-regulatory networks (including most realizator genes),revealed that a major group of genes responsive to Hox input did indeed code for realizators. Therefore, our results constitute the first experimental support of the concept postulated by Garcia-Bellido more than three decades ago. We could furthermore show that a substantial part of the Hox output is directly transferred to the realizator level, suggesting that intermediate regulators might play a smaller role than previously thought. One possible explanation why so few Hox realizators had been identified before is that most realizators will be required for general functions in many cells. Consequently, mutations in realizator genes are likely to result either in early embryonic lethality or in pleiotropic effects, making it difficult to correlate their phenotypes to those found in Hox mutants. In addition, it seems likely that realizators act redundantly or have very subtle effects,making their identification in forward genetic screens extremely difficult. Similarly, individual mutations in all known guidance factors for border cell migration in Drosophila produce either no, or only mild, defects and thus they could be identified only by expression profiling studies(Wang et al., 2006). In this study we have quantitatively identified Hox realizator genes by a comparative microarray analysis, which now can serve as a resource to study the mechanisms of segmental morphogenesis. Focusing on the differentiation of the maxillary segment, we were able to functionally correlate all major morphological defects observed in Dfd mutants with newly identified Dfd downstream genes, many of which code for realizators, demonstrating the validity of this approach.

Another important discovery of our analysis is the enormous specificity of Hox protein action on the transcriptome in vivo, which sharply contrasts with the low DNA-binding specificity in vitro. Hox proteins have been shown to bind to very similar, relatively simple, DNA sequences containing a TAAT core sequence in vitro (Biggin and McGinnis,1997; Carr and Biggin,1999; Ekker et al.,1994; Walter and Biggin,1996), whereas many of the identified Hox downstream genes are uniquely regulated by only a single Hox protein. This contrast may be explained by our observation that the majority of genes are primarily regulated at only one of the two stages, implicating that Hox proteins excessively interact with the regulatory environment in which they are embedded. Support for the notion that co-factor interactions have a major impact on Hox output also comes from a recent study, which has provided direct evidence that Hox proteins gain the ability to regulate their target genes in a context-specific manner by interaction with known cell- and/or tissue-specific transcription factors in vivo(Gebelein et al., 2004). In addition, this study also suggests that a large number of transcription factors might function as Hox co-factors, which could dictate the outcome of Hox gene action. Along these lines, we found that ubiquitous overexpression of Hox proteins never caused ubiquitous activation of downstream genes, but that ectopic expression was always locally restricted, suggesting that regional transcription factors are essential for Hox output. This is also reflected in our finding that conserved clusters of Hox binding sites in the regulatory regions of direct targets frequently contain binding sites for unrelated transcription factors. Taken together, these results support the hypothesis put forward by Michael Akam in 1998, that “we should think of the Hox genes with their short and relatively non-specific target sequences, as cofactors that modify the actions of other more specific transcription factors, rather than proteins in need of cofactors themselvess”(Akam, 1998).

We thank Jan U. Lohmann for help with the microarray analysis and quantitative real-time PCR, Benjamin Berman for receiving the stand-alone binaries and support for Cis-analyst, Tancred Frickey for help with gene annotations, Jürgen Berger for advice with the SEMs, Stephan Ossowski for help with databases, Petra Stöbe for assistance with gelshifts and Jan U. Lohmann and Detlef Weigel for discussion and for critically reading the manuscript. Supported by the Deutsche Forschungsgemeinschaft (I.L.; LO 844/3-1) and by the Max Planck Society.