Gerbera aurantiaca Sch.Bip. (Asteraceae: Mutisieae) is a long-lived perennial herb endemic to the mistbelt grasslands of KwaZulu-Natal and Mpumalanga in eastern South Africa (Johnson et al., 2005). The species is known from approximately 25 scattered, isolated populations ranging in size from a few to hundreds of clones, and is listed as endangered due to habitat transformation resulting from commercial forestry, agriculture and urban development (Raimondo et al., 2009). The populations occur in grasslands scattered along a roughly north-south axis between –26° and –30° S latitude along the eastern region of SA at elevations of between 900 and 1,700 m (Figure 1). G. aurantiaca populations are found on or near rocky outcrops in organic–rich topsoil overlying well-drained, acidic, nutrient-deficient sub-soils associated with a high degree of leaching typical of the high rainfall and cool temperatures of mistbelt grasslands (Fey, 2010).

Reproduction of G. aurantiaca is both clonal, by means of underground stems, and sexual with flowering taking place in austral spring (September to November). The radiate capitula contain about 200 florets, with approximately 100 female ray florets in three outer whorls, and with the central disc florets hermaphrodite but functionally male. Capitula are not nyctinastic and remain open throughout flowering. The species is largely self-incompatible and highly dependent on pollinators for fruit set (Johnson et al., 2004). Populations are typically predominantly red-flowered, but both completely yellow-flowered and color polymorphic populations with yellow, orange and red-flowered inflorescences also occur (Johnson et al., 2005). The chemical basis of capitulum coloration of G. aurantiaca has not been analyzed as far as we know but the flower color of G. jamesonii a sister species, and the related commercially important ornamental Gerbera hybrids, is derived from carotenoids which cause yellow and orange coloration, and the flavonoids, pelargonidin and cyanidin, which are responsible for the red coloration (Valadon and Mummery, 1967; Tyrach and Horn, 1997). A feature of G. aurantiaca capitula is the dark center during the initial female stages of flowering when the disc florets are covered by the overarching dark purple pappus hairs and dark anther caps. Pollen is the primary floral reward as the florets produce little or no nectar and have no discernible scent.

The hairy monkey beetle Eriesthis vulpina Burmeister (Coleoptera: Scarabaeidae: Hopliini), appears to be the most important pollinator of G. aurantiaca (Johnson et al., 2004). Its use of color as a cue for locating flowers was evidenced by experiments in which large numbers of these beetles were captured in red plastic dishes placed in populations (Johnson et al., 2004). Eriesthis beetles almost completely ignored blue plastic dishes placed alongside red dishes in these populations (c. 100 beetles attracted to red dishes vs. one beetle attracted to blue dishes; SD Johnson, unpublished data), suggesting that they strongly prefer red over blue colors. The honeybee, Apis mellifera scutellata Lepeletier (Hymenoptera: Apidae) is also a common visitor in the yellow-flowered population.

We recorded the proportions of capitulum color morphs (as apparent in human vision) along transects across each of 23 populations ranging from southern KwaZulu-Natal to Mpumalanga (Figure 1 and Supplementary Table 1). At least one hundred clones were counted in larger populations and all clones were counted in the smaller populations. To minimize subjective bias in color allocation, one person (IMJ) carried out all counts. Ray floret samples were matched to a flower color chart (RHS, 2007) and broadly assigned to orange, red or yellow. Geographic position (south latitude and east longitude, WGS 84) and elevation (m.a.s.l) were measured using a handheld Garmin Etrex GPS. A generalized linear model (GLM) with a binomial error distribution, logit link function, and correction for overdispersion was applied to model the proportion of clones with red capitula in each population against latitude, longitude and elevation and we used likelihood ratios to test significance of fixed effects. All GLMs in this study were implemented in SPSS ver. 27 (IBM corp.).

We selected six populations representative of the color variation across the distribution range (Figure 1) to investigate pollinator color preferences; of these three (BY, VC and UH) were predominantly red-flowered, one (NG) almost entirely yellow-flowered and two (RH and LK) polymorphic with a mixture of red, orange and yellow-flowered plants (Table 1). Here we measured spectral reflectance (300–700 nm) of the ray florets using an Ocean Optics S2,000 spectrometer (Ocean Optics Inc., Dunedin, Fla.), Ocean Optics DT-mini deuterium tungsten halogen light source and fibre optic reflection probe (QR-400-7-UVVIS; 400 lm) held at 45° to the object surface in a probe holder (RPH-1). An Ocean Optics WS-1 diffuse reflectance standard was used to calibrate the spectrometer (Johnson and Andersson, 2002). Spectral reflectance readings were averaged from the mid-adaxial surface of three ray florets from at least 20 different plants from each population. All spectral reflectance curves are available from the Floral Reflectance Database (Arnold et al., 2010). Since we were interested in color as a floral signal and insect flower visitors most readily perceive rapid changes in spectral reflectance (Chittka and Menzel, 1992; Dyer et al., 2012) we used inflection or marker points which identify the wavelengths where change in spectral reflectance is maximal. We used the online Spectral MP (Dorin et al., 2020) to class capitulum color as yellow (inflection points from 520 to 530 nm), orange (inflection points from 530 to 600 nm) or red (inflection points from 600 to 625 nm). We recorded the positions of the different color morphs in the polymorphic RH population with a handheld Garmin Etrex GPS to investigate spatial structuring.

To test whether other floral traits that might influence pollinator attraction were associated with ray floret color we measured the capitulum diameter and ray floret length and ray floret inflection points for inflorescences from 50 clones collected randomly in the polymorphic RH population and tested for trait correlations using Pearson tests. We used finite mixture analysis allowing for unequal variance and with BIC criteria implemented in the r package mclust 5.4.6 (Scrucca et al., 2016) to assess whether the distribution of inflection points in populations were best explained by a single Gaussian distribution or by two or more Gaussian distributions. Both equal variance and unequal variance models were tested. We also tested for significant deviations from unimodality of inflection points in each population using Hartigan’s diptest implemented in the r package diptest (Maechler, 2016).

We sourced the bioclimatic variables temperature (MAT: mean annual temperature °C) and precipitation (MAP: mean annual precipitation mm.) from www.worldclim.org/bioclim (Hijmans et al., 2005) for 23 populations across the distribution range of G. aurantiaca (Supplementary Table 1). A generalized linear model (GLM) with a binomial error distribution, logit link function, and correction for overdispersion was applied to model the proportion of clones with red capitula in each population against temperature and precipitation.

Since flower color in polymorphic species may be associated with differences in edaphic factors (soil characteristics) (Horovitz, 1976; Rajakaruna and Bohm, 1999) we collected soil samples from 14 population locations of G. aurantiaca. For each population 15 augered subsamples of the top 15 cm layer taken randomly across the population patch were combined, air dried, stored at room temperature and analyzed for pH (KCL), exchange acidity, total cations, acid saturation, organic carbon (C), calcium (Ca), copper (Cu), potassium (K), magnesium (Mg), manganese (Mn), nitrogen (N₂), phosphorus (P), zinc (Zn), and clay content (Supplementary Table 2) by the Soil Fertility and Analytical Services Laboratory, KwaZulu-Natal Dept. of Agriculture (Manson and Roberts, 2000). We used Principal Components Analysis (Braak and Smilauer, 2002) to assess the relationship between the 14 edaphic variables (standardized to avoid using different scales for comparative purposes). Population scores were plotted onto the principal components axes and correlation vectors were used to examine relationships between the population samples and the first two PCA axes of the soil ordination. Linear regression was carried out using the first four regression factor scores and the proportion of red clones in the 14 populations.

To test whether the capitulum color of G. aurantiaca plants from different color morph populations changes when they are grown together under identical climatic and edaphic conditions, we cultivated both ramets and seeds collected from red and yellow flowered populations in a common garden at the KwaZulu-Natal National Botanical Garden in Pietermaritzburg where conditions (soil characteristics, temperature and precipitation) were identical. We did not include orange morphs since we were not aware of their extent in the polymorphic populations at the start of this experiment. Sample sizes were limited due to permit restrictions relating to the threatened status of the species. We collected flowering ramets from each of five widely spaced clones growing in red (BY) and yellow (NG) source populations and transplanted these into the common garden nursery bed. Ramets were used as they provided a baseline spectral reflectance and because they take less time to flower in cultivation than seed raised plants. In addition, seeds were collected from clones at least five meters apart in the red-flowered LM and the yellow-flowered NG populations and marked seedlings from each population were transplanted randomly into the common garden nursery bed. Reflectance spectra of ray florets from each individual were measured at collection (in the case of ramets) and at intervals once flowering had occurred for both groups over a 5 year period, We compared mean inflection points of the spectral reflectance curves (Dorin et al., 2020) of ray florets for both ramets and seed grown plants using a GLM with normal distribution with source population and year as independent variables.

Pan traps have been widely used for passive sampling of flower-visiting insects and may give an indication of their color preferences as well as their abundance (Picker and Midgley, 1996; Leong and Thorp, 1999; Shrestha et al., 2019). Here we tested color discrimination of potential pollinators using colored pan trap arrays at the six study populations. Red and yellow traps were chosen as they represent the extremes of capitulum colors that occur naturally in the range of G. aurantiaca. We did not test for orange due the unavailability of traps with suitable reflectance spectra. Spectral reflection curves were used to calculate the inflection points (Dorin et al., 2020) of the traps and compare these to those of G. aurantiaca capitula. We used twenty sets of red and yellow plastic traps (11 cm in diameter and 8.5 cm deep) filled with 150 ml of water and placed randomly in pairs 15 cm apart amongst flowering G. aurantiaca plants. We carried out these experiments during peak flowering on sunny days between 07:00 and 16:00 h in late October when insects were most active. Captured insects were identified at least to family level and counted. Beetles appeared to be unaffected by being trapped and most were released after recording on account of their role as important pollinators of this threatened study species. Pan trap catches of the dominant insect visitors, E. vulpina beetles and honeybees, were analyzed using a GLM with negative binomial distribution and log link with trap color and population site as independent variables. For the analysis of catches, one value at sites where no insects were caught in traps had to be adjusted from a choice for yellow to one for red in order for the model to converge (Zuur et al., 2009).

Since hymenopteran vision is well-researched we used the bee color hexagon (Chittka and Menzel, 1992; Peitsch et al., 1992) to assess how honeybees would perceive the color of capitula of G. aurantiaca and the pan traps. The mean reflectance spectrum calculated from three ray florets from each of 20 individuals of each color form and those from the red and yellow pan traps was plotted in the bee color hexagon. Background color was calculated from the spectra of G. aurantiaca leaves. Coleopteran vision has been less well studied and due to this we did not assess their perception of capitulum color.

To determine if E. vulpina hopliine scarab beetles, the dominant insect visitors to G. aurantiaca, exhibit color preferences in a polymorphic population, we recorded the frequency of these beetles visiting capitula of different color morphs at the RH population along three 200 m transects during peak flowering in three separate years. As the intensity of flowering varied from year to year, we recorded different total inflorescence numbers per transect (2008 N = 155, 2011 N = 471, 2013 N = 291). The proportion of capitula of different colors (yellow, orange, and red) that were occupied by beetles was compared using a logistic GLM with binomial distribution and logit link function.

We measured fecundity (mean fruit set per capitulum) in relation to morph color in the polymorphic RH population. In 2011 (N = 61: yellow = 19, orange = 28, red = 24), and 2013 (N = 56: yellow = 9, orange = 13, red = 34) capitula were bagged and labeled at the end of the male phase and harvested after 3 weeks. Filled fruit, easily distinguished by their larger size, darker color and firmness were counted for each bagged capitulum. Fruit set per capitulum for each year was analyzed using a GLM with a negative binomial distribution and log link function. Color morph, year and the interaction of color morph and year were fixed factors in this model.

In human color vision 15 of the 23 G. aurantiaca populations measured had predominantly red-flowered clones, 5 were predominantly yellow-flowered and 3 were polymorphic with mixed yellow, orange and red-flowered clones. The red flowered populations occur mainly in the southern part of the distribution range with one polymorphic population (LK), the yellow cluster in the center, and northern populations are mainly color polymorphic with varying proportions of red, orange and yellow capitula (Figure 1 and Supplementary Table 1). The proportion of clones with red capitula in the 23 populations was significantly correlated with latitude (χ² = 9.460, df = 22, P = 0.002), and longitude (χ² = 6.791, df = 22, P = 0.009) but not with elevation (χ² = 1.805, df = 22, P = 0.179).

Mean inflection points of the reflectance spectra of the six representative populations were between 609 and 615 nm for the red-flowered populations (BY 615 nm, VC 611 nm, UH 609 nm, LK 600 nm) populations, 524 nm for the yellow-flowered NG population and 597 nm for the polymorphic RH population (Figures 2A–F and Supplementary Table 2). No UV reflectance (300–400 nm) was recorded (Figure 2). Although the northern color polymorphic RH and central LK populations showed a wide range of color variation we were able to assign plants to particular color classes as the distribution of inflection points in this population was broad (Figures 2D,F), with most inflection points clustered around 620 nm with a smaller peak around 540 nm for RH and 600–620 with a smaller cluster at 530–570 nm for LK. The frequency of inflection points in the NG, BY, VC and UH populations fitted a single Gaussian distribution better than two (ΔBIC range = 2.24–7.80), but the distribution of inflection points fitted two Gaussian distributions much better than one Gaussian distribution in the RH population (ΔBIC = 39.8) and LK population (ΔBIC = 23.5). Fits to three or more Gaussian distributions were not supported for any population, but in the RH and LK populations, models with fits to three Gaussian distributions fitted better than models with fits to a single distribution (ΔBIC = 16.4–33.6) but were still supported less than were models fitted to two Gaussian distributions (ΔBIC = 6.2–7.1). However, we found no significant deviations from unimodality for the distributions of inflection points in any of the populations (P > 0.095). We also found no significant correlations between inflection points and capitulum diameter (r = –0.15, P = 0.289) or ray length (r = –0.029, P = 0.845) (Supplementary Table 3).

Environmental conditions were similar across the distribution range of G. aurantiaca with climatic variables for the 23 population localities ranging from 14.0 to 17.6°C for mean annual temperature (MAT), 878-1057 mm for mean annual precipitation (MAP) (Supplementary Table 1). The proportion of clones with red capitula in populations was not significantly correlated with temperature (χ² = 0.662, df = 22, P = 0.416) or precipitation (χ² = 2.702, df = 22, P = 0.1).

Differences in soil chemistry between sites were small (pH ranged from 3.94 to 4.69, N from 0.12 to 0.6% and organic C content 2.5–6%) (Supplementary Table 4). Principal components of the edaphic variables were summarized into the first and second axes of the principal component analysis (PCA), which accounted for 64.9% of the total variation and were mainly explained by pH-acidity gradient (axis 1) and Org C/N gradient (axis 2) (Supplementary Figure 1 and Supplementary Table 5). The PCA biplot (Figure 3) shows an apparent pH-acidity gradient, and the two northernmost sites on the top left stand out because of their more weathered soils with higher Mg and lower than average organic carbon content. pH and Organic C can be used as largely independent surrogate variables to represent variation of other associated variables along these two typical soil chemistry gradients. We found no significant relationship between the proportion of red clones in a population and PCA Factor 1 (r² = 0.037, P = 0.512), Factor 2 (r² = 0.24, P = 0.596), Factor 3 (r² = 0.008, P = 0.977) or Factor 4 (r² = 0.05, P = 0.441) using linear regression. We further observed no spatial structure in the polymorphic RH population with different color morphs scattered randomly (Figure 4A) and often growing less than one meter from each other (Figure 4E).

Both ramets and seed-raised plants retained the capitulum color of the parent populations when grown in a common garden over a 5-year period. Flowering ramets collected in the red-flowered BY (mean inflection point 616.0 ± 1.7 nm, N = 5) and yellow-flowered NG (mean inflection point 525.6 ± 1.8, N = 5) populations showed no change from their original ray floret morph color (χ² = 3397.82, df = 1, P < 0.0001) with no significant change in morph color recorded over the 5-year observation period (χ² = 1.319, df = 1, P = 0.517) and no significant interaction between population and year (χ² = 1.530, df = 3, P = 0.465) (Figure 5A).

Seed-raised plants flowered 4 years after germination and displayed the predominant capitulum color of the source populations: LM (red with mean inflection point 606.80 ± 1.8 nm, N = 18) and NG (yellow with mean inflection point 523.4 ± 1.9 nm, N = 18). There was a significant difference between the mean inflection points of the two seed populations (χ² = 2012.457, df = 1, P < 0.0001) but there was no significant effect of time (χ² = 0.51, df = 1, P = 0.821) or interaction between seed population and year (χ² = 0.003, df = 2, P = 0.957) (Figure 5B).

A total of 682 E. vulpina beetles, and 39 honeybees, as well as some insects from other taxonomic groups were caught in the red and yellow pan traps. At all populations where beetles were present (UH, LK, NG, and RH) significantly more individuals of E. vulpina were caught in yellow than red traps (Figure 6A and Supplementary Table 6). The highest number of beetles overall was trapped at the red-flowered UH population, followed by RH (flower color polymorphic), NG (yellow), and LK (predominantly red), while no individuals of E. vulpina were trapped at the two southernmost sites, BY and VC (red). We recorded significant effects of population (χ² = 97.951, df = 5, P = 0.001), and of trap color (χ² = 14.904, df = 1, P = 0.001) on the number of E. vulpina beetles caught, but the interaction between population and trap color was not significant (χ² = 7.182, df = 6, P = 0.207) suggesting that there were no marked differences in beetle color preference between populations.

Very few honeybees were captured in pan traps, but despite this low statistical power we detected a slight significant effect of trap color (χ² = 5.172, df = 1, P = 0.023), but not of population (χ² = 9.464, df = 5, P = 0.092). There was a significant interaction between population and trap color (χ² = 15.352, df = 6, P = 0.009) with more bees trapped in yellow than red traps at four of the populations (BY, VC, UH, and LK), equal numbers at NG and none trapped at RH although they were frequently observed on other plants at the site (Figure 6B and Supplementary Table 6).

When plotted in the model of hymenopteran color perception, loci of G. aurantiaca ray floret reflectance spectra were positioned primarily in the green and UV green segments of the hexagon with the florets of the red color forms clustered more closely to the origin of the hexagon than the yellow forms, suggesting that the yellow capitula are more visible to bees than the red. The loci of orange forms were intermediate between red and yellow in bee color space (Figure 6C). The spectral reflectances of the pan traps (Supplementary Figure 2) were similar to those of the red and yellow ray florets in terms of hymenopteran color perception (Figure 6C) and inflection points (red trap = 601 nm, floret = 608 nm; yellow trap = 519 nm, floret = 524).

Analysis of the incidence of E. vulpina beetles present on orange, red and yellow colored capitula (Figures 4B–D) at the mixed color site (RH) during transect surveys in three separate years indicated that the beetles did not discriminate among color morphs, with no significant differences found between the proportions of orange, red and yellow capitula that were occupied by beetles (χ² = 0.410, df = 2, P = 0.814). The overall proportion of capitula occupied by beetles varied significantly among years (χ² = 24.90, df = 2, P < 0.0001) with a higher incidence of beetles on capitula in 2013 than in 2008 or 2011 (Figure 7A). There was no interaction between year and color morph (χ² = 0.643, df = 4, P = 0.958), implying that responses to colors by beetles did not vary among years.

We did not detect any significant differences in mean fruit set among color morphs in 2011 or 2013 (χ² = 0.349, df = 2, P = 0.840) or year (χ² = 292, df = 1, P = 0.589). There was no significant interaction between year and color morph (χ² = 0.216, df = 3, P = 0.898) (Figure 7B).