INTRODUCTION
Resource polymorphism is a phenomenon whereby a single species exhibits two or more morphs (morphotypes) showing differences in morphology, behaviour, coloration, or life history characteristics (e.g., growth characters) (Smith & Skulason, 1996). These differences are considered to be an adaptation to habitat heterogeneity through the differentiation of feeding biology and habitat utilization. Resource polymorphism is not uncommon in freshwater fish species; examples include Gasterosteus (McPhail, 1992),Coregonus (Ostbye et al., 2005; Ostbye et al., 2006)(Ostbye et al. 2005; Ostbye et al. 2006), and Salvelinus species (Jonsson & Jonsson, 2001; Savvaitova, 1995). For example, two trophic morphs ofGasterosteus aculeatus coexist in several postglacial lakes, where one morph utilizes the littoral area and feeds on benthic invertebrates while the other occupies the open-water area and feeds on zooplankton (Schluter & Mcphail, 1992). The morphological adaptations of such morphs appear to be similar among cases; that is, the planktivorous morph often has longer and denser gill rakers, an elongated head and a longer lower jaw, while the benthivorous morph is characterized by shorter and sparser gill rakers, a blunt, round head and a shorter lower jaw (Fraser, Adams, & Huntingford, 1998; Jonsson & Jonsson, 2001; Smith & Skulason, 1996). In addition to the typical benthivorous-planktivorous pairs, three or more morphs occur in some lake systems. For instance, five and four sympatric morphs have been found in Coregonus lavaretusin Fennoscandia and Salvelinus alpinus in Thingvallavatn, respectively (Hudson, Vonlanthen, & Seehausen, 2010; Jonsson & Jonsson, 2001; Praebel et al., 2013).
Phenotypic plasticity, a process by which individuals express alternate phenotypes in response to different environments, is viewed as a mechanism of resource polymorphism (Seehausen et al., 2014; Skulason et al., 2019). Praebel et al. (2013) found that European whitefish (C. lavaretus ) formed three sympatric morphs that differed in gill raker counts due to rapid adaptive radiation into the littoral, pelagic and profundal lacustrine environments in three northern Fennoscandian lakes. The differences in size, diet, and jaw features between lean and huronicus morphs of lake charr (Salvelinus namaycush ) are typical examples of adaptation to shallow- and deep-water environmental conditions (Chavarie et al., 2016). The morphs of the African barb (Labeobarbus gananensis complex), which differ in mouth morphology, gill rakers, diet and gut length, display a typical example of adaptation involving distinct feeding strategies (Levin et al., 2018). Muir et al. (2016) evaluated four morphs ofS. namaycush in many small post-glacial lakes throughout the Holarctic. Although there were differences in the number of morphs among lakes, the observed morphological variation resulted from adaptation to the separation of food resources (piscivorous and invertivorous feeding strategies), which was driven by lacustrine environmental diversity.
Cyprinids, the largest family of vertebrates, also exhibit resource polymorphism, with multiple morphs occurring in a single freshwater system. In the Genale River (Ethiopian highlands, East Africa), a barb (L. gananensis ) complex was found to be composed of six forms, five of which were related to mouth morphology, which represents a typical form of adaptive radiation in response to different resources (Levin et al., 2018). Schizothoracins, the largest and most diverse group of the Qinghai-Tibetan Plateau ichthyofauna (Chen & Cao, 2000), also show morphological diversity within species or even sympatric speciation within a single lake. For instance, four morphs ofSchizopygopsis stoliczkai were described in Lake Yashilkul, Pamir (Savvaitova et al., 1987). Furthermore, Zhao et al. (2009) demonstrated sympatric speciation between Gymnocypris eckloni eckloni andG. e. scoliostomus in Lake Sunmcuo, a small glacial lake on the Tibetan Plateau.
Two morphs were found in Schizopygopsis thermalis Herzenstein 1891 (Cyprinidae: Schizothoracinae) in Lake Amdo Tsonak Co during our team field investigations in 2014-2018. The two morphs correspond to a resource axis in the lake: one form (planktivorous) predominately feeds on plankton and inhabits pelagic lake habitats, and the other form (benthivorous) mainly feeds on benthos and dwells in the benthic zone and a tributary (Nagchu River) of the lake. The former morph possesses a normal lower jaw, a terminal mouth, and moderately or highly dense gill rakers, while the latter morph is characterized by a shortened lower jaw, an inferior mouth with a sharpened horny edge on the lower jaw and short and sparse gill rakers.
However, the two morphs of S. thermalis have not been verified via morphological analysis. In addition, their biological characteristics (e.g., growth, feeding habit, and reproductive traits) are still unclear. To address these gaps in knowledge and better understand the ecological mechanisms of polymorphism in S. thermalis , the specific objectives of this study were to (1) characterize the morphological variation of the two significantly distinct morphs coexisting in Lake Amdo Tsonak Co and quantitatively analyse their morphological characters by a combination of morphometric and traditional linear measures; (2) define the biological characters (e.g., growth, feeding habit, and reproductive traits) of the two morphs of S. thermalis ; and (3) elucidate the potential ecological mechanism of resource polymorphism inS. thermali s.
MATERIAL AND METHODS
2.1 Study area
Lake Amdo Tsonak Co (31.55-32.08°N, 91.25-91.33°E) is a resource-poor and freshwater headwater lake of the Salween River system on the Tibetan Plateau, with a surface area of 182 km2, an elevation of 4,587 m above sea level, and a maximum depth of over 20 m. The annual mean air temperature is -3-0 °C, the annual total rainfall is 350-420 mm, the mean conductivity is 581.5 μs/cm (534.2-563.6 μs/cm), the mean salinity is 0.29‰, the average pH is 8.61 (range: 8.5-8.7), the average oxygen concentration is 6.53 mg/L (range: 5.63-7.21 mg/L), the mean phosphorus concentration is 0.028 mg/L, and the mean nitrogen concentration is 0.59 mg/L. The water chemistry values were measured and assessed by our colleagues.
2.2 Field sampling
Fishes were collected from Lake Amdo Tsonak Co and its tributary (Nagchu River) (Figure 1) in May and September 2017 and April and July 2018 with gill nets (mesh size: 30 mm), cast nets and hand nets. We labelled and photographed the left lateral side of each captured fish in the field. The standard length (SL, 0.1 mm), total length (TL, 0.1 mm), total weight (W ), sex, and stage of maturity of each specimen were recorded in the field. Only adult individuals were considered in the subsequent analyses due to the less obvious characteristics of the juveniles. We also measured gill raker number and length and pharyngeal tooth row number and recorded the presence or absence of parasites and the degree of mucus cavity development. In addition, we also measured the width of the horny edge of the lower jaw of some samples with the method shown in Figure 2 and described their lower jaw traits.
2.3 Morphologicalanalysis
With respect to the analysis of both morphometric and linear traits, 154 specimens were analysed in this study. The left lateral side of the specimens was photographed using a Canon PowerShot G12 camera under natural light conditions with a ruler for scale. To obtain meaningful landmarks (homologous points), specimens were placed with the fins straightened and mouth closed in an orthogonal position. Specifically, the camera was fixed on a tripod with the lens parallel to the surface of the samples. Data on body shape and linear measurements for each fish were collected using digital landmarks on the photographs. A total of 24 landmarks were marked on each photograph (Figure 3) with tpsUtil and tpsDig2 (Rohlf, 2010). Body shapes were estimated by extracting the individual centroid size in MorphoJ 1.06 (Klingenberg, 2008). Simultaneously, new coordinates (xy ) for each fish were extracted using generalized Procrustes superimposition for subsequent analysis (Markevich, Esin, & Anisimova, 2018). To increase the accuracy of body shape estimates, some linear head traits (e.g., snout length, upper jaw length, lower jaw length, head length, head height and distance between the anterior termini of the upper and lower jaws) were measured. Allometry is common in fish, indicating that morphology and body size are typically related (Zelditch, Swiderski, & Sheets, 2004). Thus, a multivariate regression of body shape (Procrustes coordinates) on centroid size was used to correct for allometric effects, and regression residuals were used in geometric morphometric analyses (Elmer et al., 2014). Principal component analysis (PCA) was conducted via MorphoJ 1.06 to assess body shape variation (a geometric morphometric trait) among individuals without a priori grouping and to capture the maximum amount of variation with the smallest number of variables. The abbreviations and standardization of the linear head traits were as follows: HeadL = head length, SnoutL = snout length, Upper2 = upper jaw length, LowerL = lower jaw length and JawD = the distance between the anterior termini of the upper and lower jaws, which were normalized by SL , and HeadH = head height and UpperL1 = upper jaw length, which were standardized by head length and lower jaw length, respectively. Morphs were initially identified with unweighted pair-group method with arithmetic mean (UPGMA) cluster analysis using Past 3.2.6 (http://folk.uio.no/ohammer/past/) based on 7 linear traits and scores from the first two principal components (PCs) of body shape from geometric morphometric analysis (Chavarie, Howland, & Tonn, 2013). To quantify the importance of each variable for the ordination axes and thus to summarize the variation in the morphs identified with the cluster analysis, PCA of body shape PCs and linear traits was then conducted. This procedure was performed with the ggbiplot package in R (version 0.55; Vu, 2011). Discriminant function analysis (DFA) and posterior jackknife cross-validation were performed with SPSS 22.0 on morphs defined by cluster analysis to determine whether the morphs were significantly different. The efficacy of the DFA was evaluated with Wilks’ λ, which varies between 0 and 1, with zero indicating perfect identification. Finally, we performed multivariate analysis of variance (MANOVA) with Tukey’s honestly significant difference (HSD) post hoc comparisons (SPSS 22.0) on body shape PCs and linear traits of samples between morphs to test the validity of the discriminant analysis results. Morphological variation between the sexes was estimated by implementing t-tests of body shape PCs and linear traits. To visualize body shape differences between morphs, we reconstructed body shapes using landmark coordinates based on DFA of only geometric morphometric data (Procrustes coordinates) (Elmer et al., 2014; Jakubavičiūtė, De Blick, Dainys, Ložys, & Olsson, 2018). For countable variables (e.g., gill raker number and pharyngeal tooth row number) and the measurable variable (e.g., gill raker length, which was standardized by SL), the Kolmogorov-Smirnov test and Levene’s test were performed. Analysis of variance (ANOVA) was conducted for variables that were normally distributed with variance homogeneity, and nonparametric tests were implemented for variables (after logx transformation) that showed a non-normal distribution or variance heterogeneity. To test whether the presence or absence of parasites affected body shape, analysis of covariance (ANCOVA) was performed using parasites as a covariate, morph as a factor and linear traits and PCs as dependent variables.
2.4 Diet analysis
To evaluate the feeding habits of the fish, we inspected the gut contents (since cyprinid fish do not have an obvious stomach, we chose the anterior intestine: the first bend from the pharynx to the intestine, where the food had not yet been digested) of 25 individuals of one morph (planktivorous) and 30 individuals of the other morph (benthivorous) and dissected them under a dissecting stereomicroscope before other gelatinous substances were excluded, such as gastric juices. Prey items were identified at the phylum or genus level. Then, we divided all prey items into six categories: zooplankton, small fishes (including fishes and their remains), hydrophilic insects, periphytic algae, zoobenthos and others (hydrophyte debris, organic debris and small grains of sand). First, the diet compositions were estimated by occurrence rate and wet weight percentage.
F % = (Ni /Ntotal ) × 100%
F % is the percentage of occurrence of prey i ,Ni is the frequency of occurrence of preyi, andNtotalrepresents the total number of gut samples with food.
W % = (Wc /Wtotal ) ×100%
W % is the wet weight percentage of prey category c (one of the six food categories), Wc represents the wet weight of prey category c , and Wtotalrepresents the total weight of all prey items in each sample. Schoener’s index (Dxy ) (Schoener, 1970) of proportional diet overlap was calculated and used to evaluate the difference in food composition between the two morphs.
Cxy =1﹣0.5{∑|PxcPyc |}
Cxy represents the diet overlap index between the two forms (x and y ). Pxc andPyc represent the shared food category cof form x and y (W % ), respectively. Values range between 0 (no diet overlap) and 1 (complete diet overlap), and values greater than 0.6 generally indicate biologically significant overlap (Wallace, 1981). Finally, the nonparametric Kruskal-Wallis test (Zar, 1999) was implemented to estimate the food composition (W % ) difference between the two morphs.
2.5 Growth analysis
In total, 140 specimens (lapillus otoliths) were used to study the growth characteristics of the two morphs (planktivorous, N = 66; benthivorous, N = 74). Preparation of otolith sections and age determination were performed by an experienced worker. Each otolith was interpreted three times, and otoliths without at least 2 identical interpretations were excluded from the analysis. Photographs of all otolith sections were captured using MicroPublisher (5.0 RTV) under a light microscope (BH2; Olympus Optical, Tokyo, Japan). Otolith radius and ring diameter were measured with Autook (Image Analyser 2.0).
The relationship of SL with otolith radius was described by Frase-Lee regression, and the back-calculated SL of all ages was obtained using the modified Frase-Lee function (Johnson & Noltie, 1997):
logeLi = a + (logeLc -a )(logeOi /logeOc )
where Lc is the SL of the specimens,a is the intercept of the SL -otolith radius regression,Oc is the radius of the otolith at capture,Oi is the radius of the otolith at age iand Li is the back-calculated SL at agei .
Back-calculated SL was used to fit a von Bertalanffy growth function (VBGF) (von Bertanlanffy, 1938) and obtain the growth parameters of each morph:
Lt = L (1 - exp(-k ×(t-t0 )))
where Lt is the back-calculated SL at aget , L is the asymptotic SL ,k is the growth coefficient, t is age, andt0 is the theoretical age at zero length.
To compare growth parameters, the growth performance index (φ ) was calculated according to the equation of Munro & Pauly (1983):
φ = log10K +2log10L
where K is the growth coefficient and Lis the asymptotic SL .
Because of the nonlinear formulation of the VBGFs, a general linear model could not be used for ANCOVA. Instead, an analysis of residuals of the sum of squares (ARSS) was performed to compare the VBGFs between sexes and morphs (Chen, Jackson, & Harvey, 1992), and the degree of fit was denoted by the correlation coefficient and coefficient of determination (R2 ).