Symmetric and Asymmetric Components of Shape Variation in the Diatom Genus Frustulia (Bacillariophyta)

: Irregularities in cell division can produce asymmetry in symmetric structures, such as outlines of diatom cells, which can reﬂect genetic, environmental, or random variability in developmental processes. This study examined 12 phylogenetic lineages of the diatom genus Frustulia using landmark-based geometric morphometrics to assess the variation between cell segments separated by apical and transapical axes. Although asymmetric variation within cells di ﬀ ered in some lineages, these irregularities most likely did not reﬂect the evolutionary history of the lineages. The intraclonal phenotypic plasticity of diatom frustules was induced rather by nongenetic factors, i.e., inherited valve abnormalities, constraints of siliceous cell walls, and random developmental instability during morphogenesis. The positive correlations between the symmetric and asymmetric components of shape variation indicated that the morphogenesis of diatom cell walls a ﬀ ected irregularities within cells and variability among the symmetrized cells to a similar extent.


Introduction
Diatoms are eukaryotic microorganisms with bipartite siliceous cell walls that form shells (frustules) resembling a Petri dish [1]. Each part (the larger epitheca and smaller hypotheca) consists of a large end piece (valve) and a number of smaller strips (girdle bands). The common occurrence of diatoms in heterogenous habitats and their morphological variability has led to the description of more than 8000 species [2]. Nevertheless, the vast majority of the species diversity remains most likely undiscovered or hidden within morphologically defined species (i.e., morphospecies), as has been shown by DNA sequencing, breeding experiments, ecophysiological responses of populations, and morphometric techniques (for references, see [3]). Diatom classification relies principally on the overall shape, symmetry of the cell, and wall ornamentation. Similar shapes or symmetries may, however, undergo convergent and parallel evolution. Phylogenetic reconstructions [4][5][6], as well as ontogenetic studies [7][8][9], have been shown to be useful for the identification of homologies and for understanding the evolution of diatom morphology.
Diatoms are unique in their allometric shape changes during ontogeny [10][11][12]. Considerable miniaturization of cells usually occurs during the vegetative phase, and the maximum cell size is restored by sexual reproduction. Allometric shape changes are usually caused by the mechanism of cell wall formation: each of the two daughter cells inherits one theca from the parent cell, which becomes an epitheca, and a hypotheca is newly formed within the confines of the slightly larger parental theca. Maximum size is usually generated from a zygote (auxospore) de novo [1,12]. Although shape changes during ontogeny may obscure interspecific differences, species-specific allometric trajectories may be recognized and quantified using morphometric methods [13][14][15]. Cell division can produce deviations in symmetric structures such as valve outlines that may reflect genetic, environmental, or random variability in developmental processes [16][17][18]. In diatoms, pronounced intraspecific phenotypic plasticity may arise by the transmission and subsequent amplification of tiny deviations in the course of the vegetative part of the life cycle [19]. Toxic substances, e.g., heavy metals, may induce apparent valve deformities across species that have been used as a proxy for bioindications of environmental conditions [20][21][22].
The genus Frustulia Rabenhorst occurs worldwide in freshwater or brackish benthic habitats such as mires, littoral zones of lakes, and running waters [23,24]. A characteristic morphological feature of this raphid pennate diatom is a specific structure running centrally along the apical axis (longitudinal rib) and a linear-lanceolate-rhomboidal valve shape with rounded, rostrate, or capitate apices ( Figure 1, [1,23]). Morphometric and phylogenetic studies showed that not all Frustulia morphospecies can be unambiguously identified: some morphospecies show overlapping morphological variability, and some phylogenetic lineages do not share unique morphology [25][26][27]. Despite the symmetric appearance of the valves, the internal valve structures of several Frustulia lineages can be asymmetric (e.g., longitudinal rib in Figure 1a,b). We aimed to investigate whether tiny intraclonal irregularities in valve outlines reflect evolutionary changes in morphogenesis. Thus, we analyzed asymmetric and symmetric components of shape variation of 12 Frustulia lineages (which most likely represent evolutionary species) to examine the relationships between phylogeny and (a)symmetry. We aimed to investigate (1) whether clades or lineages differ in patterns of (a)symmetry; and (2) whether patterns of qualitative (a)symmetric changes differ among clades or lineages. In our study, asymmetric shape components of variation refer to intraclonal irregularities within cells, i.e., deviations of cell parts (separated by apical and/or transapical axes of symmetry) from a symmetric shape. Symmetric components of shape variation refer to intraclonal differences among cells that have been symmetrized along both axes of symmetry prior to quantification, i.e., allometry can be especially important in this component of variation. lineages. In our study, asymmetric shape components of variation refer to intraclonal irregularities within cells, i.e., deviations of cell parts (separated by apical and/or transapical axes of symmetry) from a symmetric shape. Symmetric components of shape variation refer to intraclonal differences among cells that have been symmetrized along both axes of symmetry prior to quantification, i.e., allometry can be especially important in this component of variation.

Characteristics of Strains and Data Acquisition
We examined 15 monoclonal strains representing 12 phylogenetic lineages of the genus Frustulia (three lineages were represented by two strains and nine by one strain). Since not all strains were represented by a unique morphotype, their names are composed of two species names (F. crassinervia-saxonica: F. crassinervia (Brebisson) Lange-Bertalot & Krammer, and F. saxonica Rabenhorts) with lineage attribute (III-VII) sensu [25,26]. All of the strains, which were mounted onto permanent slides using a Naphrax (Brunel Microscopes Ltd., Wiltshire, UK), originated from our previous research [25,28]. The treatment of strains was similar: slightly acidic, oligotrophic liquid medium (OGM) enriched with sodium metasilicate; 18 • C temperature and continuous lighting; harvesting of cells in the exponential phase of the growth; and cleaning of frustules by the incineration method (for details of the culture conditions, see [25,26]). We selected lineages on the basis of three criteria: a monophyletic origin [27][28][29], availability of more than 50 frustules or valves per strain, and a cell size above the minimum of their known size range (Table S1). We did not include F. cf. magaliesmontana Cholnoky in our study because two halves of the frustules were overlapping and it was not possible to distinguish separate valve outlines. We photographed strains using an Olympus BX51 (UPlan FLN 100×/1.30 oil objective, differential interference contrast imaging) with Olympus Z5060 equipment (Tokyo, Japan). We unified the orientation of cells and cut off the redundant background using Adobe Photoshop 2020 (Adobe Systems Incorporated, San Jose, CA, USA).
We digitized 50 images per strain using tpsUtil ver. 1.6 [30]. Unfortunately, an automatic tool such as SHERPA (SHapE Recognition, Processing and Analysis) software [31] could not be used for the detection of cell outlines as valves that remained within a frustule could not be focused in one plane. Thus, an automatic tool would most likely detect as an outline a mixture of shadows, overlapping valves, and parts of girdles. We registered the shapes by 96 landmarks using tpsDig ver. 2.3 [30] in clockwise (CW) and counter-clockwise (CCW) orientations. We placed two landmarks at the apices and 47 semi-landmarks along each half of the valve ( Figure 1). Originally, we planned to put four landmarks in fixed positions designating symmetric quadrants (divided by two axes of reflection symmetry). However, the determination of landmarks at the intersection of the cell outline with the transapical axis could not be unambiguously made, even with the help of the central constriction of the longitudinal rib. Although some species (e.g., F. erifuga Lange-Bertalot & Krammer, F. paulii Kilroy & Urbánková) have asymmetric internal structures, and one lineage has curved apices (F. curvata Kulichová & Urbánková), we considered all strains to be biradially symmetric in order to focus on interstrain comparisons.

Geometric Morphometric Analyses
We followed a process of landmark-based morphometrics published elsewhere [16,18,32]. Briefly, CW and CCW files of individual strains were used to test for measurement errors due to imaging and digitizing. In subsequent analyses, averaged CW and CCW coordinates were used for both the calculation of centroid sizes (CS) [33] and for the general Procrustes superimposition [34]. During Procrustes analysis, the positions of semi-landmarks were not allowed to slide along segment outlines, but their positions were determined by equal Euclidean distances within each specimen [35,36]. Symmetric and asymmetric components of shape variation were decomposed on the basis of tangent Procrustes distances among the perfectly symmetrical cell shape and reflected/relabeled configurations of a cell along the apical (horizontal), transapical (vertical), and transversal axes [37,38]. These configurations were used also for relative warp analyses (RWAs) that enabled the reconstruction of (a)symmetric patterns in shape changes. Asymmetric variation between pairs of quadrants (cymbelloid, gomphonemoid, and sigmoid symmetry) was quantified within each cell. We express levels of absolute asymmetry as total asymmetry (resulting from the sum of asymmetries) and median asymmetry. We named asymmetric components according to the symmetry of the diatom genera: cymbelloid (Cymbella Agardh, asymmetric along the apical axis and symmetric along the transapical axis, dorsiventral shape), gomphonemoid (Gomphonema Ehrenberg, symmetric along the apical axis and asymmetric along the transapical axis, heteropolar shape), and sigmoid (Pleurosigma W. Smith or Gyrosigma Hassall, asymmetric along both axes, S-like shape). The morphometric data were analyzed using tpsRegr ver. 1.46 [30] and package geomorph ver. 3.0.4 [39] implemented in R ver. 3.4.3 [40].
We determined descriptive statistics and tested significance between (a)symmetric components of shape variation by Spearman's rank correlation tests and Mann-Whitney pairwise comparisons (both with Bonferroni corrections) using PAST ver. 3.20 [41]. We also compared the relationships between (a)symmetric components of shape variation and CS.

Reconstruction of Phylogeny and Trait Mapping
The phylogeny of the genus Frustulia was computed based on four molecular markers: nuclear small subunit (SSU) ribosomal DNA (rDNA), nuclear large subunit (LSU) rDNA, the chloroplast-encoded large subunit of Rubisco gene (rbcL) and the chloroplast-encoded photosystem II protein D1 gene (psbA). We used three species as an outgroup for the genus, i.e., Amphipleura pellucida (Kützing) Kützing, Berkeleya rutilans (Trentepohl ex Roth) Grunow, and Climaconeis riddlae A.K.S.A. Prasad. All available sequences were downloaded from the GenBank database. Sequences of psbA were sequenced following reference [26] using primers psbA-F and psbA-R1 [42] with the MyTaq™ polymerase and the following PCR settings: initial denaturation at 94 • C for 1 min; 35 cycles of 94 • C for 15 s, 55 • C for 15 s, and 74 • C for 10 s; and final extension at 72 • C for 7 min. To reduce missing data in the alignment, we imputed missing sequences from strains which belonged to the same species based on the more divergent markers (rbcL or LSU). The sequence details are presented in Table S2. The newly obtained DNA sequences have the following GenBank accession numbers: MW039451-MW039479.
Trait evolution was examined in R. Species that were not used in the morphometric analysis were pruned from the tree (package ape, [51]). We used different approaches to examine phylogenetic signals in the shape components. First, we computed Pagel's λ (package phytools, [52]). The ability to estimate Pagel's λ from the dataset was assessed using the phylogenetic Monte Carlo (pmc) package [53]. Analyses were done with 1000 replicates. Second, two models of trait evolution were compared using the pmc package: a white noise (WN) model represented the null model, i.e., random distribution of the traits, and a Brownian motion (BM) model represented the phylogenetic signal in our data.

Quantification of Components of Shape Variation
The largest differences within strains were affected by symmetric differences between cells (total: 0.008-0.027) whereas asymmetric components of shape variation within cells were much less important (Table 1, Figure S1). Measurement errors due to imaging and digitizing were within the range 3-12% for all strains except F. septentrionalis Lange-Bertalot & Metzeltin, for which they reached 21% (Table 1). Among the asymmetric components, more pronounced appeared to be the cymbelloid shape (total: 0.007-0.012; median: 0.004; Figure 2). The gomphonemoid and sigmoid shapes had comparable absolute values (total: 0.004-0.011; median: 0.003). Significantly higher values of cymbelloid shape components in comparison with gomphonemoid and/or sigmoid components were found only within five lineages (F. crassinervia-saxonica lineages III, IV, VII; F. gondwana Lange-Bertalot & Beier, and F. septentrionalis; Table S3). The values quantifying (a)symmetric components of shape variation showed significant positive intercorrelations (Figure 3), except for two insignificant relationships for median values of the sigmoid shape (Table 2). Cell sizes were negatively correlated with values of both total and median symmetry; the correlations between CS and asymmetric components were not significant ( Table 2).

Reconstruction of Shape Changes
RWAs, which were based on one symmetric and three reflected/relabeled original configurations, described in the first 10 axes 94-98% of shape variability in strains. The first axis (46-76% of shape variability) was associated with symmetric changes causing widening or narrowing of cells along the whole cell outline except for the apical part (Figure 4a,b; Figure S1). The second symmetric axis described antagonistic changes in the central part of cell outlines and parts below the apices, i.e., causing an oval shape or lanceolate shape with/without rostrate apices (Figure 4c,d). The first axis associated with cymbelloid shape described both the deflection of the apices in the same direction and the differentiation of the outline to the dorsal (more convex) and ventral (less convex) sides along the apical axis (Figure 4e). The first axis associated with gomphonemoid shape was consistent among all strains and caused antagonistic changes along the transapical axis that do not affect the central part of the outlines and apices (Figure 4f). The first sigmoid axis described transversal shape changes that affect the whole cell outline (Figure 4g, Figure S1).  Figure S2).

Reconstruction of Shape Changes
RWAs, which were based on one symmetric and three reflected/relabeled original configurations, described in the first 10 axes 94-98% of shape variability in strains. The first axis (46-76% of shape variability) was associated with symmetric changes causing widening or narrowing of cells along the whole cell outline except for the apical part (Figure 4a,b; Figure S1). The second symmetric axis described antagonistic changes in the central part of cell outlines and parts below the apices, i.e., causing an oval shape or lanceolate shape with/without rostrate apices (Figure 4c,d). The first axis associated with cymbelloid shape described both the deflection of the apices in the same direction and the differentiation of the outline to the dorsal (more convex) and ventral (less convex) sides along the apical axis (Figure 4e). The first axis associated with gomphonemoid shape was consistent among all strains and caused antagonistic changes along the transapical axis that do not affect the central part of the outlines and apices (Figure 4f). The first sigmoid axis described transversal shape changes that affect the whole cell outline (Figure 4g, Figure S1).

Reconstruction of Shape Changes
RWAs, which were based on one symmetric and three reflected/relabeled original configurations, described in the first 10 axes 94-98% of shape variability in strains. The first axis (46-76% of shape variability) was associated with symmetric changes causing widening or narrowing of cells along the whole cell outline except for the apical part (Figure 4a,b; Figure S1). The second symmetric axis described antagonistic changes in the central part of cell outlines and parts below the apices, i.e., causing an oval shape or lanceolate shape with/without rostrate apices (Figure 4c,d). The first axis associated with cymbelloid shape described both the deflection of the apices in the same direction and the differentiation of the outline to the dorsal (more convex) and ventral (less convex) sides along the apical axis (Figure 4e). The first axis associated with gomphonemoid shape was consistent among all strains and caused antagonistic changes along the transapical axis that do not affect the central part of the outlines and apices (Figure 4f). The first sigmoid axis described transversal shape changes that affect the whole cell outline (Figure 4g, Figure S1).

Phylogeny and Shape Components
The strains examined by geometric morphometrics were separated into four supported clades: (1) F. crassinervia-saxonica lineages IV, V, VI, VII; (2) F. crassinervia-saxonica lineage III, F. gaertnerae, F. septentrionalis; (3) F. maoriana; and (4) F. curvata, F. erifuga, F. gondwana ( Figure S2). The phylogenetic position of F. paulii was not statistically supported. Computation of Pagel's λ suggested a lack of phylogenetic signal in (a)symmetric shape components (Table S4). However, permutation analysis showed that our dataset was too small to estimate this parameter. The median for all analyses was 0 and all estimates had wide confidence intervals (Table S4). Permutation analysis comparing the WN and BM models of evolution showed that our data were not sufficient for differentiation of the examined models, except for two cases where the distribution of total and median values of cymbelloid symmetry conformed to the WN, i.e., the distribution of cymbelloid symmetry was random in relation to phylogeny.

Shape Variation within Lineages
Our study showed that symmetric variation among cells within Frustulia lineages was the most pronounced component of shape variation (Table 1; Figure S1), i.e., differences between symmetric valves (symmetric along both the apical and transapical axes) were more important than variation between valve segments. In agreement with diatom studies that examined allometric shape changes during the life cycle [14,16], symmetric changes of Frustulia strains along the first RW axis were associated with widening or narrowing of frustules (Figure 4, Figure S1). Therefore, allometric variability among valves was dominant in our strains, even though the size range of strains represented only part of the life-cycle variability (Table S1).
The levels of asymmetry varied among cells (Figure 2), differing significantly between or within a few lineages (Table S3). Nongenetic variation in cell shapes could be affected artificially by imprecision in shape recognition, as valves were not always laying in one plane of focus, and valve outlines were sometimes distorted by girdle bands, or by culture conditions that may affect strain morphology differently (see [54] for references). Unfortunately, our collection of permanent slides did not allow us to purposefully minimize the effect of the abovementioned factors. We reduced the effect of ontogenetic shape changes by selecting strains that were both above their minimum size (fewer cell divisions after sexual reproduction) and exhibiting limited variability in size (Table S1). Insignificant correlations between asymmetric components and size-shape trajectories in several biradially symmetric strains of different genera indicate that at least the asymmetric variation is not profoundly influenced by ontogenetic changes [16,55]. Comparisons of strain pairs from three lineages indicated similar patterns in morphological (a)symmetry within lineages/species (Table 1, Figure 2). F. maoriana Lange-Bertalot & Beier strains shared increased values of all components of (a)symmetry, F. crassinervia-saxonica III shared increased levels of cymbelloid shape, and F. erifuga strains had consistently low levels of (a)symmetry. Furthermore, asymmetric components were not significantly different between strains representing the same lineage, except for the cymbelloid shape component within strain F350 (Table S3). Based on the comparisons of pairs of strains, we do not consider (a)symmetric variation to be strain specific.

Phylogeny and Asymmetric Variation
We expected to find increased asymmetry within a lineage and a clade with asymmetric internal structures; however, this was not confirmed by our data. The lineage Frustulia paulii ( Figure S2), which has an asymmetric hyaline area in the central part of the valves (Figure 1e; for SEM images, see [56]), showed similar variation in all three components of asymmetry ( Figure 2) that did not differ significantly from all other strains (Table S3). The clade represented by F. curvata, F. gondwana, and F. erifuga ( Figure S2), which share deflection of the longitudinal ribs in one direction (Figure 1), was not consistent in terms of the variation of (a)symmetric components within the clade (Table 1, Figure 2). Among the sequenced morphologically supported lineages, markedly curved ribs are also present in F. cf. krammeri Lange-Bertalot & Metzeltin and F. gibsonea Bouchard, Hamilton, Starr & Savoie, which belong to this clade as well [23], and in distantly related F. vulgaris (Thwaites) De Toni ( Figure S2; [27,29]). It is therefore likely that the shape of the longitudinal ribs is not a homologous character, but it may aid evolutionary diversification. The remaining two clades also showed variable ( Figure S2) patterns of shape variation. One of these clades included a lineage with the lowest (a)symmetric values (F. septentrionalis,  Figure 2). Nevertheless, we cannot strongly support the lack of an evolutionary signal in asymmetric components (Table S4), probably due to the small phylogeny [53].

Common Asymmetric Shape Variation across Lineages
We were surprised by the consistency in qualitative changes in valve outlines across lineages that generally differed only in terms of the magnitude of changes ( Figure S1). We did not find patterns unique for clades, lineages, or morphotypes along the most important (a)symmetric shape components. We expected to find differences between morphotypes that would indicate a common shape of the cell outline. For instance, different types of apical ends (Figure 1) showed similar valve changes ( Figure S1) that do not indicate differences between capitate, protracted, and constricted apices. Even Frustulia curvata with asymmetric apices curved in opposite directions (Figure 1d) did not show unique shape changes and exhibited only different magnitudes of changes along the first two symmetric axes (Figure 1h). The development of curved apices in F. curvata was most likely affected by the mutual influence of sibling cells, which were pressed by apices inside mucilaginous tube colonies [56], and not by strict genetic control. Abnormalities of frustules caused by epigenetic effects tend to the correct shape after perturbation, but the abnormal shape can be inherited and even amplified in the course of mitotic cell divisions [7,19]. It is interesting to note that shape changes along the first of the asymmetric axes resemble existing diatom species of asymmetric genera (Figure 4), whereas the less important axes of asymmetry (lower values of explained shape variability) account for odd diatom shapes.
The changes in cymbelloid shape observed in Frustulia (Figure 4) are not consistent with the directional asymmetry of Luticola species, where the most pronounced changes were observed in the middle part of the valve outline [18,55]. A possible explanation for the increased cymbelloid asymmetry in the middle part of Luticola valves is in the form of an initial cell (the first frustule of maximum size developed from the auxospore) that, in some/all Luticola species, has a markedly swollen outline in the middle part [55,57]. Indeed, the final shape of a frustule is associated with the specific development of the auxospore [7]. Nonrandom qualitative asymmetric variation of the valve shapes of Frustulia and Luticola species most likely indicates canalization of cell wall development in some particular direction [17,58].
Comparisons between asymmetric components showed that the cymbelloid shape is more pronounced within the majority of lineages (Table 1, Figure 2), though higher values were not significantly different in many cases (Table S3). An increased cymbelloid component of asymmetry may be caused by inherent asymmetry in pennate diatoms as valves are formed sequentially along the apical axis: the primary side is produced prior to the secondary side [59]. The sides of the valves are recognizable by irregularities in ornamentation called Voigt discontinuities and by curvature of the apical raphe endings. Extreme shape differences between the primary and secondary sides are notable in the Cymbellales, where the more convex margin can be either the primary or the secondary side of the valve [1,60]. Diatoms that are classified as bilaterally symmetric may also exhibit more or less visible cymbelloid symmetry. For instance, the genus Luticola exhibits cryptic [18,55] as well as apparent dorsiventrality (e.g., [61]), but also balanced levels of asymmetric shape components [16]. Cymbellales share with at least one Luticola species [57] the cis symmetry (i.e., the primary sides of both valves are on the same side [59]), whereas representatives of the genera Frustulia and Amphipleura (F. saxonica, F. vulgaris, A. pellucida; [59,62]) have both cis and trans symmetry (i.e., primary sides of the valves are on the same or on opposite sides of the frustule). Given that all dorsiventral and heteropolar genera share cis symmetry of frustules (sigmoid genera exhibit cis and trans), systematic cymbelloid or gomphonemoid symmetry within Frustulia is unlikely.

Relationship between Asymmetric and Symmetric Variation
We found positive correlations between the total symmetry and total asymmetric components of shape variation (Table 2, Figure 3). These correlations were not influenced by the sizes of strains, as the average centroid sizes were significantly correlated only with the total symmetry ( Table 2). The relationships revealed in our study are in agreement with the results of a study which examined allometric shape changes in two Luticola strains during their life cycle [55]. The study showed that with gradually decreasing size, the morphological disparity of the strains increased, but this variation was not associated with asymmetric variation within cells. Thus, we consider the effect of allometry on asymmetric components of variation within diatom valves to be unimportant.
The finding of positive correlations between the asymmetric and symmetric components may be particularly useful for ecotoxicological studies where frustules are markedly asymmetric. To determine the effect of toxic substances on frustule development, various metrics analyzing variations in valve shape can probably be used without much bias, because variability among symmetrized cells is most likely equivalent to irregularities within diatom cells. Therefore, ecotoxicological studies that do not separate shape variation into individual components and/or do not use the same methodological approach (e.g., [20][21][22]) are most likely comparable in their outcomes.

Conclusions
Our study enables assessment of the possible effects of genetic and nongenetic factors on shape variation of clonal diatom strains. We identified different components (according to [17]) affecting phenotypic plasticity: (1) genetic variation is most likely associated with extreme (a)symmetric values in some lineages and/or with (2) inherited valve abnormalities that were induced by variability in environmental factors; (3) intrinsic organismal factors are associated with constraints in the development of siliceous frustules (including allometric shape changes); and (4) asymmetric variation within cells may reflect random instability of developmental processes. The positive correlations between the symmetric and asymmetric components of shape variation indicated that the morphogenesis of diatom cell walls affected irregularities within the cells and variability among the symmetrized cells to a similar extent.
Supplementary Materials: The following are available online at http://www.mdpi.com/2073-8994/12/10/1626/s1, Figure S1: Reconstructions of shape changes in 15 strains along the first axes of (a)symmetric components. Figure S2: Maximum likelihood tree based on four molecular markers. Supplementary Tables are available online at https://data.mendeley.com/datasets/wrcydst96z/1. Table S1: Size ranges of the strains under study and sizes of morphospecies from publications. Table S2: GenBank accession numbers. Table S3: Matrix of Mann-Whitney pairwise comparisons of all components of asymmetry calculated within the 15 strains. Table S4: Results from the analysis of trait evolution.
Author Contributions: J.K. formulated the aims and hypotheses, and collected and processed material for geometric morphometrics. P.U. collected and processed material for the phylogenetic parts of the article. J.K. and P.U. wrote the article in collaboration. All authors have read and agreed to the published version of the manuscript.

Funding:
The research was financially supported by Charles University.