Skip to main content
Download PDF

Heterochrony in orthodenticle expression is associated with ommatidial size variation between Drosophila species

BMC Biology volume 23, Article number: 34 (2025) Cite this article

Abstract

Background

The compound eyes of insects exhibit extensive variation in ommatidia number and size, which affects how they see and underlies adaptations in their vision to different environments and lifestyles. However, very little is known about the genetic and developmental bases of differences in eye size. We previously showed that the larger eyes of Drosophila mauritiana compared to D. simulans are generally caused by differences in ommatidia size rather than number. Furthermore, we identified an X-linked chromosomal region in D. mauritiana that results in larger eyes when introgressed into D. simulans.

Results

Here, we used a combination of fine-scale mapping and gene expression analysis to further investigate positional candidate genes on the X chromosome. We found earlier expression of orthodenticle (otd) during ommatidial maturation in D. mauritiana than in D. simulans, and we show that this gene is required for the correct organisation and size of ommatidia in D. melanogaster. We discovered that the activity of an otd eye enhancer is consistent with the difference in the expression of this gene between species, with the D. mauritiana enhancer sequence driving earlier expression than that of D. simulans. When otd expression is driven prematurely during D. melanogaster eye development, the ommatidia grow larger, supporting a possible role for the timing of otd expression in regulating ommatidial size. We also identified potential direct targets of Otd that are differentially expressed between D. mauritiana and D. simulans during ommatidial maturation.

Conclusions

Taken together, our results suggest that differential timing of otd expression may contribute to natural variation in ommatidia size between D. mauritiana and D. simulans, which provides new insights into the mechanisms underlying the regulation and evolution of compound eye size in insects.

Background

Understanding the genetic basis of phenotypic diversity is one of the central themes of evolutionary developmental biology. While the causative genes and even mutations underlying evolutionary changes in a growing list of phenotypes have been identified (e. g. [1,2,3,4,5,6,7,8,9,10]) and see [11] for a more comprehensive list), we still know relatively little about the genetic basis for the evolution of organ size.

Insects exhibit remarkable variation in the size and shape of their compound eyes, which has allowed these animals to adapt to different environments and lifestyles [12, 13]. This variation greatly affects optical parameters and visual sensation, such as the detection of different intensities, polarisation and wavelengths of light to varying degrees of contrast sensitivity and acuity [13]. Compound eyes vary in the size and/or number of ommatidia: wider ommatidia capture more light, which can increase contrast sensitivity; however, larger interommatidial angles can lead to decreased acuity [13, 14]. Conversely, having many small ommatidia with narrow interommatidial angles can enhance acuity, but this may decrease contrast sensitivity [15,16,17,18].

Differences in ommatidia number and size, as well as trade-offs between these structural features of compound eyes, have been described for a range of different insects [19,20,21,22,23,24]. Furthermore, variation in ommatidia size across the eye within species is also widely documented [21, 25,26,27,28]. This size variation suggests areas of regional specialisation, where different visual tasks rely on different parts of the eye.

Several studies have also found extensive variation in eye size within and between closely related species of Drosophila, caused by differences in ommatidia number and/or ommatidia size [4, 19, 23, 25, 29,30,31,32,33,34]. Despite the pervasive variation in eye morphology and the detailed knowledge about eye development in D. melanogaster [12, 35,36,37], little is known about the genetic and developmental bases for variation in eye size even among Drosophila species with very few exceptions (e.g. [4, 30]). Interestingly, in these exceptional cases, changes in the timing of the expression of regulatory genes that consequently alter the timing of developmental processes, referred as heterochrony (i.e. differences in the timing or pace of developmental events compared to those in the ancestor [38, 39], have been involved. However, while these reported mechanisms explained eye size differences due to ommatidia number, whether changes in the timing of expression of key transcription factors (TFs) also may play a role in differences in ommatidia size, by modifying the gene regulatory network (GRN) and timing of developmental events controlled by these TF, remained unexplored.

We previously showed that D. mauritiana has larger eyes than D. simulans mainly due to larger ommatidia [23, 29]. Quantitative trait loci (QTL) mapping of this difference identified a large-effect QTL that explains 33% of the species difference [29]. Introgression of this X-linked region from D. mauritiana into D. simulans increased ommatidial size and overall eye size of the latter species [29].

Here, we combine high-resolution mapping of this previously characterised X-linked QTL, with transcriptomic analysis of eye-antennal imaginal discs (EADs) of D. simulans and D. mauritiana, to identify positional candidate genes that are differentially expressed in the developing ommatidia between these two species. We observed a temporal difference in the onset of expression of one of these candidates, the homeobox gene orthodenticle (otd), also known as ocelliless (oc), and we confirmed that Otd is involved in ommatidia organisation and size determination. We then carried out ATAC-seq to compare putative regulatory regions of otd that may underlie the difference in expression of this gene between D. mauritiana and D. simulans. Our results indicate that differential activity of an orthologous eye enhancer of otd between these species results in earlier expression of this homeobox gene during ommatidial maturation in D. mauritiana. We hypothesise that this heterochrony in otd expression and consequently longer exposure to this transcription factor (TF) in maturing ommatidia in D. mauritiana contributes to the development of larger ommatidia in this species.

Results

Enlarged ommatidia in D. mauritiana

We previously found that the central ommatidia of D. mauritiana eyes have wider diameters than those of D. simulans [23, 29]. To examine whether this phenotypic difference is prevalent in all ommatidia across the eye, we imaged the eyes of a female D. mauritiana TAM16 and a female D. simulans y, v, f using synchrotron radiation micro-CT (SRμCT) and measured the facet diameter of ommatidia in different regions of the eye using 3D reconstructions (Fig. 1). We corroborated that while the number of ommatidia is similar between these strains of D. mauritiana and D. simulans the former has larger facets, as observed for other strains of these two species [19]. This trend is consistent across anterior, central and posterior facets but is particularly pronounced in the antero-ventral region of the eye (Fig. 1 and Additional File 1: Table S1 [23, 29]) again as seen for other strains [19].

Fig. 1
figure 1

3D reconstruction and ommatidia size measurements from SRμCT data of female D. simulans (left) and D. mauritiana (right). Facet areas of the ommatidia highlighted in the antero-ventral (green), central (purple) and dorsal-posterior (blue) region of the eye are plotted in corresponding colours (far right). Ommatidia number is 996 for the D. simulans y, v, f and 1018 for the D. mauritiana TAM16. Scale bar is 100 μm

Full size image

Differentially expressed genes in a candidate region on the X chromosome

Previously we detected a QTL region located between 2.6 Mb and 8 Mb on the X chromosome, which is responsible for 33% of the difference in ommatidia size [29]. Furthermore, introgression of approximately 8.3 Mb of this X-linked region (between the yellow (y) and the vermillion (v) loci) from D. mauritiana TAM16 into D. simulans y, v, f significantly increased the eye size of the latter, consistent with the direction of the species difference [29]. Further analysis of recombinant males with breakpoints within the introgressed region revealed significant genotype–phenotype associations towards the distal end of the introgressed region near the marker v, providing a conservative interval of about 2 Mb wherein the causative loci is likely to reside (Fig. 2a). To map the candidate region to higher resolution we generated introgression lines with breakpoints in the 2 Mb interval and compared eye area and central ommatidia diameter of y, f male progeny (with some D. mauritiana DNA in the 2 Mb interval) to that of their y, v, f sibling males (i.e. without D. mauritiana DNA). We found that y, f males had significantly larger eye size than their y, v, f siblings in introgression lines IL9.1a (one tailed t = 2.18, df = 11, p = 0.026), IL9.1b (one tailed t = 4.54, df = 11, p < 0.001) and IL9.2 (one tailed t = 2.00, df = 11, p = 0.035) but ommatidia diameter was only significantly larger in y, f males than their y, v, f siblings in introgression lines IL9.1a (one tailed t = 2.52, df = 11, p = 0.014) and b (one tailed t = 3.12, df = 11, p = 0.005) (Fig. 2a). Ommatidia number and body size did not differ between y, f males and their respective y, v, f sibling males for any of the IL lines (Additional File 2: Table S2 and Additional File 3: Fig. S1). These data suggest that the candidate QTL is located in a maximum region of about 662 kb (ChrX: 7,725,195–8,387,618 in D. simulans). This mapped region contains 62 protein coding positional candidate genes.

Fig. 2
figure 2

Differential and spatial gene expression. a Fine-scale mapping of the X chromosome QTL. Marker-phenotype association in male recombinant progeny (between y and v) from three replicate introgression lines (IL1, 3, and 4, single-marker ANOVA analysis). Red dashed line indicates the Bonferroni corrected significance threshold of 0.05. Shaded grey area represents a conservative interval of ~ 2 Mb encompassing the X linked QTL. Recombination breakpoints of the new introgression lines (IL9.1–9.3) on the X chromosome (shown for D. simulans Flybase R2.02) define the 662 kb candidate region. White, black and grey boxes indicate DNA regions from D. simulans y, v, f, D. mauritiana TAM16 or not determined, respectively (the latter define the maximum candidate region). The table indicates the number of protein coding genes that are present in the candidate region in D. simulans and D. melanogaster. Distribution of eye area (left) and ommatidia diameter (right) measurements by genotype and introgression line. Asterisks indicate levels of significance between genotypes where p < 0.05 (Additional File 2: Table S2). b Differential expression of 49 protein coding genes located in the introgressed region from a and expressed at 72 h AEL, 96 h AEL and 120 h AEL. Genes with significantly higher expression at 120 h AEL in D. simulans are highlighted in blue. Genes significantly upregulated in D. mauritiana at 120 h AEL are shown in red. c Expression of differentially expressed genes at 120 h AEL in L3 EADs of D. simulans and D. mauritiana. The cartoon illustrates an eye-antenna disc, highlighting the eye territory (in dark blue) where the ommatidia will be determined during late L3 and pupal stages. Open arrowheads indicate the MF, a: antenna, e: eye, A: anterior, P: posterior, D: dorsal, V: ventral, in c

Full size image

Drosophila eyes develop from the eye-antennal imaginal discs (EAD), a monolayered epithelium [40]. Within it, the activity of two morphogens, Hedgehog (Hh) and the BMP2 Decapentaplegic (Dpp), drives a differentiation wave that progresses from the posterior to the anterior region of the primordium and leaves in its wake differentiating retinal cells. These cells then organise into clusters that will form individual ommatidia [35]. The wave-like differentiation implies a temporal axis, in which ommatidial rows are progressively added in the anterior front of the wave. The front of the differentiation wave is characterised by a transient indentation of the epithelium, the so-called morphogenetic furrow (MF). Therefore, the MF is a landmark behind which ommatidial differentiation takes place. To assay which of these candidates are expressed during the generation of ommatidia, we performed RNA-seq experiments on the EADs of 3rd instar larvae (L3). We extracted RNA from D. mauritiana and D. simulans EADs at three different developmental points: at 72 h after egg laying (AEL; late L2, at the onset of differentiation marked by the MF), when cells in the eye primordium are proliferating and specification of the ommatidial cells has not yet started; at 96 h AEL stage (L3) when the MF has moved about half way across the eye disc and the most posterior ommatidia are already determined, and at 120 h AEL (late L3), when most ommatidia are already determined but their size, structure and shape are not yet finalised [41, 42].

Comparison of the RNA-seq data among these three developmental timepoints showed that transcriptomes of 72 h AEL EADs were the most different in comparison to transcriptomes from both 96 h AEL and 120 h AEL for both species (Additional File 4: Fig. S2 and Additional File 5: Table S3). This reflects the distinctive processes that are occurring at these developmental stages [42]. We next focused on the expression of genes located within the mapped 0.66 Mb X-linked region at 120 h AEL because at this timepoint the posterior ommatidia begin to adopt their final size. Of the 62 genes located in this region, 49 were expressed at least at one of the RNA-seq timepoints and only eight of these genes were differentially expressed between these two species at 120 h AEL (Additional File 5: Table S3): spirit, otd and Ppt1 showed higher expression in D. mauritiana, whereas CG1632, Es2, Sptr, CG12112 and CG1885 were more highly expressed in D. simulans (Fig. 2b).

We next performed in situ hybridization experiments of these eight candidate genes to investigate if they are expressed in the eye field where the ommatidia are being assembled. These assays were carried out in both D. mauritiana and D. simulans, which allowed us to determine whether the differences in expression levels observed in the RNA-seq datasets are related to differences in spatial expression (Fig. 2c). Sptr, CG12112 and spirit had no detectable expression in the relevant region posterior to the MF (Fig. 2c). Ppt1 and CG1885 were expressed both anterior to and immediately posterior to the MF. CG1632 and Es2 were ubiquitously expressed in the eye disc, with no clear regional differences. Finally, as previously shown in D. melanogaster, otd was expressed in the ocellar region of the EAD and in the most posterior region of the eye field [43]. At 120 h AEL otd is already expressed in several rows of the most posterior ommatidia of D. mauritiana eye discs, whereas otd expression is weakly detected in a smaller posterior region of the eye discs of D. simulans (Fig. 2c, Additional File 6: Fig. S3, Additional File 7: Table S4). These results were consistent with our differential expression analysis, as the analyses of spatial expression showed qualitative differences in expression levels for most of the investigated genes. Taken together, these results showed that otd is the only differentially expressed positional candidate gene that is expressed in maturing ommatidia (Fig. 2c).

Differences in otd gene expression during eye development between D. simulans and D. mauritiana

Our results suggested that otd transcription in the maturing ommatidia initiates earlier in D. mauritiana than in D. simulans (Fig. 2c). To investigate this further, we performed additional in situ hybridizations at 110 h AEL to compare the onset of otd expression in the developing ommatidia of these two species. At this developmental stage, we found that otd is already transcribed in D. mauritiana eye discs, whereas there was no detectable expression in D. simulans discs (Fig. 3a, b, Additional File 8: Fig. S4). To confirm this heterochrony in otd expression, we performed immunostainings against Otd protein in developing eye discs (Fig. 3c, d). As an exact staging based on hours AEL may not recapitulate subtle interspecific differences in developmental timing, we counted the number of ommatidial rows that were already specified (i.e. with positive Elav staining) as a proxy of developmental stage and then which of these ommatidial rows showed Otd expression. We observed that D. mauritiana eye discs displayed more Otd-positive ommatidia than D. simulans eye discs at the same stage (Elav positive ommatidia rows, Fig. 3e, Additional File 9: Table S5, F 1,47 = 30.3, p-value = 1.48 × 10−06). Thus, cells in maturing ommatidia are exposed to the activity of Otd for longer in D. mauritiana since the expression of this protein extends into the pupal stage of both species (Additional File 10: Fig. S5).

Fig. 3
figure 3

otd expression in L3 eye imaginal discs. (a-d) otd mRNA at 110 h AEL in D. mauritiana (a) and D. simulans (b). Black arrowheads indicate the MF. Asterisks indicate expression in the ocellar region. D. mauritiana already has detectable otd mRNA at 110 h (red arrowhead). c–d’ Immunostaining showing Otd protein (magenta, c’, d’) in mature ommatidia (marked in green by Elav) and the ocellar region (asterisk) in D. mauritiana (c–c’) and D. simulans (d–d’). Staining against actin (blue) was used to mark the MF. e Plot showing the number of Otd-positive ommatidia rows (y-axis) at different developmental time points (x-axis, developmental points inferred by number of rows of ommatidia). Asterisks represent statistical significance p < 0.001 of the Analysis of Covariance (ANCOVA) (F 1,47 = 30.3, p-value = 1.48 × 10−06). Scale bars: 50 μm. Anterior regions of the EADs are on the left and posterior on the right

Full size image

Otd is required for the correct arrangement and size of ommatidia in Drosophila

It was previously shown in D. melanogaster that otd is expressed in photoreceptor cells and required during pupal stages for morphogenesis of their rhabdomeres, and subsequently rhodopsin expression, as well as for the synaptic-column and layer targeting of the photoreceptors [43,44,45]. We carried out further analysis of otd function during eye development using RNAi knockdown and by generating mitotic clones of homozygous otd mutant cells. Decreasing otd mRNA by overexpressing an otd-miRNA construct in cells posterior to the MF [46] resulted in defects in the final ommatidial organisation. These defects were rescued by adding a copy of UAS-otd (Additional File 11: Fig. S6a-c), supporting the specificity of the RNAi-mediated otd attenuation. Loss of otd in clones resulted in disorganised ommatidia with perturbed shapes and sizes—often smaller than the ommatidia of controls (Additional File 11: Fig. S6d).

To test whether an acceleration of the onset of otd expression could cause an increase in ommatidial size, we drove a UAS-otd transgene with the ato3’FL-GAL4 line, which drives the expression of otd in the MF and early photoreceptors [47], and therefore earlier than the onset of endogenous otd expression. Using one-way ANOVA, we observed that there are significative differences between the experimental genotype (“ato > otd”) and the two controls (“otd; + ”, and the parental line “ato > + ”; F7,95 p-value = 0.0015). We also performed a Tukey's HSD to identify the differences between the three groups and concluded that the ommatidium area was significatively larger than each of the controls (p = 0.0154 and p = 0.0052, respectively) when we overexpressed otd, while the area of controls did not differ to each other (p = 0.9004) (Additional File 11: Fig. S6, Additional File 12: Table S6). Overall, these results show that otd expression in the photoreceptor cells of maturing ommatidia is required for the proper regulation of ommatidial organisation and size, with earlier expression resulting in larger ommatidia.

Differences in chromatin accessibility in the otd locus during eye development between D. simulans and D. mauritiana

Our mapping and expression analyses indicate that the differences in otd expression likely contribute to differences in ommatidia size between D. simulans and D. mauritiana. Given that there is no amino acid difference in the Otd homeodomain between our focal strains of D. simulans and D. mauritiana (Additional File 13: Fig. S7), our data suggest that the causative changes reside in otd regulatory regions. Due to the microsyntenic conservation between D. melanogaster, D. simulans and D. mauritiana, we considered the regulatory landscape of otd as the region between its two flanking genes, Caf1-180 and CG12772, revealed by the presence of a topological associated domain (TAD) in the corresponding D. melanogaster region (a region of 69 kb in D. mauritiana and 70 kb D. simulans, Fig. 4a, Additional File 14: Fig. S8).

Fig. 4
figure 4

Chromatin accessibility at the otd locus. a Chromatin accessibility profiles at the otd locus at 96 h AEL and 120 h AEL in D. mauritiana and D. simulans EADs. Purple and orange boxes indicate APREs detected in D. mauritiana and D. simulans, respectively, when mapped against their own reference genome. Green boxes highlight APREs with differential accessibility, identified by re-mapping and peak calling of D. mauritiana datasets against D. simulans reference genome and D. simulans datasets against D. mauritiana reference genome, using the "liftOver" tool and MACS2. Red boxes represent the region corresponding to the D. melanogaster otduvi allele region. b D. simulans and D. mauritiana otd-APRE7-8 enhancer activity in 3rd instar larvae eye imaginal discs. Red bar: limit between MF and differentiated photoreceptors/eye. Scale bars: 50um. Grey: DAPI, red: Elav, green: GFP. c Plot showing rows with full enhancer activity (y-axis) activated by otd-APRE7-8 at different developmental time points (x-axis, developmental points inferred by number of ommatidia rows) for both species. d Violin plots showing the distance between the MF and the rows with full otd-APRE7-8 activity (first row counting from the posterior margin of the disc with more than 90% of the ommatidia positive for GFP signal). The D. mauritiana otd-APRE7-8 is active earlier during the differentiation of ommatidia

Full size image

To investigate the regulation of otd in the developing eyes of D. simulans and D. mauritiana further, we performed ATAC-seq [30, 48] on D. simulans and D. mauritiana EADs at 96 and 120 h AEL. We mapped our datasets against their reference genomes and also against the other species genome (see the “Methods” section) to detect common, differentially accessible and species-specific regulatory regions (Fig. 4a). The ATAC-seq peak calling of the four datasets (two developmental stages and two species) revealed a total of 22 APREs (Associated Putative Regulatory Regions) in the otd locus, all of which were located within alignable orthologous regions in the two species (Fig. 4a).

Five of these peaks showed significant differences in accessibility between D. mauritiana and D. simulans at 96 hAEL and/or 120 hAEL: APRE 6 (D. sim chrX: 8,100,587- 8,100,808, padj = 0.00155; D. mau chrX: 7,929,372–7,929,595, padj = 0.000877), APRE 7 (D. mau chrX: 7,929,910–7,930,230, padj = 1,85E-06) and APRE 11 (D. sim chrX: 8,107,881–8,108,402, padj = 0.00976; D. mau chrX: 7,936,681–7,937,179, padj = 0.0143) in the 3rd and 1st introns of otd, respectively, and APREs 17 (D. sim chrX: 8,125,765—8,126,100, padj = 0.0418; D. mau chrX: 7,954,305–7,954,661, padj = 0.0331) and 19 (D. sim chrX: 8,129,876—8,131,170, padj = 0.0317) located upstream of otd (Fig. 4a, Additional File 15: Table S7).

To test the activity of these APREs in L3 EADs, we generated D. melanogaster lines containing the APRE 6, 7–8, 11, 16–17, and 19 sequences from D. mauritiana followed by GAL4 and crossed them to UAS-GFP (see the “Methods” section, Fig. 4b and Additional File 16: Fig. S9). APREs 6 and 11 showed no activity in the EADs, which we confirmed using reporters with the equivalent D. simulans sequences (Additional File 16: Fig. S9). APREs 16–17 and 19 drove GFP expression anterior to the MF, in the ocellar domain and some regions of the antennal disc (Additional File 16: Fig. S9, Additional File 17: Table S8). Finally, APRE7-8 drove expression in the posterior of the eye disc, in maturing ommatidia, in a similar pattern to endogenous otd expression, and consistent with APRE 7 and APRE 8 demarcating the 2.5 kb region corresponding to the D. melanogaster otduvi allele previously identified as an otd eye enhancer using reporter assays [43] (Fig. 4a and Additional File 18: Fig. S10). We therefore compared the activity of the D. mauritiana and D. simulans APRE 7–8 sequences, and we observed differences in the onset of activity of these reporter lines. According to mRNA and protein Otd expression, APRE7-8-GAL4 from D. mauritiana activated GFP expression earlier than the reporter with the D. simulans sequence (Fig. 4b–d) consistent with the difference in endogenous otd expression between these species. These results suggest that differences in the D. mauritiana and D. simulans APRE 7–8 sequences underlie the difference in otd eye expression between these two species. Note that we also tested the activity of the D. melanogaster APRE 7–8 sequence and observed it was similar to the D. simulans sequence rather than that of D. mauritiana consistent with the derived larger ommatidia of the latter species (Additional File 18: Fig. S10).

We aligned the orthologous sequences of the 2.5 kb APRE 7–8 region from different D. mauritiana, D. simulans and D. melanogaster strains and found seventeen potentially fixed differences (thirteen SNPs and four short indels) in D. mauritiana (Additional File 19: Dataset S1, Additional File 20: Fig. S11). Interestingly, most of these mutations lie in a 1.5 kb sub-fragment of this enhancer identified by Vandendries et al. (1996) and several fall in predicted binding sites for Sine oculis, Cut and Otd itself (as well as other transcription factors) (Additional File 20: Fig. S11).

Differences in Otd targets during eye development between D. simulans and D. mauritiana

Next, we investigated whether differences in the onset of expression of otd between D. simulans and D. mauritiana promoted further changes in the expression of downstream genes. To this end, we searched for the Otd-binding motif in accessible chromatin regions of genes expressed during eye development. Based on this analysis, we found 1148 putative Otd target genes. We next examined which of these accessible chromatin regions were associated with genes that were differentially expressed in our transcriptome datasets. We found that 161 APREs of the 1330 genes that are upregulated in D. mauritiana contained Otd binding motifs, and 111 out of 1249 genes upregulated in D. simulans had associated peaks containing Otd binding motifs (Fig. 5a, c, Additional File 21: Table S9).

Fig. 5
figure 5

Otd downstream targets. a 161 genes upregulated in D. mauritiana have an associated peak that contains at least one Otd motif. b GO enrichment for those genes with Otd motifs in D. mauritiana (c) 111 genes upregulated in D. simulans have an associated peak that contains Otd motif. d GO enrichment for those genes with Otd motifs in D. simulans

Full size image

We then performed Gene Ontology (GO) term enrichment analysis for these differentially expressed genes with accessible chromatin containing Otd binding motifs. The D. mauritiana dataset exhibited enrichment in terms related to developmental processes, cell morphogenesis, axonogenesis, regulation of differentiation or growth, among others (Fig. 5b, Additional File 21: Table S9). By contrast, genes that were upregulated in D. simulans with associated Otd-peaks were enriched in terms such as the MAP kinase network, signalling by interleukins and cellular response to insulin stimulus (Fig. 5d, Additional File 21: Table S9).

Discussion

While much is known about the specification and differentiation of ommatidia, very little is known about the regulation and evolution of their size (although see [49]). To investigate the genetic basis of the difference in ommatidia size between D. mauritiana and D. simulans, we carried out high-resolution introgression mapping of a previously identified X-linked QTL that explains about 33% of the difference in eye size between these two species [29]. In this region, we identified eight positional candidate genes whose expression in the developing EADs differed between D. mauritiana and D. simulans. Our analysis of the spatial expression of these eight genes showed expression in the region of the MF (CG12112), as well as genes with broad expression in the presumptive retinal region (CG1632, Es2), which could contribute to eye size differences between species. However, we decided to focus on factors expressed in maturing photoreceptors in the posterior portion of the disc. Based on these criteria otd was the best candidate gene in this region for contributing to the difference in ommatidia size and thus overall eye size between these species.

otd/Otx genes play several important roles during eye development in both invertebrates and vertebrates [43, 50,51,52,53]. During Drosophila eye development, Otd regulates genes for cell adhesion and cytoskeletal organisation, which is essential for the correct development of the photoreceptor cells and ommatidia maturation as well as subtype specification through regulation of rhodopsin expression [44, 51, 53, 54]. Mutations in otd perturb morphogenesis of the photoreceptor cells [43, 44]. Intriguingly, the removal of photoreceptor cells changes ommatidia size [55]. We suggest that although otd is not expressed in the lens-secreting cone cells, it indirectly affects the organisation of these cells and thus ommatidia size through regulating the maturation and organisation of the underlying photoreceptor cells. We have shown that knockdown or loss of otd in D. melanogaster perturbs ommatidia size specification, but it remains to be directly tested if variation in the expression of this gene underlies larger and smaller ommatidia in D. mauritiana and D. simulans respectively and if otd contributes to the observed variation in ommatidia size in different regions of Drosophila eyes [19, 25]. However, the fact that early onset of otd expression in D. melanogaster using the GAL4/UAS system results in larger ommatidia strongly supports this hypothesis.

Changes in developmental timing, or heterochrony, have played an essential role in the evolution of morphologies in multiple taxa [38, 56, 57]. Classically, the term heterochrony has been used to refer to differences in the timing of developmental events and several examples of heterochrony have been described [58,59,60,61]. Most of these characterised cases showed that the genetic basis lies upstream of the mechanism responsible for the heterochrony, such as changes in proliferating rates, differences in the initial size of the primordium or distinct rates of protein stability and biochemistry [62,63,64,65]. Heterochronic shifts can also occur as a direct consequence of the causative genetic change, such as those that affect regulatory regions altering the timing of gene expression [4]. Although differences in gene expression of single TFs have the potential to completely modify the subsequent GRN, the relative contribution of such direct heterochrony in generating morphological diversity remains unknown. Our data indicate that otd is actually expressed earlier during ommatidial maturation in D. mauritiana compared to D. simulans, and that premature expression of otd results in larger ommatidia, when assayed in D. melanogaster. This suggests that cis-regulatory changes in otd lead to ommatidial cells being exposed to Otd for longer or earlier in D. mauritiana. Globally considered, our results support the notion that temporal variations in the onset of otd expression, deriving from the heterochronic action of species-specific enhancers, may contribute to the evolution of ommatidial size. The particular mechanism by which earlier (or longer) expression of otd leads to larger lenses (and presumably also larger ommatidia in general) is not well understood, although our results suggest that some rewiring of the Otd-controlled GRN may have happened in the intervening evolutionary period since the common ancestor of D. mauritiana and D. simulans diverged. Therefore, the evolutionary difference that we attribute to a heterochronic shift in otd expression could be coupled to changes in the wider GRN. Together with Ramaekers and colleagues (2019), our study shows how morphological diversity in closely related species may be achieved by subtly altering the temporal expression of a single TF. Importantly, in both cases, these transcription factors, Ey and Otd, act upstream in the GRN controlling the developmental process, thus changes in their expression may promote major differences in downstream effectors.

Further exploration and comparison of the regulatory landscape of otd between D. mauritiana than D. simulans allowed us to identify an intronic region in the otd locus, APRE 7–8, corresponding to the orthologous sequence of the previously known “eye enhancer” from D. melanogaster [43, 66, 67]. The D. mauritiana sequence reporter line showed earlier activity than D. simulans sequence, suggesting that this region is responsible for the differences in the timing of otd expression and, therefore, function. Further work is needed to identify the TFs that directly bind to this enhancer and the changes in D. mauritiana sequence that underlie the differential activity of this element in this species (Additional File 16: Fig. S9).

We also investigated how these changes in otd expression might alter target gene expression to change ommatidia size. We identified a set of genes that are differentially expressed between these two species when the ommatidia are acquiring their final size that may be acting downstream of Otd, as they have accessible chromatin regions containing putative Otd binding motifs. We compared this set of genes to known and putative targets of Otd which have been characterised later in eye development during pupal stages [44, 51]. This comparison showed that a subset of genes with altered expression in otd mutants are also differentially expressed between D. mauritiana and D. simulans in late L3. Among other TFs (Dve, vnd, MED10) we identified several genes involved in phototransduction (e.g. slo, Slob, ninaG, inaD, ninaA) and genes encoding cytoskeleton and adhesion proteins (Act88F) (Additional File 21: Table S9). In addition, we found several genes involved in regulation of growth (e.g. DAAM, Thor, Frizzled 2, kibra or Ankyrin 2) which could have an impact in ommatidia size. This further suggests that the network downstream of Otd varies between these two species and that ultimately, these changes in the GRN, promoted by an early expression of otd in D. mauritiana, result in differences in ommatidia size between D. mauritiana and D. simulans.

Interestingly, we recently showed that D. mauritiana has higher contrast sensitivity than D. simulans [19], while the latter species has greater spatial acuity consistent with the differences in ommatidia size between these species. The trade-off between contrast sensitivity and acuity is heavily influenced by various aspects of visual ecology, such as habitat type, circadian activity patterns and lifestyle. Thus, substantial functional consequences with strong ecological implications could be linked to changes in the expression of individual genes such as otd.

Conclusions

Our data suggest that changes in the timing of otd expression underlie differences in ommatidia size and thus overall eye size between D. mauritiana and D. simulans. Our work provides new insights into ommatidia size regulation and the evolution of eye size. Together with evidence from other studies showing that changes in the timing of ey expression contributes to differences in ommatidia number in Drosophila [4], we have provided new insights into the genetic and developmental mechanisms that underlie the large diversity in Drosophila eye size [4, 19, 23, 25, 29,30,31,32,33,34]. Moreover, this evidence suggests that changes in the temporal expression of upstream transcription factors could be a widespread mechanism for morphological evolution.

Methods

Fly stocks and clonal analysis

D. simulans strain yellow (y), vermillion (v), forked (f) was obtained from the Drosophila Species Stock Center, San Diego, California (Stock no.14021–0251.146). D. mauritiana TAM16 is a wild-type inbred strain [29]. UAS-miR-otd and UAS-otd (III) were kindly provided by Henry Sun [47]. GMR-GAL4 [68] was used to drive expression of the transgenes. To generate mitotic clones of mutant otd in developing eyes, we used otd[YH13], neoFRT19A/FM7c and RFP, neoFRT19A; ey-Flp, which were obtained from Bloomington Stock Centre (#8675 and #67173 respectively). The ato3’FL-GAL4 is described in [47] and drives expression in the MF and early photoreceptors (Additional File 11: Fig. S6). The UAS-otd line is the stock with Flybase number FBst0005541 (P(UAS-oc.F)RF1, w1118). The UAS-GFP line was like in [69]. To generate UAS-otd;; ato3’FL-GAL4 flies (“ato > otd”) we crossed UAS-otd females to ato3’FL-GAL4/TM6B males. As controls, we used UAS-otd;;TM6B siblings (“otd, + ”) and the parental ato3’FL-GAL4/TM6B strain (“ato > + ”). To determine the area of ommatidia, we used females. Heads were mounted and photographed as in [70]. We measured the area of three strips of 3X5 facets using the “polygonal tool” of ImageJ from the anterior-ventral region (adjacent to the second and third antennal segments) of the two eyes and then the values were averaged. The average area/15 was the facet area of each individual. This was carried out for the experimental genotype (UAS-otd;; ato3’FL-GAL4 flies) and the two controls (UAS-otd;;TM6B and ato3’FL-GAL4/TM6B). We performed a one-way ANOVA to test differences between the three different samples. We also carried out a Tukey's HSD to identify which specific groups differed significantly from each other.

otd mutant clones were induced in developing eyes using the Flp/FRT system [71]. Female flies of the genotype otd[YH13], neoFRT19A/FM7c were crossed with males of the genotype RFP, neoFRT19A; ey-Flp. Female F1 progeny were examined for the lack of the FM7c balancer and these flies were prepared for SEM analysis.

Synchrotron radiation microtomography

Fly heads were removed from the body and placed into fixative (2% PFA, 2.5% GA in 0.1 M sodium cacodylate buffer overnight at 4 °C. Heads were washed in water, then placed into 1% osmium tetroxide for 48 h at 4 °C, then washed and dehydrated in increasing concentrations of ethanol up to 100%. Heads were then infiltrated with increasing 812 Epon resin concentrations up to 100% over 5 days and polymerised in embedding moulds for 24 h at 70 °C.

Heads were scanned at the TOMCAT beamline of the Swiss Light Source (Paul Scherrer Institute, Switzerland [72]). Scans were performed using a 16 keV monochromatic beam with a 20 μm LuAG:Ce scintillator. Resin blocks were trimmed and mounted using soft wax and scanned using × 20 combined magnification (effective pixel size 325 nm) and a propagation distance of 25 mm. Two thousand projections were taken as the heads rotated through 180°, each with 200 ms exposure. Projections were reconstructed into 8-bit tiff stacks and Paganin filtered (delta = 1−8, beta = 2−9 [73] using custom in-house software [74]. Tiff stacks were segmented in Amira (v2019.2, Thermo Fisher) for measurements of facet area.

SEM microscopy

Fly heads were fixed in Bouin’s for 2 h. Then, 1/3 of total volume was replaced by 100% ethanol to fully immerse heads in Bouin’s and were left to fix overnight. Heads were washed and dehydrated 2 × 70% ethanol overnight, 2 × in 100% ethanol and finally critical point dried and mounted onto sticky carbon tabs on SEM stubs, gold coated and imaged in a Hitachi S-3400N SEM with secondary electrons at 5 kV.

Markers and introgression lines

Previously, we generated three replicate introgression lines (IL1, IL3 and IL4) by introgressing the region between y (X:178,985) and v (X:10,264,389, Supplementary Table S2b) from D. mauritiana TAM16 into D. simulans y, v, f [29]. Males recombinant within the introgressed region (with phenotypes: y, f or v, f) from the 3 replicate lines were collected at backross 7. These individuals were genotyped with eleven new additional markers (Additional File 2: Table S2).

Significant association between each marker and eye size was tested (F-test, type III sum of squares SS) by performing a single-marker ANOVA on the residuals of eye area regressed onto T1 tibia length for each replicate (introgression line (IL 1,3 and 4; n = 20–60, Additional File 2: Table S2). Multiple testing was corrected using Bonferroni correction. All ANOVA models were fitted in the R statistical environment (R Development Core Team 2012) using the CAR package [75].

To narrow down the 2 Mb region, the X chromosome region between y and v from D. mauritiana TAM16 was re-introgressed into D. simulans y, v, f (as in [29]) and y, f females were backcrossed from multiple replicate lines to y, v, f males for a further nine generations. At the end of the egg-laying cycle of that generation, we collected mothers and genotyped them for molecular markers located in the 2 Mb region (Additional File 2: Table S2). Four mothers with breakpoints within this region were identified. Two of them were siblings (IL9.1a and IL9.1b) and they had the same 4th great-grandmother as IL9.3 and the same 7th great-grandmother as IL9.2. Male progeny available for each of these females was collected and genotyped and phenotyped for eye area, ommatidia diameter, ommatidia number and T1 tibia length as described previously in [23] (Additional File 2: Table S2). To determine if the D. mauritiana DNA in the 2 Mb region resulted in larger eyes and larger ommatidia, y, f males (i.e. with some D. mauritiana DNA in the 2 Mb interval) were compared to that of their y, v, f sibling males (i.e. without D. mauritiana DNA) for each introgression line using one-tailed, two-sample, equal-variance t-tests.

In situ hybridisation and immunohistochemistry

In situ hybridizations were carried out using a standard protocol with DIG-labelled antisense RNA probes. EADs were dissected and fixed at 120 h AEL for 40 min in 4% formaldehyde. To be able to compare the expression patterns avoiding technical differences (i.e. probe affinity and probe concentration), we first aligned the sequences from D. mauritiana and D. simulans and designed primers to generate RNA probes within fragments with at least 95% of similarity between them (Additional File 22: Table S10). This design allowed us to perform the in situ hybridization experiments using the same specific probes for each of the candidate genes at the same concentration for both species. The nitro blue tetrazolium/5-bromo-4-chloro-3′-indolyphosphate (NBT-BCIP) reaction was stopped at the same time. Candidate gene sequences were cloned into a TOPO PCR4 (spirit, otd, Ppt1, CG1632, Es2 and CG12112) or pCRII (CG1885, Sptr) vectors (Invitrogen) using specific primer pairs (Additional File 22: Table S10), respectively, following the manufacturer’s protocol. M13 forward and reverse primers were used to linearise the DNA. According to the vector and orientation of the fragments T3, T7 or SP6 RNA polymerase were used to generate the DIG-labelled riboprobes.

Immunostainings with Rabbit anti-Otd [44] and rat anti-Elav (7E8A10, Hybridoma bank) were performed at 1:1500 and 1:100 dilutions respectively using standard protocols, followed by anti-rat-Cy3 (Jackson Immuno Research) and anti-rabbit-Alexafluor 647 (Molecular probes) secondary AB staining, at 1:200. The actin cytoskeleton was stained with Alexafluor 488-Phalloidin (Molecular Probes) at 1:40 dilution for 30 min after secondary antibody incubation. Discs were mounted in Prolong Gold antifade reagent, supplemented with DAPI (Molecular Probes), and captured with a Zeiss LSM 510 confocal microscope. Images were processed using NIH ImageJ software. We used an Analysis of Covariance (ANCOVA) to test for differences between strains for Otd-positive ommatidia while adjusting for differences in development stage by using the number of ommatidial rows as a proxy for the latter. The ANCOVA was performed using base R v4.0.2 (R Core Team, 2020).

RNA-seq

Flies were raised at 25 °C with a 12 h:12 h dark:light cycle and their eggs were collected in 2 h time periods. Freshly hatched L1 larvae were transferred into fresh vials in density-controlled conditions (30 freshly hatched L1 larvae per vial). EADs were dissected at three different developmental time points: 72, 96 and 120 h AEL and stored in RNALater (Qiagen, Venlo, Netherlands). Three biological replicates for each sample were generated. Total RNA was isolated using RNeasy Mini Kit (Qiagen). RNA quality was determined using the Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA) microfluidic electrophoresis.

Library preparation for RNA-seq was performed using the TruSeq RNA Sample Preparation Kit (Illumina, catalog ID RS-122–2002) starting from 500 ng of total RNA. Accurate quantitation of cDNA libraries was performed using the QuantiFluor™dsDNA System (Promega, Madison, Wisconsin, USA). The size range of final cDNA libraries was determined by applying the DNA 1000 chip on the Bioanalyzer 2100 from Agilent (280 bp). cDNA libraries were amplified and sequenced using cBot and HiSeq 2000 (Illumina): only 120 h EAD samples were sequenced as paired-end (PE) reads (2 × 100 bp), all the other samples were sequenced in single-end (SE) reads (1 × 50 bp). Sequence images were transformed to bcl files using the software BaseCaller (Illumina). The bcl files were demultiplexed to fastq files with CASAVA (version 1.8.2).

Quality control analysis using FastQC software (version 0.10.1, Babraham Bioinformatics) was performed. All RNAseq reads are accessible in the Short Read Archive through umbrella BioProject PRJNA666691 (containing PRJNA374838 and PRJNA666524). To allow using the same mapping parameters for all samples, PE 100 bp reads were converted into SE 50 bp by splitting the reads in half and merging right and left reads into a single file prior to read mapping.

The reciprocally re-annotated references described in [41] were used to map the species-specific reads. Bowtie2 [76] was used to map the reads to each reference (–very-sensitive-local –N 1), and the idxstats command from SAMtools v0.1.19 [77] was used to summarise the number of mapped reads. HTSFilter [78] was used with default parameters to filter out genes with very low expression in all samples. For the remaining genes in each pair-wise comparison, differential expression was calculated using DESeq2 v1.2.7. with default parameters [79].

ATAC-seq library preparation and sequencing

Samples were obtained following the same procedure as for the RNA-seq experiments: flies were raised at 25 °C with a 12 h:12 h dark:light cycle. Freshly hatched L1 larvae were transferred into vials with density-controlled conditions. EADs were dissected at 96 and 120 h AEL and maintained in ice cold PBS. Imaginal disc cells were lysed in 50 μl Lysis Buffer (10 mM Tris–HCl, pH = 7.5; 10 mM NaCl; 3 mM MgCl2; 0.1% IGEPAL). Nuclei were collected by centrifugation at 500 g for 5 min. 75,000 nuclei were suspended in 50 μl Tagmentation Mix [25 μl Buffer (20 mM Tris- CH3COO−, pH = 7.6; 10 mM MgCl2; 20% Dimethylformamide); 2.5 μl Tn5 Transposase; 22.5 μl H2O] and incubated at 37 °C for 30 min. After addition of 3 μl 2 M NaAC, pH = 5.2 DNA was purified using a QIAGEN MinElute Kit. PCR amplification for library preparation was done for 14 cycles with NEBNext High Fidelity Kit; primers were used according to [48]. Paired end 50 bp sequencing was carried out by the Transcriptome and Genome Analysis Laboratory Goettingen, Germany.

ATAC-seq peak calling and differential binding site analysis

ATAC-seq raw reads were generated from the following samples (2 replicates each): D. simulans larvae at 96 and 120 hAEL and D. mauritiana larvae at 96 and 120 hAEL. These reads were mapped to strain-specific genomes of D. mauritiana and D. simulans [41] using Bowtie2 (version 2.3.4.1) [76] with the parameter –X2000. The Samtools suite v0.1.19 [77] was used to convert *.sam to *.bam files and to further process the mapped reads. Duplicates were removed using Picard (version 2.20.2) with the parameter REMOVE_DUPLICATE = TRUE. Bam files were then converted to bed files using the Bedtools (version 2.24) bamtobed command. Reads were centred according to [80]. These reads were then converted to the D. simulans or D. mauritiana coordinate system using liftOver (1.14.0) with custom prepared chain files, one for the conversion of D. mauritiana coordinates to D. simulans coordinates and one for the conversion of D. simulans coordinates to D. mauritiana coordinates. Peaks were then called using MACS2 (version 2.1.2, [81]) with the following parameters: –shift – 100, extsize 200, -q 0.01.

We used the Diffbind package (version 2.12.0, [82]) in R (version 3.6.1.) to search for differentially accessible ATAC-seq regions. A consensus peak set of 19,872 peaks (96 h AEL) and 15,868 peaks (120 h AEL) was used for all samples and the reads were counted for each identified peak with the dba.count command. For each time point separately we used the dba.analyze command with default parameters to get differentially accessible peaks between the two species. This command uses by default the DESeq2 analysis. All plots were generated with the DiffBind package.

Reporter assays

All APREs with the exception of D. mauritiana-APRE6, D. simulans-APRE11 and D. mauritiana-APRE11 were cloned upstream of GAL4 into the pBPGUw backbone (a gift from Gerald Rubin, Addgene plasmid 17,575) by Genewiz. D. mauritiana-APRE6, D. simulans-APRE11 and D. mauritiana-APRE11 were first sub-cloned into the pENTR-TOPO-D plasmid backbone then shuttled into the pBPGUw backbone using Gateway cloning. Sequences of all APREs used can be found in Additional File 23: Table S11.

otd-APRE7-8-GAL4 lines with D. mauritiana and D. simulans sequences were crossed with UAS-nlsGFP to characterise the activity of APRE7-8. EADs were dissected and fixed according to standard protocols. Immunostaining was performed with mouse anti-Elav (Hybridoma Bank) at 1:100 and rabbit anti-GFP (Molecular probes) at 1:1000, followed by secondary antibody incubation with anti-rabbit Alexa Fluor 488 and anti-mouse Alexa Fluor 568, both at 1:400. Imaging was performed using a Leica Stellaris microscope with a 40X objective and zoom 2. Images were processed using the NIH ImageJ software.

We used Elav staining to count the number of ommatidial rows as a measure of developmental stage. We define the first row showing full enhancer activity as the first row counting from the posterior margin of the disc with more than 90% of the ommatidia positive for enhancer signal/GFP signal. For all discs, we determined the developmental time and the first row with full enhancer activity. Then, we calculated the distance between the MF and the first row of full enhancer activity as the difference between the total ommatidial rows and the first ommatidial row with full enhancer activity.

We then plotted the first row with full enhancer activity at different developmental times and fitted a regression line to model the observed linear trend for each species, with calculated R value for correlation and p value of fit.

Violin plots were used to show the distribution of distances between the morphogenetic furrow and the first row of full enhancer activity. We used Welch’s t-test to test whether the difference in distances between species was significant. In the case of the three-species comparison (Additional File 18: Fig. S10), we use the nonparametric Games-Howell test to make pairwise comparisons among the three species. We used an analysis of covariance (ANCOVA) to test for differences between strains for Otd-positive ommatidia rows while adjusting for differences in development stage by using the number of ommatidial rows as a proxy for the latter. The ANCOVA was performed using base R v4.0.2 (R Core Team, 2020).

Prediction of transcription factor binding sites and potential target genes

To search for TFBM of potential otd regulators, we used the universalmotif package (version 1.20.0) using its function to scan sequences. The JASPAR core and CISBP databases were used for screening DNA sequences of ATAC-seq peaks in the otd APRE7-8 sequence with all possible TFBSs from Drosophila with a threshold of p-value > 0.001. Potential fixed changes in the 2.5 kb APRE7-8 sequence of D. mauritiana were inferred from aligning the sequences of this region from strains of D. mauritiana (Red3, mau12, mav2 and TAM16), D. simulans (m3 and w501) and the D. melanogaster sequence using Clustal Omega with default parameters [77].

To define a list of potential Otd target genes, we used an Otd-motif (Dmelanogaster-FlyFactorSurvey-Oc_Cell_FBgn0004102) from the MotifDB package (version 1.16.1), which provides a collection of available transcription factors in R (version 3.3.3). We searched for Otd binding sites in accessible chromatin regions with the findMotifsGenome.pl command implemented in the HOMER (version V4.10.4, [83]) in all samples. All peaks with a predicted Otd motif were annotated to an associated gene using the annotatePeaks.pl command by HOMER and combined all time points and both species into one file. We then looked for the number of genes with an annotated Otd motif and found 1,148 unique genes, which we overlapped with our RNA-seq dataset to find out which of these target genes were differentially expressed between the two species. GO term enrichment analysis of putative Otd target genes was performed using the online tool Metascape [84].

We used the online STRING database that integrates all known and predicted associations between proteins based on evidence from a variety of sources [85], to construct networks of DEG encoded proteins. To visualize the network and map genes/prot with Otd motifs, we applied the Cytoscape software [86].

Data availability

All RNAseq and ATACseq reads are accessible in the Short Read Archive through umbrella BioProject PRJNA666691 [30] (https://www.ncbi.nlm.nih.gov/bioproject/?term=PRJNA666691, containing PRJNA374838 and PRJNA666524).

References

  1. Arnoult L, Su KFY, Manoel D, Minervino C, Magriña J, Gompel N, et al. Emergence and diversification of fly pigmentation through evolution of a gene regulatory module. Science. 2013;339(6126):1423–6.

    Article  CAS  PubMed  Google Scholar 

  2. Gompel N, Prud’homme B, Wittkopp PJ, Kassner VA, Carroll SB. Chance caught on the wing: cis-regulatory evolution and the origin of pigment patterns in Drosophila. Nature. 2005;433(7025):481–7.

    Article  CAS  PubMed  Google Scholar 

  3. Klaassen H, Wang Y, Adamski K, Rohner N, Kowalko JE. CRISPR mutagenesis confirms the role of oca2 in melanin pigmentation in Astyanax mexicanus. Dev Biol. 2018;441(2):313–8.

    Article  CAS  PubMed  Google Scholar 

  4. Ramaekers A, Claeys A, Kapun M, Mouchel-Vielh E, Potier D, Weinberger S, et al. Altering the temporal regulation of one transcription factor drives evolutionary trade-offs between head sensory organs. Dev Cell. 2019;50(6):780-792.e7.

    Article  CAS  PubMed  Google Scholar 

  5. Ridgway AM, Hood EJ, Jimenez JF, Nunes MDS, McGregor AP. Sox21b underlies the rapid diversification of a novel male genital structure between Drosophila species. Curr Biol CB. 2024;34(5):1114-1121.e7.

    Article  CAS  PubMed  Google Scholar 

  6. Santos ME, Le Bouquin A, Crumière AJJ, Khila A. Taxon-restricted genes at the origin of a novel trait allowing access to a new environment. Science. 2017;358:386–90.

    Article  CAS  PubMed  Google Scholar 

  7. Stern DL, Frankel N. The structure and evolution of cis -regulatory regions: the shavenbaby story. Philos Trans R Soc B Biol Sci. 2013;368(1632):20130028.

    Article  Google Scholar 

  8. Wucherpfennig JI, Howes TR, Au JN, Au EH, Roberts Kingman GA, Brady SD, et al. Evolution of stickleback spines through independent cis-regulatory changes at HOXDB. Nat Ecol Evol. 2022;6(10):1537–52.

    Article  PubMed  PubMed Central  Google Scholar 

  9. Xia B, Zhang W, Zhao G, Zhang X, Bai J, Brosh R, et al. On the genetic basis of tail-loss evolution in humans and apes. Nature. 2024;626(8001):1042–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Arif S, Murat S, Almudi I, Nunes MDS, Bortolamiol-Becet D, McGregor NS, et al. Evolution of mir-92a underlies natural morphological variation in Drosophila melanogaster. Curr Biol CB. 2013;23(6):523–8.

    Article  CAS  PubMed  Google Scholar 

  11. Courtier-Orgogozo V, Arnoult L, Prigent SR, Wiltgen S, Martin A. Gephebase, a database of genotype-phenotype relationships for natural and domesticated variation in Eukaryotes. Nucleic Acids Res. 2020;48(D1):D696-703.

    Article  CAS  PubMed  Google Scholar 

  12. Kittelmann M, McGregor AP. Looking across the gap: Understanding the evolution of eyes and vision among insects. BioEssays News Rev Mol Cell Dev Biol. 2024;46(5):e2300240.

    Article  Google Scholar 

  13. Land, Nilsson. Animal eyes. 2nd ed. Oxford: Oxford university press; 2012. (Oxford animal biology series).

  14. Land MF. Visual acuity in insects. Annu Rev Entomol. 1997;42:147–77.

    Article  CAS  PubMed  Google Scholar 

  15. Currea JP, Smith JL, Theobald JC. Small fruit flies sacrifice temporal acuity to maintain contrast sensitivity. Vision Res. 2018;149:1–8.

    Article  PubMed  PubMed Central  Google Scholar 

  16. Palavalli-Nettimi R, Theobald JC. Small eyes in dim light: implications to spatio-temporal visual abilities in Drosophila melanogaster. Vision Res. 2020;169:33–40.

    Article  PubMed  Google Scholar 

  17. Warrant E, Nilsson DE, editors. Invertebrate vision. Cambridge, UK ; New York: Cambridge University Press; 2006. 547 p.

  18. Warrant EJ. Seeing better at night: life style, eye design and the optimum strategy of spatial and temporal summation. Vision Res. 1999;39(9):1611–30.

    Article  CAS  PubMed  Google Scholar 

  19. Buffry AD, Currea JP, Franke-Gerth FA, Palavalli-Nettimi R, Bodey AJ, Rau C, et al. Evolution of compound eye morphology underlies differences in vision between closely related Drosophila species. BMC Biol. 2024;22(1):67.

    Article  PubMed  PubMed Central  Google Scholar 

  20. Duncan AB, Salazar BA, Garcia SR, Brandley NC. A sexual dimorphism in the spatial vision of North American band-winged grasshoppers. Integr Org Biol Oxf Engl. 2021;3(1):obab008.

    Article  CAS  Google Scholar 

  21. Gonzalez-Bellido PT, Wardill TJ, Juusola M. Compound eyes and retinal information processing in miniature dipteran species match their specific ecological demands. Proc Natl Acad Sci U S A. 2011;108(10):4224–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Horridge GA. The compound eye of insects. Sci Am. 1977;237(1):108–20.

    Article  Google Scholar 

  23. Posnien N, Hopfen C, Hilbrant M, Ramos-Womack M, Murat S, Schönauer A, et al. Evolution of eye morphology and rhodopsin expression in the Drosophila melanogaster species subgroup. PLoS ONE. 2012;7(5):e37346.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Wakakuwa M, Stavenga DG, Arikawa K. Spectral organization of ommatidia in flower-visiting insects. Photochem Photobiol. 2007;83(1):27–34.

    Article  CAS  PubMed  Google Scholar 

  25. Gaspar P, Arif S, Sumner-Rooney L, Kittelmann M, Bodey AJ, Stern DL, et al. Characterization of the genetic architecture underlying eye size variation within Drosophila melanogaster and Drosophila simulans. G3 Bethesda Md. 2020;10(3):1005–18.

    Article  CAS  PubMed  Google Scholar 

  26. Land MF. Variations in the structure and design of compound eyes. In: Stavenga DG, Hardie RC, editors. Facets of Vision. Springer: Berlin Heidelberg; 1989. p. 90–111.

    Chapter  Google Scholar 

  27. Perl CD, Niven JE. Differential scaling within an insect compound eye. Biol Lett. 2016;12(3):20160042.

    Article  PubMed  PubMed Central  Google Scholar 

  28. Streinzer M, Brockmann A, Nagaraja N, Spaethe J. Sex and caste-specific variation in compound eye morphology of five honeybee species. PLoS ONE. 2013;8(2):e57702.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Arif S, Hilbrant M, Hopfen C, Almudi I, Nunes MDS, Posnien N, et al. Genetic and developmental analysis of differences in eye and face morphology between Drosophila simulans and Drosophila mauritiana. Evol Dev. 2013;15(4):257–67.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Buchberger E, Bilen A, Ayaz S, Salamanca D, Matas de Las Heras C, Niksic A, et al. Variation in pleiotropic hub gene expression is associated with interspecific differences in head shape and eye size in Drosophila. Mol Biol Evol. 2021;38(5):1924–42.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Hilbrant M, Almudi I, Leite DJ, Kuncheria L, Posnien N, Nunes MDS, et al. Sexual dimorphism and natural variation within and among species in the Drosophila retinal mosaic. BMC Evol Biol. 2014;14:240.

    Article  PubMed  PubMed Central  Google Scholar 

  32. Keesey IW, Grabe V, Gruber L, Koerte S, Obiero GF, Bolton G, et al. Inverse resource allocation between vision and olfaction across the genus Drosophila. Nat Commun. 2019;10(1):1162.

    Article  PubMed  PubMed Central  Google Scholar 

  33. Norry FM, Gomez FH. Quantitative trait loci and antagonistic associations for two developmentally related traits in the Drosophila head. J Insect Sci. 2017;17(1):19.

    Article  PubMed  PubMed Central  Google Scholar 

  34. Reis M, Wiegleb G, Claude J, Lata R, Horchler B, Ha NT, et al. Multiple loci linked to inversions are associated with eye size variation in species of the Drosophila virilis phylad. Sci Rep. 2020;10(1):12832.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Casares F, Almudi I. Fast and Furious 800. The retinal determination gene network in Drosophila. In: Castelli-Gair Hombría J, Bovolenta P, editors. Organogenetic Gene Networks. Cham: Springer International Publishing; 2016. p. 95–124.

  36. Casares F, McGregor AP. The evolution and development of eye size in flies. Wiley Interdiscip Rev Dev Biol. 2021;10(2):e380.

    Article  PubMed  Google Scholar 

  37. Kumar JP. The fly eye: Through the looking glass. Dev Dyn Off Publ Am Assoc Anat. 2018;247(1):111–23.

    CAS  Google Scholar 

  38. Alberch P, Gould SJ, Oster GF, Wake DB. Size and shape in ontogeny and phylogeny. Paleobiology. 1979;5:296–317.

    Article  Google Scholar 

  39. McNamara KJ. Heterochrony: the evolution of development. Evol Educ Outreach. 2012;5(2):203–18.

    Article  Google Scholar 

  40. Haynie JL, Bryant PJ. Development of the eye-antenna imaginal disc and morphogenesis of the adult head in Drosophila melanogaster. J Exp Zool. 1986;237(3):293–308.

    Article  CAS  PubMed  Google Scholar 

  41. Torres-Oliva M, Almudi I, McGregor AP, Posnien N. A robust (re-)annotation approach to generate unbiased mapping references for RNA-seq-based analyses of differential expression across closely related species. BMC Genomics. 2016;17:392.

    Article  PubMed  PubMed Central  Google Scholar 

  42. Torres-Oliva M, Schneider J, Wiegleb G, Kaufholz F, Posnien N. Dynamic genome wide expression profiling of Drosophila head development reveals a novel role of Hunchback in retinal glia cell development and blood-brain barrier integrity. PLoS Genet. 2018;14(1):e1007180.

    Article  PubMed  PubMed Central  Google Scholar 

  43. Vandendries ER, Johnson D, Reinke R. orthodenticle is required for photoreceptor cell development in the Drosophila eye. Dev Biol. 1996;173(1):243–55.

    Article  CAS  PubMed  Google Scholar 

  44. Fichelson P, Brigui A, Pichaud F. Orthodenticle and Kruppel homolog 1 regulate Drosophila photoreceptor maturation. Proc Natl Acad Sci U S A. 2012;109(20):7893–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Mencarelli C, Pichaud F. Orthodenticle is required for the expression of principal recognition molecules that control axon targeting in the Drosophila retina. PLoS Genet. 2015;11(6):e1005303.

    Article  PubMed  PubMed Central  Google Scholar 

  46. Wang LH, Huang YT, Tsai YC, Sun YH. The role of eyg Pax gene in the development of the head vertex in Drosophila. Dev Biol. 2010;337(2):246–58.

    Article  CAS  PubMed  Google Scholar 

  47. Yu L, Zhou Q, Pignoni F. ato-Gal4 fly lines for gene function analysis: Eya is required in late progenitors for eye morphogenesis. Genesis. 2015;53(6):347–55.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Buenrostro JD, Giresi PG, Zaba LC, Chang HY, Greenleaf WJ. Transposition of native chromatin for fast and sensitive epigenomic profiling of open chromatin, DNA-binding proteins and nucleosome position. Nat Methods. 2013;10(12):1213–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Weinkove D, Neufeld TP, Twardzik T, Waterfield MD, Leevers SJ. Regulation of imaginal disc cell size, cell number and organ size by Drosophila class I(A) phosphoinositide 3-kinase and its adaptor. Curr Biol CB. 1999;9(18):1019–29.

    Article  CAS  PubMed  Google Scholar 

  50. Ragge NK, Brown AG, Poloschek CM, Lorenz B, Henderson RA, Clarke MP, et al. Heterozygous mutations of OTX2 cause severe ocular malformations. Am J Hum Genet. 2005;76(6):1008–22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  51. Ranade SS, Yang-Zhou D, Kong SW, McDonald EC, Cook TA, Pignoni F. Analysis of the Otd-dependent transcriptome supports the evolutionary conservation of CRX/OTX/OTD functions in flies and vertebrates. Dev Biol. 2008;315(2):521–34.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Sen S, Reichert H, VijayRaghavan K. Conserved roles of ems/Emx and otd/Otx genes in olfactory and visual system development in Drosophila and mouse. Open Biol. 2013;3(5):120177.

    Article  PubMed  PubMed Central  Google Scholar 

  53. Tahayato A, Sonneville R, Pichaud F, Wernet MF, Papatsenko D, Beaufils P, et al. Otd/Crx, a dual regulator for the specification of ommatidia subtypes in the Drosophila retina. Dev Cell. 2003;5(3):391–402.

    Article  CAS  PubMed  Google Scholar 

  54. Johnston RJ, Otake Y, Sood P, Vogt N, Behnia R, Vasiliauskas D, et al. Interlocked feedforward loops control cell-type-specific Rhodopsin expression in the Drosophila eye. Cell. 2011;145(6):956–68.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  55. Miller DT, Cagan RL. Local induction of patterning and programmed cell death in the developing Drosophila retina. Dev Camb Engl. 1998;125(12):2327–35.

    CAS  Google Scholar 

  56. Alberch P, Alberch J. Heterochronic mechanisms of morphological diversification and evolutionary change in the neotropical salamander, Bolitoglossa occidentalis (Amphibia: Plethodontidae). J Morphol. 1981;167(2):249–64.

    Article  PubMed  Google Scholar 

  57. McKinney ML, editor. Heterochrony in evolution: a multidisciplinary approach. Boston, MA: Springer US; 1988. (Stehli FG, Jones DS, editors. Topics in Geobiology; vol. 7).

  58. Briscoe J, Small S. Morphogen rules: design principles of gradient-mediated embryo patterning. Dev Camb Engl. 2015;142(23):3996–4009.

    CAS  Google Scholar 

  59. Ebisuya M, Briscoe J. What does time mean in development? Dev Camb Engl. 2018;145(12):dev164368.

    Google Scholar 

  60. Farnworth MS, Eckermann KN, Bucher G. Sequence heterochrony led to a gain of functionality in an immature stage of the central complex: a fly–beetle insight. PLOS Biol. 2020;18(10):e3000881.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Keyte AL, Smith KK. Heterochrony and developmental timing mechanisms: changing ontogenies in evolution. Semin Cell Dev Biol. 2014;34:99–107.

    Article  PubMed  Google Scholar 

  62. Gomez C, Ozbudak EM, Wunderlich J, Baumann D, Lewis J, Pourquié O. Control of segment number in vertebrate embryos. Nature. 2008;454(7202):335–9.

    Article  CAS  PubMed  Google Scholar 

  63. Kicheva A, Bollenbach T, Ribeiro A, Valle HP, Lovell-Badge R, Episkopou V, et al. Coordination of progenitor specification and growth in mouse and chick spinal cord. Science. 2014;345(6204):1254927.

    Article  PubMed  PubMed Central  Google Scholar 

  64. Matsuda M, Hayashi H, Garcia-Ojalvo J, Yoshioka-Kobayashi K, Kageyama R, Yamanaka Y, et al. Species-specific segmentation clock periods are due to differential biochemical reaction speeds. Science. 2020;369(6510):1450–5.

    Article  CAS  PubMed  Google Scholar 

  65. Rayon T, Stamataki D, Perez-Carrasco R, Garcia-Perez L, Barrington C, Melchionda M, et al. Species-specific pace of development is associated with differences in protein stability. Science. 2020;369(6510):eaba7667.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  66. Bernardo-Garcia FJ, Fritsch C, Sprecher SG. The transcription factor Glass links eye field specification with photoreceptor differentiation in Drosophila. Dev Camb Engl. 2016;143(8):1413–23.

    CAS  Google Scholar 

  67. Hauck B, Gehring WJ, Walldorf U. Functional analysis of an eye specific enhancer of the eyeless gene in Drosophila. Proc Natl Acad Sci U S A. 1999;96(2):564–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  68. Hay BA, Wolff T, Rubin GM. Expression of baculovirus P35 prevents cell death in Drosophila. Dev Camb Engl. 1994;120(8):2121–9.

    CAS  Google Scholar 

  69. Bessa J, Casares F. Restricted teashirt expression confers eye-specific responsiveness to Dpp and Wg signals during eye specification in Drosophila. Development. 2005;132(22):5011–20.

    Article  CAS  PubMed  Google Scholar 

  70. Navarro T, Iannini A, Neto M, Campoy-Lopez A, Muñoz-García J, Pereira PS, et al. Feedback control of organ size precision is mediated by BMP2-regulated apoptosis in the Drosophila eye. PLOS Biol. 2024;22(1):e3002450.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Struhl G, Basler K. Organizing activity of wingless protein in Drosophila. Cell. 1993;72(4):527–40.

    Article  CAS  PubMed  Google Scholar 

  72. Stampanoni M, Groso A, Isenegger A, Mikuljan G, Chen Q, Bertrand A, et al. Trends in synchrotron-based tomographic imaging: the SLS experience. In: Bonse U, editor. San Diego, California, USA; 2006. p. 63180M

  73. Paganin D, Mayo SC, Gureyev TE, Miller PR, Wilkins SW. Simultaneous phase and amplitude extraction from a single defocused image of a homogeneous object. J Microsc. 2002;206(Pt 1):33–40.

    Article  CAS  PubMed  Google Scholar 

  74. Marone F, Stampanoni M. Regridding reconstruction algorithm for real-time tomographic imaging. J Synchrotron Radiat. 2012;19(Pt 6):1029–37.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  75. Fox J, Weisberg S. An R companion to applied regression. Third edition. Los Angeles: SAGE; 2019. 577 p.

  76. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9(4):357–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  77. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, et al. The sequence alignment/map format and SAMtools. Bioinformatics. 2009;25(16):2078–9.

    Article  PubMed  PubMed Central  Google Scholar 

  78. Rau A, Gallopin M, Celeux G, Jaffrézic F. Data-based filtering for replicated high-throughput transcriptome sequencing experiments. Bioinforma Oxf Engl. 2013;29(17):2146–52.

    Article  CAS  Google Scholar 

  79. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.

    Article  PubMed  PubMed Central  Google Scholar 

  80. Buenrostro JD, Wu B, Chang HY, Greenleaf WJ. ATAC-seq: a method for assaying chromatin accessibility genome-wide. Curr Protoc Mol Biol. 2015;109(1):21–9.

    Article  PubMed Central  Google Scholar 

  81. Zhang Y, Liu T, Meyer CA, Eeckhoute J, Johnson DS, Bernstein BE, et al. Model-based analysis of ChIP-Seq (MACS). Genome Biol. 2008;9(9):R137.

    Article  PubMed  PubMed Central  Google Scholar 

  82. Ross-Innes CS, Stark R, Teschendorff AE, Holmes KA, Ali HR, Dunning MJ, et al. Differential oestrogen receptor binding is associated with clinical outcome in breast cancer. Nature. 2012;481(7381):389–93.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  83. Heinz S, Benner C, Spann N, Bertolino E, Lin YC, Laslo P, et al. Simple combinations of lineage-determining transcription factors prime cis-regulatory elements required for macrophage and B cell identities. Mol Cell. 2010;38(4):576–89.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  84. Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, et al. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun. 2019;10(1):1523.

    Article  PubMed  PubMed Central  Google Scholar 

  85. Szklarczyk D, Gable AL, Nastou KC, Lyon D, Kirsch R, Pyysalo S, et al. Correction to ‘The STRING database in 2021: customizable protein-protein networks, and functional characterization of user-uploaded gene/measurement sets.’ Nucleic Acids Res. 2021;49(18):10800.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  86. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13(11):2498–504.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We acknowledge the Paul Scherrer Institut, Villigen, Switzerland for provision of synchrotron radiation beamtime at the TOMCAT beamline X02DA of the SLS, and the invaluable assistance of Dr Christian Schlepütz. We thank the Transcriptome and Genome Core Unit, University Medical Center Göttingen (UMG) for support with next generation sequencing, the Advanced Light Microscopy and Image Analysis (ALMIA) platform, at CABD, for imaging support, and Franck Pichaud (UCL, UK) for access to otd microarray data. We also thank Madeleine Aase-Remedios, Antonella Iannini, Tomás Navarro and Sebastian Kittelmann for assistance with sequence analysis and statistics.

Funding

This work was funded by grants from the ERC (242553) to A.P.M. and BBSRC (BB/M020967/1, BB/T000317/1) to A.P.M. and I.A. and A.P.M. and M.K respectively. M.T.-O. was supported by a “la Caixa-DAAD” fellowship and and GG by MICINN contract PRE2019-090291. N.P. and E.B. were funded by the Deutsche Forschungsgesellschaft (DFG, Grant No. PO 1648/3–1 to N.P.).

Author information

Author notes

  1. Montserrat Torres-Oliva, Elisa Buchberger and Alexandra D. Buffry contributed equally to this work.

Authors and Affiliations

  1. Department of Developmental Biology, University Göttingen, Justus-Von-Liebig-Weg 11, 37077, Göttingen, Germany

    Montserrat Torres-Oliva, Elisa Buchberger, Georg C. Bullinger & Nico Posnien

  2. Present address: Institute of Clinical Molecular Biology, Christian-Albrechts-University of Kiel, University Hospital Schleswig-Holstein, Kiel, Germany

    Montserrat Torres-Oliva

  3. Department of Biological and Medical Sciences, Oxford Brookes University, Oxford, OX3 0BP, UK

    Alexandra D. Buffry, Maike Kittelmann, Pedro Gaspar, Maria D. S. Nunes & Isabel Almudi

  4. Andalusian Centre for Developmental Biology (CABD), CSIC/Universidad Pablo de Olavide/JA. Ctra. de Utrera Km 1, 41013, Seville, Spain

    Genoveva Guerrero & Fernando Casares

  5. Museum Für Naturkunde, Leibniz-Institut Für Evolutions- Und Biodiversitätsforschung, Invalidenstraße 43, 10115, Berlin, Germany

    Lauren Sumner-Rooney

  6. MRC London Institute of Medical Sciences, London, W12 0NN, UK

    Pedro Gaspar

  7. Centre for Functional Genomics & Department of Biological Sciences, Oxford Brookes University, Oxford, OX3 0BP, UK

    Javier Figueras Jimenez & Saad Arif

  8. Department of Biosciences, Durham University, Durham, DH1 3LE, UK

    Alistair P. McGregor

  9. Department of Genetics, Microbiology and Statistics, Universitat de Barcelona, Diagonal 643, 08028, Barcelona, Spain

    Isabel Almudi

  10. Institut de Recerca de La Biodiversitat (IRBio), Universitat de Barcelona, Diagonal 643, Barcelona, 08028, Spain

    Isabel Almudi

Authors
  1. Montserrat Torres-Oliva
    View author publications

    You can also search for this author inPubMed Google Scholar

  2. Elisa Buchberger
    View author publications

    You can also search for this author inPubMed Google Scholar

  3. Alexandra D. Buffry
    View author publications

    You can also search for this author inPubMed Google Scholar

  4. Maike Kittelmann
    View author publications

    You can also search for this author inPubMed Google Scholar

  5. Genoveva Guerrero
    View author publications

    You can also search for this author inPubMed Google Scholar

  6. Lauren Sumner-Rooney
    View author publications

    You can also search for this author inPubMed Google Scholar

  7. Pedro Gaspar
    View author publications

    You can also search for this author inPubMed Google Scholar

  8. Georg C. Bullinger
    View author publications

    You can also search for this author inPubMed Google Scholar

  9. Javier Figueras Jimenez
    View author publications

    You can also search for this author inPubMed Google Scholar

  10. Fernando Casares
    View author publications

    You can also search for this author inPubMed Google Scholar

  11. Saad Arif
    View author publications

    You can also search for this author inPubMed Google Scholar

  12. Nico Posnien
    View author publications

    You can also search for this author inPubMed Google Scholar

  13. Maria D. S. Nunes
    View author publications

    You can also search for this author inPubMed Google Scholar

  14. Alistair P. McGregor
    View author publications

    You can also search for this author inPubMed Google Scholar

  15. Isabel Almudi
    View author publications

    You can also search for this author inPubMed Google Scholar

Contributions

M.T.-O. performed and analysed RNA-seq experiments, interpreted the results and prepared figures. E.B. Bioinformatics analysis and figures. A.D.B. performed functional tests of otd perturbations and reporter studies, interpreted the results. M.K. Synchrotron radiation microtomography and interpretation of the data. G.G. reporter studies and interpretation. L.S.-R. Synchrotron radiation microtomography. P.G. gene expression experiments G.C.B. gene expression experiments. J.F.J. enhancer sequence analysis. F.C. enhancer data interpretation and supervision. S.A. fine-mapping and introgression. M.D.S.N. fine-mapping and introgression. N.P. supervision, funding and data analysis interpretation. A.P.M. and I.A. work conceptualization, supervision, design, data analysis and interpretation, figures, funding and manuscript writing with the help of all other authors.

Corresponding authors

Correspondence to Saad Arif, Nico Posnien, Maria D. S. Nunes, Alistair P. McGregor or Isabel Almudi.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Table S1. Measurements of ommatidia size of D. mauritiana and D. simulans eyes.

Additional file 2: Table S2. Mapping and introgression data.

12915_2025_2136_MOESM3_ESM.pdf

Additional file 3: Fig. S1. Ommatidia number and body size of the IL lines. Ommatidia number (left) and T1 Tibia length (right) was not significantly different between y, f males and their y, v, f sibling males for each introgression line (two-tailed, two-sample, equal-variance t-tests).

Additional file 4: Fig. S2. RNA-seq datasets. (a) Heat-map of all RNA-seq samples. (b) PCA plot of all RNA-seq samples.

12915_2025_2136_MOESM5_ESM.xlsx

Additional file 5: Table S3. Pair-wise differential gene expression analysis (D.sim vs. D.mau) results for each time point (72 hAEL, 96 hAEL and 120 hAEL).

12915_2025_2136_MOESM6_ESM.pdf

Additional file 6: Fig. S3. otd expression in 120 hAEL eye imaginal discs. (a) D. mauritiana EAD at 120h. Intensity plots show otd expression intensity measured in the red square area with ImageJ plot intensity tool. (b) D. simulans EAD at 120h. (c) Table and plot showing the area under the curve for each of the discs in which the intensity of otd signal was measured. Black arrowhead: Morphogenetic furrow; oc: ocellar region; dev eye: eye region.

12915_2025_2136_MOESM7_ESM.xlsx

Additional file 7: Table S4. Measurements of plot intensities of otd in situ hybridization in D. mauritiana and D. simulans eye discs.

12915_2025_2136_MOESM8_ESM.pdf

Additional file 8: Fig. S4. otd expression in 110 hAEL eye imaginal discs. (a) D. mauritiana EAD at 110h. (b) D. simulans EAD at 110 hAEL. Red arrowheads highlight expression of otd in some D. mauritiana discs. Black arrowhead: Morphogenetic furrow; oc: ocellar region; dev eye: eye region.

Additional file 9: Table S5. Measurements of Otd-positive ommatidia in D. mauritiana and D. simulans eye discs

12915_2025_2136_MOESM10_ESM.pdf

Additional file 10: Fig. S5. otd expression in pupal eyes. (a) D. mauritiana pupal eye (48 h APF) stained with Phalloidin (Actin, a’), anti-Elav marks photoreceptors (a”) and Otd is expressed in all ommatidia (a”’). (b) Otd is present in all ommatidia in D. simulans pupal eyes, (b’) Actin (Phalloidin) highlights ommatidia area, (b”) anti-Elav marks photoreceptors and Otd protein is shown in b”’.

12915_2025_2136_MOESM11_ESM.pdf

Additional file 11: Fig. S6. Effects of loss and gain of otd expression in ommatidial structure and facet size in D. melanogaster. Loss of otd causes defects in ommatidia structure. (a-c) Knockdown of otd by expressing UAS-miR-otd with GMR-GAL4 driver in D. melanogaster eye. (a-a’) GMR-GAL. (b-b’) GMR-GAL/ UAS-miR-otd shows a rough eye phenotype due to defects in ommatidia arrangements. (c-c’) GMR-GAL/ UAS-miR-otd phenotype is rescued by co-expressing UAS-otd. (d-d’) otd imutant mitotic clones also show defects in ommatidia size and organisation. (e) Comparison of facet area from eyes of UAS-otd (“+, otd”), ato3’FL-GAL4 (“ato>+) and UAS-otd; ato3’FL-GAL4 (“ato>otd”) female adults. Facet area is significantly larger in the experimental genotype (“ato>otd”) than in the two control genotypes. See main text for details.

Additional file 12: Table S6. ato-Gal4> UAS-otd ommatidia size.

Additional file 13: Fig. S7. Alignment of D. mauritiana TAM16 and D. simulans y, v, f Otd protein sequences.

12915_2025_2136_MOESM14_ESM.pdf

Additional file 14: Fig. S8. Topological Domain Associated with otd locus in D. melanogaster from http://chorogenome.ie-freiburg.mpg.de/.

Additional file 15: Table S7. Differential peak calling for D. mau and D. sim 96 hAEL and 120 hAEL ATAC-seq data.

12915_2025_2136_MOESM16_ESM.pdf

Additional file 16: Fig. S9. Activity of tested APREs in the otd locus of D. simulans and D. mauritiana. L3 instar eye-antenna imaginal discs from D. melanogaster lines containing D. mauritiana-APRE4-GAL4, D. mauritiana-APRE6-GAL4, D. mauritiana-APRE7-8-GAL4, D. simulans-APRE11-GAL4, D. mauritiana-APRE13-14-GAL4, D. mauritiana-APRE16-17-GAL4, D. mauritiana-APRE19-GAL4 crossed to UAS-GFP and stained with DAPI (grey). Scale bar = 50μm.

12915_2025_2136_MOESM17_ESM.xlsx

Additional file 17: Table S8. TFBM in differential open chromatin regions in otd locus of D. simulans and D. mauritiana.

12915_2025_2136_MOESM18_ESM.pdf

Additional file 18: Fig. S10. D. melanogaster, D. simulans and D. mauritiana otd-APRE7-8 enhancer activity in 3rd instar larvae eye imaginal discs. (a) D. melanogaster, D. simulans and D. mauritiana otd-APRE7-8 3rd instar larvae eye imaginal discs immunostained with anti-Elav antibody to visualise the progression of ommatidia maturation as proxy for the developmental time. (b) Plot showing the number of GFP-positive ommatidia rows (x-axis) activated by otd-APRE7-8 at different developmental time points (y-axis, developmental points inferred by number of ommatidia rows) for both species. (c) Violin plots showing the distance between MF and the first row of otd-APRE7-8 activity. D. mauritiana otd-APRE7-8 is active earlier during the differentiation of ommatidia.

12915_2025_2136_MOESM19_ESM.pdf

Additional file 19: Dataset S1. Clustal Omega alignment of the 2.5 kb APRE7-8 region from selected strains of D. mauritiana, D. simulans and D. melanogaster including the sequences used in enhancer reporter constructs for each species (otd-APRE7-8).

12915_2025_2136_MOESM20_ESM.pdf

Additional file 20: Fig. S11. Schematic showing the position of the 2.5 kb APRE7-8 within the otd locus on the X chromosome of D. melanogaster. A 1.5 kb region shown to drive similar expression to the 2.5 kb region is indicated by a purple bar (37). Below is shown the ATAC-seq profile for this region and the positions of potential fixed mutations specific to D. mauritiana (arrows indicate SNPs and rectangles indicate short indels) with predicted binding sites for Otd, So and Cut in D. melanogaster indicated by arrows. Binding sites indicated with an asterisk contain mutations in D. mauritiana.

12915_2025_2136_MOESM21_ESM.xlsx

Additional file 21: Table S9. Differentially expressed genes with predicted Otd TFBMs in their associated open chromatin regions.

Additional file 22: Table S10. Primers used to generate probes for in situ hybridization.

Additional file 23: Table S11. Sequences used to generate reporter lines.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Torres-Oliva, M., Buchberger, E., Buffry, A.D. et al. Heterochrony in orthodenticle expression is associated with ommatidial size variation between Drosophila species. BMC Biol 23, 34 (2025). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12915-025-02136-8

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12915-025-02136-8

Keywords

BMC Biology

ISSN: 1741-7007

Contact us