当前位置: 首页 > 医学版 > 期刊论文 > 基础医学 > 分子生物学进展 > 2004年 > 第2期 > 正文
编号:11259336
Evolution of Developmental Genes: Molecular Microevolution of Enhancer Sequences at the Ubx Locus in Drosophila and Its Impact on Developmen
     Institute of Genetics, University of Nottingham, Queens Medical Centre, Nottingham, United Kingdom

    E-mail: john.brookfield@nottingham.ac.uk.

    Abstract

    Homeotic genes, which function to specify segment identity along the anterior-posterior axis of embryos, are controlled by extensive batteries of enhancer sequences. We have investigated patterns of interspecific and intraspecific molecular variation in three enhancers of the Ultrabithorax (Ubx) locus, which are bx-32.8, pbx32.7, and bxd4.1, from the Drosophila melanogaster species group. These enhancer sequences control Ubx expression by binding to multiple transcription factors encoded by gap, pair-rule, and dorsoventrally expressed genes. Sequence comparisons reveal purifying selection acting on all three enhancers, both in bases binding transcription factors and in bases whose functions are as yet unknown. Neutrality tests largely fail to reject a neutral evolution model. However, using a matrix similarity value to reflect the binding affinity of the protein-binding sites, interspecific and intraspecific variation that may have potential to affect the binding affinity of the sequences homologous to those binding transcription factors in D. melanogaster are discovered, suggesting evolutionary flexibility in the way in which these sequences function in the control of development. As a means of measuring the impact of intraspecific variation on observable phenotypes, we have induced Ubx mutant phenocopies with embryonic ether treatment, and find strong and highly significant variation between D. melanogaster strains in their phenocopy frequencies. This variation shows no significant correlation with the strengths of the mutant phenotypes when the strains are heterozygous with a Ubx null mutation. Estimated phylogenetic trees have been constructed for the three enhancer regions investigated. Neither of the two phenotypic traits investigated shows any significant associations with the phylogeny of any of the three enhancers.

    Key Words: Drosophila ? Ultrabithorax ? phenocopy ? developmental evolution ? molecular evolution ? enhancer sequences

    Introduction

    Animal development is controlled by the creation of regular and repeatable patterns of gene expression, which varies between tissues and between stages of development. The process of differential gene expression is brought about in large part by the interaction of transcription factors with their target sequences, which make up enhancers, creating cis-acting control of transcription (Arnone and Davidson 1997). The evolution of development therefore must involve changes in transcription factors, in their cis-acting targets, or in both (Tautz 2000). However, a number of studies have shown that the amino acid sequences of many developmentally important transcription factors have been highly conserved among moderately diverged species (Kassis et al. 1986; Treier, Pfeifle, and Tautz 1989; Pan, Valentine, and Courey 1994; Sackerson 1995; Terol, Perez-Alonso, and de Frutos 1995). This is in accord with what one might theoretically expect, because a change in the functional properties of a transcription factor would be expected to have an impact in its interactions with multiple target genes, and the overall effect would be very unlikely to be advantageous. Changes in cis-acting control sequences, however, are likely to create a more specific and restricted effect on gene expression, which would have a higher probability of being advantageous to the organism. The control of development is likely to evolve over the short term mainly by changes in cis-acting enhancer sequences (Warren et al. 1994; Averof and Patel 1997; Weatherbee et al. 1999; Stern 2000).

    Typically, regulatory sequences can evolve unexpectedly rapidly between diverged species so that only sequences from closely related species can be reliably aligned. The comparison of the regulatory sequences between moderately diverged species often shows an interspersion pattern of small conserved sequences, most of which are identified as binding sites for transcription factors (Ludwig, Patel, and Kreitman 1998; Wilde and Akam 1987; Kassis et al. 1989; Pfeifle, and Tautz 1989; Maier, Preiss, and Powell 1990; Langeland and Carroll 1993; Lukowicz et al. 1994; Pan, Valentine, and Courey 1994; Sackerson 1995; Kreitman 1998; Wolff et al. 1999; Kim 2001; Dellino, Tatout, and Pirrotta 2002). It is interesting that, while the sequences of regulatory elements have rapidly evolved, their functions in controlling restricted gene expression have been found to be highly conserved between species. In transgenic experiments, the regulatory sequences from one species can drive the expression of a reporter gene in the same way as do the native regulatory sequences (Langeland and Carroll 1993; Lukowicz et al. 1994; Pan, Valentine, and Courey 1994; Ludwig, Patel, and Kreitman 1998; Wolff et al. 1999). Furthermore, a transformation of the coding and regulatory sequences from one species into another can often rescue the mutant phenotype in the same way as can the native gene (Maier, Preiss, and Powell 1990; Lukowicz et al. 1994). The fact that nucleotide changes are sometimes observed in protein-binding sites that have been shown to be functionally important in one species suggests that regulatory sequences may evolve under stabilizing selection (Ludwig et al. 2000). The functional redundancy of regulatory sequences may allow slightly deleterious mutations to be fixed by genetic drift, and then adaptive changes would be selected to compensate for these slightly deleterious mutations. Evidence for stabilizing selection acting on regulatory sequences has been shown in the case of the Drosophila even-skipped stripe 2 enhancer (Ludwig et al. 2000). The even-skipped stripe 2 enhancer of D. melanogaster and D. pseudoobscura can individually control the correct expression of even-skipped in transgenic D. melanogaster embryos. However, the chimeras constructed by combining halves of the enhancer from the two species no longer drive the correct even-skipped expression. Based on this finding, the authors suggested that stabilizing selection has maintained the overall function of the even-skipped stripe 2 enhancer by selecting for compensatory mutations, and this may be a general evolutionary mode for many other regulatory sequences. Such a pattern of compensatory changes regulated by stabilizing selection has allowed enhancer sequences to become remarkably diverged between moderately diverged species and yet still be able to create the same gene expression pattern.

    That enhancer sequences can evolve unexpectedly quickly between species restricts the phylogenetic domain of inference about developmental control that can be derived from model systems such as Drosophila. However, while this reduces the capacity of Drosophila to act as a model system for phylogenetically ancient patterns of expression, our great understanding of the evolutionary process in this genus suggests that it will be a useful tool to study the microevolutionary processes that operate on this kind of sequence. Investigating the ways in which the functionally important parts of enhancer sequences evolve within the genus will give us some insight into the ways in which an evolving array of transcription factor binding sites can nevertheless maintain a gene expression pattern. Evolutionary change can potentially be brought about in two different ways—either by neutral or nearly-neutral changes fixed by genetic drift or by changes driven by natural selection. In determining the relative role of these two processes in sequence change, some insight can sometimes be gained from knowledge of the function of the sequences, but another very powerful set of investigative tools comes from population genetics. Interspecific variation may not, in itself, reveal whether changes observed are neutral or selected, but these two models have differing predictions for the patterns of variability within populations, and so the simultaneous consideration of within-species and between-species divergence can often detect evidence of natural selection.

    Ultrabithorax (Ubx) is one of the homeotic genes in Drosophila (Lewis 1978). It functions to control segment identity of the parasegments 5 and 6. Temporal and spatial expression of Ubx is regulated by interaction between transcription factors and regulatory elements located in the intron and the 5' flanking region. These regulatory elements are classified into four clusters, which are anterobithorax (abx), bithorax (bx), bithoraxoid (bxd), and postbithorax (pbx). Mutations in these regulatory elements cause transformation of the third thoracic segment toward the corresponding part of the second thoracic segment. The amino acid sequences of Ubx are highly conserved, showing 94% overall amino acid identity between D. melanogaster and D. virilis (Bomze and López 1994), and 78% identity between D. melanogaster and Musca (Wilde and Akam 1987). In addition, the process of Ubx alternative splicing and the time-specific and tissue-specific expression of each UBX isoform have been highly conserved among Drosophila species (Bomze and López 1994), suggesting that the mechanisms of Ubx regulation are also highly conserved. However, it is found that the subtle morphological differences in the trichome patterns on the posterior femur of the second leg between D. melanogaster and its sibling species D. simulans are due to the differences in the levels of Ubx expression in the posterior femur during late pupal development (Stern 1998). These difference have been shown to be driven by the divergence of cis-regulatory sequences at the Ubx locus between the two species, because the rescue of a Ubx null allele by its own wild-type allele cannot be duplicated by the wild-type allele from another sibling species, even a species showing an identical Ubx amino acid sequence. In addition, it is evident that polymorphisms at the Ubx locus partly play a role in controlling the frequency of bithorax phenocopies induced by ether vapor (Gibson and Hogness 1996).

    Here, we study the process of enhancer sequence evolution by investigating the patterns of microevolutionary change in three Ubx enhancer sequences, which are named bx-32.8, bxd4.1, and pbx32.7, both within and between species. The locations of these three enhancer sequences in the Ubx gene map are shown in figure 1. The bx-32.8 enhancer is 730 bp in length and is located approximately 32.8 kb downstream from the transcription start site. This enhancer is equivalent to the BRE fragment, which functions to regulate Ubx expression in epidermal cells in the even parasegments from PS6 through PS12 (Qian, Capovilla, and Pirotta 1991). The bxd4.1 enhancer is 549 bp in length and is located 4.1 kb upstream from the transcription start site. This enhancer contains no HUNCHBACK binding sites, and does not regulate Ubx expression in a parasegmental pattern like the bx32.8 enhancer. Its specific function has not been investigated. However, it is located adjacent to the regulatory sequence controlling Ubx expression in visceral mesoderm (the sequence between 1.7 to 3.1 kb upstream of the Ubx transcription start site) (Bienz et al. 1988). The pbx32.7 enhancer is 1,036 bp in length and is aligned with the U31961 sequence at base positions 209118–210153. It is equivalent to the PBX (or 2206Hf) fragment, and it functions to regulate Ubx expression in epidermis and mesoderm in the even parasegments from PS6 through PS12 (Müller and Bienz 1992; Pirrotta et al. 1995). The initial expression of Ubx is regulated by products of gap and pair-rule genes. The gap gene product HUNCHBACK binds to the bx and the pbx elements and represses expression in the anterior half of the embryo, whereas repression by TAILLESS defines the posterior limit of the Ubx expression pattern (Qian, Capovilla and Pirrotta 1991, 1993; Muller and Bienz 1992). FTZ binding to sites in both bx and pbx stimulates expression in the even-numbered parasegments posteriorly from PS6, whereas the subsequent strengthening of expression in the anterior of each parasegment is thought to result from EN activation. TWIST stimulates expression from bx in the mesoderm (Qian, Capovilla, and Pirrotta 1993).

    FIG. 1. The Ultrabithorax gene map, showing the positions of the enhancer sequences studied. The line indicates the DNA sequence, and numbers represent the DNA walk-in kilobases. The transcription start site is shown as a bent arrow. The gene contains four exons, which are indicated by horizontal blocks (those shadowed in gray are the coding sequences). Beneath the DNA line are some Ubx mutations. Triangles indicate insertions, and brackets indicate deletions. Above the DNA line are DNA fragments investigated for their function in regulating Ubx expression in embryonic tissues (dark bars) and in imaginal discs (gray bars). The open bars indicate Polycomb response elements (PREs). The references are (1) Simon, Peifer and O'Connor (1990); (2) Irvine, Helfand, and Hogness (1991); (3) Müller and Bienz (1992); (4) Qian, Capovilla and Pirrotta (1993); (5) Pirrotta et al. (1995); (6) Chiang et al. (1995). The enhancers studied, namely bx-32.8 (730 bp), bxd4.1 (549 bp), and pbx32.7 (1,036 bp), are indicated by black blocks beneath the DNA line

    In this study, a matrix similarity value, based on Quandt et al. (1995), was used to reflect the binding affinity of the protein-binding sites and thus to assess dynamic changes in the protein-binding sites over the microevolutionary scale. We also investigate the evidence that wild-type strains of D. melanogaster vary genetically in Ubx function as measured by the phenocopy frequency after ether treatment and by the ability of the strains to complement a null allele. The variation in the two phenotypes among the strains was compared to Ubx phylogenetic trees.

    Materials and Methods

    Fly Species and Strains

    The Drosophila species studied were the eight closely related species of the melanogaster species subgroup. The D. melanogaster strains were collected from natural locations and kindly provided as follows. Canton S, Draveil, Kishinev, Loua, and Petit Bourg were collected from USA, France, Moldova, Zaire, and Guadeloupe (West Indies), respectively. These strains were provided by P. Capy (Gif). Harwich was collected from USA and provided by M. Ashburner (Cambridge). NW12B was collected from North Wootton (UK) and provided by P. Corish (Leicester). Samarkand, Tahbilk Cellar, Tahiti, and Tex were collected from Uzbekistan, Australia, South Pacific, and Texas (USA), respectively. These strains were provided by B. Clarke (Nottingham). The Monty strains were collected from Montpellier (France), and the RD3F55 and RD8F55 strains were collected from Nottingham (UK). Each D. melanogaster strain was inbred in the laboratory in a long history of single brother-sister mating.

    The extraction of chromosomes II and III was carried out in some D. melanogaster strains by crossing an individual fly from a D. melanogaster wild-type stock to C23a [In(2R) Cy sp2/In(2LR) Pm(bwV1) ds33k b dp; In(3LR) D cxF ru h/Sb], then crossing an individual F1 female or male carrying the Pm and D mutations to C23a again, and collecting F2 progeny carrying the Pm and D mutations and letting them cross to each other, finally collecting wild-type F3 progeny, which are now homozygous for chromosomes II and III isolated from their original stocks. The chromosome extraction stocks are given new names by adding a number after their original strain's name, as follows: Draveil(2), Draveil(3), Kishinev(7), Loua(4), Monty3(2), Monty5(6), Monty8(1), Monty9(3), Monty9(15), Monty11(4), Monty20(2), Petit Bourg(3), Tahiti (7), Tex(10).

    Other isofemale Drosophila species and strains were obtained from the Species Resource Centre. The seven D. simulans strains are Sim161 (14021-0251.161), Sim162 (14021-0251.162), Sim163 (14021-0251.163), Sim165 (14021-0251.165), Sim166 (14021-0251.166), Sim167 (14021-0251.167), and Sim175 (14021-0251.175). D. sechellia is Fly base stock number 14021-0248.2. The two D. mauritiana strains are line 4 (14021-0241.4) and line 5 (14021-0241.5). Single brother-sister inbreeding from the original D. teissieri (14021-0257.0) stock was carried out for five generations in the laboratory, deriving two D. teissieri strains, 1 and 5. D. yakuba, D. orena, and D. erecta are Fly base stock numbers 14021-0261.0, 14021-0245.0, and 14021-0224.0, respectively.

    DNA Extraction

    Genomic DNA extraction used commercial DNA extraction kits from Qiagen (QIAamp Tissue Kit, QIAGEN DNA Mini Kit, and QIAGEN DNA Blood Mini Kit) and Gentra system, Flowgen (Puregene DNA Isolation Kit). The methods of DNA extraction followed the protocols given with these commercial kits. However, in the last step of the extraction, the DNA was eluted with sterile distilled water preheated to 70°C. The DNA solution was stored at –20°C. All these commercial kits provided suitable DNA for polymerase chain reactions (PCR).

    PCR and DNA Sequencing

    The following primers were used to amplify the bx-32.8, bxd4.1, and pbx32.7 enhancer sequences. Two pairs of primers were used to amplify two overlapped sequences of the pbx32.7 enhancer, giving a total length of 1,036 bp (table 1).

    Table 1 Primers Used to Amplify the bx-32.8, bx4.1, and pbx32.7 Enhancer Sequences.

    The PCR was conducted in a 50 μl reaction volume (5 μl 10 x PCR buffer, 4 μl 25 mM MgCl2, 1 μl 10 mM dNTPs, 2 μl of each 10 pmol Primer, 0.4 μl Red Hot DNA Polymerase/Advanced Biotechnologies, 2 μl 40 ng template DNA, and 33.6 μl sterile distilled water). Reactions were carried out for 1 cycle of 2 min at 94°C, followed by 30 cycles of 1 min at 94°C for denaturation, 30 s at the annealing temperature (depending on which primers were being used), and 1 min at 72°C for extension, then followed by one cycle of 10 min at 72°C and 2 min at 4°C.

    The PCR product was purified using the low melting point (LMP) agarose DNA purification/freeze-and-thaw method. The PCR product was separated on a 1% LMP agarose gel, and then the band was cut out. The gel slice containing PCR product was melted at 65°C in 400 μl TE buffer for 10 min, and snap frozen in liquid Nitrogen for 2 min. The solution was thawed at 37°C for 30 min. The LMP gel was separated from the solution by spinning at 13,000 rpm in a bench centrifuge for 5 min. Supernatant containing DNA was removed and saved in a new tube. DNA was precipitated with 100 μl 0.3 M sodium acetate and 900 μl 100% ethanol at –20°C overnight, and then spun down at 13,000 rpm for 15 min. The pellet was washed with 1 ml 70% ethanol, dried, and dissolved in 30 μl sterile distilled water. The purified PCR products were kept at –20°C until used.

    The purified PCR products were sequenced using the same primers that were used to perform the PCR reactions. DNA sequencing was performed by using ABI Prism dRhodamine Terminator Cycle Sequencing Ready Reaction Kit (PE Applied Biosystems) and run on an ABI 377 automated DNA sequencer. All sequences were investigated for both strands.

    Sequence Analysis

    The bx-32.8 sequence aligned to the complementary strand of the GenBank DNA sequence, accession number U31961, at base positions 275750–276480. (The G at base 275754 in this sequence is absent in all other sequences and is probably an error.) The bxd4.1 and pbx32.7 sequences aligned to the U31961 sequence at base positions 238196–238744 and 209118–210153, respectively. Pairwise alignments using the "align" program (Myers and Miller 1988) were performed to examine the sequencing accuracy, to edit, and to complete the sequences. Multiple alignments of the completed sequences from all species and strains were produced using the ClustalW package (Thompson, Higgins, and Gibson 1994; Higgins and Sharp 1998), and manually edited. Neighbor-Joining trees with 1,000 bootstrap replications were constructed using the method of Saitou and Nei (1987) and plotted with the TreeView program (Page 1996).

    From the multiple alignments, the numbers of fixed changes between the species and the numbers of polymorphisms were counted. The minimum number of mutational steps separating the sequences observed was estimated from the branches in the estimated phylogenetic trees and from the true species phylogeny (Ashburner 1989). Deletion or insertion of any number of contiguous bases is counted as one change. However, if changes, created either by slippages or by nucleotide substitutions, were seen in the sequence included in any of the indels, then the number of such changes was also counted and included in the analysis. To calculate the numbers of fixed changes per nucleotide, the length of the D. melanogaster sequence from GenBank, accession number U31961, was used as the standard length. The DnaSP version 3.51 program was used to calculate intraspecific variation ( and ) and to perform neutrality tests (Tajima's D, and Fu and Li's D, and the HKA test) (Rozas and Rozas 1999). In this program, gaps in the alignment were excluded from the analyses.

    Among the three enhancers studied, only the bx-32.8 and pbx32.7 enhancers have been previously investigated for their protein-binding sites in D. melanogaster using a DNA footprinting technique (Pirrotta et al. 1995; Qian, Capovilla, and Pirrotta 1993). Sequence analyses comparing protein-binding and non-protein-binding sites were therefore carried out on these two enhancers. The numbers of fixed changes and polymorphisms were used for a neutrality test equivalent to the McDonald-Kreitman test (McDonald and Kreitman 1991) when applied to these non-coding sequences (Jenkins, Ortori, and Brookfield 1995).

    Estimation of Matrix Similarity of Protein Binding Site

    The calculation of a matrix similarity value is based on Quandt et al. (1995). This is involved in the construction of a binding consensus matrix for a transcription factor and the calculation of a consensus index value (Ci) for each base position within the consensus matrix. The nucleotide distribution matrixes for the binding consensuses of the FUSHI TARAZU (Ftz), HUNCHBACK (Hb), ENGRAILED (En), TAILLESS (Tll), and TWIST (Twi) proteins were created by the multiple alignments of 52 Ftz, 67 Hb, 13 En, 27 Tll, and 13 Twi binding sites, respectively. The references for these protein-binding sites are available on the Brookfield Laboratory Web site (http://evol.nolt.ac.uk/~jb_lab/index.htm); and the references for the Hb binding sites were kindly provided by C. Bergman. Note that the binding sites in the bx-32.8 and pbx32.7 enhancers are not included in the alignments. The consensus index value (Ci) represents conservation of each base position within the consensus matrix. It is calculated as follows:

    where P(i, b) is the frequency of nucleotide b at position i.

    The maximum value of Ci is 100. This happens when the base position in the consensus shows total conservation of base identity. The minimum value of Ci is 0, which happens when the position shows an equal distribution of all four nucleotides and gaps. A matrix similarity value was calculated as follows:

    where Ci(j) is the consensus index value of position j, n is the length of the consensus matrix, score (b, j) is the matrix value for base b at position j, and max_score(j) is max{score(b, j)}; b A, C, G, T, gap.

    The maximum value of a matrix similarity is 1. This happens in the case that the candidate sequence is identical to the most conserved nucleotide at every position of the consensus matrix. Because the matrix value (score) is multiplied by the Ci value, mismatches at highly conserved positions will affect the matrix similarity value more than do mismatches at less conserved positions.

    The core sequence within each protein-binding site of D. melanogaster was identified by its having the highest matrix similarity value. To achieve this, the matrix similarity value was calculated for every candidate sequence by successively shifting one base along the protein-binding sequence, as well as its complementary strand. This was performed by our program written in the Perl language. To evaluate the potential effects of evolutionary changes on the function of these protein-binding sites in other Drosophila species, the matrix similarity values were calculated for the sequences homologous to the D. melanogaster core sequences. Changes in the matrix similarity values compared to those of D. melanogaster suggest what might be the functional effect of nucleotide changes on the binding affinity of the sites.

    Measurement of Interstrain Variation in Bithorax Phenocopy Frequency

    The bithorax phenocopies are first observed in embryos treated with ether at syncytial blastoderm (1:00–2:00 h), are at the highest frequency in embryos treated at cellular blastoderm to gastrulation (2:30–3:30 h), and are not observed with treatments after germ band extension (5:30 h) (Bownes and Seiler 1977). In this study, a 10-min ether treatment was conducted on 3.0 ± 0.5-h-old embryos from different strains of D. melanogaster. Eggs were collected over a 1-h period, and placed on a 1.25% agar slice. After 2.5 h, the agar slice was transferred into a glass vial containing ether-saturated cotton wool. The eggs were exposed to ether vapor for 10 min and left in the fume cupboard for 30 min. After that, they were transferred to a new vial of cornmeal medium containing 20–40 eggs per vial. Two weeks later, the emerging adults were scored. The bithorax phenocopies show a transformation of metathoracic to mesothoracic segments, imitating the Ultrabithorax (Ubx) mutant phenotype (fig. 2). The phenocopies were transferred in 70% ethanol, dehydrated with 100% ethanol, and kept in 90% glycerol. The detailed characteristics of phenocopies were observed under a microscope.

    FIG. 2. Bithorax phenocopies showing variation in the extent of transformation of the metathoracic segment into mesothoracic-like structures. (a) Wild-type. (b and c) Phenocopies with different degrees of mesothoracic transformation, showing partial parts of transformed scutum, scutellum, and anterior compartment of wing-like structures. (d) Phenocopy dissected from a non-emerging pupa, showing a completely transformed mesothorax

    Measurement of Ubx Complementation Ability

    Single-pair crosses between a virgin female from different chromosome extraction strains of D. melanogaster and a male from the y(1) strain were carried out. The y(1) strain [(1)y; TM1: In(3LR) Me1 kniri–1 sbd1/TM2: In(3LR) emc2 Ubx130 es] was from Bloomington Stock Center (Indiana, USA), stock number 1628. The Ubx130 heterozygous progenies were identified by the presence of atypical bristles on the capitellums of the halteres (fig. 3). Halteres of the heterozygotes were removed with small needles, transferred to a drop of water on a glass microscope slide, and covered with a coverslip. The numbers of atypical bristles on the halteres on each side were scored under a light microscope.

    FIG. 3. Effect of the heterozygous Ubx130 mutation on haltere morphology. (a) Wild-type haltere. (b–d) Ubx130/Ubx+ halteres showing atypical bristles on the anterior margin of the capitellum, and an enlargement in size compared with the wild-type haltere. The numbers of atypical bristles vary among the Ubx130/Ubx+ halteres

    Results

    Evolutionary Changes in the Ubx Enhancer Sequences

    Neighbor-Joining trees based on the sequences of the bx-32.8, bxd4.1, and pbx32.7 enhancers are shown in figure 4(a–c). The overall picture shows that the topology of the trees confirms the widely accepted phylogeny of the melanogaster species subgroup, with D. erecta as the outgroup species (Ashburner 1989). The exception is the bxd4.1 tree, in which D. melanogaster is found on the same branch as D. simulans. However, the sequences of D. melanogaster are still very well grouped together: D. mauritiana and D. sechellia are sometimes found within the D. simulans branch, for example, as is seen in the bx-32.8 and bxd4.1 trees. This phylogenetic pattern is a result of these three species being closely related species that share recent coancestry. The variation in the positions of the D. melanogaster strains among the three trees shows evidence of recombination within the Ubx locus during evolutionary time. This variation is also confirmed by assessment of linkage disequilibrium. In the D. melanogaster sequences, there are four polymorphic sites in the bx-32.8 enhancer, seven in the pbx32.7 enhancer, and nine in the bxd4.1 enhancer, where the rarer variant is found more than once in the sample. Six pairwise linkage disequilibrium tests in bx-32.8 revealed one linkage disequilibrium that was significant at the 1% level and three others that were significant at the 5% level. For pbx32.7, all 21 pairwise linkage disequilibria were significant at the 1% level. Of 36 tests in bxd4.1, nine were significant at the 1% level and five others were significant at the 5% level. However, we found no significant evidence of linkage disequilibrium between variants in different enhancers.

    FIG. 4. Neighbor-Joining trees based on the bx-32.8 (a), bxd4.1 (b), and pbx32.7 (c) sequences. Each tree includes 23 sequences of D. melanogaster; 7 sequences of D. simulans, D. sechellia; 2 sequences of D. mauritiana; 2 sequences of D. teissieri, D. yakuba, D. orena, and D. erecta. The numbers represent percentages of 1,000 bootstrap estimates that support the branch

    Total numbers of fixed changes and of polymorphisms in the three enhancers are shown in Table 2. The numbers of fixed changes per nucleotide between D. melanogaster and D. simulans for the three enhancers (0.015 for bx-32.8, 0.009 for bxd4.1, and 0.032 for pbx32.7) are smaller than the average divergence rate of synonymous sites between the two species (0.072) (Sharp and Li 1989). This reveals that sequences of the three enhancers have evolved under purifying selection. The total numbers of fixed changes per nucleotide between all eight Drosophila species, for the three enhancers, are found to be significantly different among the enhancers (2 = 115.30, 2 df, P < 0.001). This suggests that the enhancers have evolved at different rates, with bx-32.8 the most conserved enhancer of the three, followed by bxd4.1 and pbx32.7.

    Table 2 Total Numbers of Nucleotide Variations in Ubx Enhancers.

    Nucleotide Polymorphisms and Neutrality Tests

    Table 3 shows a summary of nucleotide diversities within the D. melanogaster and the D. simulans samples, as well as the results from the neutrality tests. The D. simulans samples reveal higher nucleotide variation measured by and than do those from D. melanogaster, for all three enhancers. This is consistent with the study of Moriyama and Powell (1996), which showed that, on average, D. melanogaster has a lower degree of DNA polymorphism than D. simulans. In addition, the levels of nucleotide diversity among the three enhancers in the D. melanogaster samples are similar to the average level of nucleotide diversity relative to their local recombination rate (Begun and Aquadro 1992; also see Aquadro and Begun 1993).

    Table 3 Nucleotide Diversities and Neutrality Tests.

    Reflecting the numbers of fixed changes per nucleotide, bx-32.8 shows the lowest nucleotide diversity, followed by bxd4.1 and pbx32.7. Hudson-Kreitman-Aguade (HKA) tests for pairwise comparison between the three enhancers fail to reject the null hypothesis of neutral evolution (data not shown). The neutrality tests, Tajima's D and Fu and Li's D with an outgroup species, were performed to test for a non-neutral distribution of observed allele frequencies within the population (Fu and Li 1993; Tajima 1989). A significantly negative Tajima's D value is found for the bx-32.8 enhancer in the D. melanogaster samples. Correspondingly, significantly negative Fu and Li's D values are found for the bx-32.8 enhancer both for the D. melanogaster and the D. simulans samples. This suggests that there may have been a selective sweep operating on the bx-32.8 enhancer in D. melanogaster, and possibly a second in D. simulans. However, negative D values can be created by population expansion, which, without selection, results in high numbers of variable sites relative to the gene diversity. D. melanogaster and D. simulans are fruit fly species which have been commensal with humans for many thousands of years. As the human population has increased, populations of D. melanogaster and D. simulans have also expanded, which could have created the significantly negative D values of the bx-32.8 enhancer in the D. melanogaster and the D. simulans samples.

    Evolutionary Comparison of Protein-Binding and Non-Protein-Binding Sites

    Sequences of around 499 bp of the bx-32.8 enhancer and around 579 bp of the pbx32.7 enhancer have been investigated with a DNA footprinting technique to identify the DNA binding sites for some of the transcription factors encoded by gap, pair-rule, and dorsoventrally expressed genes (Pirrotta et al. 1995; Qian, Capovilla, and Pirrotta 1993). The variability within these two sequences is shown in figure 5a and 5b. Based on the results from the DNA footprinting, bases in the bx-32.8 and pbx32.7 enhancers can be classified into two types; protein-binding and non-protein-binding sites. Table 4 shows a summary of numbers of fixed changes between species and of polymorphisms in the protein-binding and non-protein-binding sites. The question is whether sequences of the protein-binding sites are more conserved than those of the non-protein-binding sites, because one might expect sequences of the protein-binding sites to evolve under a functional constraint. Interestingly, the results from table 4 show that numbers of fixed changes per nucleotide of the two base types are not significantly different from each other for either the bx-32.8 or pbx32.7 enhancers (2 = 0.852, 1 df, P = 0.356 and 2 = 0.004, 1 df, P = 0.949, respectively). This suggests that the sequences of both protein-binding and non-protein-binding sites in the bx-32.8 and pbx32.7 enhancers have evolved at the same rate.

    FIG. 5. Nucleotide variations in 499 bp of the bx-32.8 enhancer (a) and 579 bp of the pbx32.7 enhancer (b). A dot and dash indicate identity and a gap, respectively, in the aligned sequences. Above the sequences are shown the nucleotide positions. Numbers of fixed changes and polymorphisms at each base position are shown below the sequences. An insertion or a deletion of contiguous bases is counted as one change. Shaded boxes indicate bases in the protein-binding sites, which are referred from Qian, Capovilla, and Pirrotta (1993) and Pirrotta et al. (1995). Dark-shaded boxes indicate bases in the core sequences of the protein-binding sites. Changes in the core sequences which affect the matrix similarity values are highlighted in white

    Table 4 Evolutionary Comparisons Between Protein-Binding and Non-Protein-Binding Sites.

    Using data from intraspecific and interspecific variation, a neutrality test equivalent to the McDonald-Kreitman test can be applied to these noncoding sequences (Jenkins, Ortori, and Brookfield 1995). Ratios of numbers of fixed changes to polymorphisms between protein-binding and non-protein-binding sites are compared. The results show a significant difference for the bx-32.8 enhancer (P = 0.027) but no significant difference for the pbx32.7 enhancer (P = 0.093). This suggests that there might be some forms of positive selection operating on the bx-32.8 enhancer. It is noteworthy that this significant result in the bx-32.8 enhancer is due to the higher than expected numbers of polymorphisms found in the protein-binding sites and fixed changes seen in the non-protein-binding sites. It is therefore possible that either there is selection to maintain polymorphisms in the protein-binding sites or the sequences of the non-protein-binding sites have evolved adaptively. This observation is hard to explain, however, because little is known about the functions of the sequences of the non-protein-binding sites. It is possible that some of the "non-protein-binding" sites, in reality, bind unknown proteins. If as few as three of the 31 non-protein-binding sites in the bx-32.8 sequence that show fixed changes between species were reclassified as protein-binding, the result for bx-32.8 would no longer be significant.

    Changes in Protein-Binding Sites

    Both bx-32.8 and pbx32.7 enhancers contain DNA binding sites for the same transcription factors, which are Ftz, Hb, En, Tll, and Twi. Nucleotide variability in each type of transcription factor binding site in the bx-32.8 and pbx32.7 enhancers is shown in Figure 6a and 6b. The numbers of nucleotide changes of each protein-binding type, derived from the total numbers of changes in the bx-32.8 and pbx32.7 enhancers, are shown in table 5. One can ask whether the DNA binding sites of each transcription factor have evolved at the same rate and whether there is positive selection operating on any particular types of protein-binding sites. It is found that the numbers of fixed changes per nucleotide are not significantly different between the different types of proteins (2 = 1.518, 4 df, P = 0.823). This suggests that the sequences binding each type of protein have evolved at the same rate. Furthermore, the ratios of number of fixed changes to polymorphisms are also not significantly different between the types of protein bound (2 = 0.437, 4 df, P = 0.979), suggesting that there is no positive selection on any particular protein-binding type.

    FIG. 6. Changes to the protein-binding sites of the bx-32.8 (a) and pbx32.7 (b) enhancers. Shaded boxes indicate the core sequences, which give the highest matrix similarity value. The matrix similarity value of each interacting site is shown beside the sequence alignment. Nucleotide changes which cause changes in the matrix similarity values are indicated by dark boxes. The binding affinity of each protein-binding site in the bx-32.8 enhancer previously investigated by the DNA footprinting assay (Qian, Capovilla, and Pirrotta 1993) is described, but there is no information for the protein-binding sites in the pbx32.7 enhancer

    Table 5 Evolutionary Comparisons Between Bases Binding to Different Types of Transcription Factors.

    The core sequence within each protein-binding site was identified (fig. 6a and 6b). It is found that some of the large protein-binding sites contain more than one core sequence—for example, the Ftz3, Ftz4, Hb3, and Tll3 sites of the bx-32.8 enhancer. It is possible that those core sequences are individual protein-binding sites that are located close to each other. The matrix similarity value calculated for each core sequence seems to correlate with the binding affinity of the binding site investigated by a DNA footprinting assay. This can be seen from the comparison of the matrix similarity values of the binding sites in the bx-32.8 enhancer with the footprinting results of Qian, Capovilla, and Pirrotta (1993). For example, the weak Ftz1 and Ftz2 binding sites show low matrix similarity values (0.33 and 0.52, respectively) whereas the strong Tll1, Tll2, and Twi1 binding sites show high matrix similarity values of 0.95, 0.95 and 0.99. However, exceptions can be seen for some of the binding sites. For instance, the Twi3, Twi4, and En2 sites were shown to be weak binding sites, and their matrix similarity values are as high as 0.76, 0.76, and 0.73. This may be because the binding affinity of the site also depends on the structure of the DNA strand. In general, however, the matrix similarity can be used to reflect the binding affinity of the DNA target site to the transcription factors. A high matrix similarity value reflects a strong protein-binding site, whereas a low matrix similarity reflects a weak protein-binding site. It is interesting that not all of the conserved binding sites show high matrix similarity values. This implies that some of the weak binding sites are functionally important and have therefore been conserved among Drosophila species.

    Figure 6a and b shows some examples of specific changes in the protein-binding sites which might disrupt the ability of the binding sites to interact with transcription factors. This disruption is predicted by the nucleotide changes being in the core sequence of the binding sites and causing a great change in the matrix similarity values. For example, the fixed substitutions in Twi2 binding sites of the bx-32.8 enhancer in D. teissieri, D. yakuba, D. orena, and D. erecta cause a reduction of the matrix similarity value of the sequences lower than that of D. melanogaster (fig. 6a). It is possible that this bx-32.8 Twi2 binding site is no longer a protein-binding site in those Drosophila species, or it may be only a weak binding site. Another example is in the Hb2 and Twi1 binding sites of the pbx32.7 enhancer (fig. 6b). The fixed changes in the pbx32.7 Hb2 of D. yakuba, D. orena, and D. erecta, and those in the pbx32.7 Twi1 of D. teissieri and D. yakuba, also cause a great reduction of the matrix similarity value of the sequences. Many of the polymorphisms are also found to greatly affect the matrix similarity values, as can be seen in the bx-32.8 Ftz1, Ftz2, Twi2, Twi4 binding sites and the pbx32.7 Ftz4, and Hb2 binding sites. Interestingly, some of those polymorphisms are found within the species while the binding sites are highly conserved among Drosophila species. For instance, the "T" polymorphism in the Twi4 binding site of the bx-32.8 enhancer in D. melanogaster causes a reduction of the matrix similarity value from 0.76 to 0.53, and the "G" polymorphism in the bx-32.8 Ftz2 binding site of the D. simulans species complex causes a reduction of the matrix similarity value from 0.52 to 0.27.

    Numbers of fixed changes per nucleotide between bases in the core sequences and those outside the core sequences, including the non-protein-binding sites, have been compared (data not shown). Interestingly, although in this case only small regions of the protein-binding sites are considered, the numbers of fixed changes per nucleotide in the core sequences and those outside the core sequences are not significantly different from each other for either the bx-32.8 or pbx32.7 enhancers (2 = 1.540, 1 df, P = 0.215 and 2 = 0.809, 1 df, P = 0.368, respectively). This finding shows that the sequences outside the core sequences have evolved at the same rate as those in the core sequences. In addition, comparison of the fixed changes to polymorphism ratios between nucleotide changes which affect and do not affect the matrix similarity values yields an insignificant result (P = 0.537), suggesting that the pattern of nucleotide variation in the core sequences is consistent with a neutral model (data not shown).

    Variation in Ubx Homeotic Stability

    Table 6 shows the rates of induction of bithorax phenocopies for each of the D. melanogaster strains, along with the number of atypical bristles on halteres of strains heterozygous for Ubx130. The strains show remarkable variation in the phenocopy frequencies, with Monty8(1) showing the highest frequency of phenocopies. The transformation of the metathoracic to the mesothoracic segment in phenocopies varies in size (fig. 2). We found that the higher the phenocopy frequency, the more severe the transformation appears to be. The expressivity of bithorax phenocopies observed is consistent with the expressivity reported in previous studies (Bownes and Seiler 1977; Capdevila and García-Bellido 1978; Ho, Bolton, and Saunders 1983) in many respects, such that

    The numbers of individuals showing transformations on both sides, and the numbers of individuals showing transformations on both haltere and leg imaginal discs of the same side, are higher than the expected numbers under an assumption of independence.

    The effects of ether vapor are seen more in the haltere imaginal disc than in the leg imaginal disc.

    Phenocopies show transformation of the anterior compartment more frequently than the posterior compartment.

    Table 6 Means and Standard Errors of The Frequencies of Phenocopies and the Numbers of Atypical Bristles for Each of the D. melanogaster Strains.

    As with the phenocopy induction rates, the numbers of atypical bristles on the halteres of flies heterozygous for Ubx130 also vary significantly among the strains. The mean number of bristles varies between 1–9 bristle(s) per haltere. Heterozygous Loua(4) shows the highest number of atypical bristles, suggesting that this strain has the lowest ability to rescue the Ubx phenotype. We even found that a few heterozygous Loua(4) flies showed transformation of the thoracic compartment, or parts of the haltere which had a wing-like appearance (data not shown). We also found, based on information from the Monty5(6), Draveil(3), Petit Bourg(3), Draveil(2), and Monty8(1) strains, that the frequency of phenocopies is not significantly correlated across strains with the number of atypical bristles (r = –0.271, P = 0.660). The Monty8(1) strain, which has a very high frequency of phenocopies, shows a low number of atypical bristles when heterozygous with Ubx130.

    Discussion

    TOP

    Abstract

    Introduction

    Materials and Methods

    Results

    Discussion

    Acknowledgements

    Literature Cited

    Conservation and Divergence of Ubx Enhancer Sequences

    The sequences of the bx-32.8, bxd4.1, and pbx32.7 enhancers are found to be highly conserved. This reveals a strong functional constraint of the enhancers in the Drosophila melanogaster subgroup species. The differences in conservation levels among the three enhancers suggest that they have evolved subject to different degrees of constraint. The bx-32.8 enhancer shows the highest degree of constraint among the three enhancers studied. Unlike the others, the bx-32.8 enhancer is located in an intron. The requirement for transcription through this sequence, and the requirement for accurate splicing of this sequence from the transcript, might potentially impose further selective constraint on an enhancer sequence in an intron. Such constraint would add to the selection arising from its impact on the initiation of transcription. Using information about the protein-binding sites in the bx-32.8 and pbx32.7 enhancers, we have found that the sequences of the target sites of each transcription factor in the enhancers have evolved at the same rate. The sequences of both protein-binding and non-protein-binding sites reveal the same rate of divergence for either the bx-32.8 or the pbx32.7 enhancers. This observation can be seen, even in the comparison of the core sequences with the others, suggesting an equal degree of functional constraint for the two base types. It is possible that the sequences of the non-protein-binding sites function as binding sites for other transcription factors, or they may function to maintain the secondary structure of the enhancers. Both bx-32.8 and pbx32.7 enhancers drive Ubx expression during early embryogenesis (Qian, Capovilla, and Pirrotta 1991; Pirrotta et al. 1995). It is possible that other gap and pair-rule genes or dorsoventrally expressed genes might bind to the bx-32.8 and pbx32.7 enhancers, but they have not been experimentally investigated.

    Neutrality tests largely fail to reject a neutral evolution model. Although significantly negative Tajima's D, and Fu and Li's D for D. melanogaster and D. simulans are found for the bx-32.8 enhancer, this is probably a result of population expansion rather than an effect of selection. The evolutionary comparisons of nucleotide variation between types of transcription-factor-bound protein-binding, and non-protein-binding sites, and the changes affecting and not affecting the matrix similarity values of the core sequences, show that the pattern of nucleotide variation in the bx-32.8 and pbx32.7 enhancers is consistent with a neutral evolution model. One exception is in the comparison of nucleotide variation in the protein-binding and non-protein-binding sites of the bx-32.8 enhancer, which shows a significant deviation from neutrality. The data show a significant excess of protein-binding polymorphisms and/or non-protein-binding fixed changes in the bx-32.8 enhancer (table 4). In principle, if this significance is not simply the outcome of multiple testing without Bonferroni correction, two types of explanation are possible: either there is selection to maintain diversity in the protein-binding sites in this bx region, or there has been adaptive change in the non-protein-binding sites.

    Sequences of the protein-binding sites vary among the Drosophila species. The fact that the nucleotide changes are not shared by all the protein-binding sites for a given transcription factor implies that the changes observed do not result from any change in the binding site that is the target for binding by an altered transcription factor in another species. In this study, the matrix similarity value was used to reveal dynamic changes of the protein-binding sites between species or strains. The change in matrix similarity value resulting from a substitution depends on the Ci value of its position. Changes in a highly conserved position (with a high Ci value) may greatly change the matrix similarity values. Because the only way a base change can fail to change the matrix similarity value is if the scores of the two bases at the position are the same, most changes will affect the matrix similarity values. Although there is no evidence to measure quantitatively the degree to which a change in the matrix similarity value causes a difference in the biochemical function of a site, the matrix similarity value correlates with the binding affinity of the binding site investigated by DNA footprinting to the extent that a strong or a weak protein-binding site can be revealed by its having a high or a low matrix similarity value (see Results). Important changes in the protein-binding sites are thus predicted if the binding site homologs show a great difference in their matrix similarity values. In this case, we found 4 of the 17 binding sites of the bx-32.8 enhancer (Ftz1, Ftz2, Twi2, and Twi4) and 3 of the 22 binding sites in the pbx32.7 enhancer showing a great difference in the matrix similarity values.

    Although these data suggest dynamic changes in some protein-binding sites in the bx-32.8 and pbx32.7 enhancers over the microevolutionary scale of Drosophila evolution, it is possible that these dynamic changes in the protein-binding sites are, at most, only slightly deleterious and do not affect the functions of the enhancers. It is speculated that stabilizing selection may also play a role in the evolution of Ubx regulatory sequences, in the same way seen in the even-skipped stripe 2 enhancer (Ludwig et al. 2000), where selection has operated to maintain the overall functions of the Ubx enhancers during evolution, while allowing a rapid turnover of their sequences to occur. This may happen through the fixation of slightly deleterious mutations being compensated by adaptive changes at other sites. In this case, failure to achieve functional compensation will result in changes in the gene-expression pattern. Alternatively, it can also be imagined that sequence turnover may occur in the way in which new protein-binding sites are created by the fixation of neutral or selected mutations, and thus, some of the former protein-binding sites that are no longer necessary can be lost. These two models of enhancer evolution also explain the outcome of the gain and loss of protein-binding sites in regulatory sequences over the macroevolutionary scale, which is observed in regulatory sequences of many regulatory genes including Ubx (Dellino, Tatout, and Pirrotta 2002; Wilde and Akam 1987).

    Impact of Intraspecific Variation at the Ubx Locus

    The failure to detect a sign of positive selection does not necessarily mean that all of the variability observed in Ubx sequences is neutral. It is possible that the neutrality tests are not sensitive enough to detect any sign of selection. For example, in the McDonald and Kreitman test, and other tests for adaptive change, if only a few adaptive changes occur, these will be insufficient to make the results of the tests significant (Brookfield and Sharp 1994). In this study, to measure functional variation of the Ubx alleles among D. melanogaster strains, we used two phenotypic tests, the induction of bithorax phenocopy and the Ubx130/Ubx+ complementation test. The results show that the D. melanogaster strains differ greatly in their frequencies of phenocopy induction, as well as in the numbers of atypical bristles on their haltere when they are heterozygous with Ubx130. However, variation in the phenocopy frequency and the atypical bristle numbers is not correlated with the positions of the strains on the Ubx trees in figure 4. Strains that show significantly different phenocopy frequencies, or that have significantly different numbers of atypical haltere bristles (when heterozygous with Ubx130), can be found to have identical Ubx regulatory sequences, whereas strains that show the same levels of phenocopy frequencies and numbers of haltere atypical bristles can be found on different branches in the trees. It is possible that the variation of these two phenotypes is caused by variation in the other regulatory regions of Ubx, in addition to bx-32.8, bxd4.1 and pbx32.7, or it may be that the variation of the phenotypes is a result of the interaction of variations in several regions of the genome in addition to Ubx.

    The differences in the phenocopy frequencies among the strains cannot be explained by their Ubx alleles having differing levels of complementation ability because there is no significant correlation between the phenocopy frequency and the strength of the Ubx130/Ubx+ mutant phenotype (P = 0.660). An obvious example is that Monty8(1) shows a very high phenocopy frequency but has a low number of atypical bristles on the haltere when the strain is heterozygous with Ubx130. It is possible that these two phenotypes are the result of different mechanisms. The differences in the numbers of atypical bristles on the haltere among Ubx130 heterozygous strains may be caused by differences in the Ubx alleles of each strain, or it may be the result of differences in the response of downstream genes to one dose of UBX protein. However, the differences in the ether susceptibility between the strains may be caused by differences in the temporal and spatial regulation of Ubx expression, which may include the regulation of both the initial activation/repression and the maintenance of activating/repressing states. Because the sensitive period for ether induction of phenocopy is very specific, if the different strains were to be laying eggs at differing stages of development, they would differ systematically in the stage of development at which their eggs were being exposed to ether. If so, this could potentially create consistent differences between strains in their frequencies of phenocopies. However, other data (Phinchongsakuldit and Brookfield, in preparation) imply that at least some of the variation in phenocopy frequency maps to Ubx, which would not be expected under this model.

    A QTL (quantitative trait loci) study for haltere/wing margin transformation in Ubx1 heterozygotes showed that the variation in genes mapping to the left arm of chromosome III, at bands 62–65, is responsible for the variation in the atypical bristle numbers on the Ubx1/Ubx+ halteres (Gibson, Wemple, and van Helden 1999). This result gives an example of the genetic variation outside the Ubx locus that affects the variation of Ubx homeotic stability. Both Ubx130 and Ubx1 are null mutations, but, whereas Ubx130 is associated with many inversions on the third chromosome, Ubx1 is a point mutation caused by a doc element insertion in the first exon. However, the two mutations have the same phenotypic effects of creating an enlarged haltere with atypical bristles in heterozygotes. It is possible that the locus discovered by Gibson, Wemple, and van Helden may also affect the phenotype of individuals heterozygous for the Ubx130 allele. It was shown that the variation in the ether susceptibility is at least partly caused by polymorphism in the Ubx regulatory sequences (Gibson and Hogness 1996). Other genetic backgrounds contributing to the variation in ether susceptibility have not been identified.

    Acknowledgements

    We thank Keith Bradman and Roy Chaudhuri for their assistance and suggestions with the sequence analysis computer programs. We also thank Michael Akam for his comments and suggestions and Casey Bergman for supplying us with binding-site information. This work was supported by the Biotechnology and Biological Sciences Research Council through a grant awarded to S.MacA., and by the Thai Government (Development and Promotion of Science and Technology Talent Project) through a grant awarded to J.P.

    Literature Cited

    Aquadro, C. F., and D. J. Begun. 1993. Evidence for and implications of genetic hitchhiking in the Drosophila genome. Pp 159–178 in N. Takahata and A. G. Clark, eds. Mechanism of molecular evolution. Sinauer Associates, Sunderland. Mass.

    Arnone, M. I., and E. H. Davidson. 1997. The hardwiring of development: organization and function of genomic regulatory systems. Development 124:1851-1864.

    Ashburner, M. 1989. Drosophila: a laboratory handbook. Cold Spring Harbor Laboratory Press, Cold Spring Harbor, NY.

    Averof, M., and N. H. Patel. 1997. Crustacean appendage evolution associated with changes in Hox gene expression. Nature 388:682-686.

    Begun, D. J., and C. F. Aquadro. 1992. Levels of naturally occurring DNA polymorphism correlate with recombination rate in D. melanogaster. Nature 356:519-520.

    Bienz, M., G. Saari, G. Tremml, C. Muller, J. Zust, and P. A. Lawrence. 1988. Differential regulation of Ultrabithorax in two germ layers of Drosophila. Cell 53:567-576.

    Bomze, H. M., and A. J. López. 1994. Evolutionary conservation of the structure and expression of alternatively spliced Ultrabithorax isoforms from Drosophila. Genetics 136:965-977.

    Bownes, M., and M. Seiler. 1977. Developmental effects of exposing Drosophila embryos to ether vapour. J. Exp. Zool. 199:9-24.

    Brookfield, J. F. Y., and P. M. Sharp. 1994. Neutralism and selectionism face up to DNA data. Trends Genet 10:109-111.

    Capdevila, M. P., and A. García-Bellido. 1978. Phenocopies of bithorax mutants. Wilhelm Roux's Arch. 185:105-126.

    Chiang, A., M. B. O'Connor, R. Paro, J. Simon, and W. Bender. 1995. Discrete polycomb-binding sites in each parasegmental domain of the bithorax complex. Development 121:1681-1689.

    Dellino, G. I., C. Tatout, and V. Pirrotta. 2002. Extensive conservation of sequences and chromatin structures in the bxd Polycomb response element among Drosophilid species. Int. J. Dev. Biol. 46:133-141.

    Fu, Y. X., and W. H. Li. 1993. Statistical tests of neutrality of mutations. Genetics 133:693-709.

    Gibson, G., and D. S. Hogness. 1996. Effect of polymorphism in the Drosophila regulatory gene Ultrabithorax on homeotic stability. Science 271:200-203.

    Gibson, G., M. Wemple, and S. Van Helden. 1999. Potential variance affecting homeotic Ultrabithorax and Antennapedia phenotypes in Drosophila melanogaster. Genetics 151:1081-1091.

    Higgins, D. G., and P. M. Sharp. 1988. Clustal: a package for performing multiple sequence alignment on a microcomputer. Gene 73:237-244.

    Ho, M. W., E. Bolton, and P. T. Saunders. 1983. Bithorax phenocopy and pattern formation: I. Spatiotemporal characteristics of the phenocopy response. Exp. Cell Biol. 51:282-290.

    Irvine, K. D., S. L. Helfand, and D. S. Hogness. 1991. The large upstream control region of the Drosophila homeotic gene Ultrabithorax. Development 111:407-424.

    Jenkins, D. L., C. A. Ortori, and J. F. Y. Brookfield. 1995. A test for adaptive change in DNA sequences controlling transcription. Proc. R. Soc. Lond. Ser. B 261:203-207.

    Kassis, J. A., C. Desplan, D. K. Wright, and P. H. O'Farrell. 1989. Evolutionary conservation of homeodomain-binding sites and other sequences upstream and within the major transcription unit of the Drosophila segmentation gene engrailed. Mol. Cell. Biol. 9:4304-4311.

    Kassis, J. A., S. J. Poole, D. K. Wright, and P. H. O'Farrell. 1986. Sequence conservation in the protein coding and intron regions of the engrailed transcription unit. EMBO J. 5:3583-3589.

    Kim, J. 2001. Macro-evolution of the hairy enhancer in Drosophila species. J. Exp. Zool. (Mol. Dev. Evol.) 291:175-185.

    Langeland, J. A., and S. B. Carroll. 1993. Conservation of regulatory elements controlling hairy pair-rule stripe formation. Development 117:585-596.

    Lewis, E. B. 1978. A gene complex controlling segmentation in Drosophila. Nature 276:565-570.

    Ludwig, M. Z., C. Bergman, N. H. Patel, and M. Kreitman. 2000. Evidence for stabilizing selection in a eukaryotic enhancer element. Nature 403:564-567.

    Ludwig, M. Z., N. H. Patel, and M. Kreitman. 1998. Functional analysis of eve stripe 2 enhancer evolution in Drosophila: rules governing conservation and change. Development 125:949-958.

    Lukowicz, W., C. Schr?der, G. Glaser, M. Hülskamp, and D. Tautz. 1994. Regulatory and coding regions of the segmentation gene hunchback are functionally conserved between Drosophila virilis and Drosophila melanogaster. Mech. Dev. 45:105-115.

    Maier, D., A. Preiss, and J. R. Powell. 1990. Regulation of the segmentation gene fushi tarazu has been functionally conserved in Drosophila. EMBO J. 9:3957-3966.

    McDonald, J. H., and M. Kreitman. 1991. Adaptive protein evolution at the Adh locus in Drosophila. Nature 351:652-654.

    Moriyama, E. N., and J. R. Powell. 1996. Intraspecific nuclear DNA variation in Drosophila. Mol. Biol. Evol. 13:261-277.

    Müller, J., and M. Bienz. 1992. Sharp anterior boundary of homeotic gene expression conferred by the fushi tarazu protein. EMBO J. 11:3653-3661.

    Myers, E. W., and W. Miller. 1988. Optimal alignment in linear space. Comp. Appl. Bios. 4:11-17.

    Page, R. D. M. 1996. TreeView: an application to display phylogenetic trees on personal computers. Comp. Appl. Bios. 12:357-358.

    Pan, D., S. A. Valentine, and A. J. Courey. 1994. The bipartite D. melanogaster twist promoter is reorganized in D. virilis. Mech. Dev. 46:41-53.

    Pirrotta, V., C. S. Chan, D. McCabe, and S. Qian. 1995. Distinct parasegmental and imaginal enhancers and the establishment of the expression pattern of the Ubx gene. Genetics 141:1439-1450.

    Qian, S., M. Capovilla, and V. Pirrotta. 1991. The bx region enhancer, a distinct cis-control element of the Drosophila Ubx gene and its regulation by hunchback and other segmentation genes. EMBO J. 10:1415-1425.

    Qian, S., M. Capovilla, and V. Pirrotta. 1993. Molecular mechanisms of pattern formation by the BRE enhancer of the Ubx gene. EMBO J. 12:3865-3877.

    Quandt, K., K. Frech, H. Karas, E. Wingender, and T. Werner. 1995. MatInd and Matlnspector: new fast and versatile tools for detection of consensus matches in nucleotide sequence data. Nucleic Acids Res. 23:4878-4884.

    Rozas, J., and R. Rozas. 1999. DnaSP version 3: an integrated program for molecular population genetics and molecular evolution analysis. Bioinformatics 15:174-175.

    Sackerson, C. 1995. Patterns of conservation and divergence at the even-skipped locus of Drosophila. Mech. Dev. 51:199-215.

    Saitou, N., and M. Nei. 1987. The Neighbor-Joining method—a new method for reconstructing phylogenetic trees. Mol. Biol. Evol. 4:406-425.

    Sharp, P. M., and W. H. Li. 1989. On the rate of DNA sequence evolution in Drosophila. J. Mol. Biol. 28:398-402.

    Simon, J., M. Peifer, W. Bender, and M. O'Connor. 1990. Regulatory elements of the bithorax complex that control expression along the anterior-posterior axis. EMBO J. 9:3945-3956.

    Stern, D. L. 1998. A role of Ultrabithorax in morphological differences between Drosophila species. Nature 396:463-466.

    Stern, D. L. 2000. Evolutionary developmental biology and the problem of variation. Evolution 54:1079-1091.

    Tajima, F. 1989. Statistical method of testing the neutral mutation hypothesis by DNA polymorphism. Genetics 123:585-595.

    Tautz, D. 2000. Evolution of transcriptional regulation. Curr. Opin. Gen. Dev. 10:575-579.

    Terol, J., M. Perez-Alonso, and R. deFrutos. 1995. Molecular characterization of the zerknüllt region of the Antennapedia complex of D. subobscura. Chromosoma 103:613-624.

    Thompson, J. D., D. G. Higgins, and T. J. Gibson. 1994. CLUSTAL—W—Improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Res. 22:4673-4680.

    Treier, M., C. Pfeifle, and D. Tautz. 1989. Comparison of the gap segmentation gene hunchback between Drosophila melanogaster and Drosophila virilis reveals novel modes of evolutionary change. EMBO J. 8:1517-1525.

    Warren, R. W., L. Nagy, J. Selegue, J. Gates, and S. B. Carroll. 1994. Evolution of homeotic gene regulation and function in flies and butterflies. Nature 372:458-461.

    Weatherbee, S. D., H. F. Nijhout, L. W. Grunert, G. Halder, R. Galant, J. Selegue, and S. B. Carroll. 1999. Ultrabithorax function in butterfly wings and the evolution of insect wing patterns. Curr. Biol. 9:109-115.

    Wilde, C. D., and M. E. Akam. 1987. Conserved sequence elements in the 5' region of the Ultrabithorax transcription unit. EMBO J. 6:1393-1401.

    Wolff, C., M. Pepling, P. Gergen, and M. Klingler. 1999. Structure and evolution of a pair-rule interaction element: runt regulatory sequences in D. melanogaster and D. virilis. Mech. Dev. 80:87-99.(Jaros Phinchongsakuldit, )