Characterization of partial Hox gene sequences in annual fish of the subfamily Cynolebiatinae ( Cyprinodontiformes , Rivulidae )

Hox genes encode a family of transcription factors implicated in conferring regional identity along the anteroposterior axis in developing animal embryos. These genes are organized in genomic clusters, expressed collinearly and highly conserved in vertebrates. Among teleost, South American annual killifishes of the Cynolebiatinae subfamily represent an excellent model in development studies because their embryos are capable of undergoing reversible developmental arrest (diapause) at three well-defined morphological stages. They are also an excellent model for evolutionary studies due to the high rates of mutation of their mitochondrial genome, their karyotypic divergence and their morphological variability. In this study, three partial homeobox sequences were isolated from different species of the Cynolebiatinae subfamily. Phylogenetic analyses and sequence comparisons revealed that they belong to the anterior Hox complex group, specifically to paralogue groups 1 and 3. This is the first time that partial Hox genes have been described in species of the Cynolebiatinae subfamily.


Introduction
The homeobox genes play important roles in the developmental processes of many multicellular organisms.The homeobox was originally described as a conserved protein-coding sequence of about 180 base pairs, in many developmental control genes of Drosophila.The encoded proteins are evolutionarily conserved transcription factors that share a 60 amino acids DNA-binding domain (the homeodomain, or HD) and regulate axial patterning segment or cell identity and proliferation (Lewis, 1978;McGinnis et al., 1984a, b).Homeodomain proteins regulate diverse developmental programs by modulating expression patterns of target genes in a temporal, spatial and tissue-specific manner (Gehring et al., 1994a).X-ray crystallographic and NMR spectroscopic analyses on several members of this family revealed that the structure consists of a tri α-helical core and an N-terminal arm that becomes partially ordered in the presence of DNA (Qian et al., 1989;Kissinger et al., 1990;Gehring et al., 1994a;Wolberger, 1996;Passner et al., 1999).The N-terminal region, just upstream of helix 1, contacts the minor groove of DNA, whereas helix 3, known as the "recognition helix", binds to the major groove.The mode of HD-DNA interaction appears to be highly conserved within HD-containing proteins (Gehring et al., 1994b;Chauvet et al., 2000).
In vertebrates, the Hox genes (a class of homeobox genes) are organized in linked chromosomal clusters and show a striking colinearity in their 5'-3' chromosomal position.They also show remarkably similar expression patterns along the anteroposterior (AP) axis, specifying the identities of different body regions during embryogenesis (McGinnis and Krumlauf, 1992).While invertebrates studied to date have a single Hox cluster, tetrapod vertebrates have four separated Hox clusters (termed clusters A to D) lying on four different chromosomes.Within each cluster, the genes have been classified into 13 paralogue groups according to sequence homology and location in the genome (Krumlauf, 1994).The expression pattern of Hox genes during development and homeodomain sequence comparisons indicate relatedness between anterior (Hox 1-3), medial (Hox 4-8) and posterior (Hox 9-13) paralogue groups (Ruddle et al., 1994).Some paralogue groups can be recognized by the homeodomain sequence alone, others only by using characteristic residues outside the homeodomain (Sharkey et al., 1997).
Although it is widely assumed that vertebrates have four Hox clusters, studies on ray-finned fish have shown that they have more Hox clusters than tetrapods.These extranumeral Hox clusters result from a genome duplication event that is specific for the fish (actinopterygian) lineage.Sets of seven Hox clusters have been described in Danio rerio (Amores et al., 1998;Prince et al., 1998), in three pufferfish species, Takifugu rubripes, Tetraodon nigroviridis and Sphenoides nephelus (Aparicio et al., 1997(Aparicio et al., , 2002;;Jaillon et al., 2004;Amores et al., 2004) and in Oryzias latipes (Naruse et al., 2000).However, preliminary data indicates four Hox clusters in the non-annual killifish Fundulus heteroclitus (Misof and Wagner, 1996).The four-cluster situation is also retained in the sarcopterygian lineage (Koh et al., 2003), in basal ray-finned fishes, such as the bichir Polypterus senegalus (Ledje et al., 2002), and in the horn shark Heterodontus francisci (Kim et al., 2000).The study of Hox clusters, therefore, provides the opportunity to understand the relationship between molecular evolution, genome organization and gene expression.
The Cynolebiatinae subfamily (Cyprinodontiformes: Rivulidae) is a specious group of annual fish distributed from northeastern Brazil to northeastern and southern Argentina, Uruguay and Paraguay (Costa, 1995).They live in temporary ponds and each generation completes a full life cycle within one year.The population survives dry seasons in the form of eggs buried in the mud.During the subsequent rainy season, the ponds refill, the eggs hatch and the larvae rapidly grow to sexual maturity and reproduce (Wourms, 1967(Wourms, , 1972a)).The developmental pattern of annual killifish embryos is unique compared to other teleosts, and is characterized by dispersion and subsequent reaggregation of blastomeres and the occurrence of embryonic diapause.Diapause is a state of developmental arrest that precedes the onset of unfavorable environmental conditions, is promoted by genetic and environmental factors and typically occurs as part of the natural developmental program.Thus, embryos may enter diapause even under conditions considered optimal for development.Annual killifish embryos can undergo reversible developmental diapause at one or all of three well-defined morphological stages (diapause I, II and III) (Wourms, 1967(Wourms, , 1972a, b;, b;Hand and Podrabsky, 2000).Studies in diapausing embryos of the annual killifish Austrofundulus limnaeus reveal that during diapause II protein synthesis is substan-tially diminished due to depressed metabolism.At this stage embryos appear to be more resistant to environmental stresses such as anoxia and dehydration (Podrabsky andHand, 1999, 2000;Podrabsky et al., 2001).
To date, early development has been studied in some Austrolebias species (subfamily Cynolebiatinae, sensu Costa 1995;1998) and available data show that members of this genus display a dispersion-reaggregation process and that they undergo facultative arrest at diapause I and II and obligate arrest at diapause III (Wourms, 1972a;Carter and Wourms, 1991;Arezo et al., 2005).
The development of annual fish embryos exhibits several characteristics that are of interest and there is a lack of useful information about their developmental genetic control.The study presented here is the first description concerning Hox genes in different teleost fish species belonging to the Cynolebiatinae.We established their cluster affiliation and investigated if amino acid changes in homeodomain primary sequence affect its ternary structure and the homeodomain-DNA interaction.Our results suggest that the sequences isolated by us are highly homologues to the anterior genes of the Hox complex of other vertebrates.

Specimens
We studied four Austrolebias species (A. bellottii, A. charrua, A. cheradophilus and A. viarius) from temporary ponds in Uruguay and A. adloffi from the southernmost Brazilian state of Rio Grande do Sul (which borders Uruguay) kindly supplied by L. Malabarba (1991).A specimen of Brevoortia aurea (Clupeiformes: Alosinae) was also examined (Table 1).Tissues and voucher specimens are deposited in the Evolutionary Genetics Section, Science Faculty, University of the Republic (Sección Genética Evolutiva, Facultad de Ciencias, Universidad de la República), Montevideo, Uruguay.

Genomic DNA extraction, amplification, cloning and sequencing
Genomic DNA was isolated from ethanol-fixed liver tissue using sodium chloride protein precipitation followed by ethanol precipitation of total DNA (modified from Medrano et al., 1990).To ascertain the quality of DNA, an aliquot of the extract DNA was run in ethidium bromide stained 1% (w/v) agarose gel using 1XTAE (Tris-Acetate-EDTA).
Every PCR fragment was subsequently cloned into the pGEM ® -T Easy vector using the reagents and protocols supplied by the manufacturer (Promega).The recombinant DNA was obtained by DNA minipreparation of individual clones by alkaline lysis (Sambrook et al., 1989).The clones were screened for correctly sized inserts by electrophoresis in ethidium bromide stained 0.8% (w/v) agarose gel using 1XTAE, followed by a PCR with the specific Sog1 and Sog2 primers.Sequencing reactions were performed on each template using the Sog2 primer and were analyzed in an automated ABI PRISM 377 DNA Sequencer.Nucleotide sequences from A. bellottii, A. cheradophilus and A. viarius were conceptually translated using the commercially available computer software Gene Runner 3.05 for Windows (Hastings Software).The amino acid sequences were compared against the National Center for Biotechnology Information (NCBI) protein database by using the BLASTp program on the Basic Local Alignment Search Tool (BLAST) network service (Altschul et al., 1990).

Southern blot analysis
To verify the existence of Hox sequences in the Cynolebiatinae subfamily we first performed a Southern blot analysis of genomic DNA from two Austrolebias sister taxa (A.viarius and A. charrua) using the A. cheradophilus homeobox amplified fragment as a probe to screen for homologous sequences.We also carried out a Southern blot analysis of genomic DNA from three other Austrolebias species (A.adloffi, A. bellottii and A. charrua) and B. aurea using the A. bellottii amplified fragment as a probe.
For Southern blotting analysis genomic DNA was isolated as describe above from different species of Austrolebias and the distantly related taxon B. aurea (Table 1).Approximately 10 µg of DNA was digested with Hinf-I (10000 U/mL) (Promega) at 37 °C for 24 h.After digestion, the DNA was separated using 0.5X TBE (Tris-Borate-EDTA) and 1.5% (w/v) agarose gel stained with ethidium bromide and transferred onto Hybond-N+ nylon membranes (Amersham Life Science) and the blot hybridized consecutively with two different probes using the ECL direct nucleic acid labeling and detection systems (catalogue no.RPN3001; Amersham Biosciences).The probes used were the vectors containing the amplified A. cheradophilus and A. bellottii fragment, both digested with 20000 U mL -1 of the Nde-I restriction enzyme (Biolabs).Probe hybridization was carried out overnight in 0.5 M NaCl at 42 °C and the hybridized bands visualized after 2 h exposure according to the instructions of the kit.

Phylogenetic analyses
Homeodomain amino acid sequences from A. bellottii, A. viarius, A. cheradophilus and partial homeodomain sequences for anterior (1-3), medial (4-8) and posterior (9-13) Hox paralogue groups from other metazoans obtained from GenBank, were aligned by using the CLUSTALX 1.8 program (Thompson et al., 1997) (Table 2).Partial homeodomain sequences from A. gymnoventris and A. luteoflammulatus, previously isolated in our laboratory (Pereiro and Garcia, unpublished data) were also included in our alignment.The MEGA 3 program (Kumar et al., 2004) was used to construct a neighborjoining (NJ) tree (Saitou and Nei, 1987) using a p-distance model and complete deletion.Bootstrap values for the nodes were determined by analyzing 1000 bootstrap replicates to estimate the strength of the groupings.

Nucleotide substitution, structural and functional analyses
For nucleotide substitution analysis we used the DnaSP 4.00 computer program (Rozas et al., 2003) to compare synonymous (Ks) and non-synonymous (Ka) sites, to calculate ratio rates (ω = Ka/Ks) and to evaluate the codon usage bias between Cynolebiatinae species and their respective paralogue groups, in order to detect positive selection in highly conserved Hox gene fragments.
For structural and functional analysis the homeodomain three-dimensional structure was predicted by comparative protein modeling methods using the automated Swiss Model program (Schwede et al., 2003).The structure was modeled on the basis of its structural similarity with the human HoxB1 homeodomain (Protein data Bank entry 1972).The degree of identity between the template and our isolated partial Hox sequences enabled a preliminary model to be generated using the SWISS-MODEL program.The Swiss-PdbViewer 3.7 program, service pack 5, (Guex 496 Hox genes of annual fish and Peitsch, 1997) was then used to analyze and visualize the structures.

PCR products and Southern blot Hox genes analyses
The PCR amplifications using the Sog1 and Sog2 specific primers detected fragments of 105 bp in A. bellottii, 115 bp in A. cheradophilus and 115 bp in A. viarius.The deduced amino acid sequences analyses (BLASTp and sequence comparison) identified them as Hox genes.The A. cheradophilus and A. viarius amino acid sequence alignments were identical, so we will not differentiate between these Hox genes in the text.The nucleotide and deduced amino acid sequences of these subfamily species have been deposited in the GenBank database (Table 2).
The Southern blot analysis detected a single band in both A. viarius and A. charrua using the A. cheradophilus homeobox amplified fragment as a probe to screen for homologous Hox sequences (Figure 1a).A single band was also detected in B. aurea, A. adloffi, A. bellottii and A. charrua using the A. bellottii amplified fragment as a probe (Figure 1b).

Synonymous vs. Non-synonymous nucleotide substitutions
Although amino acid substitution rates increase, they never exceed those of synonymous substitutions because relaxed purifying selection is acting (Hartl, 2000).However, the non-synonymous (Ka) to synonymous (Ks) substitution ratio (ω) was significantly larger than 1, if positive selection is acting.Our analysis of synonymous nucleotide substitutions reveals that the species analyzed had a Ks value of 26, whereas the other PG1 species had a Ks value of 23.For non-synonymous substitutions, we obtained a Ka value of 89 for the Cynolebiatinae subfamily species and 81 for the remaining PG1 species (Table 2).In both cases, Cynolebiatinae subfamily species presented higher values with respect to other species of the paralogue group 1.Furthermore, when we compared the ratio rates (= Ka/Ks) between the Cynolebiatinae subfamily and other PG1 species we found that they varied from 0.14 to 1.88, indicating that a weak footprint of positive selection was acting (Hartl, 2000).
We also calculated the codon bias index (CBI) to estimate the codon usage bias.When the Cynolebiatinae subfamily sequences were excluded, the CBI was 0.670 for the PG1 species, 0.736 for PG3 species and 0.846 for Cynolebiatinae subfamily species.Additionally, it is worth noting that the G+C3s value for the Austrolebias species was 0.913, whereas 0.703 was obtained for the other PG1 species.

Structural and functional analyses
To ascertain which homeodomain (HD) amino acids are conserved within species of this annual fish subfamily, we compared Hox protein sequences from species that grouped with them (Figure 3).The alignment of PG1 partial HD sequences revealed that, as compared with the human HoxB1 sequence, the sequences for the fish studied had the following substitutions: A. cheradophilus and A. viarius had  and Ser-37 substitutions; and A. luteoflammulatus had Phe-16, Ser-23 and Thr-27 substitutions (Figure 3a).Furthermore, the comparison of the primary structures of Austrolebias Hox1 sequences revealed some differences between them.While A. cheradophilus shared high sequence identity with A. gymnoventris, sequence similarity was considerably lower when compared with A. luteoflammulatus.Interestingly, these sequences shared some unexpected residues with other PG1 species.For example, the Thr-27 present in A. cheradophilus, A. gymnoventris and A. luteoflammulatus was present in almost all sequences analyzed, except in HsHoxB1 and LmHoxD1, whereas the Val-34 present in A. cheradophilus and A. gymnoventris was only present in DrHoxB1b sequence while Ser-37 was present in both HsHoxA1 and DrHoxA1a.On the other hand, Ser-23 present in the A. luteoflammulatus sequence only was also present in TrHoxB1a.However, it is worth noting that conserved and paralogue-characteristic amino acids were present in all the sequences we analyzed.The analysis of the PG3 HD sequences revealed that A. bellottii had Arg-24, Cys-27, Pro-29, Met-34, Asn-36, Leu-37, Asn-39, Thr-41 and Arg-43 substitutions compared to the human HoxB1 sequence.However, when compared to human HoxA3 we founded that the only difference was a Val-45 substitution (Figure 3b).In this sequence, conserved and paralogue-characteristic amino acids were also present.The HD 3D structure modeling of A. bellottii, A. cheradophilus, A. gymnoventris and A. luteoflammulatus Hox proteins showed that the ternary structure of the HD and the HD-DNA interaction had not been altered (Figure 4).

Discussion
Although Hox genes have been characterized in many species of teleost fish, this study is the first report about partial Hox gene isolation and characterization in annual killifishes belonging to the subfamily Cynolebiatinae.
We employed PCR amplification with specific primers to isolate partial homeodomain sequences from three Austrolebias species, A. bellottii, A. cheradophilus and A. viarius.These sequences and those previously isolated from A. gymnoventris and A. luteoflammulatus (Pereiro and García, unpublished data) represent the first homeobox fragments reported for Cynolebiatinae subfamily species.The isolated sequences show a high percentage of similarity with those of other vertebrates.Accordingly with the homeobox sequence conservation in virtually all metazoans (McGinnis et al., 1984a), our Southern blot analyses revealed that homeobox sequences are conserved not only among species of the Cynolebiatinae subfamily but also in the unrelated clupeid fish Brevoortia aurea (Figure 1).
Because of homeodomain amino acid conservation the phylogenetic tree topology (Figure 2) followed by comparative analyses (Figure 3) allowed us to predict the distribution of Austrolebias species partial Hox sequences within the anterior group of Hox genes.So far, we have identified partial homeodomain sequences of three Hox1 sequences and a Hox3 sequence in annual fish.Interestingly, even though homeodomain sequences are extraordinarily conserved throughout evolution, comparative analyses have shown unexpected amino acid shifts among Cynolebiatinae subfamily species.The alignment of the HD sequences from the PG1 group reveals that two of these residues reside in two different Austrolebias species,  2.
Val-36A in A. gymnoventris and Phe-16L in A. luteoflammulatus (Figure 3a).Note that Leu-16 is a conserved and paralogue-characteristic amino acid present in members of the PG1 (Sharkey et al., 1997).Similarly, analysis of the HD sequences from the PG3 group revealed that A. bellottii HD shared high sequence identity with members of clusters A and D, although this species presents a Val-45I residue (Figure 3b).Although some residues distinguish Hox paralogue groups from one another, it has been shown that the HD structures are remarkably similar and that the mode of HD-DNA interaction has been extraordinarily well conserved (Sharkey et al., 1997).Because of this and the fact that the homeodomain of the human HoxB1 protein and those of Cynolebiatinae subfamily species share a greater percentage of amino acid identity (Figure 3), we used the three-dimensional structure of HoxB1 as a template for homology modeling (Figure 4).We found no ternary structure changes even if primary sequences diverge.Nevertheless, in vivo assays in which the function of homologous 500 Hox genes of annual fish proteins is tested in the context of site-directed mutagenesis (Golding and Dean, 1998) may resolve the issue.
As we expected, Hox genes from Cynolebiatinae subfamily species showed the signature of adaptive amino acid replacements since non-synonymous substitutions greatly exceeded the synonymous substitution sites among distantly related taxa of the paralogue group 1.Although they were not statistically significant, the rates of both synonymous and non-synonymous nucleotide substitutions were higher in the Cynolebiatinae subfamily species than in the other PG1 species.Despite the fact that positive selection was acting, the substitution pattern detected among Cynolebiatinae subfamily and PG1 species could be indicating that for more distantly related species two or more independent mutations (multiple hits) might have occurred at the same site.Thus, in a molecular phylogenetic context, amino acid sites were identical not because of identity by descent from a common ancestor but because different types of mutation (e.g.parallel, convergent or reverse mutations) could arise at the same site.
Our analyses based on a restricted highly conserved region of the HD showed that both synonymous and nonsynonymous rates increased in many comparisons among Cynolebiatinae subfamily species and other PG1 Hox sequences.However, the highest codon usage bias value among them (CBI = 0.846) could be indicating that some synonymous substitutions were also constrained for these Hox genes.Further substitution pattern analysis among Cynolebiatinae Hox genes could be important in detecting the characteristic evolutionary signature of duplicated divergent evolution among the Hox genes, as has previously been found in the ray-finned fish HoxA cluster (Wagner et al., 2005).
Many molecular evolutionary questions remain open about the selective constraints acting on Hox genes of annual fish of the Cynolebiatinae subfamily since high levels of mitochondrial sequence divergence and high rates of amino acid replacements among coding regions have been previously detected (García et al., 2000(García et al., , 2002)).In this case, Hox genes were assigned to a series of physiological and life-history variables such as generation time, life span, age of first reproduction, rate of population increase and metabolic rate that are certainly not exclusively confined to mitochondria.
While the limited homeobox sequences described in the present paper allowed us to assign them to specific paralogue groups, the Austrolebias partial Hox sequences could not be assigned to specific Hox clusters..More sequence data, e.g., from the entire homeobox region and flanking regions, are needed to assign these genes to their correct clusters.Furthermore, it would be interesting to determine the number of Hox clusters in even one species of the Cynolebiatinae subfamily, since at least seven Hox clusters have been described for many teleosts (Amores et al., 1998(Amores et al., , 2004;;Aparicio et al., 2002;Naruse et al., 2004) and only four in the killifish Fundulus heteroclitus (Misof and Wagner 1996).Given that the development of annual killifishes shows reversible arrest at three well-defined stages, it would also be interesting to compare Hox gene expression between species of the subfamily Cynolebiatinae.

Figure 3
Figure 3 -a) Alignment of partial homeodomain amino acid sequences from species that clustered with the Austrolebias species (A.cheradophilus, A. gymnoventris, A. luteoflammulatus and A. viarius) b) Alignment of partial homeodomain amino acid sequences from species that clustered with Austrolebias bellottii.In both alignments, the homeodomain sequence of human HoxB1 (bold type) is shown at the top of the figure for reference.Identical amino acids are shaded and dots indicate amino-acid identity to human HoxB1 homeodomain.

Gutiérrez et al. 501 Figure 4
Figure 4 -Swiss-PdbViewer (ribbon diagram) representation of the homology-modeled homeodomain-DNA (HD-DNA) complex for (a) Austrolebias cheradophilus, (b) Austrolebias gymnoventris, (c) Austrolebias luteoflammulatus and (d) Austrolebias bellottii.Selected residues within the homeodomain that were identified in our alignments are indicated in each model.In all four species none of the altered residues was predicted to interact directly with the DNA.

Table 1 -
Capture sites of the Austrolebias (A) and Brevoortia (B) specimens analyzed by us.The A. adloffi specimens came from Brazil, all other specimens came from Uruguay.

Table 2 -
Abbreviations and GenBank accession numbers of species used in Figures2 and 3.