Ecology

We have analyzed the inbreeding-purging process in four captive populations of different ungulate species with effective sizes ranging 4–40 and with available pedigrees as well as survival and productivity records. This allows us to explore the role of inbreeding and purging in determining the evolution of fitness traits in a range of scenarios relevant in the context of conservation.In A. lervia (Ne ≈ 4), purging is expected only for the most severely deleterious alleles (those giving dNe  > 1, which implies d  > 0.25 as, for example, in completely recessive alleles with deleterious homozygous disadvantage s  > 0.5). Thus, it could be that purging has not been detected for this species because such severely deleterious alleles had been purged during the demographic decline in the wild, before the foundation of the captive population. This would be consistent with the low and non-significant inbreeding load estimated in this species. It is also possible that these estimates are non-significant due to the relatively small number of individuals available.G. cuvieri and N. dama have significant initial inbreeding loads that, adding up the direct and maternal components, is about 1.25 in both cases, which is on the order of other estimates published for captive populations (Ralls et al. 1988). Since in both species Ne  > 10, purging should be efficient against less severely deleterious alleles than in A. lervia (d  > 0.1). Purging is detected for both species with very low P values. This result is in agreement with Moreno et al. (2015), who suggested that purging had occurred in G. cuvieri as they found an increased juvenile survival parallel to an increased inbreeding coefficient. The relative contribution of severe and mild deleterious effects to the inbreeding load of populations is under a scientific debate with direct implications in conservation biology (Ralls et al. 2020, Kyriazis et al. 2021, Pérez-Pereira et al. 2021). The large d estimates obtained in our analysis indicate that a substantial fraction of the initial inbreeding load is being purged under modest effective population sizes, implying that such substantial fraction is due to relatively severe deleterious mutations in these two populations. As far as we are aware, these are the first estimates of this purging parameter obtained in managed, non-experimental populations. Previous estimates of d were obtained in D. melanogaster bottlenecked populations, first for egg-to-pupae viability in lines with Ne = 6 or 12 under noncompetitive conditions (d = 0.09, Bersabé and García-Dorado 2013), and second in lines with higher Ne ≈ 40–50 under more competitive conditions, giving a larger estimate of d, of the order of that estimated in these two ungulate endangered species (d ≈ 0.3, López-Cortegano et al. 2016).Regarding G. dorcas, given its larger population size, purging is expected even against alleles with mild recessive component of the deleterious effect (d  > 0.025). However, although a significant (if modest) inbreeding load was estimated, no significant purging was detected. Nevertheless, the number of equivalent complete generations by the end of the pedigree (EqG = 7) was smaller than our proposed minimum number of generations required to detect purging (tm = 10). This suggests that, due to the large size of this population, more generations are needed to detect purging.The results above support the use of tm to get an approximate idea about when a pedigree is too shallow for purging to be detected. Should the number of generations available be larger than tm, IP predictions could additionally be computed to search the d values that can be expected to produce detectable purging. Supplementary Fig. S3 shows that the true number of generations required to detect purging becomes increasingly larger than tm for alleles with smaller d values, as they suffer weaker purging each time they are exposed in homozygosis. The tm approach helps to understand the failure of many studies to detect purging. Such is the case of the extensive meta-analyses on 119 zoo populations by Boakes et al. (2007), where the median Ne value was 22.6 while the median number of generations was t = 3 meaning that, for most species, at least 5 more generations were needed before purging could be detectable. On the contrary, and in agreement with this tm approach, purging was experimentally detected in lines of D. melanogaster with Ne = 43 (i.e., tm ≈ 10) where, after an initial period of inbreeding depression, fitness experienced a substantial recovery beginning between generations 10 and 20 (López-Cortegano et al. 2016).A reason why detecting purging in captive populations is challenging is that a fitness rebound can also be due to adaptation to captive conditions or to environmental effects, such as those derived from improved husbandry (Clifford et al. 2007). In fact, this might have been the case in Speke’s gazelle breeding program, where the observed rebound of fitness was first ascribed to purging (Templeton and Read 1984, 1998), while Kalinowski et al. (2000) suggested that husbandry improvements could also be responsible for these findings. Our estimates of d and δ, however, are based on the association between the fitness trait and purged inbreeding at the individual level (Wi, gi) which, in our data, is mainly expressed within cohorts while average survival showed little variation through time. In addition, the analyses included temporal factors (YOB or POM) that should have removed confounding effects from adaptation to captivity or improved husbandry. Therefore, adaptive processes or time-dependent environmental factors are not expected to have biased our IP estimates.For productivity, the estimates of inbreeding load were high (overall inbreeding load ~5, P value  More

Insects and bamboosFive internodes (length: mean ± SD = 44.8 ± 1.1 cm, n = 5; diameter in the middle part of internodes: 21.4 ± 0.8 mm, n = 5) of five living mature culms of P. simonii bamboo were sampled at Kawaminami, Miyazaki Prefecture, Japan [32°9′ N, 131°29′ E] on 6 June, 2019. Per internode, four semi-cylindrical strips (ca. 15 × 2 cm) were made and stored at − 25 °C until use.To obtain fungus-free larvae of D. bucculenta, we sampled five beetle eggs from P. simonii bamboo obtained at Toyota, Aichi Prefecture, Japan [35°9′ N, 137°13′ E] on 9 May, 2019 in the laboratory from ovipositing females collected at Kawaminami on 10 and 11 April, 2019. The eggs were immersed in 99.5% ethanol for 10 s followed by 70% ethanol for 10 s for surface sterilization and then individually placed on potato dextrose agar (PDA) (Difco, Detroit, MI, USA) plates. The plates were incubated at 25 °C in the dark until 30 days after larval hatching to confirm the absence of the formation of yeast or other microbial colonies. Consequently, all five larvae hatched successfully and aseptically.The bamboo used in this study was morphologically identified using the literature29. This is native to the study areas and no other host bamboo species are distributed there29. Therefore, no voucher specimen of this bamboo has been deposited in a publicly available herbarium. No specific permits were required for the described field studies. The location is not privately-owned or protected in any way. The field studies did not involve endangered or protected species. All applicable international, national, and/or institutional guidelines for the care and use of animals and plants were followed. This study is reported in accordance with ARRIVE guidelines.Component analyses of bamboo tissuesFor YP and LP, the yeast W. anomalus originating from D. bucculenta in Kawaminami (strain: DBL05Kawaminami) was cultured on a 9-cm PDA plate to obtain enough biomass for further experiments. Afterwards, yeast cells were suspended in ca. 10 mL of sterilized water, and were inoculated on the inner surface of the autoclaved internode strips using an autoclaved tissue paper immersed with the yeast suspension. For LP, additionally, the fungus-free 2nd instar larvae (weight: mean ± SD = 2.4 ± 0.4 mg, n = 5) were individually placed on the yeast-inoculated strips. Each of these yeast-inoculated and yeast-and-larva-inoculated strips was then put in a sterilized test tube (3.0 cm in diameter and 20 cm tall) with moistened cotton placed at the bottom. Each of the test tubes was covered with a sterilized polypropylene cap, sealed with Parafilm Sealing Film (Pechiney Plastic Packaging, Chicago, IL, USA) on which three small holes were made using a fire-sterilized insect pin to avoid oxygen shortage, and individually put in a plastic zipper bag. These yeasts and insects were incubated at 25 °C in the dark for 47 days for YP (n = 5), and 47 (n = 4) and 73 (n = 1) days until these larvae reached adulthood for LP (adult elytral length: mean ± SD = 9.2 ± 0.4 mm, n = 5). Microbial contamination was invisible to the naked eye.For FP, YP and LP, the inner surface (up to 0.3 mm in thickness, dry weight: 336 to 935 mg) of a strip was sampled using a small U-shaped gouge. In the case of FX, first, the pith of a strip was completely removed, and then xylem tissue (up to 0.5 mm in thickness, dry weight: 729 to 872 mg) was sampled using a small U-shaped gouge. These tissues were individually sampled from five strips derived from five different internodes for each tissue type.Samples were extracted by aqueous ethanol and hydrolyzed by sulfuric acid with reference to the literature30,31,32 as follows. Four types of samples were freeze-dried and pulverized using a rotor-speed mill (Fritsch, PULVERISETTE 14, 0.2 mm mesh). About 80 mg of powdered sample was extracted using 5-mL 80% ethanol aqueous solution (aq.) at 63 °C three times. The volume of the extracts was adjusted to 25 mL, filtered, and analyzed using ion exchange chromatography measurements (extractable sugar analysis). Their extracted residues were hydrolyzed using sulfuric acid as follows: 50-mg samples were immersed in 1.64-g 72% sulfuric acid aq. at 30 °C for 2 h, boiled in 39.4-g 3% sulfuric acid aq. for 4 h, and filtered to collect sulfuric acid residues as sulfuric acid lignin fractions. The volumes of the filtrates were fixed to 100 mL, passed through a sulfuric acid-removing filter (DIONEX OnGuard IIA), and submitted to ion exchange chromatography measurements (structural sugar analysis). For the uronic acid measurements, the sulfuric acid-removing filter was not used.Ion exchange chromatography measurements were conducted using a DIONEX ICS-3000 apparatus. The measurement conditions were as follows: column, CarboPac PA-1 (2.0 mm I.D. × 250 mm L, Dionex corp.); flow rate, 0.3 mL min−1; column temperature, 30 °C; injection volume, 25 µL; eluent, H2O (solvent A), 100 mM NaOHaq. (solvent B), aqueous solution containing 100 mM NaOH and 1.0 M CH3COONa (solvent C), and aqueous solution containing 100 mM NaOH and 150 mM CH3COONa (solvent D). The gradient conditions for monomers, dimers, and uronic acids were as follows: for monomers, with a gradient of B 0.5% C 0% 45 min, C 100% 10 min, B 100% 10 min, B 0.5% C 0% 20 min; for dimers, with a gradient of B 50% C 0% 50 min, C 100% 10 min, B 100% 10 min, B 50% C 0% 15 min; for uronic acids, with a gradient of D 100% 10 min. These extraction, hydrolysis, and measurement procedures were conducted using n = 5 samples. For the structural sugars, their yield was calculated as the dehydrated state. The values of other extractives % were calculated by the subtraction of total extractable sugars % from total extractives %.Elemental analysis (carbon, hydrogen, nitrogen) was conducted by 2400 CHNS Organic Elemental Analyzer (PerkinElmer Japan, Yokohama, Japan). About 1-mg dried samples were burned completely and the produced CO2, H2O, and N2 (after reduction of NOx species) gasses were quantified by a thermal conductivity detector.Means of components of bamboo tissues were compared among tissue types using the Steel–Dwass test after the Kruskal–Wallis test. Calculations were performed using R 3.5.133.Carbon assimilation testThe yeast W. anomalus (DBL05Kawaminami) was cultured aerobically in 20 mL of yeast nitrogen base (YNB) (Difco) containing 0.5% glucose at 25 °C in the dark for 2 days with shaking at 85 rpm. The culture media were centrifuged and cell pellets were suspended in sterile water, in which the OD600 was adjusted to 0.10. Fifty μL of the cell suspension was added into a tube (2 mL) with 1 mL of each of 14 different media containing YNB and one of the following carbon sources: d-glucose, d-galactose, d-mannose, d-xylose, l-arabinose, d-fructose, d-galacturonic acid, d-glucuronic acid, sucrose, cellobiose, starch from corn, xylan from corn, carboxymethyl cellulose, and no carbon source (n = 5 to 6). The concentration of each carbon source was 0.5 g L−1, except for xylan at 1.5 g L−1. The tubes were shaken at 85 rpm and incubated at 25 °C in the dark for 7 days. Afterwards, the presence of visible pellets of yeasts and OD600 were recorded to determine the growth of the strain. The degree of assimilation was scored according to the presence of the pellets and the difference in the turbidity increase (ΔOD600) between culture media containing no and a given carbon source as follows: no growth (without a pellet, ΔOD600  More

1.Hughes, T. et al. Coral reefs in the Anthropocene. Nature 546, 82–90 (2017).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
2.Moberg, F. & Folke, C. Ecological goods and services of coral reef ecosystems. Ecol. Econ. 29, 215–233 (1999).Article 

Google Scholar 
3.Brander, L. M., Van Beukering, P. & Cesar, H. S. The recreational value of coral reefs: A meta-analysis. Ecol. Econ. 63, 209–218 (2007).Article 

Google Scholar 
4.Jackson, J. B. et al. Historical overfishing and the recent collapse of coastal ecosystems. Science 293, 629–637 (2001).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
5.Bruno, J. F., Sweatman, H., Precht, W. F., Selig, E. R. & Schutte, V. G. Assessing evidence of phase shifts from coral to macroalgal dominance on coral reefs. Ecology 90, 1478–1484 (2009).PubMed 
Article 
PubMed Central 

Google Scholar 
6.Hughes, T. Catastrophes, phase shifts, and large-scale degradation of a Caribbean Coral Reef. Science 265, 1547–1551 (1994).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
7.Hughes, T. P., Bellwood, D. R., Folke, C. S., McCook, L. J. & Pandolfi, J. M. No-take areas, herbivory and coral reef resilience. Trends Ecol. Evol. 22, 1–3 (2007).PubMed 
Article 
PubMed Central 

Google Scholar 
8.Gardner, T. A., Côté, I. M., Gill, J. A., Grant, A. & Watkinson, A. R. Long-term region-wide declines in Caribbean corals. Science 301, 958–960 (2003).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
9.Edgar, G. J. et al. Reef Life Survey: Establishing the ecological basis for conservation of shallow marine life. Biol. Conserv. 252, 108855 (2020).Article 

Google Scholar 
10.Chabanet, P., Bigot, L., Garnier, R., Tessier, E. & Moyne-Picard, M. Coral reef monitoring at Reunion island (Western Indian Ocean) using the GCRMN method. Proc. 9th Int. Coral Reef Symp. 2, 873–878 (2000).
Google Scholar 
11.Lang, J. C., Marks, K. W., Kramer, P. A., Kramer, P. R. & Ginsburg, R. N. AGRRA Protocols Version 5.4. (2010).12.Cortés, J. et al. The CARICOMP network of Caribbean Marine Laboratories (1985–2007): History, key findings, and lessons learned. Front. Mar. Sci. 5, 519 (2019).Article 

Google Scholar 
13.Dethier, M. N., Graham, E. S., Cohen, S. & Tear, L. M. Visual versus random-point percent cover estimations: ‘objective’ is not always better. Mar. Ecol. Prog. Ser. 96, 93–100 (1993).Article 

Google Scholar 
14.Beijbom, O., Edmunds, P. J., Kline, D. I., Mitchell, B. G. & Kriegman, D. Automated annotation of coral reef survey images. In 2012 IEEE Conference on Computer Vision and Pattern Recognition 1170–1177 (IEEE, 2012).15.Beijbom, O. et al. Towards automated annotation of benthic survey images: Variability of human experts and operational modes of automation. PloS One 10, e0130312 (2015).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
16.Williams, I. D. et al. Leveraging automated image analysis tools to transform our capacity to assess status and trends of coral reefs. Front. Mar. Sci. 6, 222 (2019).Article 

Google Scholar 
17.Roelfsema, C. et al. Benthic and coral reef community field data for Heron Reef, Southern Great Barrier Reef, Australia, 2002–2018. Sci. Data 8, 1–7 (2021).MathSciNet 
Article 

Google Scholar 
18.Gonzalez-Rivero, M. et al. Monitoring of coral reefs using artificial intelligence: A feasible and cost-effective approach. Remote Sens. 12, 489 (2020).Article 

Google Scholar 
19.Cairns, S. D. Stony corals (Cnidaria: Hydrozoa, Scleractinia) of Carrie Bow Cay, Belize. Smithson. Contrib. Mar. Sci. 21, 271–302 (1982).
Google Scholar 
20.Rutzler, K. & Macintyre, I. G. The Atlantic Barrier Reef Ecosystem at Carrie Bow Cay, Belize, 1: Structure and Communities (Smithsonian Institution Press, 1982). https://doi.org/10.5479/si.01960768.12.539.Book 

Google Scholar 
21.Rutzler, K. Caribbean coral reef ecosystems: Thirty-five years of smithsonian marine science in Belize. In Proceedings of the Smithsonian Marine Science Symposium (2009).22.McField, M. et al. Mesoamerican Reef Report Card. (2020).23.Cox, C. E. et al. Genetic testing reveals some mislabeling but general compliance with a ban on herbivorous fish harvesting in Belize. Conserv. Lett. 6, 132–140 (2013).Article 

Google Scholar 
24.Pebesma, E. J. & Bivand, R. S. Classes and methods for spatial data in R. R News 5, 9–13 (2005).
Google Scholar 
25.Pebesma, E. Simple features for R: Standardized support for spatial vector data. R J. 10, 439–446 (2018).Article 

Google Scholar 
26.R Core Team. A language and environment for statistical computing. (2020).27.Althaus, F. et al. A standardised vocabulary for identifying benthic biota and substrata from underwater imagery: The CATAMI classification scheme. PLoS One 10, e0141039 (2015).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
28.Goatley, C. H., Bonaldo, R. M., Fox, R. J. & Bellwood, D. R. Sediments and herbivory as sensitive indicators of coral reef degradation. Ecol. Soc. 21(1), 29 (2016).Article 

Google Scholar 
29.Connell, S., Foster, M. & Airoldi, L. What are algal turfs? Towards a better description of turfs. Mar. Ecol. Prog. Ser. 495, 299–307 (2014).Article 

Google Scholar 
30.Lozada-Misa, P., Schumacher, B. D. & Vargas-Angel, B. Analysis of benthic survey images via CoralNet: A summary of standard operating procedures and guidelines. Pacific
Islands Fish. Sci. Cent. Natl. Mar. Fish. Serv. https://doi.org/10.7289/V5%2FAR-PIFSC-H-17-02 (2017).Article 

Google Scholar 
31.Obura, D. & Grimsditch, G. Resilience Assessment of Coral Reefs: Assessment Protocol for Coral Reefs, Focusing on Coral Bleaching and Thermal Stress (Citeseer, 2009).
Google Scholar 
32.Broeke, J., Pérez, J. M. M., & Pascau, J. Image processing with ImageJ. (Packt Publishing Ltd, 2015).
Google Scholar 
33.Wood, S. N. Generalized Additive Models: An Introduction with R (CRC Press, 2017).MATH 
Book 

Google Scholar 
34.Fasiolo, M., Nedellec, R., Goude, Y. & Wood, S. N. Scalable visualization methods for modern generalized additive models. J. Comput. Graph. Stat. 29, 78–86 (2020).MathSciNet 
Article 

Google Scholar 
35.Oksanen, J. et al. Community ecology package. R Package Version 2, (2013).36.Arnold, S. N. & Steneck, R. S. Settling into an increasingly hostile world: The rapidly closing ‘“Recruitment Window”’ for corals. PLoS One. https://doi.org/10.1371/journal.pone.0028681 (2011).Article 
PubMed 
PubMed Central 

Google Scholar 
37.Mumby, P. J. & Harborne, A. R. Marine reserves enhance the recovery of corals on Caribbean reefs. PLoS One 5, e8657 (2010).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
38.Adam, T. C., Burkepile, D. E., Ruttenberg, B. I. & Paddack, M. J. Herbivory and the resilience of Caribbean coral reefs: Knowledge gaps and implications for management. Mar. Ecol. Prog. Ser. 520, 1–20 (2015).Article 

Google Scholar 
39.Suchley, A., McField, M. D. & Alvarez-Filip, L. Rapidly increasing macroalgal cover not related to herbivorous fishes on Mesoamerican reefs. PeerJ 4, e2084 (2016).PubMed 
PubMed Central 
Article 

Google Scholar 
40.Arnold, S. N., Steneck, R. S. & Mumby, P. J. Running the gauntlet: Inhibitory effects of algal turfs on the processes of coral recruitment. Mar. Ecol. Prog. Ser. 414, 91–105 (2010).Article 

Google Scholar 
41.Box, S. J. & Mumby, P. J. Effect of macroalgal competition on growth and survival of juvenile Caribbean corals. Mar. Ecol. Prog. Ser. 342, 139–149 (2007).Article 

Google Scholar 
42.Williams, I. & Polunin, N. Large-scale associations between macroalgal cover and grazer biomass on mid-depth reefs in the Caribbean. Coral Reefs 19, 358–366 (2001).Article 

Google Scholar 
43.Newman, M. J., Paredes, G. A., Sala, E. & Jackson, J. B. Structure of Caribbean coral reef communities across a large gradient of fish biomass. Ecol. Lett. 9, 1216–1227 (2006).PubMed 
Article 
PubMed Central 

Google Scholar 
44.Mumby, P. J., Steneck, R. S., Roff, G. & Paul, V. J. Marine reserves, fisheries ban, and 20 years of positive change in a coral reef ecosystem. Conserv. Biol. https://doi.org/10.1111/cobi.13738 (2021).Article 
PubMed 
PubMed Central 

Google Scholar 
45.Aronson, R., Precht, W., Toscano, M. & Koltes, K. The 1998 bleaching event and its aftermath on a coral reef in Belize. Mar. Biol. 141, 435–447 (2002).Article 

Google Scholar 
46.Green, D. H., Edmunds, P. J. & Carpenter, R. C. Increasing relative abundance of Porites astreoides on Caribbean reefs mediated by an overall decline in coral cover. Mar. Ecol. Prog. Ser. 359, 1–10 (2008).Article 

Google Scholar 
47.Roff, G., Joseph, J. & Mumby, P. J. Multi-decadal changes in structural complexity following mass coral mortality on a Caribbean reef. Biogeosciences 17, 5909–5918 (2020).Article 

Google Scholar 
48.Graham, N. & Nash, K. The importance of structural complexity in coral reef ecosystems. Coral Reefs 32, 315–326 (2013).Article 

Google Scholar 
49.Aronson, R. B. & Precht, W. F. White-band disease and the changing face of Caribbean coral reefs. Ecol. Etiol. New. Emerg. Mar. Dis. 159, 25–38 (2001).
Google Scholar 
50.Aronson, R. B., Macintyre, I. G., Precht, W. F., Murdoch, T. J. & Wapnick, C. M. The expanding scale of species turnover events on coral reefs in Belize. Ecol. Monogr. 72, 233–249 (2002).Article 

Google Scholar 
51.McField, M. et al. Status of the Mesoamerican Reef after the 2005 coral bleaching event. Status Caribb. Coral Reefs Bleach. Hurric. In 45–60 (2005).52.Arias-González, J. E. et al. A coral-algal phase shift in Mesoamerica not driven by changes in herbivorous fish abundance. PLoS One 12, e0174855 (2017).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
53.Miller, S. Automatically Annotating 175,000+ Images with the CoralNet API. CoralNet (Accessed 23 August 2021); https://coralnet.ucsd.edu/blog/automatically-annotating-175000-images-with-the-coralnet-api/ (2020).54.Muller, E. M., Sartor, C., Alcaraz, N. I. & van Woesik, R. Spatial epidemiology of the stony-coral-tissue-loss disease in Florida. Front. Mar. Sci. 7, 163 (2020).Article 

Google Scholar 
55.Alvarez-Filip, L., Estrada-Saldívar, N., Pérez-Cervantes, E., Molina-Hernández, A. & González-Barrios, F. J. A rapid spread of the stony coral tissue loss disease outbreak in the Mexican Caribbean. PeerJ 7, e8069 (2019).PubMed 
PubMed Central 
Article 

Google Scholar 
56.Weil, E. et al. Spread of the new coral disease “SCTLD” into the Caribbean: implications for Puerto Rico. Reef Encount. 34, 38–43 (2019).
Google Scholar 
57.Heres, M. M., Farmer, B. H., Elmer, F. & Hertler, H. Ecological consequences of Stony Coral Tissue Loss Disease in the Turks and Caicos Islands. Coral Reefs 40, 609–624 (2021).Article 

Google Scholar 
58.Walton, C. J., Hayes, N. K. & Gilliam, D. S. Impacts of a regional, multi-year, multi-species coral disease outbreak in Southeast Florida. Front. Mar. Sci. 5, 323 (2018).Article 

Google Scholar  More

Morphology and light-emitting behavior of the undescribed Japanese Terebellidae wormIn 2016, some of the present authors were exploring shallow coastal waters (depth less than 1 m) to observe the ecological behaviors of marine animals in the Noto Peninsula, when they discovered unknown violet-light-emitting worms. At the sampling point, the worms were living in small holes (a few centimeters in diameter) or in cracks in rocks covered by sand at the shallow sea bottom (Supplementary Fig. S1). We successfully video-recorded their emission of violet light from the whole tentacle stretching into sea water when stimulated by air bubbling at night (Fig. 1A–C; Supplementary Videos S1 and S2). The violet-light emission consisted of rapid flashes with variable duration in the order of milliseconds (Supplementary Video S3), as observed for the worm P. perplexus in response to stimulation17. From our morphological observation, we identified the violet-light-emitting worm as a member of Polycirrus on the basis of the following characteristics18: (1) a sheetlike prostomium covering the upper lip; (2) avicular unicini on some neuropodia; (3) no branchiae. The specimens also have the following characteristics: (1) neurochaetae beginning on last notochaetigerous segment, chaetiger 14; (2) uncini with a long neck and concave base; (3) notopodial pre- and post-chaetal lobes both similar shape. These characters are also found in Polycirrus disjunctus Hutchings and Glasby18; however some of the characters in parapodial lobes and chaetae have differentiation. Thus, we concluded that this species should be treated as an undescribed species. Further comparative observation is needed to describe the species. At this time we treated the Polycirrus species observed in this study as Polycirrus sp. ISK. Application of an electric pulse also caused clear light emission from the tentacles of the living worm (Fig. 1D; Supplementary Video S3), and the luminescence spectrum showed that its λmax was 444 nm or slightly longer, depending on the individual (Fig. 2A). We also found that light emission was efficiently induced by the addition of KCl solution and observed the time course of light emission with rapid fluctuations with variable duration in the order of milliseconds for up to 30 s (Fig. 2B). The flash pattern was similar to that observed in a study of P. perplexus17. In the genus Polycirrus, P. medius and P. nervosus in Japan have been described18,19. However, the morphological features of the species in the present study differed from these species on the basis of our observations described above.Figure 1Photographs of Polycirrus sp. ISK. (A) Polycirrus sp. ISK in its natural habitat with bright-field illumination. (B) Bioluminescence of Polycirrus sp. ISK in its natural habitat without bright-field illumination. The worms were stimulated by air bubbling from SCUBA gear. (C) A single worm with stretched tentacles. Tentacles are indicated by white arrows. (D) The worm with light emission at the tentacles. This worm was stimulated by an electric shock. Scale bars = 100 mm for A and B, 10 mm for C and D. Each photograph was extracted from the videos recorded with the following settings: sensitivity, ISO 51200 or 11 lx; white balance, 4300 K or 5800 K; shutter speed, 1/30 s or 1/60 s; iris, F1.8-3.5; frame rate, 29.97 fps or 60 i; frame size, 1920 × 1080 pixels. Original high quality videos are available at https://youtu.be/KEsU0kWAEfg and https://youtu.be/24dxvPlBDB0Full size imageFigure 2Luminescence spectra and KCl-induced light emission of Polycirrus sp. ISK. (A) Spectrum analysis of Polycirrus sp. ISK using a living worm stimulated by an electric shock. The luminescence spectra were obtained from two different individuals. The λmax represented in closed circles and open circles were 444 nm and 446 nm, respectively. (B) Typical light-emission signal of a living worm soaked in 667 mM KCl. The black line indicates luminescence intensity after adding KCl solution, and the gray line indicates luminescence intensity before adding KCl solution.Full size imageJapanese Polycirrus spp. have not been described as luminous worms according to our review of the literature and web pages. In addition, the number of reports for new Polycirrus spp. from all over the world has been increasing, but a limited number of species are known to emit light13,17,18. Our finding of KCl-induced light emission from Polycirrus sp. ISK suggested that we can easily test the light-emitting ability of Polycirrus spp. by luminescence measurement just after adding KCl solution. A spectrum pattern has been reported for only one species, P. perplexus collected in California17, and it would be necessary for further understanding of these species to examine the light-emitting abilities and to compare light-emitting behaviors and spectrum patterns. The color of bioluminescence is often related to habitat, and light in the blue range is typical for pelagic species20. Thus, one of the points to be focused on is the ecological function of the violet-light emission of this worm inhabiting in a shallow coastal water environment. In P. perplexus, deterring predation is a possible function of luminescence based on that species’ habitat and its violet-light emission17,21. As shown in Supplementary Videos S1, S2, which are the first video records of in situ light emission of a Polycirrus species, the air bubble-stimulated luminescence of Polycirrus sp. ISK in its natural habitat also seemed to deter predation, but this explanation is still speculative.Differentially expressing genes between the tentacles and the rest of bodyA few years after discovering this worm, we found it difficult to collect enough of them to conduct common biochemical and chemical analyses because we did not find a place densely inhabited by hundreds of the worms whose wet weight was a few tens of milligrams (e.g. 16.5, 29.8, or 31.8 mg). Next, we conducted RNA-Seq analysis. In luminous animals with strong light emission, such as firefly or syllid polychaetes (Syllidae), luciferase expression is high especially at the luminous organ or in the whole body22,23. On the other hand, the light emission of Polycirrus sp. ISK was not so strong compared to that of fireflies, and the light-emitting area was limited to the tentacles. In addition, the genetic information related to the tentacles responsible for various ecological functions is still limited. Thus, in the present study we decided to purify RNA from the tentacles and the rest of body separately (Fig. 1C) and performed RNA-Seq analysis followed by a computational analysis using the MASER pipeline24. By de novo assembly, 110,775 contigs were predicted; 26.1% of them showed more than twice the expression level in the tentacles than in the rest of body, whereas 20.8% showed more than twice the expression in the rest of body than in the tentacles. When we performed a blastX search to the NCBI nr database for the contigs longer than 300 bp, 35.6% showed significant homology with registered genes with e-values of less than 1e−10. The average length for these contigs was 1384 bp, and half of them were in the range of 463–1863 bp (Supplementary Fig. S2). In the assembled sequence, we found the cytochrome oxidase subunit I (COI) gene and tried to construct a phylogenetic tree. However, the obtained phylogenetic tree was unreliable due to the low bootstrap values as shown in Supplementary Fig. S3.To focus on the tissue-specific genes, we first picked up genes with high expression levels based on high fpkm values (over 1000) and then ranked these genes based on the tissue-specificity judged by the comparison of fpkm values in tentacles and the rest of body. In tentacle-specific genes, we found that some genes coding for lectin(-like) domains were ranked in the top eight as shown in Supplementary Table S1. Of the top eight genes in the rest of body-specific genes (Supplementary Table S2), seven exhibited no similarity to any genes, and the remaining gene exhibited significant similarity to a hypothetical protein of Capitella teleta, which is a Polychaetes species with whole-genome information available25. Recently, TPM is preferably used to normalize expression level, and the value is used for statistical differential expression analysis26, and we also calculated TPM for tissue-specific genes (Supplementary Table S3).As we were unable to conduct statistical differential expression analysis due to no biological/technical replication resulted from difficulties in the sample collection, we simply compared TPM value between the tentacle and the rest of body samples. The ratio of TPM (tentacle/rest of body) was calculated, and then top 100 genes (Fig. 3A), which were highly expressing in the tentacle, were selected. Similarly, top 100 genes highly expressing in the rest of body were selected using the ratio of TPM (rest of body/tentacle) (Fig. 3B). These gene lists were annotated by gene ontology (GO) terms and analyzed using WEGO program27. WEGO results showed different expression patterns for the tentacle and the rest of body. In the tentacle, GO terms including cell adhesion, biological adhesion, small molecular binding, positive regulation of biological process, regulation of response to stimulus, carbohydrate binding, and immune response were significantly higher (Fig. 3C, D). In the rest of body, GO terms including hydrolase activity, catalytic activity, localization, and establishment of localization were significantly higher. In the top 100 genes highly expressing in the tentacle, we found 21 genes annotated as a gene coding for fucolectin by blast search (Supplementary Table S4). Fucolectin is a fucose-binding lectin involved in the innate immunity of diverse invertebrate species28. However, its function in invertebrates remains unclear, and no information is available for Terebellidae, including sequence information. Fucolectin was first identified in eel with mRNA distribution mainly in liver and gill28. In sea cucumber, expression of the fucolectin gene is confirmed in respiratory trees, muscle, and tentacle29. We were not able to see whether this gene was expressed in the respiratory organ of Polycirrus sp. ISK because a characteristic of the genus Polycirrus is the absence of branchiae18. Nevertheless, the tentacle-specific expression of fucolectin was consistent with the observation in sea cucumber, and the high expression of such proteins involved in innate immunity seemed reasonable because tentacles stretching out of their bodies can be damaged by attack of predators and thus are threatened by infectious bacteria and other pathogens11, as is the respiratory organ. In addition, localization of antimicrobial compounds in Terebellidae worms is suggested to be of antiseptic importance in damage by predation14. This study would provide indispensable information about the ecological meaning of Polycirrus sp. ISK’s life in future genetic studies.Figure 3WEGO analysis of highly expressing genes in the tentacle and the rest of body. (A) Box plot graph for the distribution of TPM value for top 100 genes highly expressing in the tentacle. Corresponding genes in each part are colored in the same gradation color according to the TPM value (red to blue form higher to lower value). (B) Box plot graph for the distribution of TPM value for top 100 genes highly expressing in the rest of body. Each gene is colored as in (A). (C) WEGO analysis of top 100 genes highly expressing in the tentacle (orange bar) and the rest of the body (blue bar). (D) P-values from Chi-square tests obtained by WEGO analysis. CC cellular component, MF molecular function, BP biological process.Full size imageTranscripts coding for luciferase-like genes in the wormTo find genes similar to the known luciferase, which is an enzyme oxidizing a specific compound called luciferin to emit light, from related species in polychaetes, we performed a blastX analysis against the Odontosyllis luciferase sequence using our RNA-Seq data. We found a gene coding for a protein that exhibited similarity to Odontosyllis luciferase, but the e-value was more than 1e−10 (Supplementary Fig. S4). In addition, the top hit for this gene analyzed by blastX was annotated to code an uncharacterized protein of Saccoglossus kowalevskii (Hemichordata), and its specific function was not predicted. Other hits were for genes from Chordata, Mollusca, and other phyla but there was no hit from Annelida. This result would suggest that the light-emission system of Polycirrus sp. ISK differs from that of the genus Odontosyllis, although further experiments using high purity Odontosyllis luciferase and the substrate will be necessary to confirm this. In further blastX analyses of representative luciferases, photoproteins, and a putative luciferase [luciferases from the ostracod Cypridina noctiluca (Accession number: BAD08210.1), the copepod Gaussia princeps (AAG54095.1), the deep-sea shrimp Oplophorus gracilirostris (BAB13775.1 and BAB13776.1), the firefly Photinus pyralis (AAA29795.1), the sea pansy Renilla reniformis (AAA29804.1); photoproteins from the hydrozoan jellyfish Aequorea victoria (AAA27720.1), the hydroid Clytia gregaria (CAA49754.1), the hydroid Obelia geniculate (AAL86372.1); a putative luciferase from the tunicate Pyrosoma atlanticum30 sequences using our RNA-Seq data], we found some tissue-nonspecific genes whose sequences exhibited similarity to firefly luciferase (FLuc) or Renilla luciferase-like protein (RLuc-like) sequences with an e-values of less than 1e−10 and percent identity of more than 50%. FLuc is a member of the acyl-adenylate-forming superfamily of enzymes responsible for firefly luciferin-dependent bioluminescence, which is found in terrestrial luminous beetles emitting light ranging from green to red31. Previously, a putative acyl-CoA synthetase protein was found in the luminous organ of firefly squid emitting blue light32, but there is no clear biochemical evidence that such protein is responsible for firefly squid’s bioluminescence. On the other hand, RLuc is responsible for coelenterazine-dependent bioluminescence, which is found in marine luminous organisms belonging to various taxa. An RLuc-like protein is found to be localized in luminous organs of the brittle star Amphiura filiformis, as revealed by taking advantage of the cross reactivity of anti-RLuc antibody to A. filiformis RLuc-like protein33. A recent study reported that recombinant RLuc-like protein found in P. atlanticum exhibited luciferase activity to coelenterazine30. However, an RLuc-like protein from sea urchin Strongylocentroutus purpuratus is confirmed to exhibit dehalogenase activity to various substrates but no luciferase activity to coelenterazine34. Therefore, it is suspected that Polycirrus sp. ISK possesses a luminescence system using an RLuc-like enzyme.Coelenterazine content in the wormTo investigate whether Polycirrus sp. ISK possesses not only a Renilla luciferase homologous gene but also coelenterazine, we analyzed an ethanolic extract of Polycirrus sp. ISK by UPLC with a UV–visible detector (Fig. 4). The obtained UPLC chromatogram did not show a peak corresponding to that of authentic coelenterazine. When further checking the chromatogram, we found the peak at a retention time similar to those of authentic coelenteramide and coelenteramine, which can be formed from coelenterazine. However, the absorption spectrum obtained by UPLC analysis and the mass spectrum obtained by MS/MS analysis were not identical to those of authentic coelenteramide or coelenteramine (Fig. 4 and Supplementary Figs. S5 and S6). In addition, when the worm extract was mixed with a recombinant RLuc, we did not detect luminescence using a luminometer. These results suggested that the luminescence system in the worm was independent of coelenterazine, although a RLuc homologous gene was found. Similarly, the existence of an RLuc homologous gene was reported in P. atlanticum, which has been suggested to use a coelenterazine-independent luminescence system relying on bacterial bioluminescent symbionts30,35. We also mixed the worm extract with a recombinant cypridinid luciferase, but we did not detect luminescence using a luminometer. This result was consistent with Harvey’s observation for P. caliendrum16. To examine whether the luminescence system is based on luciferin–luciferase reaction, which is found in various luminous animals including some syllid Odontosyllis spp.23,36,37,38,39, we prepared two different extracts of the worm using 100 mM HEPES–NaOH buffer (pH 7.4) and methanol, and subsequently subjected a mixture of the two to luminescent measurement. As a result, no light emission was detected from the mixture of the buffer and methanolic extracts of the worm. This result was also consistent with Harvey’s observation for P. caliendrum16. However, there is still a possibility that the light emission is based on luciferin–luciferase reaction, because luciferin–luciferase reaction found in fireflies or luminous mushrooms requires a cofactor such as ATP or NADPH, and we did not test all possible conditions due to the limitation of the number of collected specimens. In addition, extraction of luciferin and luciferase in the active form is sometimes difficult, as shown in previous studies37. Further studies using hundreds or more of the specimens must be performed to elucidate the mechanism underlying the violet-light emission.Figure 4Comparison of the ethanolic extract of Polycirrus sp. ISK with CTZ, CTMD, and CTM. (A) UPLC analysis of (a) the extract, (b) authentic CTZ, (c) authentic CTMD, and (d) authentic CTM using a multiwavelength detector. The black solid line indicates detection at 333 nm, and the blue solid line indicates detection at 435 nm. The compound between the red vertical dashed lines was collected for MS/MS analysis. (B) Absorption spectra of the compound from the extract, CTZ, CTMD, and CTM obtained at retention times of (a) 9.65, (b) 10.89, (c) 9.47, and (d) 9.27 shown in (A). CTZ coelenterazine, CTMD coelenteramide, CTM coelenteramine. These chemical structures are shown in Supplementary Fig. S5.Full size image More

1.Hansen, L. J. & Hoffman, J. R. Climate Savvy: Adapting Conservation and Resource Management to a Changing World, Austral Ecology (Island Press, 2011).Book 

Google Scholar 
2.Huang, Q., Fleming, C., Robb, B., Lothspeich, A. & Songer, M. How different are species distribution model predictions? Application of a new measure of dissimilarity and level of significance to giant panda Ailuropoda melanoleuca. Ecol. Inform. 46, 114–124 (2018).Article 

Google Scholar 
3.Giorgi, F. & Mearns, L. O. Probability of regional climate change based on the Reliability Ensemble Averaging (REA) method. Geophys. Res. Lett. 30, 1–4 (2003).Article 

Google Scholar 
4.Palmer, T. N., Doblas-Reyes, F. J., Hagedorn, R. & Weisheimer, A. Probabilistic prediction of climate using multi-model ensembles: From basics to applications. Philos. Trans. R. Soc. London B Biol. Sci. 360, 1991–1998 (2005).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
5.Dormann, C. F. et al. Prediction uncertainty of environmental change effects on temperate European biodiversity. Ecol. Lett. 11, 235–244 (2008).PubMed 
Article 

Google Scholar 
6.Spence, M. A. et al. A general framework for combining ecosystem models. Fish Fish. 19, 1031–1042 (2018).Article 

Google Scholar 
7.Janssen, P. H. M., Petersen, A. C., van der Sluijs, J. P., Risbey, J. S. & Ravetz, J. R. A guidance for assessing and communicating uncertainties. Water Sci. Technol. 52, 125–131 (2005).CAS 
PubMed 
Article 

Google Scholar 
8.Frewer, L. The public and effective risk communication. Toxicol. Lett. 149, 391–397 (2004).CAS 
PubMed 
Article 

Google Scholar 
9.Leiserowitz, A. A., Malbach, E. W., Roser-Renough, C., Smith, N. & Dawson, E. Climategate, public opinion, and the loss of trust. Am. Behav. Sci. 57, 818–837 (2012).Article 

Google Scholar 
10.Hyder, K. et al. Making modelling count: Increasing the contribution of shelf-seas community and ecosystem models to policy development and management. Mar. Policy 61, 291–302 (2015).Article 

Google Scholar 
11.Cartwright, S. J. et al. Communicating complex ecological models to non-scientist end users. Ecol. Model. 338, 51–59 (2016).Article 

Google Scholar 
12.Doyle, E. E. H., Johnston, D. M., Smith, R. & Paton, D. Communicating model uncertainty for natural hazards: A qualitative systematic thematic review. Int. J. Disaster Risk Reduct. 33, 449–476 (2019).Article 

Google Scholar 
13.Kamal, A. et al. Recent advances and challenges in uncertainty visualization: A survey. J. Vis. https://doi.org/10.1007/s12650-021-00755-1 (2021).Article 

Google Scholar 
14.Spiegelhalter, D. J. & Riesch, H. Don’t know, can’t know: Embracing deeper uncertainties when analysing risks. Philos. Trans. R. Soc. A 369, 4730–4750 (2011).ADS 
MathSciNet 
Article 
MATH 

Google Scholar 
15.Brodlie, K., Osorio, R. A. & Lopes, A. A review of uncertainty in data visualization. In Expanding the Frontiers of Visual Analytics and Visualization (eds Dill, J. et al.) 81–109 (Springer, 2012).Chapter 

Google Scholar 
16.MacEachren, A. M. et al. Visualizing geospatial information uncertainty: what we know and what we need to know. Cartogr. Geogr. Inf. Sci. 32, 139–160 (2005).Article 

Google Scholar 
17.Ibrekk, H. & Morgan, M. G. Graphical communication of uncertain quantities to nontechnical people. Risk Anal. 7, 519–529 (1987).Article 

Google Scholar 
18.Hawkins, E. The cascade of uncertainty in climate projections. https://www.climate-lab-book.ac.uk/2014/cascade-of-uncertainty/. 6 Feb 2014.19.Wilby, R. L. & Dessai, S. Robust adaptation to climate change. Weather 65, 180–185 (2010).ADS 
Article 

Google Scholar 
20.Daron, J., Lorenz, S., Taylor, A. & Dessai, S. Communicating future climate projections of precipitation change. Clim. Change. 166, 23 (2021).ADS 
Article 

Google Scholar 
21.Kinkeldey, C., MacEachren, A. M. & Schiewe, J. How to assess visual communication of uncertainty? A systematic review of geospatial uncertainty visualisation user studies. Cartogr. J. 51, 372–386 (2014).Article 

Google Scholar 
22.van Vuuren, D. P. et al. The representative concentration pathways: An overview. Clim. Change 109, 5–31 (2011).ADS 
Article 

Google Scholar 
23.Breslow, N. E. & Clayton, D. G. Approximate inference in generalized linear mixed models. J. Am. Stat. Assoc. 88, 9–25 (1993).MATH 

Google Scholar 
24.Bates, D., Maechler, M., Bolker, B. & Walker, S. Fitting linear mixed-effects models using lme4. J. Stat. Softw. 67, 1–48 (2015).Article 

Google Scholar 
25.Christensen, R.H.B. ordinal—Regression models for ordinal data. R package version 2018.4-19. https://CRAN.R-project.org/package=ordinal Deposited 19 Apr 2018.26.Peterson, B. & Harrell, F. E. Partial proportional odds models for ordinal response variables. Appl. Stat. 39, 205–217 (1990).Article 
MATH 

Google Scholar 
27.Momeni, A., Pincus, M. & Libien, J. Introduction to Statistical Methods in Pathology (Springer, 2018).Book 

Google Scholar 
28.Turner, H. & Firth, D. Bradley-Terry models in R: The BradleyTerry2 package. J. Stat. Softw. 48, 1–21 (2012).Article 

Google Scholar 
29.Akaike, H. Information theory and an extension of maximum likelihood principle. In Proceedings of the 2nd International Symposium on Information Theory. Tsahkadsor, Armenia, USSR, 2–8 September (1973).30.Rubinstein, R. The cross-entropy method for combinatorial and continuous optimization. Methodol. Comput. Appl. Probab. 1, 127–190 (1999).MathSciNet 
Article 

Google Scholar 
31.Kendall, M. G. A new measure of rank correlation. Biometrika 30, 81–93 (1938).Article 
MATH 

Google Scholar 
32.Venables, W. N. & Ripley, B. D. Modern Applied Statistics with S (Springer, 2002).Book 
MATH 

Google Scholar 
33.Firth, D., qvcalc: Quasi variances for factor effects in statistical statistical models. R package version 0.9-1. https://CRAN.R-project.org/package=qvcalc. Deposited 19 Sept 2017.34.Pihur, V., Datta, S. & Datta, S. RankAggreg: Weighted rank aggregation (2018). R package version 0.6-4. https://CRAN.R-project.org/package=RankAggreg. Deposited 18 Mar 2018.35.Belia, S., Fidler, F., Williams, J. & Cumming, G. Researchers misunderstand confidence intervals and standard error bars. Psychol. Methods 10, 389 (2005).PubMed 
Article 

Google Scholar 
36.Cumming, G., Fidler, F. & Vaux, D. L. Error bars in experimental biology. J. Cell Biol. 177, 7–11 (2007).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
37.Hullman, J., Resnick, P. & Adar, E. Hypothetical outcome plots outperform error bars and violin plots for inferences about reliability of variable ordering. PLoS ONE 10, e0142444 (2015).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
38.Correll, M. & Gleicher, M. Error bars considered harmful: Exploring alternate encodings for mean and error. IEEE Trans. Vis. Comput. Graph. 20, 2142–2151 (2014).PubMed 
PubMed Central 
Article 

Google Scholar 
39.Lipkus, I. M. & Hollands, J. G. The visual communication of risk. JNCI Monogr. 25, 149–163 (1999).Article 

Google Scholar 
40.Spiegelhalter, D., Pearson, M. & Short, I. Visualizing uncertainty about the future. Science 333, 1393–1400 (2011).ADS 
CAS 
PubMed 
Article 

Google Scholar 
41.Few, S. Solutions to the problem of over-plotting in graphs. Visual Business Intelligence Newsletter (2008).42.Daron, J. D., Lorenz, S., Wolski, P., Blamey, R. C. & Jack, C. Interpreting climate data visualisations to inform adaptation decisions. Clim. Risk Manag. 10, 17–26 (2015).Article 

Google Scholar 
43.Heer, J. & Agrawala, M. Multi-scale banking to 450. IEEE Trans. Vis. Comput. Graph. 12, 701–708 (2006).PubMed 
Article 

Google Scholar 
44.Lorenz, S., Dessai, S., Paavola, J. & Forster, P. M. The communication of physical science uncertainty in European National Adaptation Strategies. Clim. Change 132, 143–155 (2015).ADS 
CAS 
PubMed 
Article 

Google Scholar 
45.Quispel, A., Maes, A. & Schilperoord, J. Graph and chart aesthetics for experts and laymen in design: The role of familiarity and perceived ease of use. Inf. Vis. 15, 238–252 (2016).Article 

Google Scholar 
46.Saary, M. J. Radar plots: A useful way for presenting multivariate health care data. J. Clin. Epidemiol. 61, 311–317 (2008).PubMed 
Article 

Google Scholar 
47.Vaughan, N. E. & Gough, C. Expert assessment concludes negative emissions scenarios may not deliver. Environ. Res. Lett. 11, 95003 (2016).Article 
CAS 

Google Scholar 
48.Peltier, J. Excel Charting Dos and Don’ts (Peltier Technical Services, 2013).
Google Scholar 
49.Lohse, G. L. The role of working memory on graphical information processing. Behav. Inf. Technol. 16, 297–308 (1997).Article 

Google Scholar 
50.Grainger, S., Mao, F. & Buytaert, W. Environmental data visualisation for non-scientific contexts: Literature review and design framework. Environ. Model. Softw. 85, 299–318 (2016).Article 

Google Scholar 
51.Few, S. Heatmaps: to bin or not to bin? Visual Business Intelligence Newsletter (2017).52.Moreland, K. Diverging color maps for scientific visualization. In Advances in Visual Computing (eds Bebis, G. et al.) 92–103 (Springer, 2009).Chapter 

Google Scholar 
53.Alrehiely, M., Eslambolchilar, P. & Borgo, R. Evaluating different visualization designs for personal health data. In Proceedings of the 32nd International BCS Human Computer Interaction Conference, Belfast, UK, 4–6 July (2018).54.Saket, B., Endert, A. & Demiralp, C. Task-based effectiveness of basic visualizations. IEEE Trans. Visual Comput. Graph. 25, 2505–2512 (2019).Article 

Google Scholar 
55.Dahshan, M., Polys, N. F., Jayne, R. S. & Pollyea, R. M. Making sense of scientific simulation ensembles with semantic interaction. Comput. Graph. Forum. 39, 325–343 (2020).Article 

Google Scholar  More

Water useMany potato physiological features (photosynthesis intensity, leaf water potential) morphological and agronomic features as the Soil Plant Analysis Development (SPAD) and dry matter content can be used as indicators of potato water stress. In this result water consumption and the average daily amount of water used for irrigation differed over the growing season, but differences also occurred between varieties and the humidity level (Table 1). When irrigating the Julinka variety at all stages of the growing season, regardless of the established pF values, water consumption per pot was higher. The average dose of water supplied per pot was 9.7%, 30.7% and 26.6% greater than for the Denar variety, at humidity levels 1, 2 and 3, respectively. The highest water consumption was observed during the potato growth period from BBCH 40/400 to 69/609 and ranged from 0.39 l/pot /day (level 1) to 0.99 l/ pot/day (level 3).Table 1 Water consumption per pot within potato growing stages (in liters) and average consumption of water per pot (in brackets).Full size tableThe highest water consumption in both potato varieties occurred in July (11–18 July). Analyzing the remaining two months of the irrigation period, it can be seen that in June the plants used less water than in July. Seasonal irrigation doses in mid-early potato of studies of Rolbiecki et al. (2015)9 ranged from 40 to 170 mm, and the highest daily values of field water consumption (over 3 mm) occurred in July, similar to the results in this research.Depending on the irrigation system, water consumption efficiency in potato varies from 5.4 to 12 kg m−316,24. Drip irrigation is one of the most effective methods and ranged from 6.3 to 8.6 kg m−3 (Sharma 2007)25. Different values for average WUE index’ in potato cultivation were obtained by Ati et al. (2012)26, and indicated value ranged from 5.9 to 12.2 kg m−3. In present research, average WUE index’ for the Denar variety was from 0.00 l day−1 in the 1st period to 0.79 l day−1 in the 5th harvest period, while for the Julinka it was from 0.49 to 0.92 l day−1, respectively.In the research by Zin El Abedin et al. (2019)27 the amount of water used for irrigating potato amounted to 1505 mm and 1062 mm for FI (full irrigation) and PRD (partial root zone drying) variants, respectively. The use of 50% of water consumption in the PRD reduced water productivity (WP), as compared to water stress in the form of excess FI and deficit irrigation (DI). A large amount of water in conditions of water deficit causes losses due to evaporation and leads to degradation of the soil environment. In turn, in this research the highest water consumption in both varieties was found at level 3, 39.60 l for the Denar variety and 50.15 l for the Julinka variety.Pszczółkowski et al. (2009)28 showed that early potato varieties water requirements in the period from May 1 to August 31 amounted to 336.4 mm, with greatest requirements in July (108—119.6 mm). In our research, the amount of water used depended on the assumed humidity level and amounted from 19.60 × 103 to 39.60 × 103 cm3 for the Denar variety and between 21.50×103 to 50.15 × 103 cm3 for the Julinka (Table 1).Total potato and tuber massThe total weight of plants aboveground—(stems with leaves) and underground (tubers, stolons and roots) was greater in water humidity level 1 than in humidity levels 2 and 3 (Table 2). Administration of increased amounts of water in the later stages of potato growth resulted in inhibition of biomass growth, mainly for the Julinka variety. At the 5th harvest time, at humidity level 3, the total weight and the weight of tubers were 59.2% and 54.7% lower than those obtained at level 1, respectively. At the same time, the difference for Denar was 11.9% and 18.8%, respectively. Begum et al. (2015, 2018)16,22 and Reyes-Cabrera et al. (2016)5 showed that the production of total and commercial tuber yield was strongly dependent on the total biomass production and its structure.Table 2 Potato total biomass and tuber increase depending on water humidity level (g per plant).Full size tableA three-factor analysis of variance showed that the total weight as well as the weight of potato tubers differed significantly by the humidity level and the variety. A significant effect was found for humidity level on the total weight and tuber weight for the Denar variety and tuber weight for the Julinka variety (Table 3).Table 3 Variance analysis for total biomass and tuber of potato depending on factors (significance verified by the Fisher test).Full size tableAnalysis of variance showed a significant impact of the variety on potato plant weight, while it did not show significant interaction of weight and weight of tubers between measurement dates. No significant effect was obtained for interaction between the factors studied (Table 3).Wang et al. (2009)29, concluded that the use of irrigation significantly contributed to an increase total and commercial tubers of medium-early Folva variety yield and its quality. Ossowski et al. (2013)30, shown that irrigation had a significant effect on medium-early potato varieties: Barycz, Mors, Triada tuber yield. When using drip irrigation, yield increased by 26%. In turn, Mazurczyk et al. (2007)31 showed that drip irrigation increased the tuber yield from 29.4–37.5 to 45.1–54.4 t·ha−1.Over the period from the 1st to the 5th harvest date, the total plant biomass increased from 3.5-fold (Julinka—level 3) to 7.2-fold (Julinka—level 1). On the first harvest, Denar did not produce tubers at levels 1 and 3, and for level 2 its weight was the lowest (6 g from a pot). The increase in tuber weight to the last harvest date was the highest for level 2: 23.9- and 22.9-fold, in Denar and Julinka varieties, respectively. At level 3, the growth dynamics of tubers was the lowest: 11.7 times for the Julinka and 9.1 times for the Denar variety (measured from the second harvest date). The highest total biomass increases and tuber weight was found between the 3rd, 4th and 5th dates when humidity was at levels 1 and 2, and between the 3rd and 4th dates at level 3.Kumari et al. (2011, 2018)1,2 concluded that drip irrigation significantly contributed to an increase in potato tuber yield 18% greater than with other irrigation methods. Xu et al. (2010)32 achieved higher yields using the same irrigation system (40–48 t ha−1), and potato tuber weight was reduced under the slight water stress. Potato reacts to stress when soil water tension exceeds 20 kPa24. In a study by Amer et al. (2016)33 potato tuber yield also decreased with the application of excessive irrigation, resulting in greater stress, increased vegetative growth and potential leaching of nutrients from the root zone.Changes also occur in the quality of potato tubers, such as the shape, skin smoothness and chemical composition34.In the research carried out by Zin El-Abedin et al. (2019)27 differences were found in potato tuber yield depending on the irrigation variant. At FI, the highest tuber yields of 31.77–35.91 Mg ha−1 were obtained. Water deficiency reduced tuber yield, in DI variants, by 53.24–65.15% as compared to the FI. Similar results were obtained by Kumari et al. (2011)1. In the present research, the tuber weight of the Denar variety in the fifth term in level 1, increased by 26% compared to the irrigation at level 2 and was a 24% increase for the Julinka variety under similar conditions. At humidity level 3 there was a decrease in total biomass by 12% and 59% (for Denar and Julinka, respectively) in comparison obtained at level 1. In the research Liu et al. (2006)35 the aboveground biomass reached the highest values in excess water conditions.Potato varieties react differently to the humidity of the soil. Mahmood et al. (2016)36 response of potato varieties diversity to soil water deficit, also Hassanapanah (2010)17 showed the reaction of potato varieties to stress conditions. In our study, a higher total and tubers weight was found for the Julinka variety than for the Denar variety.Regardless of the humidity level and variety, the trends in the biomass yield structure were similar (Fig. 3). A downward trend from the 1st to 5th harvest period was shown for roots and stolons. This varied from 5 to 18% at the beginning of the study to 2–5% by the 5th period. It should be noted that under level 3, especially for the Denar variety, the percentage of roots and stolons was at a constant, low level. The percentage of stems with leaves decreased from 68–90% at the first harvest time to 40–55% at 5th. The dynamics of the decline in the share of stems and leaves was highest at humidity level 3. The tuber percentage was from 0 to 20% for the 1st period to 40–60% for the 5th period.Figure 3Potato biomass structure changes depending on humidity level and tuber harvest term (percentage).Full size imageThe Denar variety, regardless of the humidity level, was characterized by a greater share of stems and leaves. For the Julinka, the tuber percentage at the last harvest was at the same or higher than in the case of stems and leaves. At humidity levels 2 and 3, tubers accounted up to 60% of the harvested biomass.The growth of stem and stolon biomass was noticeable at all stages of potato development (Table 2); greater dynamics were found in the growth of tuber mass (Fig. 3). Under level 3, the growth of the biomass of stems with leaves and stolons was slower than in level 2 of water was used.Water use efficiencyAverage daily doses of water used for the Denar and Julinka varieties in potato harvesting periods are shown in Fig. 4. The volume of water was determined each time for the corresponding level of humidity (1, 2 and 3). Based on the data obtained, a proportional increase in water consumption was found for both potato varieties. The most intensive increase in water consumption was noted at humidity level 3. The W index corresponding to the average daily dose of water calculated for the Denar variety varied from 0.40 l day−1 in the 1st period (O1) to 0.79 l day−1 in the 5th harvest period (O5), whereas for the Julinka it was from 0.49 l day−1 (O1) to 0.92 l day−1 (O5). The W values for the level 3 changed for the Denar variety from 0.23 l day−1 in (O1) to 0.38 l day−1 (O5), while for the Julinka from 0.28 l day-1 (O1) to 0.28 l day-1, respectively (O5). The difference in the intensity of water consumption increase for humidity levels was expressed by varying the values of simple directional coefficients approximating empirical data. The highest values of these coefficients were obtained for the humidity level 1. The directional coefficient for the Denar was 0.0077 day−1, and for the Julinka variety 0.009 day−1. For humidity level 3, these values are 4 and 6 times lower: 0.002 day−1 (Denar) and 0.0014 day−1 (Julinka), respectively.Figure 4Average daily water consumption for potato varieties, at three soil humidity levels (1, 2, 3) and in each of five growing stages (O1), (O2), (O3), (O4), (O5).Full size imageThe average daily water consumption throughout the growing season calculated from potato planting is shown in Fig. 5. The average daily water use was the highest for both varieties at humidity level 3. Index W1 for the Denar was 0.53 l day−1, while for Julinka was higher—0.70 l day−1. The water consumption for the humidity level 1 was about 2 times lower: for the Denar—0.27 l day−1 and for Julinka—0.29 l day−1.Figure 5Average daily water consumption for potato varieties, at three soil humidity levels (1, 2, 3), cumulative calculation from potato planting.Full size imageAhmadi et al. (2017)37 used various irrigation schedule strategies for water demand measurements at evapotranspiration. Water demand has been fully or partially satisfied in static and dynamic modes. The research presents dynamics of vapor pressure deficit (VPD) throughout the growing season. The value of VPD in the first days after planting the potato was about 0.5 kPa while in 70 days maximum value was noted (2.5 kPa), and at the end of the growing season (after 150 days) about 1.5 kPa. Due to the shorter potato growing season in present research, no decrease in water demand was noticed up to about 70 days and, as in the results of the research presented by Ahmadi et al. (2017)37, a steady increase in water demand was noted. Similar results were obtained by King et al. (2020)38 and the largest water deficit was found in the middle of vegetation, after 70–80 days after planting35,39.Values for average daily increase in potato tuber weight (index W2) in individual vegetation periods are presented below (Fig. 6). No approximation of functional models to empirical data is possible; hence, the conclusions are based on a description. In the 1st period, i.e. until day 24 (O1), tuber weight gains were smaller than in the other periods. Depending on the humidity level, these amounted to 2.0 to 3.5 g day−1 for the Denar variety, and 2.7 to 3.9 g day−1 for the Julinka. The differences for Denar were 1.5 g day−1 and for Julinka 1.2 g day−1. In the 2nd irrigation period (O2), average daily increase in potato tuber weight was the highest, from 5.9 g day−1 for level 2 to 7.9 g day−1 for level 3. Average daily tuber weight gain was 13% higher for level 1 than for level 2.Figure 6Average daily potato varieties tuber increase, at three soil humidity levels (1, 2, 3), in each of five potato growing stages (O1), (O2), (O3), (O4), (O5).Full size imageThe average daily weight gain of tubers of potato varieties (W3), calculated incrementally from the beginning of the experiment (Fig. 7). For the entire growing season, this indicator for the Denar variety was the highest for the humidity level 1st (5.7 g day−1), at the level 3rd (5.1 g day−1) and the lowest at the level 2nd (4.3 g day−1). The average daily weight gain of potato tubers of the Julinka was definitely highest for the first humidity level (8.1 g day−1).Figure 7Average daily potato varieties tuber increase, at three soil humidity levels (1, 2, 3), cumulative calculation from potato planting.Full size imageThe ratio of the average daily water consumption to the average weight gain of potato tuber (W4) for individual periods is given in Fig. 8. For humidity level 1 for Denar and Julinka varieties, the values decreased with the growing period of vegetation. In the period (O1), 0.079 l of water was used for the Denar variety and 0.075 l for the Julinka for an increase in potato tuber weight of 1 g. In the next stages of the growing season, this index ranged from 0.35 to 0.45 l g−1 for the Denar variety, for the Julinka it was definitely smaller and range from 0.25 to 0.34 l g−1. At humidity level 1, Julinka used less water than Denar to produce the same weight of tubers. At humidity level 2, the volume of water used at the beginning of growth was also the largest for the Denar variety (0.159 l g−1). This amount was two times higher than the volume at level 1. In subsequent periods, the indicator changed and ranged from 0.059 to 0.105 l g−1. For the Julinka variety, water consumption varied in individual periods from 0.085 to 0.113 l g−1 and showed no trend. At humidity level 3, Denar used the greatest amount of water, as compared to levels 1 and 2, and showing no trend. The Julinka variety used even more water at the same humidity level. This amount ranged from 0.164 to 0.298 l g−1 and, unlike in previous cases, it showed an upward trend with plant development.Figure 8Ratio of average daily water consumption to average daily tuber mass increase dependent on three soil humidity levels (1), (2), (3), in each of five potato growing stages (O1), (O2), (O3), (O4), (O5).Full size imageJovanovic et al. (2010)40 divided the potato growing season into five stages related to growth phases. There were no increases in the weight of leaves and stems, while the tuber weight, regardless of the irrigation method (PRD and FI), increased steadily. The weight of tubers in the last harvest, as compared to the first, increased five-fold. A similar relationship was obtained in the work of Shahnazari et al. (2007)41. This research also took account of different levels of humidity using the strategies of PRD and FI, also considering soil retention characteristics (pF curve). The research showed a clear steady increase in potato tuber weight in each harvest.The ratio of the average daily water consumption to the average weight gain of potato tuber varieties calculated cumulatively from the planting (Fig. 9). The W5 value (0.114 l g−1) for the Denar variety at the end of the growing season was the highest for the 3rd humidity level and was about two times higher than at level 1. Water consumption efficiency for the Denar variety was the highest at humidity level 1. The sequence of W5 values is similar for the Julinka, with the difference that for the 3rd level it was 0.205 l g−1; i.e. six times higher than the indicator for level 1. Water consumption efficiency for the Julinka variety was definitely highest at humidity level 1.Figure 9Ratio of average daily water consumption to average daily tuber mass increase dependent on three soil humidity levels (1), (2), (3), cumulative calculation from potato planting.Full size imageBadr et al. (2010)42 analyzed the tuber yield, using two irrigation systems: surface and subsurface drip line. The total volume of water applied during the growing season was the differentiating factor. Results showed that as the volume of water applied during the growing season increased, the yield increased. When the subsurface line was used, applying 75 mm of water during the growing season, the total yield was approx. 27.5 t ha−1, and 32.5 t ha−1 for 325 mm. The effect of water amount on increase in yield was greater for the surface drip line. After applying 75 mm, the yield was 17.5 t ha−1, and 40 t ha−1 (for 325 mm). Similar results were obtained in the work of Linker et al. (2016)43. Regardless of the frequency, amount and total size of irrigation treatments, a proportional increase in the size of crops was observed with increasing doses of water.Shahnazari et al. (2007)41 planned several harvest dates (H0–H4) throughout the entire growing season, analyzing the irrigation efficiency indicator (average WUE index’). Regardless of the irrigation technique, and taking into account, above all, the amount of water administered, the value of the average WUE index’ indicator was the highest in the period H2–H3, similar results were found in our own research. More

1.Parmesan, C. & Yohe, G. Nature 421, 37–42 (2003).CAS 
Article 

Google Scholar 
2.Charmantier, A. & Gienapp, P. Evol. Appl. 7, 15–28 (2014).Article 

Google Scholar 
3.Merilä, J. & Hendry, A. P. Evol. Appl. 7, 1–14 (2014).Article 

Google Scholar 
4.Gienapp, P., Teplitsky, C., Alho, J. S., Mills, J. A. & Merilä, J. Mol. Ecol. 17, 167–178 (2008).CAS 
Article 

Google Scholar 
5.Bonamour, S., Chevin, L.-M., Charmantier, A. & Teplitsky, C. Philos. Trans. R. Soc. Lond. B Biol. Sci. 374, 20180178 (2019).Article 

Google Scholar 
6.Visser, M. E., Caro, S. P., van Oers, K., Schaper, S. V. & Helm, B. Philos. Trans. R. Soc. Lond. B Biol. Sci. 365, 3113–3127 (2010).CAS 
Article 

Google Scholar 
7.Caro, S. P., Schaper, S. V., Hut, R. A., Ball, G. F. & Visser, M. E. PLoS Biol. 11, e1001517 (2013).CAS 
Article 

Google Scholar 
8.Visser, M. E. & Gienapp, P. Nat. Ecol. Evol. 3, 879–885 (2019).Article 

Google Scholar 
9.Cole, E. F., Regan, C. E. & Sheldon, B. C. Nat. Clim. Change https://doi.org/10.1038/s41558-021-01140-4 (2021).Article 

Google Scholar 
10.Radchuk, V. et al. Nat. Commun. 10, 3109 (2019).Article 

Google Scholar 
11.Charmantier, A. et al. Science 320, 800–803 (2008).CAS 
Article 

Google Scholar 
12.Dewitt, T. J., Sih, A. & Wilson, D. S. Trends Ecol. Evol. 13, 77–81 (1998).CAS 
Article 

Google Scholar 
13.Díaz, S. et al. (eds) Summary for Policymakers of the Global Assessment Report of the Intergovernmental Science-Policy Platform on Biodiversity and Ecosystem Services (IPBES Secretariat, 2019). More

All experimental procedures were permitted by the Institutional Animal Care and Usage Committees of the Grassland Science College of Gansu Agricultural University (GSC-IACUC-2015-0011). Our experiments were conducted according to their guidelines, which are in accordance with the Guide for the Care and Use of Laboratory Animals (the Constitution of Experimental Animal Ethics Committee of Gansu Agricultural University). All experiments were performed in accordance with ARRIVE guidelines.Animals and laboratory conditionsAdults of both sexes (three males and three females) were captured in April 2015. Specifically, the animals were captured in Mayin Tan (37°12′N, 102°46′E; Tianzhu Tibetan Autonomous County, China) using live traps25 set at fresh surface mounds. The individuals were then transported to the laboratory and housed in an acrylic box with a pipeline covered with soil. The box and pipeline were covered by black cloth to simulate the dark environment of plateau zokors. The temperature in the room was maintained between 20 and 25 °C. Food was supplied daily and consisted of potatoes, lettuce and carrot. After three days of acclimatization to the laboratory the animals were used in the different experiments. Our laboratory is located 2 km away from the field site. At the end of the experiments all animals were returned to the capture site in good health.Laboratory testing arenaThe experimental setup for the laboratory experiments was as follows (Fig. 1): A transparent Perspex tube (8 cm × 8 cm × 80 cm) was joined to the side of the dark acrylic box (40 cm × 40 cm × 40 cm). A rubber stopper was inserted into one end of the tube to avoid effects from the external environment. Treatment apparatus was placed into the rubber stopper (see “Laboratory treatment apparatus” section, below), and, to avoid the apparatus being damaged by the animals, wire mesh (8 cm × 8 cm × 0.5 cm) was placed about 15 cm from one end of the tube. A mercury thermometer was inserted into the tube in the middle to monitor the tube’s temperature. Between experiments with different animals, the box and tube were wiped with 95% alcohol and then with distilled water.Figure 1Schematic drawing of the setup used to test burrow-sealing behavior in plateau zokors in the laboratory. (1) acrylic box covered with soil 30 cm in depth; (2) experimental animal; (3) mercury thermometer; (4) transparent Perspex tube; (5) the pipe’s support clip; (6) wire mesh (8 cm × 8 cm × 0.5 cm); (7) rubber stopper (8 cm × 8 cm × 5 cm).Full size imageLaboratory treatment apparatusA rubber stopper with seven holes was used for plugging one end of the tube (Fig. 2). The oxygen concentration, light, temperature, sound and gas flow were considered in this design.Figure 2Schematic diagram of the rubber stopper used to simulate the entrance plug of the burrow. (1) power supply; (2) light bulb switch; (3) electric wire switch; (4) oxygen cylinder; (5) in situ three-parameter soil gas analyzer; (6) voice recorder; (7) negative pressure drainage device; (8) rubber plug; (9) LED bulb; (10) the iron rod; (11) heating cord; (12) AVOXIVY speaker with 5 cm diameter.Full size imageOxygen treatmentTo avoid the oxygen that was delivered into the tube causing the gas to flow too strongly, become drier, and create a sound, a steel oxygen cylinder and thin hose (0.3 cm in diameter) were selected, and one end of the hose was connected directly to the oxygen cylinder with a humidifier bottle, while the opposite end was inserted into the rubber stopper (Fig. 2). Before beginning the experiment, we allowed the oxygen cylinder to sit for two hours at laboratory temperature to remove any temperature effects. A three-parameter soil gas analyzer (13.05.03Pro, Shanghai SAFE Biotech Co., Ltd, China) was used to monitor the oxygen concentration in the tube (Fig. 2).Light treatmentThe average light intensity—that is, 360 Lux from 8:00 am to 8:00 pm—was measured in the field. One end of a wire was connected to an LED light (1 Watt), and the other end to the power supply (Fig. 2).Temperature treatmentThe temperature in the burrow entrance in the field was about 3 °C warmer than that at a tunnel depth of 10 cm. As such, one end of a wire was connected to a heater strip and the other end to the electrical power supply (Fig. 2). A thermometer was inserted into the tube to monitor the temperature inside the tube (Fig.1). During the experiment period in the laboratory, we switched on or off to make sure the relatively constant temperature inside the tube. The temperature range inside the tube was 3.2 ± 0.27 °C .Sound treatmentWhen a burrow is opened, wind whistle can be produced around the burrow entrance. Accordingly, a voice recorder (PCM-D50, frequency response 50 Hz–40 kHz, Sony, Japan) was placed at the burrow entrance in the field to record the burrow-entrance sound, the duration of which was 30 min. In the laboratory, the two ends of a wire were connected to an AVOXIVY loudspeaker (diameter: 5 cm; impedance: 4 Ω; 50 Hz–20 kHz) and a voice recorder, respectively (Fig. 2). The recorded sound was played back with a 60 dB sound pressure level, as measured at the burrow entrance in the field (XL2 sound level meter, Nti Audion, Switzerland). The sound was repeatedly played within one hour.Gas flow treatmentTo avoid ambient atmosphere entering the tube, a negative pressure drainage ball with plastic tube (12 cm long, 2 cm in diameter) connected the tube through a rubber stopper (Fig. 2). The tunnel gas was inhaled by the ball, then we pinched the ball to blow the gas into the tunnel as gas flow treatment.Field treatment apparatusFor the field experiment, the apparatus consisted of a tube (40 cm long, 8 cm in diameter) and an alarm device. The alarm device was made up of a loudspeaker, two slide rails (15 cm long), two metal plates (approximately 7 cm in length and 3 cm in width), and three coiled metal springs (5 cm long, 2 cm in diameter). The three springs were joined to one of the metal plates, while the other metal plate was fixed on the slide rails. The two metal plates were touched by the plateau zokor when it was plugging, which triggered the alarm device, thus enabling us to know whether or not burrow-sealing behavior was occurring (Fig. 3). The aluminum tube with an oxygen device was embedded into the burrow. The soil covering the tube served as an excellent insulator, buffering the tube from the aboveground temperature (Fig. 4A). A steel oxygen cylinder and thin hose (0.3 cm in diameter) were applied by connecting one end of the hose directly to the oxygen cylinder with a humidifier bottle, and then the opposite end of the hose was inserted into the tube (Fig. 4A). A three-parameter soil gas analyzer (13.05.03Pro, Shanghai SAFE Biotech Co., Ltd, China) was used to monitor the oxygen concentration in the tube (Fig. 4A). Allowing sunlight to enter the burrow, a glass bottle, open at one end but closed at the other, was embedded into the burrow. We also used soil to cover the bottle, and there was a 5 cm gap at the surface (Fig. 4B). The aluminum tube with high thermal conductivity was embedded into the burrow. Again, we used soil to cover the bottle and retained a 20 cm gap (Fig. 4B).Figure 3Schematic drawing of the apparatus used to test the burrow-sealing behavior of plateau zokor in the field. (1) tube; (2) loudspeaker; (3) slide rail; (4) metal plate; (5) coiled metal springs.Full size imageFigure 4(A) Schematic drawing of the apparatus used in the oxygen treatment placed in the tunnel of the plateau zokor. (B) Schematic drawing of the apparatus used for the temperature and light treatments placed in the tunnel of the plateau zokor. (1) tunnel of the plateau zokor; (2) oxygen cylinder; (3) three-parameter soil gas analyzer; (4) plateau zokor.Full size imageProcedureIn the laboratory experiment, we tested three males and three females for their responses to each treatment. To avoid generating stress and habituation to treatments, zokors were tested for 12 h each day and there was one hour interval between treatments, and five days interval between round of testing for the same individual (Table 1). We performed a control experiment in which a rod was inserted into the burrow but no further treatment was applied, which allowed us to evaluate whether it was the treatment that was causing the burrow-sealing behavior. Before beginning treatment experiment, each zokor was tested 24 times (12 h × 2 days) under the control experiment. We determined the rod movement as occurrence of burrow-sealing behavior.Table 1 Times of the experiments for each treatment in the laboratory simulation.Full size tableIn the field experiment, we tested three zokors (one male, two females), and six zokors were caught in the cold season and warm season (three males and three females, respectively). We then fastened radio collars (Ag357, Biotrack, Ltd., UK) to each captured individual to allow us to track the position in foraging tunnels of each zokor. Each zokor was used three times in the experiments under each treatment, and, after finishing each experiment, we changed the position of the foraging tunnel to ensure the test tunnel was not an abandoned tunnel. According to radio-tracking data, the straight-line distance between the test tunnel and the nest for each treatment was about 5 m. We conducted a control experiment that whether plateau zokor move to the test tunnel or not during the time between treatments. In the cold season, from 4 October 2015 to 2 November 2015, the burrow-sealing behavior of each zokor was tested under different treatments during their active time (12:00–18:00) and inactive time (09:00–11:00) for a total of 27 days (Table 2). The same was done in the warm season but for a total of 18 days from 15 May 2016 to 5 June 2016, in which the active time was 14:00–20:00 and the inactive time was 08:00–13:00 (Table 2).Table 2 Times of the experiments for each treatment in the active and inactive periods of plateau zokors during the warm and cold season.Full size tableData analysisThe occurrence of burrow-sealing was recorded as “1”, and non-sealing was recorded as “0”. The frequency of burrow-sealing was the number of times the burrow was sealed divided by the total number of experiments for each treatment26, and we considered the frequency for each individual as a replicate. The latency to reseal the burrow was the period from the start of the treatment to the sealing of the burrow, and we considered each instance of latency to reseal the burrow as a repeat. The latency to reseal the burrow for non-sealing under each treatment was unavailable data and was therefore removed. The presence of a normal distribution in the initial data was determined using the Kolmogorov–Smirnov test. All data followed a normal distribution. A comparison of males and females in their frequency of sealing the burrow and in their latency to reseal the burrow under each treatment was performed with an independent-samples T-test. Multiple comparisons were made for the frequency of burrow-sealing and the latency to reseal the burrow under different treatments by using the least significant difference method at the significance level of P = 0.05. In the field experiment, the number of replicates was fewer than three for frequency and the latency to reseal the burrow, we did not conduct multiple comparisons.Preliminary statistical analysis of the data was performed using Excel 2013 and SPSS 19.0. All the figures and tables were produced in GraphPad Prism 8.0 and Excel 2013. More