Spatiotemporal niche overlap, asymmetric reproductive interference, and population genetics between the sympatric species, Rhododendron diversipilosum and Rhododendron subarcticum, in alpine fellfield habitat

ABSTRACT

Reproductive interference between sympatric-related species often causes adverse impacts on rare species, which increases the risk of local extinction, particularly in small and isolated populations. To evaluate the congeneric interactions in alpine plants, we compared the ecological and genetic properties and assessed the reproductive interference between tetraploid Rhododendron diversipilosum (widespread species) and diploid Rhododendron subarcticum (rare species in alpine fellfields) in northern Japan. In alpine fellfields, R. diversipilosum is commonly distributed close to shrubby patches, whereas R. subarcticum tends to grow in more exposed places, although they are sometimes mixed. R. subarcticum initiated flowering one week earlier; however, the flowering periods overlapped between species, indicating incomplete phenological isolation. The pollination experiment showed that both species were self-incompatible. Furthermore, heterospecific pollination occurred only for R. subarcticum; however, hybrid seeds seldom germinated. These results indicate that reproductive interference is asymmetric between species, where only the rare R. subarcticum may suffer from heterospecific pollination. Genetic analysis showed that R. subarcticum populations had lower genetic diversity but higher divergence among populations than R. diversipilosum. Although spatiotemporal niche separations may mitigate reproductive interference, the risk of local extinction could be higher in R. subarcticum populations because of isolated distributions with low genetic diversity.

KEYWORDS:

• Alpine plant
• competition
• congeneric species
• interspecific interaction
• reproductive interference

Introduction

Alpine ecosystems are characterized by highly heterogeneous environments at a fine scale because the complex terrain and geological aspects of mountains create spatial variations in microclimate, snow distribution, hydrology, and edaphic conditions (Körner 2003; Nagy and Grabherr 2009). Reflecting the combinations of diverse habitat types, alpine plant communities are often distributed as small mosaic patches, where diverse plant species grow adjacently within a local area, including both widespread and narrow-range species (Körner 1995; McEwen and Vamosi 2010). The patchy distribution of plant communities causes diverse interspecific interactions, such as competition, facilitation, and reproductive interference between interacting species at the microscale. The competitive and facilitative relationships between sympatric alpine species vary greatly depending on the growing conditions, such as environmental harshness, stability, and productivity of habitats (Choler, Michalet, and Callaway 2001; Onipchenko et al. 2009).

Reproductive interference caused by heterospecific pollen transfer often causes asymmetric fitness loss between species sharing pollinators (Ashman and Arceo-Gómez 2013). Reproductive interference strongly affects the reproduction of closely related species (Kyogoku 2015), resulting in asymmetric hybridization and niche separation between the interacting species (Tiffin, Olson, and Moyle 2001; Zhang, Montgomery, and Huang 2016). For example, the extent of reproductive interference between congeneric alpine shrub species changes drastically along the gradient of snowmelt time, reflecting the seasonal dynamics of pollinator activity and floral abundance of interacting species (Kameyama, Kasagi, and Kudo 2008). Recent studies stress the importance of reproductive interference as a driver promoting niche separation between related species (Runquist and Stanton 2013; Toll and Wiilis 2018; Christie and Strauss 2020).

The current geographic distribution of alpine plants in mid-latitude mountain regions closely reflects the migration history during the Pleistocene, where several species originating from arctic regions are distributed only in specific mountain regions as relict species (Birks 2008; Testolin et al. 2021). It should be an ecologically important issue how these species have survived and preserved their populations in isolated habitats since the last glacial period, particularly when competitive species exist. Generally, spatial and temporal niche separation is effective for the coexistence of sympatric-related species (Coyne and Orr 2004; Ramírez-Aguirre et al. 2019). Habitat use separation at the microscale and phenological isolation of resource use timing can mitigate exploitation competition for limited resources (Shigesada, Kawasaki, and Teramoto 1979; Snyder 2008) or pollinator acquisition (Waser 1978). In addition, spatial and temporal separation can reduce the risk of reproductive interference caused by heterospecific pollen transfer (Paudel et al. 2018; Christie and Strauss 2020). Because the risk of heterospecific pollination is commonly frequency dependent (Kyogoku 2015; Christie and Strauss 2020) and increases with closer phylogenetic relatedness (Streher et al. 2020), rare species are more susceptible to the interference competition with widespread congeneric species. However, little is known regarding the intensity of interspecific competition and the effectiveness of niche separation between related alpine plants.

In addition to the negative interactions, the risk of local extinction in relict alpine plants may increase when isolated small populations experience strong genetic drift (Ellstrand and Elam 1993). Though genetic drift provides a possibility of local speciation in alpine plants (Gavini, Ezcurra, and Aizen 2019), it may cause various deleterious effects on rare species through the loss of genetic diversity. Alpine plant populations in marginal regions of the geographic distribution range often have low genetic diversity (Hirao et al. 2015; Vásquez et al. 2016). This is possibly because gene flow via pollen or seed dispersal from other populations is limited by a geographic barrier (Ægisdóttir, Kuss, and Stöcklin 2009; Chen et al. 2019). Because small genetic variations may amplify the deleterious effects of inbreeding depression (Carr and Dudash 2003) and vulnerability to stochastic environmental change (Bouzat 2010), the genetic diversity of local populations could be a crucial issue for the conservation of relict alpine plants.

Congeneric ericaceous shrubs Rhododendron diversipilosum (Nakai) Harmaja and Rhododendron subarcticum Harmaja (synonym, Ledum palustre var. diversipilosum and Ledum palustre var. decumbens) are ideal species for testing interspecific interactions in alpine plants. Rhododendron diversipilosum is distributed in mid-latitudinal Far East Asia (Kuril Islands, Sakhalin), including northern Japan (Yamazaki 1993; Harmaja 1999). In contrast, R. subarcticum is widely distributed in the arctic tundra regions across Kamchatka, Alaska, northern Canada, and Greenland (Hultén 1968) but only sporadically in the mid-latitudinal alpine regions of eastern Asia. In Hokkaido of northern Japan, both species grow, but their major growing habitats are different; R. diversipilosum grows in various environments, such as low-elevation wetlands, understory of subalpine forests, and alpine shrubby zones, whereas R. subarcticum only grows in the exposed alpine fellfield habitat of the central mountain range (the Taisetsu Mountains). Thus, R. subarcticum is a rare alpine species, whereas R. diversipilosum is a major widespread species growing from the lowland to the alpine zone in Hokkaido. Although both species are sometimes mixed in the fellfield habitat of alpine zones, little is known regarding the ecological properties, mating system, genetic structure, and interactions between the species.

To assess the interactions between R. diversipilosum (widespread species) and R. subarcticum (rare species), we compared the ecological and genetic characteristics of these species in alpine fellfields of the Taisetsu Mountains. This study focused on the following questions: (1) Are there spatial and temporal niche separations between species? If so, what are the ecological characteristics related to these patterns? (2) Is there any reproductive interference between species? If so, does reproductive interference occur symmetrically or asymmetrically between species? (3) To what extent do population genetic properties vary between species in terms of genetic diversity, genetic divergence, and isolation by distance among populations?

Materials and methods

Material plants

Rhododendron diversipilosum is an ericaceous evergreen shrub that grows in moorlands, understory of subalpine forests, and alpine regions of northern Japan. This species is 30 to 70 cm in height with elliptical leaves and corymb inflorescences composed of many small white flowers (Figure 1a). Flowering in the alpine zone usually occurs during late June and mid-July; flowers are mainly visited by syrphid and other flies (Kudo 2022), and seeds are dispersed by wind in September. This species is reported to be weakly self-compatible, and its ploidy level has been identified as tetraploid (Wakui 2021). Rhododendron subarcticum is widely distributed in the arctic tundra of North America and Siberia, and its southern distribution margin is located in the alpine zone of Hokkaido. This species is morphologically similar to R. diversipilosum but smaller in size (commonly <20 cm in height) with narrower leaves (Figure 1b). Dipteran insects are main pollinators also in this species (Kudo 2022). As for the ploidy levels of this species, Löve and Löve (1982) reported it as tetraploid (2n = 52), whereas Lantai and Kihlman (1995) reported it as diploid (2n = 26). However, the mating system of this species remains unclear.

Figure 1. Pictures of two related Rhododendron species, (a) R. diversipilosum and (b) R. subarcticum, growing at alpine fellfields in Hokkaido of northern Japan.

Study sites

Field surveys were conducted at four mountain sites in the Taisetsu Mountains of Hokkaido (Mt. Hira, Hisago, Muka, and Tengu) in 2018 and 2019 (Figure 2). Two mountain sites (Mt. Teshio and Ishikari) were used for genetic comparisons. The location and elevation of each site are listed in Table 1. All research sites were located on exposed plateaus or ridges with little snow cover during winter; that is, the fellfield habitat. Dwarf shrubs (such as Loiseleuria procumbence, Diapensia lapponica, and Empetrum nigrum) and lichens (Cladonia and Cetraria spp.) dominated in the study sites, mixed with sparse patches of alpine dwarf pine (Pinus pumila). The height of P. pumila patches was generally 50 to 100 cm in the fellfield habitat. Rhododendron diversipilosum often grows around P. pumila patches, whereas R. subarcticum tends to grow in exposed places, although both species sometimes grow together. This study focused on P. pumila patches that can nurse other species in the alpine fellfields, as observed for R. diversipilosum. Pinus pumila patches act as shelters of other alpine plants against intense radiation and strong winds during the growing period in the exposed habitat. Furthermore, P. pumila patches trap snow (Okitsu and Ito 1989) by which neighboring small plants are protected from freezing temperatures under snow cover during winter.

Figure 2. (a) Location of the Taisetsu Mountains in Hokkaido and (b) study sites of this study.

Table 1. Geographical locations of individual study sites and observation year.

Site	Longitude (N)	Latitude (E)	Elevation (m)	Observation year
Mt. Hira	43.7629	143.0089	1,737	2018–2019
Mt. Hisago	43.5480	142.8667	1,678	2018–2019
Mt. Ishikari^a	43.5483	142.8733	1,876	2020
Mt. Muka	43.7065	143.1717	1,725	2018–2019
Mt. Tengu	43.8175	143.0635	1,472	2018–2019
Mt. Teshio^a	43.9641	142.8865	1,516	2018

^aThese sites were used only for genetic analysis.

Morphological and performance traits

To compare morphological and performance traits, the following seven traits were measured for approximately twenty individuals per species at each site: leaf mass per area (LMA), leaf area, leaf longevity, annual leaf production, annual shoot growth, canopy height, and flower number per shoot. Leaf traits were measured for current-year leaves. Leaf area was measured for two mature leaves per plant using a CanoScan LiDE 60 scanner and the free software ImageJ (Abràmoff, Magalhâes, and Ram 2004). Leaf mass was measured after desiccation at 60°C for 48 hours, and the LMA was calculated. Leaf longevity was estimated by making a life table of the leaf cohort based on the number of living leaves and leaf scars remaining in each annual segment on a stem (Xiao 2003). Annual leaf production was measured by multiplying the number of annual leaves by the mean mass of a single leaf. Annual shoot growth was expressed as the length of a one-year-old branch.

Generalized linear mixed models (GLMMs) were used to compare individual traits (leaf area, LMA, leaf longevity, annual leaf production, annual shoot growth, canopy height, and flower number per shoot) between species across sites, where species was set as an explanatory variable and site as a random factor. GLMM postulating a Poisson error distribution was applied for the analysis of flower number, and GLMM postulating a gamma error distribution was used for other traits. R v3.6.1 (R Core Team 2019) was used for all statistical analyses.

Spatial niche separation

To evaluate the dependency of each Rhododendron species to the patches of P. pumila, we conducted a belt transect survey. Five 1 × 1 m quadrats were set in a row from the vicinity of developed P. pumila patches (with 70–100 percent canopy cover) toward an exposed direction (see Figure 3a), and canopy projection areas of P. pumila, R. diversipilosum, and R. subarcticum were recorded in each quadrat. The projection area of each species was recorded with 10 percent accuracy and expressed as a rank from 0 (< 5 percent) to 10 (95–100 percent; hereafter abundance). The number of belt transects was eight for the Mt. Tengu site and ten for the other four sites.

Figure 3. (a) An example of belt transect setting. Five quadrates (white squares) were arranged from a marginal part of the Pinus pumila patch toward the exposed direction. (b) Relationship between the distance from P. pumila patches and the abundance of each species, expressed as a rank of plant cover at 10 percent accuracy. Black solid line indicates P. pumila, black broken line indicates Rhododendron diversipilosum, and grayish line indicates R. subarcticum. Vertical axes show the mean value of the coverage of each species. Error bars indicate standard error.

To evaluate the microscale habitat use pattern between the species, we tested the response of R. subarcticum to the abundance of R. diversipilosum and P. pumila using a GLMM postulating a Poisson error distribution. In the GLMM, the response variable was the abundance of R. subarcticum, and the explanatory variables were the abundance of R. diversipilosum and P. pumila, including their interactions, and the belt transect nested by site was set as a random factor. Thereafter, the best-fit model was selected based on the model selection using Akaike’s information criterion.

Phenological separation

The flowering overlap between R. diversipilosum and R. subarcticum was measured at four sites. This survey was conducted three times during the flowering season in 2019: early (late June), middle (early July), and late (mid-July) flowering periods of each species. We set six 1 × 1 m quadrats arbitrary on the flowering patches of these species in each flowering period, and the number of flowering inflorescences was recorded. However, at the Mt. Tengu site, we failed to observe this in the early flowering period because of earlier flowering progress than we expected. The flowering pattern of R. diversipilosum was compared with R. subarcticum using a GLMM, postulating a Poisson error distribution. In the GLMM, the number of inflorescences at flowering within a quadrat was a response variable, species (R. diversipilosum or R. subarcticum), day of the year in 2019 (both linear and squared terms), and their interactions were set as explanatory variables, and the site was set as a random factor.

Pollination experiment

To clarify the mating system of each species and the reproductive interference between the species, we conducted a hand pollination experiment in mid-July 2018 at Mt. Hisago. Four pollination treatments were conducted for each species: simple bagging, self-pollination, outcross pollination, and heterospecific pollination. We arbitrarily selected twenty inflorescences per treatment at 5-m intervals to avoid multiple sampling within the same clonal patch. All inflorescences were covered with fine-meshed bags prior to flowering. Once flowering had occurred, hand pollination was performed. After hand pollination, the treated inflorescences were bagged again. For heterospecific pollination, we deposited pollen grains of R. diversipilosum (or R. subarcticum) onto the stigmas of R. subarcticum (or R. diversipilosum) flowers by hand using forceps. Both species have porous dehiscence anthers as well as other Rhododendron species in which tetrad pollen is contained within anther locules (King and Buchmann 1995). Once anthers mature, pollen grains covered with viscin strands are exposed from the rounded apical pores when pollinators touch the anthers. Because anthers mature after opening flowers and anthers and stigma are separately arranged within a flower (herkogamy), automatic self-pollination rarely occurs in bagged inflorescences. Thus, emasculation was not conducted before hand pollination in the pollination experiment. For the hand pollination treatment, we confirmed that the stigmas of the recipient flowers were fresh (whitish-yellow color) using a hand loupe. After the pollination treatments, some inflorescences were damaged by frost events, resulting in a decreasing sample size (six to eighteen inflorescences per treatment). All fruits were harvested from the manipulated inflorescences in mid-September, before seed dispersal. The sampled fruits were dried at room temperature. One to three developed fruits were selected from each inflorescence, and the number of developed seeds and undeveloped ovules was recorded in each fruit under a microscope.

Seed set rates (seed-to-ovule ratios) of manipulated flowers were analyzed using a GLM postulating a negative binomial error distribution in each species due to excessive zero values in seed number. In the GLM, the response variable was the number of developed seeds, the explanatory variable was pollination treatment, and the total number of ovules (sum of developed and undeveloped seeds) was set as an offset term. Furthermore, Tukey’s post hoc tests were conducted to evaluate which combination of pollination treatments had a significant difference.

Germination experiment

We conducted a germination experiment using seeds obtained from the pollination experiment. All seeds were stored at 4°C under dry conditions for two months until the germination experiment. Pretreatment of cold stratification under moist conditions was not conducted because a previous study reported that seeds of a congeneric species had high germinability without cold stratification under warm conditions (Karlin and Bliss 1983). To determine the best thermal conditions for germination in each species, germination rates were compared under three thermal conditions: (1) constant 25°C, (2) constant 15°C, and (3) 12 hours at 25°C/12 hours at 10°C. For all treatments, thirty seeds were disseminated on wet filter paper in a 5-cm Petri dish under a 12-hour light period, with four replicates in each treatment. Distilled water was added to each Petri dish daily to maintain moist conditions. We conducted another germination experiment on seeds obtained from heterospecific pollination. Because only R. subarcticum produced seeds by heterospecific pollination (see results), this experiment was conducted only for seeds obtained from R. subarcticum. Germinability of the hybridized seeds was tested under a constant temperature of 25°C with a 12-hour light period. To compare germination rates, a GLMM postulating a binomial error distribution was conducted, where the explanatory variable was the thermal condition in each species. Furthermore, the comparison between the species under each thermal condition was conducted using a GLMM in which replicates were selected as a random factor.

Genetic diversity and differentiation among sites

To compare genetic differentiation among sites between the species, we collected leaf samples of R. diversipilosum from the Mt. Hira, Hisago, Muka, Tengu, and Teshio sites and R. subarcticum from the Mt. Hira, Hisago, Muka, Tengu, and Ishikari sites. Leaf sampling was conducted at 2- to 20-m intervals, depending on the population size and plant density. Leaf samples were stored in silica gel at room temperature. The number of samples per site ranged from twenty to twenty-four.

Genomic DNA was extracted from the leaf samples using the Cetyltrimethylammonium bromid method (Doyle and Doyle 1987). Eight microsatellite markers developed for R. diversipilosum (Wakui 2021; Table S3) were used for both Rhododendron species. Each marker was amplified by two-step polymerase chain reaction (PCR). Reaction mix for first-step PCR was 5 µL in total, including 2.5 µL of 2 × Type-it PCR Multiplex Master Mix, 1.0 µL of diluted DNA, 0.5 µL of forward primer mix (0.1 µM for each locus), 0.5 µL of reverse primer mix (0.1 µM for each locus), and 0.5 µL of sterile water. The PCR cycle conditions for the first step were as follows: initial 5 minutes at 95°C and nineteen cycles of 30 seconds at 95°C, 90 seconds at 60°C, and 30 seconds at 72°C. After the first step amplification, 0.5 µL of 2 × Type-it PCR Multiplex Master Mix and 0.5 µL mix of fluorescently labeled universal primers (2.4 µM for each color) and reverse primers (2.4 µM for each locus) were added. We conducted the second PCR step under the following conditions: fifteen cycles of 30 seconds at 95°C, 90 seconds at 60°C, and 30 seconds at 72°C. Microsatellite fragments were analyzed using an Applied Biosystems 3730 Genetic Analyzer with the Gene Scan 500 LIZ Size Standard, and genotypes were determined using Peak Scanner software v1.0 (Applied Biosystems, USA). After genotyping, samples with the same genotypes at all loci were treated as the same genet, and the genet number across all populations in each species was calculated.

Allelic richness (A_r), expected heterozygosity (H_e), and observed heterozygosity (H_o) were analyzed as indices of genetic diversity using SPAGeDi 1.5 (Hardy and Vekemans 2002). We conducted a Welch’s t test for A_r, H_e, and H_o between R. diversipilosum and R. subarcticum. To assess genetic differentiation among populations in each species, global F-statistics and pairwise F_st (Weir and Cockerham 1984) were calculated using SPAGeDi 1.5 (Hardy and Vekemans 2002). We inspected the geographic patterns of genetic differentiation—that is, isolation by distance (IBD)—for each species. To reveal the degree of IBD, the correlation between geographic distance and pairwise F_st value was analyzed using the Mantel test with 9,999 permutations. We calculated the geographic distance between sites using Google Earth Pro 7.3.4.8248 (https://www.google.com/earth/versions/#earth-pro).

Results

Morphological and reproductive traits

In comparing morphological and performance traits between R. diversipilosum and R. subarcticum, leaf area, annual leaf production, and plant height were significantly larger in R. diversipilosum, whereas LMA and leaf longevity were significantly larger in R. subarcticum (Figure 4, Table S1). There were no significant differences in annual shoot growth and flower number per inflorescence between species.

Figure 4. Comparisons of morphological and performance traits between Rhododendron diversipilosum and Rhododendron subarcticum at the fellfield habitat: (a) individual leaf area, (b) leaf mass per area (LMA), (c) leaf longevity, (d) annual leaf production in dry weight, (e) annual shoot elongation, (f) plant height, and (g) flower number per inflorescence. ***p < .001, n.s. = not significant.

Spatial niche separation

As expected, the abundance of R. diversipilosum was larger close to P. pumila patches, whereas R. subarcticum stayed away from P. pumila cover (Figure 3b). In the GLMM, the abundance of R. subarcticum was negatively correlated with the abundance of P. pumila (z = −2.39, p = .017), whereas the abundance of R. diversipilosum (z = −1.29, p = .20) and the interaction between P. pumila and R. diversipilosum (z = 0.59, p = .55) did not affect the abundance of R. subarcticum. In the best-fit model, the abundance of R. diversipilosum and the interaction term were removed, and only the abundance of P. pumila was negatively correlated with the abundance of R. subarcticum (z = −3.91, p < .0001). The Akaike’s information criterion of the full model was 553, whereas that of the best-fit model was 552. Thus, R. subarcticum distribution is exclusively dependent of P. pumila patches but is independent of R. diversipilosum.

Phenological difference

In the phenological survey, the flowering of R. subarcticum occurred approximately one week earlier than that of R. diversipilosum (Figure 5). A significant difference in flowering phenology was detected between the species using the GLMM (z = 6.189, p < .0001). However, the flowering periods of both species overlapped at every site to some extent, indicating incomplete phenological isolation between species. The major flowering period of R. subarcticum occurred from late June to early July, whereas that of R. diversipilosum occurred from early July to late July, although flowering progress at the Mt. Tengu site occurred one week earlier than at other sites due to the lower elevation (Table 1).

Figure 5. Comparisons of flowering periods of Rhododendron diversipilosum and R. subarcticum in each study site. (a) Mt. Hira, (b) Mt. Muka, (c) Mt. Hisago, and (d) Mt. Tengu. The number of flowering inflorescences per plot (n = 6) was measured in early, middle, and late flowering seasons of two Rhododendron species. Open column: R. diversipilosum, gray column: R. subarcticum. Error bar indicates standard error.

Pollination experiment

Among the manipulated plants in the pollination experiment (twenty inflorescences per treatment), two to ten inflorescences of R. diversipilosum and two to fourteen inflorescences of R. subarcticum were damaged by frost events during the flowering period. These inflorescences were excluded from analysis. Under natural conditions, the seed set rate of R. diversipilosum (37.2 ± 19.6 percent; mean ± SD) was significantly higher than that of R. subarcticum (11.8 ± 16.4 percent; p < .0001). In a pollination experiment conducted for R. diversipilosum, the bagging, self-pollination, and heterospecific pollination treatments resulted in very low seed set rates. In contrast, there was no significant difference in seed set rates between the outcross pollination treatment and control groups (Figure 6a, Table S2). This indicates that R. diversipilosum is physiologically self-incompatible, but there is no pollen limitation under natural conditions.

Figure 6. Results of the hand-pollination experiment conducted in 2018 at the Mt. Hisago site for (a) Rhododendron diversipilosum and (b) R. subarcticum. Seed set percentage in each treatment (Co: control, Oc: outcross pollination, He: heterospecific pollination, Ba: simple bagging, Se: self-pollination) is shown. Sample size of each treatment is shown within parentheses. Different characters indicate significant difference in each species among treatments (p < .05 by Tukey’s test).

In the pollination experiment conducted for R. subarcticum, the seed set rates of the bagging and self-pollination treatments were close to zero, but the seed set rates of the outcross and heterospecific pollination treatments were significantly higher than those of the control (Figure 6b, Table S2). These results indicate that R. subarcticum is physiologically self-incompatible but compatible with the heterospecific pollen of R. diversipilosum. Furthermore, seed production under natural conditions is limited by insufficient pollen receipt. The results of the pollination experiment indicate that heterospecific compatibility is asymmetric between species.

Germination experiment

Seeds of both R. diversipilosum and R. subarcticum showed high germination ability without cold–wet stratification under warm conditions. In both species, most seeds germinated within 5 to 15 days at constant 25°C (Figure 7). Germination rates of outcrossed seeds reached >80 percent in all temperature treatments for R. diversipilosum, but the 25/10°C treatment showed the maximum germinability (93.3 ± 3.8 percent; mean ± SD; Table 2). For R. subarcticum, the germination rate of outcrossed seeds was significantly lower at constant 15°C (25.0 ± 7.4 percent) than at constant 25°C (75.0 ± 13.9 percent) and 25/10°C (76.7 ± 16.1 percent; Table 2). Thus, R. subarcticum has a high thermal requirement for germination. Furthermore, hybrid seeds produced in R. subarcticum rarely germinated (5.0 ± 4.1 percent), even at 25°C (Figure 7), indicating very low germinability of hybrid seeds.

Figure 7. Germination progress after seed sowing under 25°C condition. Open circles connected with solid line indicate Rhododendron diversipilosum seeds, gray circles with solid line indicate R. subarcticum seeds, and black circles with dashed line indicate hybrid seeds. Mean values of four replicates with standard errors are shown.

Table 2. GLMM results of the seed germination experiment in each species (Rhododendron diversipilosum and Rhododendron subarcticum) under three thermal conditions: constantly high (25°C), constantly low (15°C), and fluctuation (25/10°C).

Variable	R. div.	z	P	R. sub.	z	P
Intercept (high)	1.74 ± 0.26	6.78	***	1.35 ± 0.50	2.69	***
Low	−0.18 ± 0.35	−0.52	ns	−2.57 ± 0.70	−3.66	***
Fluctuation	0.91 ± 0.45	2.03	*	−0.02 ± 0.70	−0.03	ns

Estimated coefficient ± standard error, z value, and P value are shown.

***p < .001, *p < .05, ns = not significant.

Genetic diversity and differentiation

Two of the eight primer sets developed for R. diversipilosum (LP_SSR_034 and 003) did not amplify the DNA fragments successfully in R. subarcticum (Table S3). Thus, we used the remaining six primer sets for genetic analysis. A total of 205 leaf samples from the six sites were PCR-amplified and sequenced. In genotyping procedures, R. subarcticum showed at most two microsatellite peaks, indicating a diploid species as described by Lantai and Kihlman (1995). In contrast, R. diversipilosum often showed four microsatellite peaks, indicating a tetraploid species as described by Wakui (2021). Despite careful sampling in the field, some ramets were genetically identified as the same genets, indicating the existence of large clone patches, particularly for R. subarcticum. Therefore, we conducted genetic analyses at the genet level (Table 3).

Table 3. Genetic characteristics of Rhodoendron diversipilosum and Rhodoendron subarcticum populations at each site.

Population	Species	N	GN	Ar	He	Ho
Mt. Tengu	R. div.	20	14	3.91	0.82	0.87
	R. sub.	20	8	2.76	0.54	0.38
Hisago	R. div.	20	15	3.63	0.78	0.67
	R. sub.	21	5	2.83	0.58	0.33
Mt. Muka	R. div.	20	11	3.66	0.79	0.7
	R. sub.	20	16	3.11	0.7	0.47
Mt. Hira	R. div.	20	3	3.59	0.8	0.79
	R. sub.	20	4	1.65	0.33	0.42
Mt. Teshio	R. div.	20	8	4.12	0.85	0.77
Mt. Ishikari	R. sub.	24	23	3.56	0.74	0.35

Sampling was conducted only for R. diversipilosum at Mt. Teshio and R. subarcticum at Mt. Ishikari.

N: number of samples, GN: genet number, A_r: allelic richness, H_e: expected heterozygosity, H_o: observed heterozygosity.

The indices of genetic diversity—that is, allelic richness (A_r), expected heterozygosity (H_e), and observed heterozygosity (H_o)—were larger in R. diversipilosum populations than in R. subarcticum populations (t = 3.01, 3.14, 8.55, p < .05, .05, .001, respectively). The global F_st value among the five sites was lower in R. diversipilosum than in R. subarcticum (0.063 vs. 0.114). The correspondence of geographic distance and pairwise F_st among sites for both species is shown in Table S4. A Mantel test revealed that significant IBD was not detected in either R. diversipilosum (r = 0.272, p = .31) or R. subarcticum (r = 0.164, p = .32; Figure 8).

Figure 8. The correlation between the geographic distance and pairwise F_st value in (a) Rhododendron diversipilosum and (b) R. subarcticum. No significant correlation was found by Mantel test in each species.

Discussion

The present study demonstrated the differences in habitat use pattern and flowering phenology between R. diversipilosum and R. subarcticum, but the spatiotemporal isolation mechanisms were not perfect. Furthermore, asymmetric reproductive interference existed between the species, where the negative effect of heterospecific pollination through inviable hybrid seed production existed only for R. subarcticum. Finally, R. subarcticum showed lower genetic diversity within populations but higher genetic divergence among populations than R. diversipilosum.

Morphological differences and spatiotemporal separations between the species

Rhododendron subarcticum has smaller leaves, shorter plant height, smaller annual leaf production, larger LMA, and longer leaf life span than R. diversipilosum. Strong wind and solar radiation in the fellfield habitat often increase the risk of physical damage and drought stress, resulting in slower growth, lower productivity, and higher risk of shoot mortality in alpine plants (Wardle 1985; Jones 1992; Körner 2003). Under such harsh environments, prostrate growth form, smaller leaf size with higher toughness, and extension of leaf longevity are adaptive responses, particularly in evergreen species (Kudo 1995; Körner 2003). Because LMA is positively related to leaf thickness, toughness, and longevity in evergreen plants (Givnish 1988; Körner et al. 1989), the larger LMA in R. subarcticum, originating from arctic tundra regions, might have preadapted to the fellfield habitat in comparison with R. diversipilosum, originating from boreal regions.

The belt transect survey demonstrated a negative correlation in the abundance between R. diversipilosum and R. subarcticum. This result indicates that the habitat use patterns of these species are different at the microscale level, although spatial segregation was not perfect. Their distributions reflected the existence of P. pumila cover, where R. diversipilosum was positively correlated, and R. subarcticum was negatively correlated. In the fellfield habitat, dwarf P. pumila patches create heterogeneous environments at the microscale. The presence of shrubs in tundra and alpine environments often modifies the growth conditions of other plants by shading, litter accumulation, snow trapping, windbreaks, moisture conditions, and nutrient cycling (Myers-Smith et al. 2011). In particular, snow trapping and windbreak functions may facilitate the growth of alpine plants under fellfield conditions. Rhododendron diversipilosum might take advantage of the facilitation effect of P. pumila in the fellfield habitat, whereas R. subarcticum was distributed to escape from the shading by P. pumila, reflecting the potentially high tolerance to the exposed alpine environment, as mentioned before.

Observation of flowering time in the fellfield habitat indicated that the flowering of R. subarcticum was initiated and completed approximately one week earlier than that of R. diversipilosum, indicating a partial separation of anthesis between the species. Several studies have reported phenological isolation between sympatric-related species in alpine zones (Yang, Gituru, and Guo 2007; Paudel et al. 2018). Because phenological isolation between the species is not perfect, as well as the habitat use pattern, heterospecific pollen transfer is likely to occur under natural conditions. Although the quantitative survey of pollinators was not conducted in this study, both species are frequently visited by syrphid flies and other flies. Assessment of the heterospecific pollination frequency under natural conditions is needed to clarify the effectiveness of the partial phenological differences as a strategy to prevent reproductive interference.

Asymmetric reproductive interference between the species

Reproductive interference is assumed to have various evolutionary and ecological consequences (Kyogoku 2015). Our pollination experiments revealed that both R. diversipilosum and R. subarcticum were self-incompatible, but the reproductive interference between the species was asymmetric. Although R. diversipilosum had a strong heterospecific incompatibility, R. subarcticum showed a high heterospecific compatibility equivalent to the level of conspecific pollination. This means no physiological barrier in R. subarcticum that prevents crossing with related species. Because hybrid seeds seldom germinated even under warm conditions, there is no reproductive advantage for R. subarcticum producing hybrid seeds. Therefore, the reproductive interference between the species causes a fitness loss only for R. subarcticum.

Asymmetric reproductive interferences between related species have been studied in several plant species. The direction of heterospecific pollen flows varies depending on the flowering order and floral abundance of competing species for the service of pollinators (Runquist and Stanton 2013; Zhang, Montgomery, and Huang 2016). For instance, unilateral pollen flow was common from later flowering species (Silene yunnanensis) to earlier flowering species (Silene asclepiadea) due to protandrous flowering (Zhang, Montgomery, and Huang 2016). Because we did not measure pollen flows between R. diversipilosum and R. subarcticum in this study, the existence of asymmetric reproductive interference during the pollination process is unclear. As the mechanisms of crossing barriers after heterospecific pollination, differences in pollen germinability on stigma, stigma/style lengths, physiological compatibility, and differential fruit abortion between conspecific and heterospecific crossings are reported (reviewed in Tiffin, Olson, and Moyle 2001). Unfortunately, any physiological or genetic mechanism of the asymmetric seed production between R. subarcticum and R. diversipilosum is not clear in this study. One possible reason for the low germinability of hybrid seeds is a triploid block (reproducing nonviable progeny of triploid; Köhler, Scheid, and Erilova 2010) caused by a crossing between diploid R. subarcticum and tetraploid R. diversipilosum. Further studies on the mechanisms of the asymmetric reproductive interference in these species are needed.

The reproductive interference between R. subarcticum and R. diversipilosum might be reduced by spatial segregation at the microscale and phenological isolation of flowering periods. Although spatial segregation is effective to reduce reproductive interference (Kuno 1992; Christie and Strauss 2020), the acceleration of phenological isolation by earlier flowering might decrease the seed set success of R. subarcticum. The pollination experiment indicated pollen limitation only in R. subarcticum. Because the flowering onset of R. subarcticum started approximately one week earlier than that of R. diversipilosum, R. subarcticum might suffer from low pollinator activity due to cool weather conditions in the early season (Kudo and Suzuki 2002). Furthermore, the risk of frost damage to flowers increases during early flowering (Kudo 2021). Actually, the frequency of frost damage was higher in R. subarcticum in 2018. Thus, there is a trade-off between phenological isolation to reduce harmful heterospecific pollen flow and the risk of frost damage and/or pollen limitation in R. subarcticum.

Comparisons of genetic properties between the species

Rhododendron subarcticum populations had lower genetic diversity than R. diversipilosum. This genetic property may reflect isolated distribution of R. subarcticum populations. Although R. subarcticum grows only in the exposed fellfield habitat of restricted mountains, R. diversipilosum is widely distributed across various environments and elevations in Hokkaido. Because the population size of R. subarcticum is much smaller than that of R. diversipilosum, R. subarcticum may have experienced strong genetic drift in each mountain. Small genetic variation and clonal growth in R. subarcticum may cause not only deleterious effects of harmful recessive genes but also mating failure. Because R. subarcticum is self-incompatible, loss of genetic variation may lead to a reduction in seed production. It is known that self-incompatible plants tend to decrease sexual seed production with decreasing genetic diversity of the populations (Bushakra et al. 1999; Uesugi et al. 2020). Another possibility for the low genetic diversity in R. subarcticum is the difference in the ploidy levels between the species. Diploids often maintain lower genetic diversity than polyploids in isolated populations (Rosche et al. 2016) because of tetraploid polysomic inheritance (Soltis and Soltis 2000), which maintains a higher genetic diversity within the polyploid population. If the effects of genetic drift are more intense in diploid species, R. subarcticum may be at a higher risk of local extinction.

Large genetic differentiation among R. subarcticum populations may reflect the isolated distributions. In contrast, small genetic differentiation in R. diversipilosum populations may be due to the continuous distribution across wider elevations and habitat types, by which stepwise gene flows are possible between mountain populations. In addition, significant IBD was not detected in R. subarcticum although gene flow is expected to occur depending on the geographic distance (Ellstrand and Elam 1993). Severe genetic drift with random allele fixation can cause high genetic differentiation (Cai et al. 2020). In R. subarcticum, different alleles might have been fixed randomly in each isolated population because of geographic barriers after migration in the last glacial period.

Conclusions

This study demonstrated asymmetric reproductive interference between R. diversipilosum and R. subarcticum (widespread vs. rare species) in the alpine fellfield habitat, where only R. subarcticum suffers from heterospecific pollination and wasted seed production by hybridization. Habitat use differentiation at the microscale and partial phenological isolation may mitigate reproductive interference, resulting in the coexistence of both species in alpine fellfields. However, low genetic diversity and unique genetic composition in R. subarcticum populations suggest a higher risk of local extinction at the southern distribution edges under global warming. Shifts in the distribution range and flowering phenology of individual species may change the relationships between interacting species. Further studies focusing on the impact of climate change on the interspecific interactions for mating properties of alpine plants are crucial.

Acknowledgments

We are grateful to A. Shibata, A. Wakui, N. Oshima, J. Shigeeda, and Y. Ito for their cooperation in field research.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Supplementary material

Supplemental data for this article can be accessed online at https://doi.org/10.1080/15230430.2023.2188716.

Additional information

Funding

This study was supported by the Sumitomo Foundation [203058], a Shirataki Geopark Research Grant 2019, and JSPS KAKENHI [22H02695].