Conserved and Divergent Aspects of Plasticity and Sexual Dimorphism in Wing Size and Shape in Three Diptera

The ability of powered flight in insects facilitated their great evolutionary success allowing them to occupy various ecological niches. Beyond this primary task, wings are often involved in various premating behaviors, such as the generation of courtship songs and the initiation of mating in flight. These specific functions imply special adaptations of wing morphology, as well as sex-specific wing morphologies. Although wing morphology has been extensively studied in Drosophila melanogaster (Meigen, 1830), a comprehensive understanding of developmental plasticity and the impact of sex on wing size and shape plasticity is missing for other Diptera. Therefore, we raised flies of the three Diptera species Drosophila melanogaster, Ceratitis capitata (Wiedemann, 1824) and Musca domestica (Linnaeus, 1758) at different environmental conditions and applied geometric morphometrics to analyze wing shape. Our data showed extensive interspecific differences in wing shape, as well as a clear sexual wing shape dimorphism in all three species. We revealed an impact of different rearing temperatures on wing shape in all three species, which was mostly explained by plasticity in wing size in D. melanogaster. Rearing densities had significant effects on allometric wing shape in D. melanogaster, while no obvious effects were observed for the other two species. Additionally, we did not find evidence for sex-specific response to different rearing conditions in D. melanogaster and C. capitata, while a male-specific impact of different rearing conditions was observed on non-allometric wing shape in M. domestica. Overall, our data strongly suggests that many aspects of wing morphology underly species-specific adaptations and we discuss potential developmental and functional implications of our results.

The ability of powered flight in insects facilitated their great evolutionary success allowing them to occupy various ecological niches. Beyond this primary task, wings are often involved in various premating behaviors, such as the generation of courtship songs and the initiation of mating in flight. These specific functions imply special adaptations of wing morphology, as well as sex-specific wing morphologies. Although wing morphology has been extensively studied in Drosophila melanogaster (Meigen, 1830), a comprehensive understanding of developmental plasticity and the impact of sex on wing size and shape plasticity is missing for other Diptera. Therefore, we raised flies of the three Diptera species Drosophila melanogaster, Ceratitis capitata (Wiedemann, 1824) and Musca domestica (Linnaeus, 1758) at different environmental conditions and applied geometric morphometrics to analyze wing shape. Our data showed extensive interspecific differences in wing shape, as well as a clear sexual wing shape dimorphism in all three species. We revealed an impact of different rearing temperatures on wing shape in all three species, which was mostly explained by plasticity in wing size in D. melanogaster. Rearing densities had significant effects on allometric wing shape in D. melanogaster, while no obvious effects were observed for the other two species. Additionally, we did not find evidence for sex-specific response to different rearing conditions in D. melanogaster and C. capitata, while a male-specific impact of different rearing conditions was observed on non-allometric wing shape in M. domestica. Overall, our data strongly suggests that many aspects of wing morphology underly speciesspecific adaptations and we discuss potential developmental and functional implications of our results.

INTRODUCTION
Natural variation in size and shape of animal body parts is often the result of adaptation to an ever-changing environment. Insect wings are an excellent model to study organ size and shape because insects are the only group of arthropods that developed the ability of a powered flight. This adaptation allowed them to occupy various ecological niches, including air, and contributed to the morphological diversity and great ecological success of the entire class. Flying helps insects to surmount long distances in a relatively short time, facilitating basic tasks, such as finding mating partners and food resources. In line with these crucial functions, wing size has been shown to be highly variable between (Houle et al., 2003;Salcedo et al., 2019) and within species (Coyne and Beecham, 1987;Imasheva et al., 1994;James et al., 1995James et al., , 1997van't Land et al., 1999;Zwaan et al., 2000;Alves and Bélo, 2002;Hoffmann and Shirriffs, 2002;Bai et al., 2016;Rohner et al., 2019). Most of our current understanding of intraspecific variation in wing size is based on extensive work in Drosophila melanogaster (Meigen, 1830). For instance, longitudinal (Coyne and Beecham, 1987;Imasheva et al., 1994;James et al., 1997;Zwaan et al., 2000) and altitudinal (Pitchers et al., 2013;Klepsatel et al., 2014) clines in wing size have been observed on all continents studied to date with positive correlations between wing size and latitude and altitude, respectively. The analysis of latitudinal clines in conjunction with laboratory experiments showed that temperature is a major environmental factor affecting wing size (Cavicchi et al., 1985;Partridge et al., 1994a;James et al., 1995;van't Land et al., 1999;Zwaan et al., 2000). The effect of temperature on wing size can be further modulated by other environmental factors. Wing length, for example, reacts more to changes in rearing temperature under limited food conditions (de Moed et al., 1997). Overall, this data obtained from D. melanogaster revealed that natural variation in wing size is caused by a combination of the genetic and environmental factors.
Variation in organ size must be accompanied by shape changes to retain fully functional properties. As shown for wing size, latitudinal and altitudinal clines for wing shape have been observed in D. melanogaster and environmental factors, such as temperature, have profound effects on wing shape (Gilchrist et al., 2000;Debat et al., 2003Debat et al., , 2009Pitchers et al., 2013). Advances in geometric morphometrics methodology to quantify shape variation nowadays allow analyzing size and shape independently to study the allometric effect on wing shape (Bookstein, 1996;Bitner-Mathé and Klaczko, 1999;Debat et al., 2003;Mitteroecker and Gunz, 2009;Klingenberg, 2016). Such analyses revealed that natural variation in D. melanogaster wing shape is predominantly caused by differences in wing size (i.e., strong allometric component of shape variation) (Gilchrist et al., 2000;Debat et al., 2003). Strong allometric effects on wing shape as well as similar clinal patterns and plastic responses to different temperature regimes have been observed in other insects, such as sepsid flies (Rohner et al., 2019) and grasshoppers (Bai et al., 2016). These findings suggests that wing size and shape must be well-coordinated during development. Indeed, in D. melanogaster both aspects of wing morphology are regulated by similar patterning and differentiation processes during larval and pupal wing imaginal disc development (Day and Lawrence, 2000;Matamoro-Vidal et al., 2015;Testa and Dworkin, 2016). However, despite this tight connection between wing size and shape, laboratory crosses in D. melanogaster showed that wing size variation is shaped by directional selection, while wing shape shows signs of optimizing selection (Gilchrist and Partridge, 2001). Similarly, experimental evolution of D. subobscura (Collin, 1936) at three different temperature regimes for 81 generations revealed clear effects of temperature adaptation on wing shape, while wing size was not affected (Santos et al., 2006). These results suggest that wing size and shape are under different selection pressures and may be less functionally connected than expected. Support for this notion also comes from studies in other insects, such as the Muscidae Polietina orbitalis (Stein, 1904;Alves et al., 2016) and Australian populations of D. serrata (Malloch, 1927;Hoffmann and Shirriffs, 2002) where wing shape varies independently from wing size. Data from different insects, thus, revealed contrasting results with respect to the impact of wing size on wing shape in different environmental conditions. To test whether allometric aspects of wing shape in different rearing conditions are indeed species-specific or rather conserved among insects, interspecific comparisons with comparable experimental designs are needed. We have previously shown that the size of different body parts, including the wings, vary in response to different rearing temperatures and densities in the three Diptera species D. melanogaster, the Mediterranean fruit fly Ceratitis capitata (Wiedemann, 1824) and the housefly Musca domestica (Linnaeus, 1758;Siomava et al., 2016). Therefore, we used these three species to reveal conserved and species-specific effects of different rearing temperatures and densities on allometric aspects of wing shape.
In modern insects, wings became engaged in other essential processes, such as mating. For instance, numerous fruit flies, such as D. melanogaster and C. capitata use their wings to perform courtship songs, reflecting the size and vigor of males to help females choosing the right mating partner (Keiser et al., 1973;Burk and Webb, 1983;Webb et al., 1983;Partridge et al., 1987). Some insects mate in the air during flight (Wilkinson and Johns, 2005), while others, such as M. domestica initiate the mating process in flight but always land prior to copulation (Murvosh et al., 1964). These different behaviors in combination with intersexual food competition and various reproductive roles of the wing pose a constant selective pressure resulting in pronounced sexual dimorphism in wing size and shape in insects (Gidaszewski et al., 2009;de Camargo et al., 2015). In addition to sexual wing size dimorphism (Siomava et al., 2016), D. melanogaster exhibits a clear wing shape dimorphism with male wings being broader than female wings (Bitner-Mathé and Klaczko, 1999;Gilchrist et al., 2000;Gidaszewski et al., 2009). Given that males of C. capitata produce courtship songs like D. melanogaster, it is not surprising that similar trends in wing size and shape have been observed in this species (Siomava et al., 2016;Pieterse et al., 2017;Lemic et al., 2020). Wing size (Siomava et al., 2016;Cortés-Suarez et al., 2021) and shape (Rohner, 2020;Cortés-Suarez et al., 2021) are also sexually dimorphic in M. domestica. Since this species does not produce songs, but relies on the initiation of mating in flight, it remains to be tested, whether sexual differences in wing shape reflect this mating behavior in comparison to song-producing Diptera, such as D. melanogaster and C. capitata. In addition to sexual dimorphisms in wing morphology, data from D. melanogaster showed that plastic responses of wing shape to different rearing temperatures is sex-specific because opposite shape changes were observed in males and females (Aytekin et al., 2009). However, in the other two species it is still unknown whether plasticity in wing shape is sexually dimorphic. Therefore, we studied plasticity in wing morphology in different Diptera species to reveal conserved and divergent features of sexual size and shape dimorphisms.
In this study, we applied geometric morphometrics to comprehensively describe and compare plasticity in wing shape in C. capitata, D. melanogaster, and M. domestica. Furthermore, we asked whether plasticity in allometric and non-allometric wing shape were sexually dimorphic in the three species. We found clear evidence for sexual shape dimorphism in all three species with major contribution of wing size on shape differences in D. melanogaster. Different rearing temperatures had a strong effect on total and non-allometric wing shape in all three species, while density effects were most pronounced in D. melanogaster. Eventually, we did not find strong arguments for sexual dimorphism in response to the different rearing conditions, although M. domestica males showed slightly stronger differences than females. We identified highly variable regions, the radio-medial (r-m) crossvein, the R2 + 3 radial vein and the basal-medial-cubital (bm-cu) crossvein, which changed similarly among the three species in response to various larval densities and temperature regimes. Overall, we reveal mostly speciesspecific aspects of plasticity and sexual dimorphism in wing morphology, and we discuss developmental and functional implications of our results.

Fly Species
All experiments were performed using three different fly species ( Figure 1A; Siomava et al., 2016). We used the highly inbred laboratory strain Drosophila melanogaster w1118, which was kept at 18 • C on standard food (400 g of malt extract, 400 g of corn flour, 50 g of soy flour, 110 g of sugar beet syrup, 51 g of agar, 90 g of yeast extract, 31.5 ml of propionic acid, and 7.5 g of Nipagin dissolved in 40 ml of Ethanol, water up to 5 l). The other two flies were Musca domestica wild type ITA1 collected in Italy, Altavilla Silentia in 2013 (Y. Wu and L. Beukeboom, GELIFES, Netherlands) and Ceratitis capitata wild type Egypt II (IAEA). The M. domestica strain was reared at room temperature (RT) (22 ± 2 • C) on food composed by 500 g of wheat bran, 75 g of wheat flour, 60 g of milk powder, 25 g of yeast extract, 872 ml of water, and 18.85 ml of Nipagin (2.86 g of Nipagin in 10 ml of Ethanol). Adult M. domestica flies were provided with sugar water. C. capitata were kept at 28 • C, 55 ± 5% RH on a diet composed by 52.5 g of yeast extract, 52.5 g of carrot powder, 2 g of Sodium benzoate, 1.75 g of agar, 2.25 ml of 32% HCl, and 5 ml of Nipagin (2.86 g of Nipagin in 10 ml of Ethanol), water up to 500 ml for larvae. For adult flies, we used a 1:3 mixture of sugar and yeast extract.

Treatment of Experimental Groups
To generate a range of sizes for each species, we varied two environmental factors known to influence overall body size -temperature and density. Prior to the experiment, D. melanogaster flies were placed at 25 • C for two days. On the third day, flies were moved from vials into egg-collection chambers and provided with apple-agar plates. After several hours, we started egg collection by removing apple-agar plates with laid eggs once per hour. Collected plates were kept at 25 • C for 24 h to allow embryonic development to complete. Freshly hatched first-instar larvae were transferred into 50 ml vials with 15 ml of fly food. Three vials containing 25 larvae each (low density) and three vials with 300 larvae each (high density) were moved to 18 • C; the second set of six vials with the same densities was left at 25 • C.
Ceratitis capitata flies were kept at 28 • C and allowed to lay eggs through a net into water. Every hour, eggs were collected and transferred to the larval food. After 22 h, first-instar larvae were transferred into small Petri dishes (diameter 55 mm) with 15 ml of the larval food in three densities: 25 (low density), 100 (medium density), or 300 (high density) larvae per plate. Two plates of each density were moved to 18 • C. The second set of six plates was left at 28 • C for further development.
Eggs of M. domestica were collected in the wet larval food at RT and after 24 h, all hatched larvae were removed from food. Only larvae hatched within the next hour were transferred into 50 ml vials with 5 g of food. Collection of larvae was repeated several times to obtain two experimental sets with three replicates of three experimental densities 10 (low density), 20 (medium density), or 40 (high density) larvae. One set of nine vials was moved to 18 • C, the other was left at RT.
After pupation, individuals of C. capitata and M. domestica were collected from the food and kept until eclosion in vials with a wet sponge, which was refreshed every other day.
The experimental temperature regimes were chosen for the following reasons. D. melanogaster is known to survive in the range 10 to 33 • C, but flies remain fertile at 12 to 30 • C with the optimum at 25 • C (Hoffmann, 2010). Reproduction temperatures in C. capitata range from 14 • C to 30 • C with the optimum at 28 • C (Duyck and Quilici, 2002;Navarro-Campos et al., 2011). M. domestica flies survive at 10 to 35 • C (Hewitt, 1914) with the optimum between 24 and 27 • C (Hafez, 1948;Chun-Hsung, 2012). The low temperature for our experiment was chose as the one above the survival and fertile minimums for all three species -18 • C. The warm temperature was aimed to be optimal for each species.

Data Collection
For each combination of temperature and density, we randomly picked 8 to 24 (median: 12.5) flies for C. capitata, 20 to 31 (median: 23) flies for D. melanogaster and 5 to 15 (median: 8) flies for M. domestica. Both wings were dissected, embedded in Roti R -Histokitt II (Roth, Buchenau) on a microscope slide, and photographed under the Leica MZ16 FA stereo microscope with the Q Imaging Micro Publisher 5.0 RTV Camera.
Wing shape was analyzed using landmark-based geometric morphometric methods (Rohlf, 1990;Bookstein, 1991). We digitized 10 anatomically homologous landmarks on wings of the three species ( Figure 1A). The landmarks were the following [nomenclature is given after (Colless and McAlpine, 1991)]: 1, branching point of veins R 1 and R S (base of R 2+3 and R 4+5 ); 2, branching point of veins R 2+3 and R 4+5 ; 3, intersection of FIGURE 1 | Schematic overview of the phylogenetic relationships Muscomorpha and interspecific wing shape variation D. melanogaster, C. capitata, and M. domestica. (A) Homologous landmarks 1-10 in C. capitata, D. melanogaster, and M. domestica are shown as black points with the respective number in red. A, anal vein; bm-cu, basal-madial-cubital crossvein; CuA, anterior cubital vein; C, costal vein; h, humeral crossvein; M, medial vein; R, radial veins; r-m, radio-medial crossvein; SC, subcosta. Branch lengths of the tree do not indicate evolutionary time or distance. (B) Wing shape variation between D. melanogaster, C. capitata, and M. domestica. Principal Component Analysis of shape (PC1 and PC2). The red wireframes depict shape changes along the main axes of variation (PC1: -0.2, 0.2; PC2: -0.2, 0.1) relative to the gray wireframes (0 for both axes). The amount of variation explained by each PC is shown in brackets and the number of individuals analyzed for each sex is provided next to each species.

Procrustes Superimposition and Interspecies Comparison
Wing images were digitized using tpsUtil and tpsDig2 (Rohlf, 2015) to obtain raw x and y landmark coordinates. Using superimposition methods, it is possible to register landmarks of a sample to a common coordinate system in three steps: translating all landmark configurations to the same centroid, scaling all configurations to the same centroid size, and rotating all configurations until the summed squared distances between the landmarks and their corresponding sample average is a minimum scaling (Slice, 2005;Mitteroecker et al., 2013). To follow these three steps, we applied the generalized Procrustes analysis (GPA) (Dryden and Mardia, 1998;Slice, 2005) as implemented in MorphoJ (version 1.06d) (Klingenberg, 2011). The wings were aligned by principal axes, the mean configuration of landmarks was computed, and each wing was projected to a linear shape tangent space. The coordinates of the aligned wings were the Procrustes coordinates. To evaluate wing asymmetry, we used the R package Geomorph (v. 3.3.1) (Adams et al., 2021) to partition shape variation by individuals, side (directional asymmetry) and the interaction between individual and side (fluctuating asymmetry). Procrustes ANOVA was used to evaluate statistical significance. Since we obtained no evidence for fluctuating asymmetry for all species (Supplementary Table 1), we averaged the coordinates for the right and left wings for further analyses. Next, we used Type III Procrustes ANOVA (with Randomization of null model residuals and 1,000 permutations) as implemented in Geomorph (v. 3.3.1) (Adams et al., 2021) to determine the effect of the species, sex, temperature, and density as well as potential interactions between these explanatory variables and shape ( Table 1). Temperature and density values were re-labeled as low, high, or intermediate to allow the comparison between species. To visualize the main components of shape variation, we performed a principal component analysis (PCA) and generated wireframes using MorphoJ (version 1.06d) ( Figure 1B).

Intraspecific Sexual Dimorphism and Effects of Rearing Conditions
Since most of the variation in shape was explained by differences between species (see section "Results"), we split the analysis to further evaluate the effects of sex, rearing temperature, and density as well as potential interactions on wing shape within each species using Procrustes ANOVA in Geomorph (v. Table 2). Magnitudes of sexual shape dimorphism were estimated using the discriminant function analysis (DFA) and expressed in units of Procrustes distance using MorphoJ (version 1.06d). DFA identifies shape features that differ the most between groups relative to within groups and it can only be applied to contrast two experimental groups. Therefore, we used this method to define sexual shape dimorphism (males and females), as well as effects of the rearing temperature (high and low) and density (high and low) in each species. Wing shape changes were visualized using wireframe graphs. To test for the significance of the observed differences, we ran a permutation test with 1,000 random permutations (Good, 1994) for each test using MorphoJ (version 1.06d).

Estimation of the Non-allometric Component of Shape
To evaluate the impact of wing size on shape variation, we estimated wing centroid size (WCS) that was computed from raw data of landmarks and measured as the square root of the sum of squared deviations of landmarks around their centroid (Bookstein, 1996;Dryden and Mardia, 1998;Slice, 2005) in MorphoJ (version 1.06d). The wing images were taken at different magnifications to account for the different body sizes of the three species. C. capitata wings were takes at a 1,280 × 960 resolution, M. domestica at a 2,560 × 1,920 resolution and for D. melanogaster images with both resolutions were obtained. Therefore, the WCS values were corrected for differences in magnification and resolution among photos before being used as covariate. The effects of sex, temperature, density, and their interactions on WCS was tested using Type III ANOVA for each species.
Since a Type III Procrustes ANOVA for all three species including WCS as covariate revealed significant interaction between species and WCS (Supplementary Table 3, sheet "All Species_total"), we estimated the non-allometric shape component for each species independently. The non-allometric Procrustes coordinates were obtained as the residuals of the multiple regression of WCS onto the Procrustes coordinates pooling the samples by sex, temperature, and density, as implemented in MorphoJ (version 1.06d) (Klingenberg, 2011(Klingenberg, , 2016). To evaluate the success of the size correction, we ran a Type II ANOVA testing for additive effects of explanatory variables prior to and after the size correction. This analysis revealed that WCS had no impact on wing shape after size correction (Supplementary Table 3, compare sheets " * _total" to " * _non-allometric"), suggesting that the approach accounted for most of the impact of size on shape. Note that species-specific Type III Procrustes ANOVAs including WCS as covariate revealed significant interactions for C. capitata and D. melanogaster (Supplementary Table 3, sheets "C. capitata_total, " "D. melanogaster_total"), suggesting complex relationships between wing size and shape. Since the regression approach to estimate allometry across all explanatory variables (i.e., sex, temperature, and density) assumes similar slopes of within-group regressions (Gidaszewski et al., 2009;Klingenberg, 2016), we analyzed those regressions in more detail. Although the directions of within-group regressions (i.e., regression scores of the regression of Procrustes coordinates onto WCS pooled by sex, temperature, and density) were comparable (Figure 3), type III ANOVA testing for interactions of explanatory variables partially resulted in significant interactions in C. capitata and D. melanogaster (Supplementary Table 4), suggesting some significant differences in slopes. Therefore, the distinction between allometric and non-allometric shape differences may not be fully accurate. However, due to a lack of a better method, we decided to proceed with the analysis and to interpret the results cautiously. The amount of shape variation explained by wing centroid size was estimated by obtaining R 2 for the linear model describing the pooled regression (i.e., the unique curve after pooling the data by groups: centered_reg_scores ∼centered_WCS).
Magnitudes of non-allometric shape differences for each explanatory variable and for sex-specific effects of different temperature and density conditions were estimated using DFA as described before. Since wing size and shape of intermediate densities for C. capitata and M. domestica overlapped with the high-and low-density conditions, respectively (Figure 3), we restricted the detailed shape analysis to the extreme values.

Wing Shape Variation in Dipteran Species
To characterize wing shape variation among Ceratitis capitata, Drosophila melanogaster, and Musca domestica, we applied geometric morphometrics based on 10 landmarks ( Figure 1A). Procrustes ANOVA revealed that most of the shape variance was caused by differences among species (Table 1). This result was also reflected by the principal component analysis (PCA) which clearly distinguished the three species with the first two principal components (PCs) accounting for almost 98% of the variation (Figure 1B). The main shape difference along PC1 (77.3% of the variation) reflected the ratio between the proximal and distal parts of the wing as well as wing length. The wireframe graphs showed that C. capitata wings were broad in the proximal part (landmarks 1-5) and narrow in the distal part (landmarks 6-10), while D. melanogaster wings showed the opposite pattern with the proximal part being heavily compressed along the anteriorposterior axis. M. domestica had an intermediate morphology being, however, more similar to C. capitata ( Figure 1B, wireframe graphs along PC1). PC2 explained 20.5% of the variation mainly accounting for the ratio between the length and width of the whole wing. C. capitata and D. melanogaster showed wider and shorter wings, while M. domestica showed more elongated wings ( Figure 1B, wireframe graphs along PC2). Together, we found evidence for extensive wing shape variation among C. capitata, D. melanogaster, and M. domestica.

Sexual Dimorphism in Wing Shape
Our interspecific shape analysis revealed a clear sexual shape dimorphism which was most pronounced for C. capitata ( Figure 1B and Table 1). The extreme sexual shape dimorphism in C. capitata and the observation that males and females occupy different relative positions along PC1 compared to D. melanogaster and M. domestica likely explain the significant interaction between species and sex ( Table 1). Species-specific Procrustes ANOVA for all three species revealed that withinspecies shape variation was affected independently by sex and rearing conditions (i.e., temperature and density, see below) (Supplementary Table 2). Therefore, we split the analysis by species and ran a DFA to determine the total shape differences caused by sex (Figures 2A-C). For all three species we found a highly significant sexual shape dimorphism (Figures 2A-C). Male wings were broader than those of females in C. capitata (Figure 2A) and D. melanogaster (Figure 2B), while the opposite trend was observed in M. domestica (Figure 2C). Wings of C. capitata females were slightly longer than male wings (Figure 2A) and in D. melanogaster and M. domestica male wings were longer (Figures 2B,C). The radio-medial crossvein (r-m) defined by landmarks 4 and 5 was different between males and females in C. capitata (Figure 2A). In D. melanogaster we observed clear sexual differences in the radial vein R2 + 3 defined by landmarks 2 and 9 and in the basal-madial-cubital crossvein (bm-cu) defined by landmarks 6 and 7 ( Figure 2B). In summary, FIGURE 2 | Total and non-allometric sexual shape dimorphism. Summary of discriminant function analyses for total (A-C) and non-allometric (D-F) shape differences between females and males. The wireframes represent differences between female (red) and male (blue) average wing shapes. The scale factor is provided next to the wireframes. The magnitude of sexual shape dimorphism is indicated in units of Procrustes distance with the corresponding p-values based on 1,000 random permutations (***P ≤ 0.001). Histograms with the distributions of the discriminant scores show shape separation into two distinct groups for each species. See also Supplementary Table 5 for details.
Frontiers in Ecology and Evolution | www.frontiersin.org FIGURE 3 | Relationship of the shape component that is most associated with wing size and wing centroid size. Scatterplots showing the relationship of the regression score (i.e., regression of wing centroid size onto Procrustes coordinates pooled by sex, temperature, and density) and wing centroid size (WCS). Linear regression lines were fitted for each category to illustrate the allometric relationships in detail. In the upper row, the data is color coded by the different rearing temperatures and in the lower row, the data is color coded by the different rearing temperatures. Equations and corresponding R 2 values in the plots of the lower row represent pooled regressions for all data combined for each species. f, females; m, male.
we found a clear total sexual shape dimorphism in all three studied species.

Influence of Sexual Size Dimorphisms on Wing Shape
Since sex had a significant influence on wing size in all three species (Table 2, see also Figure 3) (see Siomava et al., 2016 for independent wing size measures) we next asked how much of the shape differences between sexes was explained by differences in wing size. Wing size explained 16% and 14% of wing shape variation in C. capitata and M. domestica, respectively (Figure 3). Wing shape was more clearly associated with differences in wing size in D. melanogaster (34%; Figure 3). To analyze the non-allometric shape differences between sexes in more detail, we accounted for wing size (Klingenberg, 2016; see section "Materials and Methods" for details). A species-specific DFA using the non-allometric shape component revealed significant sexual shape dimorphism for all three species (Figures 2D-F). This observation was confirmed by a Procrustes ANOVA (Supplementary Table 2). In line with the minor impact of wing size on wing shape, the extant of the total compared to the non-allometric sexual shape dimorphism was very similar in C. capitata and M. domestica (Figures 2D,F). Accordingly, the differences in shape decreased after size correction in D. melanogaster ( Figure 2E). Specifically, differences in the length and the width of the wings were largely reduced in the nonallometric shape component. These results suggest that major difference between sexes observed in C. capitata and M. domestica can be explained by the non-allometric shape component, while wing size differences contribute to the sexual shape dimorphism in D. melanogaster.

Effect of Different Rearing Conditions on Wing Shape
To evaluate the effect of wing size on shape in the three species in more detail, we raised flies at different temperatures and densities. The Procrustes ANOVA for all three species revealed significant interactions between species and density and temperature, respectively ( Table 1), suggesting speciesspecific effects of the rearing conditions on wing shape. Density had significant effects on wing size in all three species, while temperature affected only D. melanogaster and C. capitata ( Table 2). Since temperature and density contributed additively to wing shape differences (Supplementary Table 2), we therefore split the analyses by species and performed species-specific DFA to study the total shape differences caused by different rearing temperatures (Figures 4A-C) and densities (Figure 5A), respectively. Additionally, we calculated the non-allometric component of shape (see section "Materials and Methods" for details) and analyzed the differences using DFA (Figures 4D-F for temperature and Figure 5B for density).
In accordance with the species-specific Procrustes ANOVA testing for effects of rearing conditions on total wing shape variation (Supplementary Table 2), the DFA clearly assigned wings to one of the two rearing temperatures (Figures 4A-C). The most obvious effect of different rearing temperatures on wing shape was observed in D. melanogaster ( Figure 4B). Wings of D. melanogaster flies raised at higher temperatures were wider in the distal-central region defined by landmarks 7-9 and distally shortened (i.e., displacement of landmark 10) ( Figure 4B). The distal contraction remained after size correction, while the width was much less affected ( Figure 4E). This observation suggests that temperature-dependent plasticity in wing width is predominantly caused by differences in wing size. Different rearing temperatures also affected distal-central wing width in C. capitata ( Figure 4A) and M. domestica ( Figure 4C). Narrower wings were observed at higher temperatures in M. domestica, while wings were narrower at lower temperatures in C. capitata (Figures 4A,C). M. domestica wings were longer at higher rearing temperatures (i.e., displacement of landmark 10) ( Figure 4C). Only minor changes were present between total and nonallometric shape differences in C. capitata, suggesting that temperature-dependent size differences had small effects on wing shape in this species. In line with non-significant effects of temperature on wing size ( Table 2), we observed no differences between the allometric and non-allometric shape in M. domestica (Figures 4D,F), In all three species, we observed temperaturedependent plasticity in the placement of the radio-medial (r-m) crossveins (defined by landmarks 4 and 5), the basal-medialcubital (bm-cu) crossveins (defined by landmarks 6 and 7), the radial vein R2 + 3 (defined by landmarks 2 and 9) and the anterior cubital (CuA1) veins (defined by landmarks 7 and 8) (Figure 4). Table 2), the DFA showed that different rearing densities only significantly affected total wing shape in D. melanogaster ( Figure 5A and Supplementary Table 5). Wings of flies raised at high densities were wider in the central-distal region (defined by landmarks 7-9) and elongated (i.e., displaced landmark 10) ( Figure 5A). We observed plasticity in the placement of the basalmedial-cubital (bm-cu) crossveins (defined by landmarks 6 and 7), the anterior cubital (CuA1) veins (defined by landmarks 7 and 8), the radial vein R2 + 3 (defined by landmarks 2 and 9) and the costal (C) vein (defined by landmarks 3 and 9) ( Figure 5A). All these differences were gone after size correction (Figure 5B), suggesting that wing shape differences in response to rearing densities were predominantly associated with variation in size. The displacement of landmark 3 was consistent in both analyses (compare Figures 5A,B), implying that this region of the wing may be affected by different rearing densities independent of wing size. Despite a weak effect of density on C. capitata wing shape in our Procrustes ANOVA (Supplementary Table 2), the DFA did not reveal significant shape differences for flies raised at high and low densities (Supplementary Table 5). The same was true for M. domestica, supporting the Procrustes ANOVA results (Supplementary Table 2). Note that we observed significant differences in wing shape at different rearing densities after size correction in C. capitata and M. domestica in our Procrustes ANOVA (Supplementary Table 2) and in the DFA (Supplementary Table 5).

As suggested by the Procrustes ANOVA (Supplementary
In summary, we showed that different rearing temperatures had a profound effect on wing shape D. melanogaster and minor effects were observed in C. capitata and M. domestica. Variation in rearing densities only clearly affected wing shape in D. melanogaster.

Sexual Dimorphism in Response to Rearing Conditions
Since the non-allometric component of shape was consistently affected by different rearing conditions (Supplementary  Tables 2, 5), we next asked for each species, whether a sexual dimorphism in response to different rearing temperatures and densities exists. In all three species, we did not observe significant interactions between sex and rearing conditions in our Procrustes ANOVA of size corrected shape variables (Supplementary Table 2). Accordingly, in C. capitata and D. melanogaster a DFA revealed very similar non-allometric shape differences between high and low temperatures (Supplementary  Table 2; " * _non-allometric"-sheets). Interestingly, the DFA revealed a significant effect of temperature (Figures 6A,B FIGURE 4 | Total and non-allometric shape differences due to different rearing temperatures. Summary of discriminant function analyses for total (A-C) and non-allometric (D-F) shape differences for flies raised at different temperatures. The wireframes represent differences between the low temperature (turquoise) and high temperature (orange) average wing shapes. The scale factor is provided next to the wireframes. The magnitude of shape variation is indicated in units of Procrustes distance with the corresponding p-values based on 1,000 random permutations (*P ≤ 0.05 and ***P ≤ 0.001). Histograms with the distributions of the discriminant scores show shape separation into two distinct groups for each species. See also Supplementary Table 5 for details. FIGURE 5 | Total and non-allometric shape differences in D. melanogaster due to different rearing densities. Summary of discriminant function analyses for total (A) and non-allometric (B) shape differences for flies raised at different densities. The wireframes represent differences between the low density (green) and high density (purple) average wing shapes. The scale factor is provided next to the wireframes. The magnitude of shape variation is indicated in units of Procrustes distance with the corresponding p-values based on 1,000 random permutations (***P ≤ 0.001). Histograms with the distributions of the discriminant scores show shape separation into two distinct groups for each species. See also Supplementary  Table 5) on non-allometric shape differences in males, while females were unaffected.
In summary, we found no evidence for a sexual dimorphism in the response to different rearing conditions in C. capitata and D. melanogaster. In contrast, we found indications for a sexual dimorphism in the non-allometric shape differences between extreme rearing temperatures and densities in M. domestica.

DISCUSSION
Wing morphology is thought to be a highly adaptive trait because wings facilitate basic tasks, such as finding food and mating partners (e.g., Wootton, 1992;Le Roy et al., 2019). In some insects, wings are engaged in mating behaviors. For instance, males of D. melanogaster and C. capitata generate species-specific courtship songs by fast and repetitive wing movements and females judge the size and vigor of a potential mating partner by the intensity of this buzzing (Burk and Webb, 1983;Sivinski et al., 1984;Webb et al., 1984;Churchill-Stanland et al., 1986;Partridge et al., 1987). In other fly species, such as M. domestica, the mating process is initiated during flight by an attack of a male against the female (i.e., mating strike). A successful mating strike usually results in the immediate landing and start of copulation (Murvosh et al., 1964). In line with these different behaviors, our interspecific PCA separated C. capitata and D. melanogaster (courtship songs) from M. domestica (mating strike) along PC2 that explained about 20% of the shape variation and captured differences in the ratio of wing width and wing length. Hence, the longer and narrower wings of M. domestica might be under selection for better flight performance (Alves and Bélo, 2002;Shyy et al., 2013), while the shorter and rounder wings of D. melanogaster and C. capitata seem to be better suited to displace more air for repeated buzzing (Burk and Webb, 1983;Webb et al., 1983;Wheeler et al., 1988;Talyn and Dowse, 2004;de Souza et al., 2015). Despite this potential link between interspecific wing shape differences and mating behavior, our data revealed that all three species were equally separated along PC1 that captured wing width along the proximal-distal axis and explained about 77% of the variation. Therefore, different mating behaviors (i.e., courtship songs vs. mating strike) are not the predominant driver for interspecific wing shape differences and many other species-specific adaptations must be at play. For instance, interspecific differences in absolute body size and shape and general flying behaviors influence selection on wing size and shape.
In the following we discuss our major findings on plasticity and sexual dimorphism in wing size and wing shape across three Diptera species with an emphasis on conserved and divergent aspects of variation in wing morphology.

Species-Specific and Conserved Aspects of Plastic Response to Different Rearing Conditions
We raised flies at different temperature and density regimes and analyzed the impact on wing size, wing shape and the relationship of both. Overall, we revealed a few common, but mostly speciesspecific patterns of plastic responses.

Plasticity in Response to Rearing Temperature
The most consistent results were obtained for different rearing temperatures that affected wing size in D. melanogaster and C. capitata and wing shape in all three species. The temperature impact on wing size is consistent with previous observations FIGURE 6 | Sex-specific non-allometric shape differences due to different rearing temperatures and densities in M. domestica. Summary of discriminant function analyses for non-allometric shape differences for flies raised at different temperatures (A,B) or densities (C,D). The wireframes represent differences between the low temperature (turquoise) and high temperature (orange) (A,B) and the low density (green) and high density (purple) (C,D) average wing shapes, respectively. The scale factor is provided next to the wireframes. The magnitude of shape variation is indicated in units of Procrustes distance with the corresponding p-values based on 1,000 random permutations ( ns P > 0.05; **P ≤ 0.01; and ***P ≤ 0.001). Histograms with the distributions of the discriminant scores show shape separation into two distinct groups for each species. See also Supplementary Table 5 for details. in these two species (Cavicchi et al., 1985;Partridge et al., 1994a;James et al., 1995;van't Land et al., 1999;Zwaan et al., 2000;Navarro-Campos et al., 2011) and the general trend that insects develop larger body sizes and wings in at lower temperatures (Bergmann, 1847;Ray, 1960). A previous study of natural M. domestica populations in South America found a positive correlation between latitude and wing size in females (Alves and Bélo, 2002). Since such latitudinal clines are often related to variable temperatures, the complete lack of association between wing size and rearing temperature in our work is a bit surprising. A potential explanation could be that the higher temperature (i.e., room temperature of 22-25 • C) in our experiment could be lower than optimal for the Italian M. domestica strain used here (Siomava et al., 2016) and we may not have covered the entire range of optimal rearing temperatures. However, the discrepancy could also be a real biological pattern because two very different M. domestica populations are being compared. Support for this notion comes from an analysis of thermal plasticity in the two sepsid fly species Sepsis fulgens (Meigen, 1826) and Sepsis punctum (Fabricius, 1794) which share comparable geographic distributions and ecological niches, but still exhibited opposite clinal patterns of wing size (Rohner et al., 2019). Thermal plasticity in wing shape was most obvious in wing width since wings of flies raised at higher temperatures were generally wider in D. melanogaster and C. capitata and narrower in M. domestica. We detected a slightly stronger plasticity in proximal landmarks in D. melanogaster, which is in contrast to previous data that revealed a strong response of the distal wing region to different temperatures (Debat et al., 2003). This discrepancy could partially be due to the more moderate temperature range used in our study (18-25 • C) because the lower temperature of 12-14 • C used by Debat et al. (2003) is already stressful as for instance wing size variation increases at 13 • C (Bubliy and Loeschcke, 2002).

Plasticity in Response to Rearing Density
Besides temperature, environmental factors related to nutrition levels and/or quality influence organ size and shape. For instance, in the wasp Ephedrus persicae (Froggatt, 1904) that parasitizes aphids, wing size and shape depends on the type of host (Bogdanović et al., 2009). Different density regimes, inducing competition for nutrition levels, clearly affected wing size in all three studied species, a common trend previously observed in D. melanogaster (Bitner-Mathé and Klaczko, 1999;Siomava et al., 2016) and in C. capitata and M. domestica (Siomava et al., 2016). Data from C. capitata natural strains strongly suggests that these laboratory experiments elicit ecologically relevant plasticity because body and wing size in C. capitata depends on the sugar content of the host fruits on which larvae develop (Krainacker et al., 1987;Navarro-Campos et al., 2011). An association of wing shape variation with the host fruit has also been observed in this species (Pieterse et al., 2017). Since different rearing densities in our experiment reflect variable nutrient quantities with constant composition and we did not observe density related wing shape effects in C. capitata, the shape variation induced by different host fruits may be the result of nutrient quality rather than quantity. In contrast, the plasticity in wing shape in D. melanogaster observed in this work is most likely a result of nutrient quantity. It will be interesting in future controlled laboratory experiments to test how nutrient quality and quantity affect wing size and shape.

Species-Specific Scaling Relationships
We observed the most obvious differences among species with respect to allometric shape differences. Overall, we could assign about 34% of shape variation in D. melanogaster to variation in wing size. This trend went down to 16% and 13% in C. capitata and M. domestica, respectively. Intriguingly, our data showed that scaling relationships were uniform in D. melanogaster across all groups (i.e., sex, temperature and density), while very complex scaling relationships, partially with significant interaction terms, were observed for the other two species (see Figure 3). These complex scaling relationships may have influenced our regression approach to partition total shape variation into an allometric and a non-allometric component (Gidaszewski et al., 2009;Klingenberg, 2016). For example, in some cases we found significant differences in the nonallometric shape component despite insignificant total shape differences (e.g., compare Figures 4A,B for C. capitata). While this rather counterintuitive finding could be an artifact, we confirmed that the size correction indeed removed the effect of centroid size on wing shape variation. Moreover, a study in the Leporinus cylindriformis (Borodin, 1929) group of rayfinned fishes revealed that the distinction of species based on geometric morphometrics data was facilitated by the inclusion of size correction (Sidlauskas et al., 2011), suggesting that indeed relevant aspects of shape variation could be recovered. To clarify this aspect further, more accurate estimates of withingroup scaling relationships are needed because the number of flies for each subgroup (i.e., a certain sex at a specific rearing temperature and density, respectively) was rather low in our nested experimental design.
The different scaling relationships among species could be caused by complex interactions between temperature and density regimes applied during our experiments as artificial thermal selection experiments in D. melanogaster revealed effects on resource intake and management and thus the ability to accept different resource levels (Bochdanovits and de Jong, 2003). For instance, while cold adapted flies process food resources more efficiently resulting in better larval survival when raised in low density regimes, they show lower larval survival at higher larval density (Partridge et al., 1994b. Moreover, the rearing regimes applied for C. capitata and M. domestica may have been more in the stressful range compared to D. melanogaster. In the latter species, increased variation in wing size has been observed at higher density (Imasheva and Bubliy, 2003) and at stressful low temperatures (13 • C) (Bubliy and Loeschcke, 2002), showing that environmental factors do not only affect morphology of the individual, but also the level of interindividual phenotypic variability. Genetic factors may also play a role in defining the species-specific scaling relationships observed in our data as the interindividual genetic diversity could be different between the three species studied here. We used a highly inbred laboratory strain of D. melanogaster (>20 years inbreeding) and more recently collected wildtype strains of C. capitata and M. domestica, respectively. For instance, the Italian M. domestica strain used for our work was established as laboratory culture only 1-2 years prior to the experiments (Siomava et al., 2016). Since mutations affecting thermal plasticity in wing development have been identified in D. melanogaster (Debat et al., 2009), it is likely that more such mutations segregate in the more recently collected and thus less inbred species. However, more "plasticity loci" segregating in C. capitata and M. domestica cannot fully explain our observations, since this hypothesis implies that these two species should show the most pronounced plastic responses. In contrast, we found the most obvious responses of wing size and shape to different rearing conditions in D. melanogaster. Therefore, the species-specific scaling relationships must be caused by additional processes. While it is generally accepted that variation in organ size entails shape changes (Debat et al., 2003;Klingenberg, 2016), the level of this coupling could be variable among species. Our D. melanogaster data confirms previous findings in this species that processes regulating wing size and shape are highly coordinated during development (Day and Lawrence, 2000;Matamoro-Vidal et al., 2015;Testa and Dworkin, 2016), while partially different processes could be at play in the other two species.
An interesting approach to learn more about the developmental basis of coordinated size and shape regulation may be to study the cellular basis of the observed plasticity. Data in D. melanogaster (Zwaan et al., 2000;Klepsatel et al., 2014) and D. subobscura (Calboli et al., 2003) revealed that different clinal patterns of wing size variation can be caused by differences in cell number or cell size, while wing size variation induced by different temperature conditions seems to be predominantly driven by variation in cell size (reviewed in Arendt, 2007). Artificial selection for different wing shapes (i.e., rounded vs. elongated wings) in D. melanogaster resulted in flies with different wing shapes, but comparable wing sizes (Menezes et al., 2013;Torquato et al., 2014). Rounder wings have more cells than elongated wings, while the cell size does not seem to contribute to the differences (Torquato et al., 2014). Intriguingly, complex interactions between environmental factors have been observed in D. melanogaster since wing length and cell size reacts more to changes in rearing temperature under limited food conditions, while cell number seems to be unaffected (de Moed et al., 1997). Therefore, the cellular basis of size and shape differences can be differentially shaped by various environmental factors. Overall, based on our species comparison we propose that the general capacity for plastic development and the processes affecting scaling relationships can be uncoupled and evolve independently. Future experiments testing whether plasticity in wing morphology in different insects is caused by variation in cell size or cell number will help formulating thorough developmental hypotheses.

Conserved Aspects of Plastic Shape Variation
Our shape analysis revealed a high level of variation in the positioning of the r-m (landmarks 4 and 5), R2 + 3 (landmarks 2 and 9), CuA1 (landmarks 7 and 8), and bm-cu (landmarks 6 and 7) veins in all treatment groups in all three species. These distal and posterior wing regions were also shown to be most variable across 25 Drosophila species (Houle et al., 2003) and between inbred strains from one D. melanogaster population (Matamoro-Vidal et al., 2015;Pitchers et al., 2019). Interestingly, our data showed that variation in these wing regions was predominantly caused by variation in wing size. Similar results were obtained for different Drosophila species as well (Cavicchi et al., 1991;Bitner-Mathé et al., 1995;Baylac and Penin, 1998;Haas and Tolley, 1998), suggesting that the development of these veins may be a hotspot for the integration of wing size variation (Cavicchi et al., 1985;Guerra et al., 1997;Pezzoli et al., 1997;Baylac and Penin, 1998). Moreover, the recurrent variation in similar wing regions implies that the development of those wing veins may be highly coordinated. It will therefore be interesting to study the cellular and molecular mechanisms underlying plasticity in wing morphology in different species to test whether the same or different developmental mechanisms are responsible for the observed patterns.
A correlation between male wing shape and mating success has been established in C. capitata because males with wider wings copulated more often successfully (de Souza et al., 2015). Interestingly, males with elongated wings have a higher mating success in D. melanogaster (Menezes et al., 2013), showing that intra-sex variation in wing shape are caused by species-specific processes. Indeed, both species differ in the frequency and quality of their courtship song (Cowling and Burnet, 1981;Burk and Webb, 1983;Briceño et al., 2002;Rybak et al., 2002). In contrast to C. capitata and D. melanogaster, we found the most obvious sexual differences in M. domestica in wing length. The slightly elongated wings in males may increase flying efficiency and could be linked to the mating strike (Murvosh et al., 1964;Alves and Bélo, 2002).
Among the three studied species we observed different effects of wing size on sex-specific wing shape. Exclusion of the allometric coefficient clearly decreased the sexual shape dimorphism in D. melanogaster, suggesting that most of the observed shape differences could be explained by differences in wing size. In contrast, sex had only a minor effect on wing size in C. capitata and M. domestica (Siomava et al., 2016) and accordingly, the impact of the allometric component on wing shape was weak. Species-specific sexual dimorphism in scaling relationships were also observed in five schizophoran Diptera (Rohner, 2020) and in Sphingidae moths, where up to 60% of the shape variation could be explained by wing size differences (de Camargo et al., 2015). An extreme example has been observed in the aphid parasitoid Ephedrus persicae, where male wings are twice as large as female wings although females have generally larger body sizes. In this species, only 5% of the sexual shape differences can be explained by wing size (Bogdanović et al., 2009), suggesting that exaggerated male wing size may be more advantageous than sex-specific wing shape. All this data implies that sexual dimorphisms in wing shape scaling may reflect species-specific adaptations rather than conserved patterns across insects.
Since sexually dimorphic organs have been shown to grow disproportionally with body size in one sex in response to environmental cues (Teder and Tammaru, 2005;Bonduriansky, 2007;Lavine et al., 2015;Siomava et al., 2016), we tested for sexually dimorphic patterns of plasticity. We did not find evidence for a sexual dimorphism in the response to different rearing temperatures and densities in all three studied species. It is important to note that our experimental design is particularly limited with respect to this question because the number of individuals in each relevant group (i.e., split by sex, density and temperature) is very low and we might not have the statistical power to detect any general trends. However, our DFA implies that non-allometric shape differences between extreme rearing temperatures and densities observed in M. domestica may be restricted to males. This trend contradicts a previous observation in this species that showed a significant positive correlation between wing size and wing width with latitude specifically in females of different populations from South America (Alves and Bélo, 2002). This discrepancy could be explained by the use of an Italian strain in our experiments. Future experiments addressing these specific questions should be performed to obtain more conclusive results.

CONCLUSION
In summary, our morphometric analysis of wing morphology in three Diptera species revealed various species-specific aspects of plasticity and sexual dimorphism in wing size and shape. Such analyses allow identifying interesting species to further study the developmental and genetic underpinnings of the observed wing morphology and link these to species-specific ecological and behavioral adaptations.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material. Raw wing images, landmark files, R scripts, and MorphoJ project files are available in an online repository (https://doi.org/10.25625/BXZHAF).

AUTHOR CONTRIBUTIONS
MR, NS, EW, and NP conceptualized the research. NS performed all fly experiments. MR and NS conducted data analysis and data visualization. NS, NP, and MR wrote the original draft. NP administered the project. All authors revised the manuscript.

FUNDING
This work has been funded by a German Academic Exchange Service (DAAD) fellowship number 91540915 to NS, Göttingen Graduate School for Neurosciences, Biophysics, and Molecular Biosciences (GGNB), and Volkswagen Foundation (project number: 85 983 to NP).