Ecology

Okogwu, O. I., Nwonumara, G. N. & Okoh, F. A. Evaluating heavy metals pollution and exposure risk through the consumption of four commercially important fish species and water from cross river ecosystem, Nigeria. Bull. Environ. Contam. Toxicol. 102, 867–872. https://doi.org/10.1007/s00128-019-02610-4 (2019).CAS 
Article 
PubMed 

Google Scholar 
Fuentes-Gandara, F., Pinedo-Hernández, J., Marrugo-Negrete, J. & Díez, S. Human health impacts of exposure to metals through extreme consumption of fish from the Colombian Caribbean Sea. Environ. Geochem. Health 40, 229–242. https://doi.org/10.1007/s10653-016-9896-z (2018).CAS 
Article 
PubMed 

Google Scholar 
Liu, X. et al. Human health risk assessment of heavy metals in soil–vegetable system: A multi-medium analysis. Sci. Total Environ. 463, 530–540 (2013).PubMed 

Google Scholar 
Lelieveld, J., Evans, J. S., Fnais, M., Giannadaki, D. & Pozzer, A. The contribution of outdoor air pollution sources to premature mortality on a global scale. Nature 525, 367–371. https://doi.org/10.1038/nature15371 (2015).CAS 
Article 
PubMed 

Google Scholar 
Huang, R.-J. et al. High secondary aerosol contribution to particulate pollution during haze events in China. Nature 514, 218–222. https://doi.org/10.1038/nature13774 (2014).CAS 
Article 
PubMed 

Google Scholar 
Rajeshkumar, S. et al. Studies on seasonal pollution of heavy metals in water, sediment, fish and oyster from the Meiliang Bay of Taihu Lake in China. Chemosphere 191, 626–638. https://doi.org/10.1016/j.chemosphere.2017.10.078 (2018).CAS 
Article 
PubMed 

Google Scholar 
Gao, X. & Chen, C.-T.A. Heavy metal pollution status in surface sediments of the coastal Bohai Bay. Water Res. 46, 1901–1911. https://doi.org/10.1016/j.watres.2012.01.007 (2012).CAS 
Article 
PubMed 

Google Scholar 
Naser, H. A. Assessment and management of heavy metal pollution in the marine environment of the Arabian Gulf: A review. Mar. Pollut. Bull. 72, 6–13. https://doi.org/10.1016/j.marpolbul.2013.04.030 (2013).CAS 
Article 
PubMed 

Google Scholar 
Zhang, Y. et al. Heavy metals in aquatic organisms of different trophic levels and their potential human health risk in Bohai Bay, China. Environ. Sci. Pollut. Res. 23, 17801–17810 (2016).CAS 

Google Scholar 
Wei, M., Yanwen, Q., Zheng, B. & Zhang, L. Heavy metal pollution in Tianjin Bohai bay, China. J. Environ. Sci. 20, 814–819 (2008).
Google Scholar 
Zhao, B. et al. Spatiotemporal variation and potential risks of seven heavy metals in seawater, sediment, and seafood in Xiangshan Bay, China (2011–2016). Chemosphere 212, 1163–1171. https://doi.org/10.1016/j.chemosphere.2018.09.020 (2018).CAS 
Article 
PubMed 

Google Scholar 
Wang, Y. & Fang, X. Analysis of the impact of heavy metal on the Chinese aquaculture and the ecological hazard. GuangDong 836, 156.152 (2016).
Google Scholar 
Pini, J., Richir, J. & Watson, G. Metal bioavailability and bioaccumulation in the polychaete Nereis (Alitta) virens (Sars): The effects of site-specific sediment characteristics. Mar. Pollut. Bull. 95, 565–575 (2015).CAS 
PubMed 

Google Scholar 
Amoozadeh, E. et al. Marine organisms as heavy metal bioindicators in the Persian Gulf and the Gulf of Oman. Environ. Sci. Pollut. Res. 21, 2386–2395 (2014).CAS 

Google Scholar 
Gu, Y.-G., Huang, H.-H., Liu, Y., Gong, X.-Y. & Liao, X.-L. Non-metric multidimensional scaling and human risks of heavy metal concentrations in wild marine organisms from the Maowei Sea, the Beibu Gulf, South China Sea. Environ. Toxicol. Pharmacol. 59, 119–124. https://doi.org/10.1016/j.etap.2018.03.002 (2018).CAS 
Article 
PubMed 

Google Scholar 
Kennedy, A., Martinez, K., Chuang, C.-C., LaPoint, K. & McIntosh, M. Saturated fatty acid-mediated inflammation and insulin resistance in adipose tissue: Mechanisms of action and implications. J. Nutr. 139, 1–4. https://doi.org/10.3945/jn.108.098269 (2008).CAS 
Article 
PubMed 

Google Scholar 
Hao, Z. et al. Heavy metal distribution and bioaccumulation ability in marine organisms from coastal regions of Hainan and Zhoushan, China. Chemosphere 226, 340–350. https://doi.org/10.1016/j.chemosphere.2019.03.132 (2019).CAS 
Article 
PubMed 

Google Scholar 
Golden, C. D. et al. Nutrition: Fall in fish catch threatens human health. Nat. News 534, 317 (2016).
Google Scholar 
Bosch, A. C., O’Neill, B., Sigge, G. O., Kerwath, S. E. & Hoffman, L. C. Heavy metals in marine fish meat and consumer health: A review. J. Sci. Food Agric. 96, 32–48 (2016).CAS 
PubMed 

Google Scholar 
Burger, J., Gochfeld, M., Jeitner, C., Pittfield, T. & Donio, M. Heavy metals in fish from the Aleutians: Interspecific and locational differences. Environ. Res. 131, 119–130. https://doi.org/10.1016/j.envres.2014.02.016 (2014).CAS 
Article 
PubMed 

Google Scholar 
Anandkumar, A., Nagarajan, R., Prabakaran, K., Chua Han, B. & Rajaram, R. Human health risk assessment and bioaccumulation of trace metals in fish species collected from the Miri coast, Sarawak, Borneo. Mar. Pollut. Bull. 133, 655–663. https://doi.org/10.1016/j.marpolbul.2018.06.033 (2018).CAS 
Article 

Google Scholar 
Murtala, B. A., Abdul, W. O. & Akinyemi, A. A. Bioaccumulation of heavy metals in fish (Hydrocynus forskahlii, Hyperopisus bebe occidentalis and Clarias gariepinus) organs in downstream Ogun coastal water, Nigeria. J. Agric. Sci. 4, 51 (2012).
Google Scholar 
Ahmed, A. S. S., Rahman, M., Sultana, S., Babu, S. M. O. F. & Sarker, M. S. I. Bioaccumulation and heavy metal concentration in tissues of some commercial fishes from the Meghna River Estuary in Bangladesh and human health implications. Mar. Pollut. Bull. 145, 436–447. https://doi.org/10.1016/j.marpolbul.2019.06.035 (2019).CAS 
Article 
PubMed 

Google Scholar 
Sun, X. et al. Source identification, geochemical normalization and influence factors of heavy metals in Yangtze River Estuary sediment. Environ. Pollut. 241, 938–949. https://doi.org/10.1016/j.envpol.2018.05.050 (2018).CAS 
Article 
PubMed 

Google Scholar 
Dadar, M., Adel, M., NasrollahzadehSaravi, H. & Fakhri, Y. Trace element concentration and its risk assessment in common kilka (Clupeonella cultriventris caspia Bordin, 1904) from southern basin of Caspian Sea. Toxin Rev. 36, 222–227 (2017).CAS 

Google Scholar 
Chakraborty, P., Raghunadh Babu, P. V., Acharyya, T. & Bandyopadhyay, D. Stress and toxicity of biologically important transition metals (Co, Ni, Cu and Zn) on phytoplankton in a tropical freshwater system: An investigation with pigment analysis by HPLC. Chemosphere 80, 548–553. https://doi.org/10.1016/j.chemosphere.2010.04.039 (2010).CAS 
Article 
PubMed 

Google Scholar 
Handy, R. Seminar Series-Society for Experimental Biology 29–60 (Cambridge University Press, 1997).
Google Scholar 
Ahmed, M. K. et al. Human health risks from heavy metals in fish of Buriganga river, Bangladesh. Springerplus 5, 1–12 (2016).
Google Scholar 
WHO. Heavy metals-environmental aspects. Environment Health Criteria. No. 85. (1989).Xu, H. et al. Long-term study of heavy metal pollution in the northern Hangzhou Bay of China: Temporal and spatial distribution, contamination evaluation, and potential ecological risk. Environ. Sci. Pollut. Res. 28, 10718–10733 (2021).CAS 

Google Scholar 
El-Moselhy, K. M., Othman, A. I., AbdEl-Azem, H. & El-Metwally, M. E. A. Bioaccumulation of heavy metals in some tissues of fish in the Red Sea, Egypt. Egypti. J. Basic Appl. Sci. 1, 97–105. https://doi.org/10.1016/j.ejbas.2014.06.001 (2014).Article 

Google Scholar 
Jezierska, B. & Witeska, M. Soil and Water Pollution Monitoring, Protection and Remediation 107–114 (Springer, 2006).
Google Scholar 
Bawuro, A. A., Voegborlo, R. B. & Adimado, A. A. Bioaccumulation of heavy metals in some tissues of fish in Lake Geriyo, Adamawa State, Nigeria. J. Environ. Public Health 2018, 1854892. https://doi.org/10.1155/2018/1854892 (2018).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
Zhuang, P., McBride, M. B., Xia, H., Li, N. & Li, Z. Health risk from heavy metals via consumption of food crops in the vicinity of Dabaoshan mine, South China. Sci. Total Environ. 407, 1551–1561. https://doi.org/10.1016/j.scitotenv.2008.10.061 (2009).CAS 
Article 
PubMed 

Google Scholar 
Hosseini, M., Nabavi, S. M. B., Nabavi, S. N. & Pour, N. A. Heavy metals (Cd Co, Cu, Ni, Pb, Fe, and Hg) content in four fish commonly consumed in Iran: Risk assessment for the consumers. Environ. Monit. Assess. 187, 237. https://doi.org/10.1007/s10661-015-4464-z (2015).CAS 
Article 
PubMed 

Google Scholar 
Prabhakaran, K., Nagarajan, R., MerlinFranco, F. & AnandKumar, A. Biomonitoring of Malaysian aquatic environments: A review of status and prospects. Ecohydrol. Hydrobiol. 17, 134–147. https://doi.org/10.1016/j.ecohyd.2017.03.001 (2017).Article 

Google Scholar 
Meche, A. et al. Determination of heavy metals by inductively coupled plasma-optical emission spectrometry in fish from the Piracicaba River in Southern Brazil. Microchem. J. 94, 171–174 (2010).CAS 

Google Scholar 
Zhang, Y. et al. Temporal and spatial changes of microbial community in an industrial effluent receiving area in Hangzhou Bay. J. Environ. Sci. 44, 57–68. https://doi.org/10.1016/j.jes.2015.11.023 (2016).CAS 
Article 

Google Scholar 
Huang, L. et al. Quantifying the spatiotemporal dynamics of industrial land uses through mining free access social datasets in the Mega Hangzhou Bay Region, China. Sustainability 10, 3463 (2018).
Google Scholar 
Pang, H.-J. et al. Contamination, distribution, and sources of heavy metals in the sediments of Andong tidal flat, Hangzhou bay, China. Continental Shelf Res. 110, 72–84. https://doi.org/10.1016/j.csr.2015.10.002 (2015).Article 

Google Scholar 
National Bureau of Statstics. Zhejiang Statistical Yearbook-2017 (China Statistics Press, 2017).
Google Scholar 
Chen, W., Zheng, Y., Chen, Y. & Mathews, C. An assessment of fishery yields from the East China Sea ecosystem. Mar. Fish. Rev. 59, 1–7 (1997).
Google Scholar 
Zhejiang Provincial Development and Reform Commission. Zhejiang Zhoushan Islands New Area Development Plan (In Chinese). (2021).Che, Y., He, Q. & Lin, W.-Q. The distributions of particulate heavy metals and its indication to the transfer of sediments in the Changjiang Estuary and Hangzhou Bay, China. Mar. Pollut. Bull. 46, 123–131 (2003).CAS 
PubMed 

Google Scholar 
Li, R. et al. Environmental health and ecological risk assessment of soil heavy metal pollution in the coastal cities of Estuarine Bay—a case study of Hangzhou Bay, China. Toxics 8, 75 (2020).CAS 
PubMed Central 

Google Scholar 
Bergami, E., Manno, C., Cappello, S., Vannuccini, M. L. & Corsi, I. Nanoplastics affect moulting and faecal pellet sinking in Antarctic krill (Euphausia superba) juveniles. Environ. Int. 143, 105999. https://doi.org/10.1016/j.envint.2020.105999 (2020).CAS 
Article 
PubMed 

Google Scholar 
Fang, H., Huang, L., Wang, J., He, G. & Reible, D. Environmental assessment of heavy metal transport and transformation in the Hangzhou Bay, China. J. Hazard. Mater. 302, 447–457 (2016).CAS 
PubMed 

Google Scholar 
Zhu, G. et al. Evaluation of ecosystem health and potential human health hazards in the Hangzhou Bay and Qiantang Estuary region through multiple assessment approaches. Environ. Pollut. 264, 114791. https://doi.org/10.1016/j.envpol.2020.114791 (2020).CAS 
Article 
PubMed 

Google Scholar 
Li, F. et al. Distribution and risk assessment of trace metals in sediments from Yangtze River estuary and Hangzhou Bay, China. Environ. Sci. Pollut. Res. 25, 855–866. https://doi.org/10.1007/s11356-017-0425-0 (2018).CAS 
Article 

Google Scholar 
Liu, L., Huang, X., Cao, W. & Yang, Y. Pollution load characteristics of the Hangzhou Bay and its surrounding areas. Ocean Dev. Manage 5, 108–112 (2012).
Google Scholar 
He, Z., Li, F., Dominech, S., Wen, X. & Yang, S. Heavy metals of surface sediments in the Changjiang (Yangtze River) Estuary: Distribution, speciation and environmental risks. J. Geochem. Explor. 198, 18–28. https://doi.org/10.1016/j.gexplo.2018.12.015 (2019).CAS 
Article 

Google Scholar 
Jin, X., Zhao, X., Meng, T. & Cui, Y. The Fishery Resources and the Environment of the Bohai Sea and Yellow Sea (Science Press, 2005).
Google Scholar 
Huang, Z. The Species and Distribution of Marine Organisms of China (Ocean Press, Beijing, 1994) (In Chinese).
Google Scholar 
Schram, F. R. Checklist of Marine Biota of China Seas. J. Crustac. Biol. 30, 339–339. https://doi.org/10.1651/09-3228.1 (2010).Article 

Google Scholar 
AQSIQ, P. in GB 17378.6–2007 (General Administration of Quality Supervision, Inspection and Quarantine of People’s Republic of China, 2007).Zhang, L. et al. Distribution and bioaccumulation of heavy metals in marine organisms in east and west Guangdong coastal regions, South China. Mar. Pollut. Bull. 101, 930–937. https://doi.org/10.1016/j.marpolbul.2015.10.041 (2015).CAS 
Article 
PubMed 

Google Scholar 
Zhong, W. et al. Health risk assessment of heavy metals in freshwater fish in the central and eastern North China. Ecotoxicol. Environ. Saf. 157, 343–349. https://doi.org/10.1016/j.ecoenv.2018.03.048 (2018).CAS 
Article 
PubMed 

Google Scholar 
Wang, Q. et al. Bioaccumulation and biomagnification of emerging bisphenol analogues in aquatic organisms from Taihu Lake, China. Sci. Total Environ. 598, 814–820. https://doi.org/10.1016/j.scitotenv.2017.04.167 (2017).CAS 
Article 
PubMed 

Google Scholar 
Arnot, J. A. & Gobas, F. A. A review of bioconcentration factor (BCF) and bioaccumulation factor (BAF) assessments for organic chemicals in aquatic organisms. Environ. Rev. 14, 257–297 (2006).CAS 

Google Scholar 
Duan, X., Zhao, X., Wang, B., Chen, Y. & Cao, S. Exposure Factors Handbook of Chinese Population (Adults) (China Environmental Science Press, 2013).
Google Scholar 
Chauhan, G. & Chauhan, U. Human health risk assessment of heavy metals via dietary intake of vegetables grown in wastewater irrigated area of Rewa, India. Int. J. Sci. Res. Publ. 4, 1–9 (2014).
Google Scholar 
USEPA. (Philadelphia PA; Washington, DC, 2007).Wang, X., Sato, T., Xing, B. & Tao, S. Health risks of heavy metals to the general public in Tianjin, China via consumption of vegetables and fish. Sci. Total Environ. 350, 28–37. https://doi.org/10.1016/j.scitotenv.2004.09.044 (2005).CAS 
Article 
PubMed 

Google Scholar 
USEPA. (2015).FAO/WHO. Wastewater Use in Agriculture. 988 (World Health Organization).Ahmed, A. S. S. et al. Bioaccumulation of heavy metals in some commercially important fishes from a tropical river estuary suggests higher potential health risk in children than adults. PLoS One 14, e0219336. https://doi.org/10.1371/journal.pone.0219336 (2019).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
Saha, N., Mollah, M. Z. I., Alam, M. F. & Safiur Rahman, M. Seasonal investigation of heavy metals in marine fishes captured from the Bay of Bengal and the implications for human health risk assessment. Food Control 70, 110–118. https://doi.org/10.1016/j.foodcont.2016.05.040 (2016).CAS 
Article 

Google Scholar 
Yin, S., Feng, C., Li, Y., Yin, L. & Shen, Z. Heavy metal pollution in the surface water of the Yangtze Estuary: A 5-year follow-up study. Chemosphere 138, 718–725. https://doi.org/10.1016/j.chemosphere.2015.07.060 (2015).CAS 
Article 
PubMed 

Google Scholar 
USEPA. Risk-based concentration table. United States Environmental Protection Agency, Washington DC, Philadelphia (2000).Hu, B. et al. Assessment of heavy metal pollution and health risks in the soil-plant-human system in the Yangtze River Delta, China. Int. J. Environ. Res. Public Health 14, 1042 (2017).PubMed Central 

Google Scholar 
USEPA. in United States Environmental Protection Agency, Washington DC, Philadelphia (2010).Kwok, C. K. et al. Bioaccumulation of heavy metals in fish and Ardeid at Pearl River Estuary, China. Ecotoxicol. Environ. Saf. 106, 62–67. https://doi.org/10.1016/j.ecoenv.2014.04.016 (2014).CAS 
Article 
PubMed 

Google Scholar 
Yu, T., Zhang, Y., Hu, X. & Meng, W. Distribution and bioaccumulation of heavy metals in aquatic organisms of different trophic levels and potential health risk assessment from Taihu lake, China. Ecotoxicol. Environ. Saf. 81, 55–64. https://doi.org/10.1016/j.ecoenv.2012.04.014 (2012).CAS 
Article 

Google Scholar 
Qiu, Y.-W., Lin, D., Liu, J.-Q. & Zeng, E. Y. Bioaccumulation of trace metals in farmed fish from South China and potential risk assessment. Ecotoxicol. Environ. Saf. 74, 284–293. https://doi.org/10.1016/j.ecoenv.2010.10.008 (2011).CAS 
Article 
PubMed 

Google Scholar 
Arulkumar, A., Paramasivam, S. & Rajaram, R. Toxic heavy metals in commercially important food fishes collected from Palk Bay, Southeastern India. Mar. Pollut. Bull. 119, 454–459. https://doi.org/10.1016/j.marpolbul.2017.03.045 (2017).CAS 
Article 
PubMed 

Google Scholar 
Jonathan, M. P. et al. Metal concentrations in demersal fish species from Santa Maria Bay, Baja California Sur, Mexico (Pacific coast). Mar. Pollut. Bull. 99, 356–361. https://doi.org/10.1016/j.marpolbul.2015.07.032 (2015).CAS 
Article 
PubMed 

Google Scholar 
Liu, H., Yang, J. & Gan, J. Trace element accumulation in bivalve mussels Anodonta woodiana from Taihu Lake, China. Arch. Environ. Contam. Toxicol. 59, 593–601. https://doi.org/10.1007/s00244-010-9521-6 (2010).CAS 
Article 
PubMed 

Google Scholar 
Wang, W. X. et al. Copper and zinc contamination in oysters: Subcellular distribution and detoxification. Environ. Toxicol. Chem. 30, 1767–1774 (2011).CAS 
PubMed 

Google Scholar 
de FreitasRebelo, M., do Amaral, M. C. R. & Pfeiffer, W. C. High Zn and Cd accumulation in the oyster Crassostrea rhizophorae, and its relevance as a sentinel species. Mar. Pollut. Bull. 46, 1354–1358 (2003).
Google Scholar 
AQSIQ, P. in GB 18421–2001 (General administration of quality supervision, inspection and quarantine of People’s Republic of China, 2001).FAO/WHO. in Fifth Session [displayed 10 February 2014]. ftp://ftp.fao.org/codex/meetings/CCCF/cccf5/cf05_INF.pdf.Nauen, C. E. Compilation of legal limits for hazardous substances in fish and fishery products. FAO Fisheries Circular (FAO). no. 764. (1983).Rajeshkumar, S. & Li, X. Bioaccumulation of heavy metals in fish species from the Meiliang Bay, Taihu Lake, China. Toxicol. Rep. 5, 288–295. https://doi.org/10.1016/j.toxrep.2018.01.007 (2018).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
Baki, M. A. et al. Concentration of heavy metals in seafood (fishes, shrimp, lobster and crabs) and human health assessment in Saint Martin Island, Bangladesh. Ecotoxicol. Environ. Saf. 159, 153–163. https://doi.org/10.1016/j.ecoenv.2018.04.035 (2018).CAS 
Article 
PubMed 

Google Scholar 
Vu, C. T., Lin, C., Yeh, G. & Villanueva, M. C. Bioaccumulation and potential sources of heavy metal contamination in fish species in Taiwan: Assessment and possible human health implications. Environ. Sci. Pollut. Res. 24, 19422–19434. https://doi.org/10.1007/s11356-017-9590-4 (2017).CAS 
Article 

Google Scholar 
Sharma, B., Singh, S. & Siddiqi, N. J. Biomedical implications of heavy metals induced imbalances in redox systems. BioMed Res. Int. 20, 14 (2014).
Google Scholar 
Feng, W., Wang, Z., Xu, H., Chen, L. & Zheng, F. Trace metal concentrations in commercial fish, crabs, and bivalves from three lagoons in the South China Sea and implications for human health. Environ. Sci. Pollut. Res. https://doi.org/10.1007/s11356-019-06712-8 (2020).Article 

Google Scholar 
Ruiz-Fernández, A. C. et al. A comparative study on metal contamination in Estero de Urias lagoon, Gulf of California, using oysters, mussels and artificial mussels: Implications on pollution monitoring and public health risk. Environ. Pollut. 243, 197–205 (2018).PubMed 

Google Scholar 
Bergstad, O. A. In Encyclopedia of Ocean Sciences (Second Edition) (ed. Steele, J. H.) 458–466 (Academic Press, 2009).
Google Scholar 
Mauchline, J. & Gordon, J. Foraging strategies of deep-sea fish. Mar. Ecol. Prog. Ser. 27, 227–238 (1986).
Google Scholar 
Li, J., He, M., Han, W. & Gu, Y. Analysis and assessment on heavy metal sources in the coastal soils developed from alluvial deposits using multivariate statistical methods. J. Hazard. Mater. 164, 976–981. https://doi.org/10.1016/j.jhazmat.2008.08.112 (2009).CAS 
Article 
PubMed 

Google Scholar 
Yu, P. Applications of hierarchical cluster analysis (CLA) and principal component analysis (PCA) in feed structure and feed molecular chemistry research, using synchrotron-based Fourier transform infrared (FTIR) microspectroscopy. J. Agric. Food Chem. 53, 7115–7127 (2005).CAS 
PubMed 

Google Scholar 
Kara, D. Evaluation of trace metal concentrations in some herbs and herbal teas by principal component analysis. Food Chem. 114, 347–354 (2009).CAS 

Google Scholar 
Chai, X. et al. Distribution, sources and assessment of heavy metals in surface sediments of the Hangzhou Bay and its adjacent areas. Acta Sci. Circum. 35, 3906–3916 (2015).CAS 

Google Scholar 
Mackay, D. & Fraser, A. Bioaccumulation of persistent organic chemicals: Mechanisms and models. Environ. Pollut. 110, 375–391. https://doi.org/10.1016/S0269-7491(00)00162-7 (2000).CAS 
Article 
PubMed 

Google Scholar 
ATSDR, T. ATSDR (Agency for toxic substances and disease registry). Prepared by Clement International Corp., under contract 205, 88–0608 (2000).Traina, A. et al. Heavy metals concentrations in some commercially key species from Sicilian coasts (Mediterranean Sea): Potential human health risk estimation. Ecotoxicol. Environ. Saf. 168, 466–478. https://doi.org/10.1016/j.ecoenv.2018.10.056 (2019).CAS 
Article 
PubMed 

Google Scholar 
Ozmen, M., Ayas, Z., Güngördü, A., Ekmekci, G. F. & Yerli, S. Ecotoxicological assessment of water pollution in Sariyar Dam Lake, Turkey. Ecotoxicol. Environ. Saf. 70, 163–173. https://doi.org/10.1016/j.ecoenv.2007.05.011 (2008).CAS 
Article 
PubMed 

Google Scholar 
Jeffrey, B. & Alison, G. Guidance for assessing chemical contaminant data for use in fish advisories. v. 1. Fish sampling and analysis-v. 4. Risk communication. (1993).Regulations, U. S. E. P. A. O. o. W. Assessing Human Health Risks from Chemically Contaminated Fish and Shellfish: A Guidance Manual. (US Environmental Protection Agency, 1989).Liu, Q., Liao, Y. & Shou, L. Concentration and potential health risk of heavy metals in seafoods collected from Sanmen Bay and its adjacent areas, China. Mar. Pollut. Bull 131, 356–364. https://doi.org/10.1016/j.marpolbul.2018.04.041 (2018).CAS 
Article 
PubMed 

Google Scholar 
Abtahi, M. et al. Heavy metals (As, Cr, Pb, Cd and Ni) concentrations in rice (Oryza sativa) from Iran and associated risk assessment: A systematic review. Toxin Rev. 36, 331–341 (2017).CAS 

Google Scholar 
WHO. WHO Technical Report Series. Evaluation of Certain Food Additives and Contaminants. Fifty-Third Report of the Joint FAO/WHO Expert Committee on Food Additives (JECFA). http://www.Who.Int/foodsafety/publications/jecfa-reports/en/ (2000).USEPA. USEPA Regional Screening Level (RSL) summary table: November 2011. (2011).Farkas, A., Salánki, J. & Specziár, A. Age-and size-specific patterns of heavy metals in the organs of freshwater fish Abramis brama L. populating a low-contaminated site. Water Res. 37, 959–964 (2003).CAS 
PubMed 

Google Scholar 
Canpolat, Ö. & Çalta, M. Heavy metals in some tissues and organs of Capoeta capoeta umbla(Heckel, 1843) fish species in relation to body size, age, sex and seasons. Fresenius Environ. Bull. 12, 961–966 (2003).CAS 

Google Scholar 
Hosseini, M., Nabavi, S. M. B., Nabavi, S. N. & Pour, N. A. Heavy metals (Cd Co, Cu, Ni, Pb, Fe, and Hg) content in four fish commonly consumed in Iran: Risk assessment for the consumers. Environ. Monit. Assess. 187, 1–7 (2015).CAS 

Google Scholar 
Jiang, X. et al. Assessment of heavy metal accumulation in freshwater fish of Dongting Lake, China: Effects of feeding habits, habitat preferences and body size. J. Environ. Sci. 112, 355–365 (2022).
Google Scholar 
Yi, Y., Tang, C., Yi, T., Yang, Z. & Zhang, S. Health risk assessment of heavy metals in fish and accumulation patterns in food web in the upper Yangtze River, China. Ecotoxicol. Environ. Saf. 145, 295–302 (2017).CAS 
PubMed 

Google Scholar 
USEPA. Assessing Human Health Risks from Chemically Contaminated Fish and Shellfish: A Guidance Manual. (US Environmental Protection Agency, 1989).Means, B. Risk-assessment guidance for superfund. Volume 1. Human health evaluation manual. Part A. Interim report (Final). (Environmental Protection Agency, Washington, DC (USA). Office of Solid Waste …, 1989).Raknuzzaman, M. et al. Trace metal contamination in commercial fish and crustaceans collected from coastal area of Bangladesh and health risk assessment. Environ. Sci. Pollut. Res. 23, 17298–17310. https://doi.org/10.1007/s11356-016-6918-4 (2016).CAS 
Article 

Google Scholar 
Kalantzi, I. et al. Metals in tissues of seabass and seabream reared in sites with oxic and anoxic substrata and risk assessment for consumers. Food Chem. 194, 659–670. https://doi.org/10.1016/j.foodchem.2015.08.072 (2016).CAS 
Article 
PubMed 

Google Scholar 
Sarkar, S., Mukherjee, S., Chattopadhyay, A. & Bhattacharya, S. Differential modulation of cellular antioxidant status in zebrafish liver and kidney exposed to low dose arsenic trioxide. Ecotoxicol. Environ. Saf. 135, 173–182. https://doi.org/10.1016/j.ecoenv.2016.09.025 (2017).CAS 
Article 
PubMed 

Google Scholar 
Mandal, B. K. & Suzuki, K. T. Arsenic round the world: A review. Talanta 58, 201–235. https://doi.org/10.1016/S0039-9140(02)00268-0 (2002).CAS 
Article 
PubMed 

Google Scholar 
Kibria, G., Hossain, M. M., Mallick, D., Lau, T. C. & Wu, R. Trace/heavy metal pollution monitoring in estuary and coastal area of Bay of Bengal, Bangladesh and implicated impacts. Mar. Pollut. Bull. 105, 393–402. https://doi.org/10.1016/j.marpolbul.2016.02.021 (2016).CAS 
Article 
PubMed 

Google Scholar 
Fang, Y. et al. Concentrations and health risks of lead, cadmium, arsenic, and mercury in rice and edible mushrooms in China. Food Chem. 147, 147–151. https://doi.org/10.1016/j.foodchem.2013.09.116 (2014).CAS 
Article 
PubMed 

Google Scholar 
Vannoort, R. & Thomson, B. New Zealand Total Diet Study—Agricultural Compound Residues (Selected Contaminant and Nutrient Elements. Ministry for Primary Industries, 2009).
Google Scholar 
Praveena, S. M., Pradhan, B. & Ismail, S. N. S. Spatial assessment of heavy metals in surface soil from Klang District (Malaysia): An example from a tropical environment. Hum. Ecol. Risk Assess. Int. J. 21, 1980–2003 (2015).CAS 

Google Scholar  More

Study sitesWe sampled foliar fungal endophytes and root fungi (root endophytes and AM fungi) in the Colorado Rockies at the Rocky Mountain Biological Laboratory, Gunnison Co., Colorado, USA (38°57’N, 106°59’W). This region has predictable decreases in air temperature (c. 0.8 °C per 100 m; [40]) and declines in soil nutrients with altitude [41], but increases in precipitation, mainly as snow [42]. The entire region is warming at rates of 0.5–1.0 °C per decade [43].To capture environmental, spatial, and grass-host specific variation in fungal guilds, we sampled 66 sites encompassing 9–13 elevations from each of six altitudinal gradients in July 2014 (Supplementary Table S1, Supplementary Fig. S1). Elevational gradients represented separate mountains in the Gunnison Basin and were located within 20 km of each other. We created a regional climate model to interpolate average number of growing degree days (GDD, base 0 °C), mean annual temperature (MAT), maximum temperature (Tmax), minimum temperature (Tmin), mean annual precipitation (MAP), and mean snow depth (MSD) for each site based on data from 29 local meteorological stations [44]. At each site, soil edaphic parameters were measured on dried soil at the UC Davis soils lab (see [24] for more details) and soil nutrients at Western Ag (Saskatoon, Canada). Soil pH was measured in a 1:1 solution with diH2O, and soil moisture was measured gravimetrically. Physical characteristics of each site (e.g., aspect, soil depth, elevation) were measured as described in Lynn et al. [44]. Environmental variation across sites was large. For example, MAT varied from 7.1 to 13.3 °C, MAP from 563 to 1171 mm, and Total N from 2 to 316 ug/g dry soil (Table S1).Host plant speciesWe focused on grasses because grasslands cover ~20% of Earth’s land surface [45] and dominate subalpine meadows of the Rocky Mountains. In addition, individual grass species spanned the entire elevational range of our study system [46], whereas tree, shrub, and forb species did not. At each location, we sampled nine adult individuals from up to 13 grass species representing five genera (Poaceae, subfamily Pooideae; Supplementary Table S1). Many sites had fewer than 13 grass species present, but all sites, except for two, had at least two grass species. Samples were composited by tissue type (leaves v. roots) and grass species within each site.Fungal compositionCollected root and leaf samples were surface sterilized (1 min in 95% ethanol, 2 min in 1% sodium hypochlorite solution, and 2 min in 70% ethanol) over ice to focus on the endophytic fungal community [34]. Following surface sterilization, samples were rinsed in purified water (Milli-Q Integral Water Purification System, EMD Millipore Corporation, Billerica, MA), stored in RNAlater, and refrigerated. All samples were then frozen in liquid nitrogen and ground using a mortar and pestle. Total DNA was extracted from ~50 mg of ground sample using QIAGEN DNeasy plant extraction kits (QIAGEN Inc., Valencia, CA).Fungal composition was characterized using barcoded primers targeting the ITS2 region for leaf and root endophytes [47], and FLR3-FLR4 primers targeting ~300 bp in the 28S region for AMF [48]. Each PCR contained 5 μL of ~1–10 ng/μL DNA template, 21.5 μL of Platinum PCR SuperMix (Thermo Fisher Scientific Inc., Waltham, MA), 1.25 μL of each primer (10 μM), 1.25 μL of 20 mg/mL BSA, and 0.44 μL of 25 mM MgCl2. For the ITS2 primers, the reactions included an initial denaturing step at 96 °C for 2 min, followed by 24 cycles of 94 °C for 30 sec, 51 °C for 40 s, and 72 °C for 2 min, with a final extension at 72 °C for 10 min. For the 28S primers, reactions started with an initial denaturing step at 93 °C for 5 min, followed by 33 cycles of 93 °C for 1 min, 55 °C for 1 min, and 72 °C for 1 min, with a final extension at 72 °C for 10 min.Three PCR replicates from each sample were pooled and then cleaned and concentrated using a ZR-96 DNA Clean & Concentrator-5 (Zymo Research Corporation, Irvine, CA). PCR was then carried out on all samples to add dual indexes and Illumina sequencing adaptors; each reaction began with an initial denaturing step at 98 °C for 30 s, followed by 7 cycles of 98 °C for 30 s, 62 °C for 30 s, and 72 °C for 30 s, with a final extension at 72 °C for 5 min. Sequencing was performed by the Genomic Sequencing and Analysis Facility at The University of Texas at Austin using paired-end 250 base Illumina MiSeq v.3 chemistry (Illumina, Inc., San Diego, CA). We aimed to obtain a minimum of 30,000 reads/sample for the ITS2 region and 20,000 reads/sample for the 28S region. All sequences are deposited in the NCBI SRA database under accession number (PRJNA639093).BioinformaticsWe processed reads to generate OTUs using commands from USEARCH (v9.2.64). Reads from previous studies [24] and this study were clustered together to improve OTU delineations for a total of 36,754,931 reads. We merged paired-end reads using the fastq_mergepairs from USEARCH with “fastq_maxdiffs” set to 20 and “fastq_maxdiffpct” set to 10 to ensure proper merging at a low error rate. The merged reads and the forward unmerged reads were trimmed at the primer sites using cutadapt with “e” set to 0.2, “m” set to 200, and untrimmed reads were discarded. Merged reads were filtered using fastq_filter from USEARCH with “fastq_maxee” set to 1.0. The forward reads were first trimmed to 230 using fastx_truncate from USEARCH with “trunclen” set to 230 and then filtered by fastq_filter from USEARCH with “fastq_maxee” set to 1.0. We then concatenated the merged and forward reads into one file and de-replicated using fastx_uniques from USEARCH with “minuniquesize” set to 2. After these steps, 11,357,274 sequences remained. We clustered these sequences to form OTUs at 97% similarity [49] using cluster_otus command from UPARSE. The reads (all reads before filtering step) of each sample were mapped to OTUs with usearch_global from USEARCH with “id” set to 0.97. We determined taxonomy for the representative OTUs using sintax from USEARCH with the database set to UNITE all eukaryotes (v. 8.2) “strand” set to both and “sintax_cutoff” set to 0.8 [50]. Representative OTUs were also blasted against Genbank with “perc_identity” set to 80 and “max_target_seqs” set to 50. All OTUs identified as “fungi” were retained, and OTUs labeled as “unknown” or “unidentified” were manually inspected based on blast results to determine retention. Our filtering criteria left between 5 and 418 OTUs per sample (Supplementary Table S2).Due to low fungal abundance in leaves [34], many leaf samples were dominated by plant sequences (average ~78% plant reads). Therefore, fungal sequence numbers in leaf samples were low, despite adequate sequencing depth to capture trends in fungal endophyte communities across sites based on prior analyses [24, 34, 35]. We included only samples that contained at least 50 fungal sequences after data processing (Leaves N = 192, Roots N = 191, AMF N = 251), and most samples had much greater sequencing depth, especially for roots (Supplementary Table S2). Nevertheless, there were no correlations between sequence read depth and richness, alpha diversity, or evenness of our samples (P  > 0.05 in all cases), and plant species did not differ in the average sequencing depth for samples (P  > 0.05). Data for each fungal OTU were transformed to the proportion of total sequence abundance to minimize any differences in sampling effort [51].Diversity and compositionWe calculated the alpha diversity metrics of richness, Shannon’s Diversity, Inverse Simpson’s Diversity, and Pielou’s Evenness. For each fungal guild, differences among plant species and elevation in alpha diversity were first determined using a general linear mixed effects model with plant species (categorical) and elevation (continuous) as fixed effects and site nested within elevation gradient (e.g., mountain identity, Supplementary Table S1, Supplementary Fig. S1) as random effects to account for the lack of statistical independence among plant species sampled at the same site and among sites located within the same mountain elevation gradient (Supplementary Fig. S1). Models were constructed using the lmer function in R package lme4 [52, 53]. To address, do fungal community patterns along environmental gradients differ among guilds: leaf endophytes, root endophytes, or arbuscular mycorrhizal fungi?, we then compared alpha diversity metrics among fungal guilds using a general linear mixed effects model with fungal guild, plant species, and elevation as fixed effects and site nested within elevation gradient as random effects. In all models, we evaluated parameter fit with analysis of deviance using Wald chi-square tests and corrected for multiple comparisons using a false discovery alpha of 0.05. Differences among grass species were determined using Tukey post-hoc tests.Because elevation is a good proxy for variation in both climate and soil parameters (Supplementary Table S1), in all community analyses, we first ran models with grass species and elevation to parse biotic versus abiotic influences on fungal OTUs, then secondly ran full variance partitioning models with all environmental covariates (Supplementary Table S1, climate, physical, soil) in addition to grass species identity and space (gradient location, Supplementary Fig. S1). Because leaf and root endophytes were sequenced using different primers than AM fungi, we could not compare composition among the three guilds directly. Instead, we compared the relative influence of biotic and abiotic drivers on fungal composition within each guild to compare patterns among guilds. To do so, we first used distance-based redundancy analysis (dbRDA) to analyze the effects of plant host species and elevation on fungal composition for general fungal communities in leaves and roots and separately for AM fungal communities in roots. All models were run on quantitative Jaccard indices of fungal composition for each guild and included site nested within elevation gradient (e.g., mountain side, Supplementary Fig. S1) as random effects. Second, to evaluate which environmental variables most strongly influenced fungal composition, we further partitioned variance in fungal composition due to grass species, climate variables (MAP, MAT, MSD, Tmax, Tmin, and GDD), soil variables (total nitrogen, total phosphorus, nitrate, ammonium, calcium, magnesium, potassium, iron, manganese, sulfur, aluminum, soil pH, soil gravimetric moisture content), physical variables (aspect degree, aspect category (e.g., cardinal direction), slope, soil depth, and elevation) and spatial variables (latitude and longitude) using the varpart function in Vegan v. 2–5.3 [54]. Plots of fungal composition by plant host were also generated using dbRDA separately for each fungal guild. Spatial variables were de-trended and tested for spatial autocorrelation using the ade4 package v. 1.7–16 [55]. When we detected significant spatial autocorrelation eigenvectors, we included these in the spatial variable matrix. To characterize how many fungal taxa occurred in multiple plant taxa and elevations, we used the VennDiagram package v. 1.6.20 [56].Turnover and rewiringTo evaluate whether fungal composition was driven by grasses associating with different fungal taxa or differing relative abundances of the same fungal taxa, we first performed a beta partitioning analysis using betapart v. 1.5.3 [57]. Each fungal guild was analyzed separately. Next, to examine turnover in the abundances of fungal functional groups (pathogens, saprotrophs, mutualists), we defined groups using the FungalTrait database, which merges previous databases into one cohesive framework of 17 functional trait types (referred to here as functional groups; [58]). We recognize that fungal functions are highly environmentally dependent and therefore these functional groups may represent potential function more than actual function. Functional group identity was ascribed to 60% of leaf endophyte and 62% of root endophyte fungal taxa. Then, cumulative abundance of proportionally transformed sequence reads in each functional group was analyzed using a general linear mixed effects model with grass species and elevation as fixed effects and site nested within elevation gradient as random effects, as above. Finally, we defined indicator species within the OTUs that comprised at least 1% of the total abundance of each fungal guild by grass host, gradient, and elevation classes (rounded to the nearest 100 m) using the indicspecies package v. 1.7.9 [59]. Functional group assignments using the FungalTrait database from above were assigned to each indicator taxon [58]. A large percentage of significant indicator taxa out of the total number of OTUs would confirm that turnover in the species identity of fungal associations is stronger than turnover in the relative abundances of the same fungal taxa.Network propertiesTo address does grass-fungal network structure track elevation?, we analyzed four properties that encompass different facets of ecological networks at the site level. First, we calculated network nestedness, or the propensity for specialists to interact with the same plant species as generalists, using the weighted NODF (Nestedness metric based on Overlap and Decreasing Fill; [60]). Second, we calculated complexity as linkage density or the average number of interactions per plant species [61]. Third, to characterize specialization, we used the H2’ Index [62]. Finally, network evenness was calculated as Alatalo’s interaction evenness [63]. In all cases, these network metrics were weighted indices to increase accuracy [64], and calculations were performed in the Bipartite package v. 2.15 [65]. To address, how much do fungal guilds differ in altitudinal variation in network structure?, we compared network-level statistics among fungal guilds using a general linear mixed effects model with fungal guild as a fixed effect, number of grass hosts as a fixed effect, and gradient as a random effect (function lmer in lme4 [52],). We compared relationships with elevation separately for each fungal guild, using general linear mixed effects models with elevation as a fixed, continuous effect, number of grass hosts within the network as a fixed, continuous effect, and gradient identity as a random effect (Supplementary Table S1, Supplementary Fig. S2). We evaluated parameter fit with analysis of deviance using Wald chi-square tests using the car package 3.0–10 in R [66].All data met model assumptions of normality of residuals and homogeneity of variance. All analyses were performed in R v. 3.5.0 [53]. More

German, C. R. & Von Damm, K. L. Treatise on Geochemistry (eds Heinrich, D. H. & Karl, K. T.) 181–222 (Pergamon, 2003).Van Dover, C. The Ecology of Deep-Sea Hydrothermal Vents (Princeton University Press, 2000).Spiess, F. N. et al. East Pacific rise: Hot springs and geophysical experiments. Science 207, 1421–1433 (1980).CAS 

Google Scholar 
Haymon, R. M. et al. Hydrothermal vent distribution along the East Pacific Rise crest 9° 09’–54’ N and its relationship to magmatic and tectonic processes on fast-spreading mid-ocean ridges. Earth Planetary Sci. Lett. 104, 513–534 (1991).
Google Scholar 
Edmonds, H. N. et al. Discovery of abundant hydrothermal venting on the ultraslow-spreading Gakkel Ridge in the Arctic Ocean. Nature 421, 252–256 (2003).CAS 

Google Scholar 
German, C. R. et al. Hydrothermal activity and seismicity at Teahitia Seamount: Reactivation of the society islands hotspot? Front. Mar. Sci 7, 73 (2020).
Google Scholar 
de Ronde, C. E. J. et al. Intra-oceanic subduction-related hydrothermal venting, Kermadec volcanic arc, New Zealand. Earth Planetary Sci. Lett. 193, 359–369 (2001).
Google Scholar 
Ishibashi, J. & Urabe, T. Backarc Basins: Tectonics and Magmatism (ed Taylor, B.) 451–495 (Springer, 1995).Fouquet, Y. et al. Hydrothermal activity and metallogenesis in the Lau back-arc basin. Nature 349, 778–781 (1991).CAS 

Google Scholar 
Boschen, R. E., Rowden, A. A., Clark, M. R. & Gardner, J. P. A. Mining of deep-sea seafloor massive sulfides: A review of the deposits, their benthic communities, impacts from mining, regulatory frameworks, and management strategies. Ocean Coastal Manage. 84, 54–67 (2013).
Google Scholar 
Lisitsyn, A. P. et al. Active hydrothermal activity at Franklin Seamount, Western Woodlark Sea (Papua New Guinea). Int. Geol. Rev. 33, 914–929 (1991).
Google Scholar 
Laurila, T. E. et al. Tectonic and magmatic controls on hydrothermal activity in the Woodlark Basin: Hydrothermalism in the Woodlark Basin. Geochem. Geophys. Geosyst. 13, Q09006 (2012).Goodliffe, A. M. et al. Synchronous reorientation of the Woodlark Basin spreading center. Earth Planetary Sci. Lett. 146, 233–242 (1997).CAS 

Google Scholar 
Martínez, F., Taylor, B. & Goodliffe, A. M. Contrasting styles of seafloor spreading in the Woodlark Basin: Indications of rift-induced secondary mantle convection. J. Geophys. Res. 104, 12909–12926 (1999).
Google Scholar 
Taylor, B., Goodliffe, A., Martinez, F. & Hey, R. Continental rifting and initial sea-floor spreading in the Woodlark Basin. Nature 374, 534–537 (1995).CAS 

Google Scholar 
Schellart, W. P., Lister, G. S. & Toy, V. G. A Late Cretaceous and Cenozoic reconstruction of the Southwest Pacific region: Tectonics controlled by subduction and slab rollback processes. Earth-Sci. Rev. 76, 191–233 (2006).
Google Scholar 
Hall, R. Cenozoic geological and plate tectonic evolution of SE Asia and the SW Pacific: Computer-based reconstructions, model and animations. J. Asian Earth Sci. 20, 353–431 (2002).
Google Scholar 
Breusing, C. et al. Allopatric and sympatric drivers of speciation in Alviniconcha hydrothermal vent snails. Mol. Biol. Evol. 37, 3469–3484 (2020).CAS 

Google Scholar 
Ondréas, H., Scalabrin, C., Fouquet, Y. & Godfroy, A. Recent high-resolution mapping of Guaymas hydrothermal fields (Southern Trough). BSGF – Earth Sci. Bull. 189, 6 (2018).
Google Scholar 
Nakamura, K. et al. Water column imaging with multibeam echo-sounding in the mid-Okinawa Trough: Implications for distribution of deep-sea hydrothermal vent sites and the cause of acoustic water column anomaly. Geochem. J. 49, 579–596 (2015).CAS 

Google Scholar 
Xu, G., Jackson, D. R. & Bemis, K. G. The relative effect of particles and turbulence on acoustic scattering from deep sea hydrothermal vent plumes revisited. J. Acoust. Soc. Am. 141, 1446–1458 (2017).
Google Scholar 
Park, S.-H. et al. Petrogenesis of basalts along the eastern Woodlark spreading center, equatorial western Pacific. Lithos 316–317, 122–136 (2018).
Google Scholar 
Chadwick, J. et al. Arc lavas on both sides of a trench: Slab window effects at the Solomon Islands triple junction, SW Pacific. Earth Planetary Sci. Lett. 279, 293–302 (2009).CAS 

Google Scholar 
Fouquet, Y. et al. Geodiversity of Hydrothermal Processes Along the Mid-Atlantic Ridge and Ultramafic-Hosted Mineralization: A New Type of Oceanic Cu-Zn-Co-Au Volcanogenic Massive Sulfide Deposit (eds Rona, P. A., Devey, C. W., Dyment, J. & Murton, B. J.) Vol. 188, 321–367 (American Geophysical Union, 2010).Von Damm, K. et al. Chemistry of submarine hydrothermal solutions at 21N, East Pacific Rise. Geochim. Cosmochim. Acta 49, 2197–2220 (1985).
Google Scholar 
Seyfried, W. E. & Bischoff, J. L. Experimental seawater-basalt interaction at 300 °C, 500 bars, chemical exchange, secondary mineral formation and implications for the transport of heavy metals. Geochim. Cosmochim. Acta 45, 135–147 (1981).CAS 

Google Scholar 
Pester, N. J., Rough, M., Ding, K. & Seyfried, W. E. A new Fe/Mn geothermometer for hydrothermal systems: Implications for high-salinity fluids at 13°N on the East Pacific Rise. Geochim. Cosmochim. Acta https://doi.org/10.1016/j.gca.2011.08.043 (2011).Podowski, E. L., Moore, T. S., Zelnio, K. A., Luther, G. W. & Fisher, C. R. Distribution of diffuse flow megafauna in two sites on the Eastern Lau Spreading Center, Tonga. Deep Sea Res. Part I: Oceanogr. Res. Papers 56, 2041–2056 (2009).CAS 

Google Scholar 
Collins, P., Kennedy, R. & Van Dover, C. A biological survey method applied to seafloor massive sulphides (SMS) with contagiously distributed hydrothermal-vent fauna. Mar. Ecol. Prog. Ser. 452, 89–107 (2012).CAS 

Google Scholar 
Desbruyères, D., Hashimoto, J. & Fabri, M.-C. Composition and biogeography of hydrothermal vent communities in Western Pacific back-arc basins. Geophys. Monogr. Ser. 166, 215–234 (2006).Reid, W. D. K. et al. Spatial differences in East scotia ridge hydrothermal vent food webs: Influences of chemistry, microbiology, and predation on trophodynamics. PLoS One 8, e65553 (2013).CAS 

Google Scholar 
Van Audenhaege, L., Fariñas-Bermejo, A., Schultz, T. & Lee Van Dover, C. An environmental baseline for food webs at deep-sea hydrothermal vents in Manus Basin (Papua New Guinea). Deep Sea Res. Part I: Oceanogr. Res. Papers https://doi.org/10.1016/j.dsr.2019.04.018 (2019).Erickson, K. L., Macko, S. A. & Van Dover, C. L. Evidence for a chemoautotrophically based food web at inactive hydrothermal vents (Manus Basin). Deep-Sea Res. Part II: Top. Stud. Oceanogr. 56, 1577–1585 (2009).CAS 

Google Scholar 
Comeault, A., Stevens, C. J. & Juniper, S. K. Mixed photosynthetic-chemosynthetic diets in vent obligate macroinvertebrates at shallow hydrothermal vents on Volcano 1, South Tonga Arc—evidence from stable isotope and fatty acid analyses. Cahiers de Biologie Marine 51, 351–359 (2010).
Google Scholar 
Bennett, S. A., Dover, C. V., Breier, J. A. & Coleman, M. Effect of depth and vent fluid composition on the carbon sources at two neighboring deep-sea hydrothermal vent fields (Mid-Cayman Rise). Deep-Sea Res. Part I: Oceanogr. Res. Papers 104, 122–133 (2015).CAS 

Google Scholar 
Levin, L. A. et al. Hydrothermal vents and methane seeps: Rethinking the sphere of influence. Front. Marine Sci. 3, 1–23 (2016).
Google Scholar 
Hügler, M. & Sievert, S. M. Beyond the Calvin cycle: Autotrophic carbon fixation in the ocean. Annu. Rev. Mar. Sci. 3, 261–289 (2011).
Google Scholar 
Wang, X., Li, C., Wang, M. & Zheng, P. Stable isotope signatures and nutritional sources of some dominant species from the PACManus hydrothermal area and the Desmos caldera. PLoS One 13, e0208887 (2018).
Google Scholar 
Tunnicliffe, V. & Southward, A. J. Growth and breeding of a primitive stalked barnacle Leucolepas longa (Cirripedia: Scalpellomorpha: Eolepadidae: Neolepadinae) inhabiting a volcanic seamount off Papua New Guinea. J. Mar. Biol. Ass. 84, 121–132 (2004).
Google Scholar 
Auzende, J. M., Pelletier, B. & Lafoy, Y. Twin active spreading ridges in the North Fiji Basin (southwest Pacific). Geology 22, 63–66 (1994).
Google Scholar 
Parson, L. M. & Wright, I. C. The Lau-Havre-Taupo back-arc basin: A southward-propagating, multi-stage evolution from rifting to spreading. Tectonophysics 263, 1–22 (1996).
Google Scholar 
Thaler, A. D. et al. Comparative population structure of two deep-sea hydrothermal-vent-associated decapods (Chorocaris sp. 2 and Munidopsis lauensis) from Southwestern Pacific back-arc basins. PLoS One 9, e101345 (2014).
Google Scholar 
Lee, W.-K., Kim, S.-J., Hou, B. K., Van Dover, C. L. & Ju, S.-J. Population genetic differentiation of the hydrothermal vent crab Austinograea alayseae (Crustacea: Bythograeidae) in the Southwest Pacific Ocean. PLoS One 14, e0215829 (2019).CAS 

Google Scholar 
Plouviez, S. et al. Amplicon sequencing of 42 nuclear loci supports directional gene flow between South Pacific populations of a hydrothermal vent limpet. Ecol. Evol. https://doi.org/10.1002/ece3.5235 (2019).Tran Lu Y, A. et al. Fine-scale genomic patterns of connectivity in the deep sea hydrothermal gastropod Ifremeria nautilei over its species range using outlier scans and demo-genetic inferences. Mol. Ecol. (In Revision).Yearsley, J. M. & Sigwart, J. D. Larval transport modeling of deep-sea invertebrates can aid the search for undiscovered populations. PLoS One 6, e23063 (2011).CAS 

Google Scholar 
Mitarai, S., Watanabe, H., Nakajima, Y., Shchepetkin, A. F. & McWilliams, J. C. Quantifying dispersal from hydrothermal vent fields in the western Pacific Ocean. Proc. Natl Acad. Sci. USA 113, 2976–2981 (2016).CAS 

Google Scholar 
Marsh, L. et al. Microdistribution of faunal assemblages at deep-sea hydrothermal vents in the southern ocean. PLoS One 7, e48348 (2012).CAS 

Google Scholar 
Jollivet, D. et al. The Biospeedo cruise: A new survey of hydrothermal vents along the south East Pacific Rise from 7°24’ S to 21°33’ S. InterRidge News 13, 20–26 (2005).Girard, F. et al. Currents and topography drive assemblage distribution on an active hydrothermal edifice. Prog. Oceanogr. 187, 102397 (2020).
Google Scholar 
Hessler, R. R. & Lonsdale, P. F. Biogeography of Mariana Trough hydrothermal vent communities. Deep Sea Res. Part A. Oceanogr. Res. Papers 38, 185–199 (1991).
Google Scholar 
Fujikura, K. Biology and earth scientific investigation by the submersible ‘Shinkai 6500’ system of deep-sea hydrothermal and lithosphere in the Mariana back-arc basin. JAMSTEC J. Deep Sea Res. 13, 1–20 (1997).
Google Scholar 
Connelly, D. P. et al. Hydrothermal vent fields and chemosynthetic biota on the world’s deepest seafloor spreading centre. Nat. Commun. 3, 620 (2012).
Google Scholar 
Cline, J. D. Spectrophotometric determination of hydrogen sulfide in natural waters. Limnol. Oceanogr. 14, 454–458 (1969).CAS 

Google Scholar 
Craddock, P. R., Rouxel, O. J., Ball, L. A. & Bach, W. Sulfur isotope measurement of sulfate and sulfide by high-resolution MC-ICP-MS. Chem. Geol. 253, 102–113 (2008).CAS 

Google Scholar 
Mateo, M. A., Serrano, O., Serrano, L. & Michener, R. H. Effects of sample preparation on stable isotope ratios of carbon and nitrogen in marine invertebrates: Implications for food web studies using stable isotopes. Oecologia 157, 105–115 (2008).
Google Scholar 
Hedges, J. I. & Stern, J. H. Carbon and nitrogen determinations of carbonate-containing solids1. Limnol. Oceanogr. 29, 657–663 (1984).CAS 

Google Scholar 
Coplen, T. B. Guidelines and recommended terms for expression of stable-isotope-ratio and gas-ratio measurement results: Guidelines and recommended terms for expressing stable isotope results. Rapid Commun. Mass Spectrom. 25, 2538–2560 (2011).CAS 

Google Scholar 
Folmer, O., Black, M., Hoeh, W., Lutz, R. & Vrijenhoek, R. DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Mol. Mar. Biol. Biotechnol. 3, 294–299 (1994).CAS 

Google Scholar 
Methou, P., Michel, L. N., Segonzac, M., Cambon-Bonavita, M.-A. & Pradillon, F. Integrative taxonomy revisits the ontogeny and trophic niches of Rimicaris vent shrimps. R. Soc. Open Sci. 7, 200837 (2020).CAS 

Google Scholar 
Leigh, J. W. & Bryant, D. Popart: Full‐feature software for haplotype network construction. Methods Ecol. Evol. 6, 1110–1116 (2015).
Google Scholar  More

Study sites and samplingWe sampled gas and soil in 29 regions throughout the A (rainy tropical), C (temperate), and D (boreal) climate types of the Köppen classification from six continents during the vegetation period between August 2011 and June 2018, following a standard protocol26. According to the protocol, the gas and soil samples were collected from locations in public domain or in previous agreement with the local community and/or property owner. The samples were exported from the origin countries and imported to Estonia, EU in cooperation with customs officers of the respective states, following the legal provisions of soil export and import, specifically exemptions for scientific purposes. To capture the full range of environmental conditions in each region, we established 76 wetland soil sites under different vegetation (mosses, sedges, grasses, herbs, trees, and bare soil) and land-use types (natural—raised bog, fen, and forest; agricultural—arable, hay field and pasture; and a peat extraction area) (Fig. 1a; Supplementary Data 1). We used a four-grade land-use intensity index to quantify the effect of land conversion: 0, no agricultural land use (natural mire, swamp, or bog forest); 1, moderate grazing or mowing (once a year or less); 2, intensive grazing or mowing (more than once a year); and 3, arable land (directly fertilized or unfertilized). The vegetation and land-use intensity types and the land-use intensity index were estimated from observations and contacts with site managers and local researchers.Within the sites, we established 1–4 stations 15–500 m apart to maximize the captured environmental variation. Each of the 196 stations were equipped with 3–5 opaque PVC 65 L truncated conical chambers 1.5–5 m apart and an observation well (perforated, 50 mm diameter PP-HT pipe wrapped in geotextile; 1 m in length). From each of the 645 chambers, N2O fluxes were measured following the static chamber method37 using PVC collars (0.5 m diameter, installed to 0.1 m depth in soil). Stabilization of 3–12 h was allowed before gas sampling to reduce the disturbance effect of inserting the collars on fluxes. The chambers were placed into water-filled rings on top of the collars. Gases were sampled from the chamber headspace into a 50 mL glass vial every 20 min during a 1-h session. The vials had been evacuated in the laboratory 2–6 days before the sampling. At least three sampling sessions per location were run within 3 days. Water-table height was recorded from the observation wells during the gas sampling at least 8 h after placement. Soil temperature was measured at the 10 and 20 cm depth.Soil samples of 150–200 g were collected from the chambers at 0–10 cm depth after the final gas sampling, and transported to laboratories in Tartu, Estonia. The homogenized samples were divided into subsamples for physical–chemical analyses and DNA extraction. The samples for chemical analyses were stored at 4 °C and microbiological samples were stored at –20 °C. DNA extraction was provided at the Tartu University environmental microbiology laboratory (see details below). Using a PP-HT plastic cylinder, intact soil cores (diameter 6.8 cm, height 6 cm) for the N2 analysis with the He–O2 method38 were collected from the topsoil (0−10 cm) inside 252 chambers at 26 sites, starting from 2014. Samples from different climates were run at respective temperatures. During transport, the soil samples were kept below the ambient soil temperature from which they were collected.Gas flux analysesThe gas samples were analyzed for N2O concentration within 2 weeks using two Shimadzu GC-2014 gas chromatographs equipped with ECD, TCD, and a Loftfield-type autosampler. The N2O fluxes were determined on linear regressions obtained from consecutive N2O concentrations taken when the chamber was closed, using p  0.05 we removed one outlier. If the regression remained insignificant but the flux value fell below the gas-chromatography measuring accuracy (regression change of N2O concentration δv  More

In the work reported here, we have examined the interaction of symbiotic partners representative of the three major diversification centers. Although P. vulgaris could establish symbiosis with diverse rhizobial lineages, Rhizobium etli seemed to predominate in nature in the bean nodules collected from the Americas8,9, while the Americas is where the origin and diversification of the host have been experimentally supported19,20. Genotypes other than R. etli that also induce nodule formation in the bean have already been taxonomically defined as species, for instance Rhizobium tropici and Rhizobium ecuadorense, both of which were isolated from areas in northwestern South America, namely Ecuador, Brazil, and Colombia.American-bean rhizobia, from soil samples retrieved by the common bean as well as isolates from nodules found in nature have possessed polymorphism in the nodC gene, disclosing three nodC genotypes namely α, (upgamma), and (updelta)9. The different nodC alleles in American strains exhibit a varying predominance in their regional distributions in correlation with the centers of bean genetic diversification. The nodC types α and (upgamma) were detected both in bean nodules and in soils from Mexico, whereas the nodC type (updelta) was clearly predominant in soil and nodules from the Southern Andes (i. e., in Bolivia and northwest Argentina9). A quantitatively balanced representation of rhizobia with nodC type α and (upgamma) was detected in soils from Ecuador, but the nodC type (upgamma) had been found to be predominantly isolated from nodules formed in nature in that area5,9,10. It should be noted that we have reported finding of equal distribution of allele nodC type α and γ among the nine R etli isolates from bean in Mexico reported by Silva et al.7,9. The occurrence of this polymorphism proved to contribute to examining rhizobial populations inhabiting the Americas and characterizing the interaction with beans in BGD centers from Mexico to the northwest of Argentina. In performing our nodC analysis, we were aware that rhizobia genes for symbiosis are carried on plasmids which might mediate horizontal transfer, however in agreement with Silva et al.7 we assumed that although genetic exchange could be important, it is not so extensive to prevent epidemic clones from arising at significant frequency. Similar findings were found in R. leguminosarum bv trifolii associated with native Trifolium species growing in nature21.Investigations in the last decade have proposed an evolutionary pathway for the host bean that provided us with a framework for examining our results on rhizobia-bean interactions and facilitated an interpretation of the results. The current model proposes the occurrence of a Mesoamerican origin from where dispersion by independent migrations over time led to the Mesoamerican and Andean gene pools and to the Ecuador-Peru wild common-bean populations2,19,20. We found a balanced competition between α and (upgamma) nodC types in beans from Mesoamerica and the southern Andes, whereas the beans from Ecuador and Peru revealed a clear affinity for nodulation with strains of nodC type α rather than with the sympatric strains nodC type (upgamma) that we assayed (R. ecuadorense, CIAT894 and Bra-5). Nevertheless, we have previously reported that native strains and isolates with respectively both nodC types α and (upgamma) were found in soils and bean nodules from Mexico9, whereas lineages harboring nodC type (upgamma) were found to be predominant in beans from the northern and central regions of Ecuador-Peru8,9. The present results, however, indicated a clear affinity of the Ecuadorean-Peruvian—i. e., AHD—beans for strains nodC type α when assessed for competition against nodC type (upgamma) (Fig. 2A). We also found that nodC type (updelta) displayed a clear predominant occupancy of nodules of the AHD beans in contrast to the scarce occupancy of nodules of the Mesoamerican and Andean beans (Fig. 2B). Taken together, these results indicate no affinity of AHD beans for sympatric rhizobial strains containing nodC type (upgamma)—despite the finding that rhizobia of nodC type (upgamma) appear to predominate in isolates of nodules formed in Ecuador9,10.We conclude that although rhizobial type nodC (upgamma) was previously found to predominate in bean nodules from Ecuador, the competitiveness of that rhizobial strain for nodulation compared to other genotypes of bean rhizobia was relatively low. A possible explanation could be that seeds might be assumed to play a key role as carriers during the dissemination of the bean throughout the American regions. Thus, we can hypothesize that at the time of bean dissemination both R. etli nodC types α and (upgamma) (R. ecuadorense and other lineages) moved in conjunction with the host from Mesoamerica to northern Ecuador-Peru, but the strains bearing nodC type (upgamma) achieved an adaptation—probably due to edaphic characteristics, environmental stresses, and other features—in such a way that that strain predominated in soils and succeeded in nodulation.Alternatively, that prevalence might arise from a host selection for a rhizobium that is more effective in nitrogen fixation in a new and different environment. A poor relationship, however, between competitiveness and efficiency was found in the pea22. In line with the concept of adaptation, the bean had been found to be preferentially nodulated by species of R. tropici in acidic soils in regions of Brazil and Africa4,23. Since the environment could also be a major influence on the host and its symbiotic interactions, the Andean area represents a cooler climate for the growth of the bean than the Mesoamerican region24,25. Furthermore, since our assays were performed in laboratory environment parameters, we do not rule out the effect -if any- by the diverse and complex soil microbial community coexisting with bean rhizobia. Within this context, two contrasting soils from Argentina which differ in geolocation and edaphic properties and the perlite substrate were assayed side by side in nodule occupancy of Negro Xamapa after inoculation with a mixture of strains nodC type α and γ (Results not shown). Our results showed that the predominance of nodC type γ in the occupation of the nodules of this variety (about 80% occupation) is not affected by the type of substrate (p = 0.5566). Yet, we assume that the performance in diverse soil and ecosystems should be further evaluated in situ. In agreement, a good coevolution of rhizobia strains with nodC type (upgamma) was detected in nodules of bean varieties from the Mesoamerica and Andean genetic pools inoculated with soil samples from Mexico, Ecuador, and Northwest of Argentina, respectively (see Table 2 in Aguilar et al., 2004) [9].With respect to the interaction in the southern Andes, we propose another interpretation that takes into consideration the bottleneck that occurred before domestication in the Andes, as was indicated by Bitocchi et al.26, which scenario enables the assumption that the adaptation and concomitant diversification involved a coevolution of the symbioses. Therefore, similar profiles of competitiveness for nodulation in Mesoamerican and Andean beans were found between nodC type (upgamma) versus nodC types α and (updelta), but a significant occupancy by the nodC type (updelta) was recorded in the Andean beans.Our work suggests that the genetics of both the host and the bacteria determine the mutual affinity and additionally indicates that symbiotic interaction is another trait of legumes sensitive to the effects of evolution and ecological adaptation to the locale environment such as the characteristics of the soil and the climate.The analysis of the genetic sequences of the bean that encode genes associated with symbiosis, revealed variation of NFR1, NFR5 and NIN over the representative accessions of the Mesoamerican, the Andean, and the AHD gene pools. It is proposed that a receptor complex composed of NFR1 and NFR5 initiates signal transduction in response to Nod-factor synthesized and released by rhizobia27. Although the variation consisted mainly in neutral-amino-acid substitutions, thus suggesting only minimal changes in the functionality, if any at all; we could cite the convincing and relevant evidence reported by Radutoiu et al.27 that the amino-acid residue 118 of the second LysM module of NFR5 is essential for the recognition of rhizobia by species of Lotus japonicus and Lotus filicaulis. Our finding that the Mesoamerican-bean NFR5 has glutamine (Q) in position 151, whereas the Andean and the AHD both have proline (P)—neither of which amino acids is neutral—would merit further investigation to evaluate if such a mutation might play a role in nodulation preference. Although this result must be considered with caution, we found that the conserved polymorphism in the NFR1 and NFR5 proteins has caused the beans representative of the genetic pool Ecuador-Peru—i. e., the AHD—to be grouped in a cluster separate from those of Mesoamerica and the Andes. What we found to be interesting was that the phylogenetic and RMSD profiles of grouping the sequences are consistent with different evolutionary pathways in beans from the AHD and the Andean areas. This observation agrees with the proposal of Randón-Anaya et al.2 that those former beans from northern Peru-Ecuador originated from an ancestral form earlier than that of Mesoamerican- and Andean-bean genotypes. In addition, by applying a suppressive subtractive hybridization approach a set of bean genes were identified in our laboratory to be expressed in early step of infection by the cognate rhizobia28. Taken these results together, we conclude that genomic regions and patterns of expression in the host appear associated with an affinity for nodulation.Within a broader context, we believe that our results on the biogeography of bean-rhizobia interactions in the region where the origin and domestication of the host plants occurred provide novel useful issues to be considered in inoculation programs, for instance those involving selection of strains and cultivars, and invite to validate these findings in follow up field trials. More

The Bayesian VMD Method we developed can classify pulsed signals with similar frequency content in poor SNR files from underwater acoustic recordings. The Method consists of two parts. The first part scans the incoming audio data as segments that potentially contain signals of interest by detecting energy peaks. It then uses the start and end of the energy peaks to isolate those areas of interest from non-signal areas of the audio file. The second part classifies the detected signals into separate categories based on their frequency content. The algorithms of our Detector and Classifier steps are self-developed, but some key components in them were inspired by previous work39,40,41.DetectorThe proposed detector uses full audio files that are 4.5 s long at a sampling rate of 100 kHz. It then finds audio file segments where potential signals of interest exist.For a given audio file, denoted by ({hat{x}}(n)), where (n=1, dots , N), and N is the total number of samples, the Laplacian Differential Operator (LDO) is applied to ({hat{x}}(n)) resulting in an enhanced version of the audio file denoted by y(n), as follows:$$begin{aligned} y(n) = frac{1}{4}frac{partial ^2 {hat{x}}}{partial n^2} end{aligned}$$
(1)
The LDO enhances the transient signals (edge detection) and filters out the low frequencies ((< 10) kHz) which are not needed for Gg and Lo pulsed signal classification. The y(n) is then transformed into a time-frequency representation using Short-time Fourier transform (STFT). The STFT was implemented on 1024 samples with 90% overlap and a 1024-point Hanning window. The magnitude of the STFT matrix s(n, f) is given as ({hat{S}}(n,f)).$$begin{aligned} {hat{S}}(n,f) = begin{bmatrix} |s_{11}| &{} dots &{} |s_{1N}|\ vdots &{} ddots \ |s_{M1}| &{} &{} |s_{MN}| end{bmatrix} end{aligned}$$ (2) Where (N) is the length of the input segment and (M) is the number of frequency bins. The dimensionality of matrix ({hat{S}}(n,f)) is reduced from 2-D to 1-D as follows:$$begin{aligned} S_{d}(n) = sum _{f=1}^{M} {hat{S}}(n,f) end{aligned}$$ (3) The resulting temporal sequence is an accumulated sum of all frequency bins from (begin{aligned} {hat{S}}(n,f) end{aligned}), so scaling is applied, as follows:$$begin{aligned} S_{d}(n) = frac{S_{d}(n)}{max{S_{d}(n)}} end{aligned}$$ (4) After finding (S_{d}(n)) from Eq. (4), the mean of (S_{d}(n)) is subtracted. Then, to determine the boundaries of the acoustic signal, an adaptive threshold is applied. The first step in developing the threshold is to vectorize the matrix ({hat{S}}(n,f)) in column order into a vector called (S_{r}(n)):$$begin{aligned} S_{r}(n) = overrightarrow{{hat{S}}(n,f)} end{aligned}$$ (5) Then, (S_r (n)) is scaled similar to (S_{d}(n)) and is sorted into ascending order, denoted by ({hat{S}}_r(n)). The changing point where the root-mean-square level of the sorted curve ({hat{S}}_r(n)) changes the most is obtained by minimizing Eq. (6)39,40,42$$begin{aligned} J(k) = sum _{i=1}^{k-1} Delta ({hat{S}}_{r,i}; chi ([{hat{S}}_{r,1} dots {hat{S}}_{r,k-1}])) + sum _{i=k}^{N} Delta ({hat{S}}_{r,i}; chi ([{hat{S}}_{r,k} dots {hat{S}}_{r,N}])) end{aligned}$$ (6) where (k) and N are the index of the changing point and the length of the sorted curve ({hat{S}}_r (n)), respectively, and$$begin{aligned} sum _{i=u}^{v} Delta ({hat{S}}_{r,i}; chi ([{hat{S}}_{r,u} dots {hat{S}}_{r,v}])) = (u-v+1)log left( frac{1}{u-v+1}sum _{n=u}^{v}{hat{S}}_{r,n},^{2}right) end{aligned}$$ (7) The threshold, (lambda), is the value of ({hat{S}}_r (k)) which equals the noise floor estimation, and can be represented as follows:$$begin{aligned} begin{aligned} {mathcal {H}}_{0}: S_d(n) < lambda \ {mathcal {H}}_{1}: S_d(n) ge lambda end{aligned} end{aligned}$$ (8) where ({mathcal {H}}_0) and ({mathcal {H}}_1) are the hypothesis that the activity was below or above the threshold, respectively. The calculated threshold can vary for each file, thus making it adaptable if ambient noise conditions change between files. The threshold (lambda) is then projected onto the temporal sequence (S_{d}(n)) to extract the boundaries of the regions of the acoustic signal that comprised the detected energy peak. The start and end points of each acoustic signal are determined as the first and last points that are greater than (lambda) in amplitude.The boundaries of the detected segments are scaled by the sampling rate to obtain start and end times which will be used to extract the audio file segments from the original data file in the classification step. Figure 4 illustrates the layout of the the proposed detector.Figure 4Block diagram of the proposed detector.Full size imageClassifierOnce segments with energy peaks were identified, they were scanned by the team’s bioacoustics expert, and any segments confirmed to contain only Gg or Lo signals were sifted out for use in testing the accuracy of the Bayesian VMD Method classifier.In this paper, the metric weight was defined for classification purposes. The weight for a parameter (varvec{theta _i}) given its measurement (varvec{y_i}) is defined as$$begin{aligned} w(varvec{theta _i} mid varvec{y_i}) = P_{varvec{Theta mid Y}}(varvec{theta _i} mid varvec{y_i}) * varvec{p_i} end{aligned}$$ (9) where (varvec{theta _i}) is the probability density function (PDF) of (varvec{y_i}), (varvec{y_i}) is one measurement in the measurement vector (varvec{y}), (P_{varvec{Theta mid Y}}(varvec{theta _i} mid varvec{y_i})) is the posterior probability of the parameter (varvec{theta _i}) given the measurement (varvec{y_i}), and (varvec{p_i}) is the scaled prominence value of (varvec{y_i}).When a detected audio file segment is fed into the Bayesian VMD classifier, the classification process starts with a feature extraction step. During this step, peak and notch frequencies and their prominence values were obtained from the VMD-Hilbert spectrum of the segment. The prominence measures how much a peak stands out due to its intrinsic height or how much a notch stands out due to its depth and its location relative to surrounding peaks or notches. In general, peaks that are taller and more isolated have a higher “prominence” (p) than peaks that are shorter or surrounded by other peaks.In the feature extraction step, VMD decomposed the input audio segment into a set of IMFs. The HHT was then applied to all IMFs to create a Hilbert spectrum with a frequency resolution of 50 Hz. The Hilbert spectrum is a matrix, (H(n,f)) that contains the instantaneous energies, (h(n,f)).$$begin{aligned} H(n,f) = begin{bmatrix} h_{11} &{}dots &{} h_{1R} \ vdots &{} ddots \ h_{Q1} &{} &{} h_{QR} end{bmatrix} end{aligned}$$ (10) where r is the length of the input segment and q is the number of frequency bins in (H).The matrix (H (n,f)) is then converted from a 2-D array to a 1-D spectral representation by summing all instantaneous energy values in each frequency bin, as follows:$$begin{aligned} H(f) = sum _{n=1}^{R} H(n,f) end{aligned}$$ (11) The energy summation sequence was converted to a base-10 logarithmic scale and then smoothed by passing through a 17-point median filter and an 11-point moving average filter for the purpose of easily extracting features. All peaks and notches in the sequence whose prominence values exceeded the threshold of 0.5 were located, and their frequency values and prominence values were then stored as extracted features from the input signal (see Fig. 5).Figure 5Example of locating peak and notch frequencies and how prominent they are compared to other peaks and notches. The wave form in (a) is the smoothed energy summation sequence from the Hilbert spectrum of the Lo signal in Fig. 1. Subplot (b) is a flipped version of the energy summation sequence for the convenience of extracting notch frequencies and their prominence values. The length of the red line represents the prominence value of a peak or notch.Full size imageFor testing the effectiveness of the VMD feature extractor, a second set of features were extracted from the FFT-based power spectrum using the same input signals with the Welch’s algorithm. The FFT-based spectrum was calculated on 2048 samples with 50% overlap and a 2048-point Hanning window with 48.82 Hz frequency resolution. The power spectral density sequence was then converted to dB and went through a 21-point median filter and a 15-point moving average filter. Feature extraction followed the same strategies as in VMD feature extractor except using a prominence threshold of 2 dB.Next, the measured features, frequencies (Hz) of the peaks and notches (henceforth referred to as “measured peaks and notches”), were matched with the probability distribution functions (PDFs) of peaks and notches (henceforth referred to as “parameter peaks and notches”) from Soldevilla et al. (2008). The matching between measured and parameter peaks and notches was done in preparation of weight calculations, and it was implemented for both Gg and Lo. There are four Gaussian PDFs for parameter peaks and three for parameter notches for each species in Soldevilla et al. (2008) (Table 2). A 95% confidence interval of a Gaussian PDF was used here as a frequency range defined as 1.96 standard deviations to the left and right of its mean value. When measured peaks and notches were matched to parameter peaks and notches, only the peak or notch that fell within a 95% confidence interval were kept. Any peaks or notches outside the 95% confidence intervals were discarded.Because there are overlaps between the 95% confidence intervals of 22.4 kHz and 25.5 kHz parameter peaks of Gg and between 33.7 kHz and 37.3 kHz parameter peaks of Lo (see Table 2), it is likely that some measured peaks will fall in the overlapping areas. In this paper, the maximum a posterior (MAP) estimation41 was used to determine which PDF results in the measured peak in an overlapping area. For a measured peak that falls into an overlapping area, two parameter peaks’ PDFs are plugged in the MAP estimation equation sequentially, and then the measured peak will be matched with the PDF that maximizes the posterior probability of it given the measured peak.After the preliminary match, if more than one measured peak or notch remains in any one PDF confidence interval, the measured peak and notch with the highest prominence value is selected as the real measured peak or notch of this PDF, and the redundant ones are discarded. Finally, all remaining peak prominence values and notch prominence values were scaled to be between 0 and 1, respectively.Once peak and notch matching and selection was finished, Bayesian weights were calculated to select the most likely species. From Bayes’s rule, the posterior probability of a parameter given its measurement is proportional to the product of the likelihood function of the measurement given the parameter and the prior probability of the parameter41, as shown in Eq. (12).$$begin{aligned} P_{varvec{Theta mid Y}}(varvec{theta _i} mid varvec{y_i}) propto f_{varvec{Y mid Theta }}(varvec{y_i} mid varvec{theta _i}) P_{varvec{Theta }}(varvec{theta _i}) end{aligned}$$ (12) therefore, substitution of the posterior probability in Eq. (9) yields$$begin{aligned} w(varvec{theta _i} mid varvec{y_i}) = f_{varvec{Y mid Theta }}(varvec{y_i} mid varvec{theta _i}) *P_{varvec{Theta }}(varvec{theta _i}) * varvec{p_i} end{aligned}$$ (13) Figure 6Example of feature matching. The top plots show a set of measured peaks and notches matched with both Gg’s PDFs (a) and Lo’s PDFs (b) parameter peaks and notches like in Fig. 5 during the feature matching and selection step. Middle plots show how closely to the parameter PDFs that the measured peaks match either Gg (c) or Lo (d) and their weight calculations. The width of each PDF represents its 95% confidence interval, and the ordinate represents the weight value. Subplots (e) and (f) show the same weight calculations for notches. The final weight value is the summation of all weight values of peaks and notches matched with Gg or Lo.Full size imageWith all PDFs and a priori probabilities from Soldevilla et al. (2008), the weight value in terms of Gg and Lo given a set of measurements, (varvec{y}), was obtained by Eqs. (13) and (14)$$begin{aligned} w(Gg mid varvec{y}) = sum _{forall i} w(varvec{theta _i} mid varvec{y_i}) qquad w(Lo mid varvec{y}) = sum _{forall j} w(varvec{theta _j} mid varvec{y_j}) end{aligned}$$ (14) where (varvec{y_i}) and (varvec{y_j}) are the remaining measured peaks and notches that were matched with Gg’s PDFs and Lo’s PDFs after the matching and matching step. The feature matching and selection results and the weight calculation process are shown in Fig. 6.The last step was a comparison between weight values in terms of Gg and Lo. If (w(Lo mid varvec{y}) > w(Gg mid varvec{y})), the signal was labeled an Lo signal; otherwise, it was labeled a Gg signal. The classifier is illustrated in Fig. 7. The weight values are significant to three digits because weights are normally smaller than 1.000 and three significant digits was sufficient for comparing all calculated weight values for these audio files. In the case that the weight comparison is equal to three significant digits (even though this never happened in these 174 signals), the Bayesian VMD algorithm will automatically classify the input as a Gg signal given that the highest precision (85.91%) by the Bayesian VMD Method was achieved on Gg.Figure 7Block diagram of the Bayesian VMD Method classifier.Full size image More

Our research determined noise to share a potentially important negative relationship with native bird richness and abundance and appears to be the most limiting factor in green space occupancy by native bird species, more so than the type and amount of vegetation present in urban green spaces, and more so than urbanization itself, represented as building height and cover surrounding green spaces. Thus, noise is potentially acting as an invisible source of habitat degradation, limiting the bird species capable of inhabiting an area, regardless of whether the appropriate vegetative conditions exist.As predicted, native urban avoiders reached their maximum abundances in PAR, which, given their high vegetation cover and large size, act as patches of natural habitat in cities. Native urban utilizers tended to be found in more suburban areas, and urban dwellers, both native and exotic, were detected in green spaces of all noise levels. All exotic bird species were urban dwellers, referring to their high tolerance to urbanization5,25, thus reaching the high abundances observed, particularly in SGS.SGS possessed higher average noise levels and greater exotic bird abundance than PAR, which presented significantly higher numbers of native bird richness and abundance. The potential influence of noise on native bird species first becomes evident when we consider that native bird abundance tended to rise above the generally high abundance of exotic birds when average noise levels in green spaces reached below 52 dB (it should be noted that, according to the Chilean Noise Norm No. 146, the maximum allowable noise levels generated by fixed sources in residential areas of Santiago is 55 dB during the day, 7 a.m.–9 p.m.). The negative relations between noise and urban avoider, urban utilizer, and urban dweller species richness and abundance further indicate how noise may be regulating the native bird species present in green spaces, affecting urban avoider richness the most and urban dweller richness the least, while influencing the abundance of all native bird species rather similarly. Meanwhile, building height surrounding green spaces negatively influenced native urban avoider and urban dweller richness and abundance, with the greatest influence on urban dweller abundance, yet all native birds were less likely to be detected in green spaces surrounded by buildings over 10 m tall on average.The importance of vegetation for native bird communities also cannot be denied, given that native birds reached higher abundances than exotic birds when vegetation cover reached an average NDVI value greater than 0.5. Results from this study thus suggest that exotic birds begin to replace native birds in terms of abundance as noise levels rise in urban green spaces, vegetation cover decreases, and building height surrounding green spaces increases, with native urban avoider species being the least tolerant to the influences of urbanization, and, consequently, the first to disappear when noise levels and building height become too great. The observed negative relationship between native bird species richness and maximum noise levels, and the positive relationship with vegetation cover, are comparable to results seen in other Neotropical cities24,26, yet our results indicate that the relationships between these variables and bird abundance are stronger. This may indicate how bird abundance fluctuates in green spaces as some birds temporarily leave during noisy events or become quieter and more cryptic under noisy conditions26, while noise also negatively influences bird species richness by filtering the species that can inhabit areas of varying noise levels.Detection probability models found native bird detectability to mostly increase with vegetation cover and tree cover in urban green spaces, except for the common diuca finch, whose detectability decreased with rising tree cover. Some of the bird species that displayed the lowest detection probabilities, such as the picui ground dove and fire-eyed diucon (Xolmis pyrope), are not frequently found in cities and possess vocalizations that are unlikely to be heard well in high-noise areas due to their low frequencies, making them more easily masked by the anthrophony, characterized by its low frequency and high intensity31. Consequently, birds whose vocalizations are similar in frequency and amplitude to the anthrophony were more commonly or exclusively found in green spaces that registered low noise levels, their detectability also decreasing with rising noise, as was the case with the fire-eyed diucon.Urban green space occupancy by native bird species was mainly influenced by average maximum noise levels recorded in green spaces. Of the modeled native species, the long-tailed meadowlark and the picui ground dove, an urban avoider and an urban utilizer species respectively, were the species most sensitive to noise, their probability of occupying green spaces with average maximum noise levels over 55 dB decreasing rapidly and approaching zero when over 65 dB. Meanwhile, the austral thrush, an urban dweller species, was by far the most tolerant to noise of the native birds, its presence probability just beginning to decrease when average maximum noise levels reached over 73 dB in green spaces. The differing tendencies of urban avoiders, urban utilizers, and urban dwellers to occupy green spaces of varying noise levels is thus evident, with native urban dweller species more likely to occupy higher noise urban green spaces than urban avoiders and utilizers, seemingly more adapted to the high noise levels that come with inhabiting a busy city. Nonetheless, although native urban dwellers displayed greater noise tolerances than urban avoiders and utilizers, their presence in city parks can also be expected to diminish if noise levels become too high, which for the most tolerant of the native birds, means reaching an average maximum level of 73 dB or more, but 55 dB or more for less tolerant species.No relation was found between vegetation cover and noise, and some of the highest noise levels were recorded in PAR. This suggests that PAR, often considered to be quiet and peaceful areas to escape the busyness of city life, can reach noise levels as high as those recorded in SGS, reducing the quality of the greatest sources of natural habitat for birds and other wildlife in cities.The results from this study regarding the influence of noise on bird communities support previous studies indicating that birds may be excluded from suitable habitats on account of the acoustic conditions of the local environment12,15. Despite abundant vegetation in PAR and some SGS, certain bird species, particularly urban avoiders and utilizers, were less likely to occupy areas that presented high noise levels. However, it is important to consider other potential influencing factors, such as predators (e.g., dogs and cats) and food availability, both of which could be linked to pedestrians and could therefore also increase noise levels in green spaces. Furthermore, in an effort to focus on the influence of anthropogenic variables on urban birds (i.e., urban morphology, noise, and vegetation type and cover), this study did not consider the size of urban green spaces as a variable in occupancy modeling, but as the results of this study suggest and others in Latin America have shown23,32, green space size is likely an influencing factor that should be considered in future studies. Another variable worth considering would be road coverage, which undoubtedly plays a role in noise levels, particularly for SGS.Measures to control the COVID-19 pandemic have significantly reduced noise levels in major cities worldwide33,34,35. Noise reduction in the San Francisco Bay Area, characterized by a Mediterranean climate like Santiago, resulted in songbirds rapidly occupying newly available acoustic niches within urban soundscapes and maximizing communication through higher performance songs35. Consequently, native bird species not commonly found in high-noise areas, mainly urban avoider and utilizer species, may now be found in greater abundance at the community level in urban green spaces where they had been scarce or non-existent during this study, conducted pre-pandemic. Furthermore, if average noise levels dropped below 52 dB in Santiago green spaces due to region-wide shut-down measures, native birds may reach higher abundances than exotic birds. The negative effects of urban noise on bird communities are extensive, yet recent research indicating birds’ rapid adaptability and improved vocal performance when noise levels are significantly lowered provides hope. Native bird species susceptible to noise may stand a chance despite growing urbanization, if noise levels in urban green spaces are regulated.Rapid urban expansion in Latin America places natural ecosystems at great risk, reducing or altogether eliminating natural habitats for native birds and other wildlife, making urban green spaces necessary for their persistence, especially in biodiversity hotspots like central Chile. As this study illustrates, noise associated with urbanization plays a significant role in influencing green space occupancy by native bird species, and, quite possibly, other animal species dependent on acoustic signaling (e.g., amphibians and mammals). Given the recreational role of urban green spaces in cities, noise regulation within these areas should be considered, while also considering how city morphology may impact bird communities. This study exemplifies how, in addition to noise, the size of urban green spaces and the vegetation cover in them, particularly tree cover, are vital aspects to consider in city planning in order to preserve native bird communities in urban systems. Large urban parks held significantly richer bird communities than small green spaces, with greater native bird richness and abundance. Therefore, it is imperative that science and city planning collaborate to develop cities with networks of large green spaces with abundant tree cover, surrounded by smaller urban morphology, where noise is regulated and maintained at tolerable levels for native birds. There is a clear need to move towards biophilic city planning to harmonize urban growth and the protection and expansion of networks of green areas that generate habitat for birds that, in turn, provide important ecosystem services to cities. More

Department of Zoology, University of Cambridge, Cambridge, UKYucheng Wang, Bianca De Sanctis, Ana Prohaska, Daniel Money & Eske WillerslevLundbeck Foundation GeoGenetics Centre, GLOBE Institute, University of Copenhagen, Copenhagen, DenmarkYucheng Wang, Mikkel Winther Pedersen, Fernando Racimo, Antonio Fernandez-Guerra, Alexandra Rouillard, Anthony H. Ruter, Hugh McColl, Nicolaj Krog Larsen, James Haile, Lasse Vinner, Thorfinn Sand Korneliussen, Jialu Cao, David J. Meltzer, Kurt H. Kjær & Eske WillerslevThe Arctic University Museum of Norway, UiT— The Arctic University of Norway, Tromsø, NorwayInger Greve Alsos, Eric Coissac, Marie Kristine Føreid Merkel, Youri Lammers & Galina GusarovaDepartment of Genetics, University of Cambridge, Cambridge, UKBianca De Sanctis & Richard DurbinUniversité Grenoble Alpes, Université Savoie Mont Blanc, CNRS, LECA, Grenoble, FranceEric CoissacCenter for Macroecology, Evolution and Climate, GLOBE Institute, University of Copenhagen, Copenhagen, DenmarkHannah Lois Owens, Carsten Rahbek & David Nogues BravoDepartment of Geosciences, UiT—The Arctic University of Norway, Tromsø, NorwayAlexandra RouillardUniversité Paris-Saclay, CEA, CNRS, Institute for Integrative Biology of the Cell (I2BC), Gif-sur-Yvette, FranceAdriana AlbertiGénomique Métabolique, Genoscope, Institut François Jacob, CEA, CNRS, Université Evry, Université Paris-Saclay, Evry, FranceAdriana Alberti, France Denoeud & Patrick WinckerInstitute of Earth Sciences, St Petersburg State University, St Petersburg, RussiaAnna A. Cherezova & Grigory B. FedorovArctic and Antarctic Research Institute, St Petersburg, RussiaAnna A. Cherezova & Grigory B. FedorovSchool of Geography and Environmental Science, University of Southampton, Southampton, UKMary E. EdwardsAlaska Quaternary Center, University of Alaska Fairbanks, Fairbanks, AK, USAMary E. EdwardsCentre d’Anthropobiologie et de Génomique de Toulouse, Université Paul Sabatier, Faculté de Médecine Purpan, Toulouse, FranceLudovic OrlandoNational Research University, Higher School of Economics, Moscow, RussiaThorfinn Sand KorneliussenDepartment of Geography and Environment, University of Hawaii, Honolulu, HI, USADavid W. BeilmanDepartment of Geosciences and Natural Resource Management, University of Copenhagen, Copenhagen, DenmarkAnders A. BjørkCarlsberg Research Laboratory, Copenhagen, DenmarkChristoph Dockter & Birgitte SkadhaugeCenter for Environmental Management of Military Lands, Colorado State University, Fort Collins, CO, USAJulie EsdaleFaculty of Biology, St Petersburg State University, St Petersburg, RussiaGalina GusarovaDepartment of Glaciology and Climate, Geological Survey of Denmark and Greenland, Copenhagen, DenmarkKristian K. KjeldsenDepartment of Earth Science, University of Bergen, Bergen, NorwayJan Mangerud & John Inge SvendsenBjerknes Centre for Climate Research, Bergen, NorwayJan Mangerud & John Inge SvendsenUS National Park Service, Gates of the Arctic National Park and Preserve, Fairbanks, AK, USAJeffrey T. RasicZoological Institute, , Russian Academy of Sciences, St Petersburg, RussiaAlexei TikhonovResource and Environmental Research Center, Chinese Academy of Fishery Sciences, Beijing, ChinaYingchun XingCollege of Plant Science, Jilin University, Changchun, ChinaYubin ZhangDepartment of Earth and Atmospheric Sciences, University of Alberta, Edmonton, Alberta, CanadaDuane G. FroeseCenter for Global Mountain Biodiversity, GLOBE Institute, University of Copenhagen, Copenhagen, DenmarkCarsten RahbekSchool of Environment, Earth and Ecosystem Sciences, The Open University, Milton Keynes, UKPhilip B. Holden & Neil R. EdwardsDepartment of Anthropology, Southern Methodist University, Dallas, TX, USADavid J. MeltzerDepartment of Geology, Quaternary Sciences, Lund University, Lund, SwedenPer MöllerWellcome Trust Sanger Institute, Wellcome Genome Campus, Cambridge, UKEske WillerslevMARUM, University of Bremen, Bremen, GermanyEske Willerslev More