key: cord-0028741-jz8xlpoj
authors: Echeverri, Alejandra; Smith, Jeffrey R.; MacArthur-Waltz, Dylan; Lauck, Katherine S.; Anderson, Christopher B.; Monge Vargas, Rafael; Alvarado Quesada, Irene; Wood, Spencer A.; Chaplin-Kramer, Rebecca; Daily, Gretchen C.
title: Biodiversity and infrastructure interact to drive tourism to and within Costa Rica
date: 2022-03-04
journal: Proc Natl Acad Sci U S A
DOI: 10.1073/pnas.2107662119
sha: bc4c27916a26b939e314da976d9c1a2bfb598a9c
doc_id: 28741
cord_uid: jz8xlpoj

Nature-based tourism has potential to sustain biodiversity and economic development, yet the degree to which biodiversity drives tourism patterns, especially relative to infrastructure, is poorly understood. Here, we examine relationships between different types of biodiversity and different types of tourism in Costa Rica to address three questions. First, what is the contribution of species richness in explaining patterns of tourism in protected areas and country-wide in Costa Rica? Second, how similar are the patterns for birdwatching tourism compared to those of overall tourism? Third, where in the country is biodiversity contributing more than other factors to birdwatching tourism and to overall tourism? We integrated environmental data and species occurrence records to build species distribution models for 66 species of amphibians, reptiles, and mammals, and for 699 bird species. We used built infrastructure variables (hotel density and distance to roads), protected area size, distance to protected areas, and distance to water as covariates to evaluate the relative importance of biodiversity in predicting birdwatching tourism (via eBird checklists) and overall tourism (via Flickr photographs) within Costa Rica. We found that while the role of infrastructure is larger than any other variable, it alone is not sufficient to explain birdwatching and tourism patterns. Including biodiversity adds predictive power and alters spatial patterns of predicted tourism. Our results suggest that investments in infrastructure must be paired with successful biodiversity conservation for tourism to generate the economic revenue that countries like Costa Rica derive from it, now and into the future.

Nature-based tourism has potential to sustain biodiversity and economic development, yet the degree to which biodiversity drives tourism patterns, especially relative to infrastructure, is poorly understood. Here, we examine relationships between different types of biodiversity and different types of tourism in Costa Rica to address three questions. First, what is the contribution of species richness in explaining patterns of tourism in protected areas and country-wide in Costa Rica? Second, how similar are the patterns for birdwatching tourism compared to those of overall tourism? Third, where in the country is biodiversity contributing more than other factors to birdwatching tourism and to overall tourism? We integrated environmental data and species occurrence records to build species distribution models for 66 species of amphibians, reptiles, and mammals, and for 699 bird species. We used built infrastructure variables (hotel density and distance to roads), protected area size, distance to protected areas, and distance to water as covariates to evaluate the relative importance of biodiversity in predicting birdwatching tourism (via eBird checklists) and overall tourism (via Flickr photographs) within Costa Rica. We found that while the role of infrastructure is larger than any other variable, it alone is not sufficient to explain birdwatching and tourism patterns. Including biodiversity adds predictive power and alters spatial patterns of predicted tourism. Our results suggest that investments in infrastructure must be paired with successful biodiversity conservation for tourism to generate the economic revenue that countries like Costa Rica derive from it, now and into the future.

conservation j earth observations j recreation j rural livelihoods j species distribution models T he tourism sector is well-poised to generate win-win approaches to biodiversity conservation and sustainable development given its nonextractive nature and its dependence on scenic beauty (1) . However, for sustainable tourism to succeed as a strategy for biodiversity conservation, the role that biodiversity plays in driving tourism patterns needs to be better understood. On one hand, wildlife and nature motivate a significant portion of global tourism (2) , and protected areas with higher species richness tend to attract more tourists and yield higher economic benefits (3) . On the other hand, tourism hotspots also tend to occur in places where more human-built infrastructure (e.g., hotels, roads, and airports) enables access (4, 5) . Studies have reached mixed conclusions on the relative importance of biodiversity and accessibility for tourism, and little is known about how they work in concert (6) (7) (8) (9) (10) (11) . Given the potential negative impacts of infrastructure on biodiversity conservation, their relative contributions to tourism deserves explicit study, particularly in developing countries where both biodiversity conservation and sustainable development are urgently needed (7, 12) .

Drivers of tourism patterns across landscapes have been explored through questionnaire surveys and structured interviews that ask tourists about their affinity for landscape features (13) (14) (15) and through spatial models that predict recreation using photographs (e.g., geographically weighted regression, MaxEnt) (16, 17) . Recently, geo-tagged photographs and species lists shared on social media platforms have become popular tools for tourism-focused research (18) (19) (20) (21) (22) . These studies, however, typically focus narrowly on the role of single taxa (22, 23) and landscape attributes without accounting for species diversity (14, 20) , or they focus only on the role of infrastructure as a driver of tourism (5) . Recent advances in satellite Earth observations make it possible to capture more of the ecosystem heterogeneity that can drive variability in species distributions, compared to more conventional modeling based on land cover (24) (25) (26) . An integrated approach is needed, linking Significance Tourism accounts for roughly 10% of global gross domestic product, with nature-based tourism its fastest-growing sector in the past 10 years. Nature-based tourism can theoretically contribute to local and sustainable development by creating attractive livelihoods that support biodiversity conservation, but whether tourists prefer to visit more biodiverse destinations is poorly understood. We examine this question in Costa Rica and find that more biodiverse places tend indeed to attract more tourists, especially where there is infrastructure that makes these places more accessible. Safeguarding terrestrial biodiversity is critical to preserving the substantial economic benefits that countries derive from tourism. Investments in both biodiversity conservation and infrastructure are needed to allow biodiverse countries to rely on tourism for their sustainable development. species richness of multiple taxa along with infrastructure variables, both modeled and mapped through high resolution Earth observations. Such an approach could be scaled up to larger regions and applied globally, helping to identify where biodiversity is playing a significant role in driving tourism, such that governments and the tourism sector can prioritize investments in biodiversity conservation.

Here, we ask three questions in the iconic case of Costa Rica. First, what is the contribution of species richness (of vertebrate taxa) in explaining patterns of tourism in protected areas and also country-wide? Second, how similar are the patterns for birdwatching tourism compared to those of overall tourism? Third, where in the country is biodiversity contributing more than other factors to birdwatching tourism and to overall tourism? We predict that vertebrate species richness is more important for driving tourism in protected areas than in the rest of the country, because nature-seeking tourists often go to protected areas to find wildlife (3) . We also expect that birdwatching tourism is predicted by richness of threatened and endemic bird species rather than total species richness, given birdwatchers' preferences for rare birds (27) . We predict a saturating relationship between species richness and tourism, because beyond a large number of species additional species are unlikely to contribute more to tourism (28) . Finally, we predict that national-level tourism is better explained by infrastructure (such as roads and hotels) and distance to water than by biodiversity, because tourists going to Costa Rica often seek activities such as surfing and relaxing in beach resorts (29, 30) . We predict nonlinear effects of proximity to roads and water, because a place is deemed inaccessible if it is further away from roads, and a beach tourist destination is also either close to water or not a destination at all. Access diminishes rapidly over a few miles (31) .

To answer our research questions, we analyze patterns of tourism at two different spatial scales. First, we analyze tourism patterns in protected areas only. Many protected areas provide reliable data on visitation rates and biodiversity tends to be higher inside than outside protected areas (29) . However, given that they are often visited by tourists who are already interested in biodiversity, protected areas are not representative of all tourism patterns (30) . Second, we evaluate the relative importance of biodiversity to all tourism across Costa Rica (excluding offshore islands). Investigating tourism across the whole nation may give a better understanding of how biodiversity contributes to tourism writ large and not only for tourists with a predisposition for finding wildlife. At both the protected area and national scales, we use a modified MaxEnt model that integrates species distribution models for 66 terrestrial vertebrate species (including amphibians, reptiles, and mammals) and 699 bird species, based on remotely sensed climate and habitat variables, with spatial patterns of infrastructure (hotel density and distance to roads) and distance to water. We measure tourism in two ways for both scales: using eBird checklists as a proxy for birdwatching (32) , and using Flickr photographs as a proxy for all international and domestic tourism (19) .

Costa Rica is an ideal country to explore these questions because tourism represents 7% of the national gross domestic product and employs 3% of the working population directly and a further 9% indirectly (33) . Approximately 70% of all international visitors to Costa Rica state that the wildlife, dramatic scenery, and opportunities for adventure sports are the main motivation for visiting the country (33, 34) . However, the importance of biodiversity as a factor that influences tourism to, or domestically within, the country has not been evaluated (apart from very local studies) (35, 36) . The Central Bank of Costa Rica is currently piloting a nature-based tourism account under the United Nations System of Environmental Economic Accounts (UN SEEA). The state of the art with this methodology is to attribute value to different ecosystems, which may vary widely in their biodiversity. Understanding the relationships between biodiversity and tourism in Costa Rica is a key step toward maintaining the vibrant ecotourism industry and can serve as an example for other biodiverse nations that often look to Costa Rica as a leader in sustainable development (37) .

Tourism inside Protected Areas and Country-wide. As expected, vertebrate biodiversity (species richness of amphibians, reptiles, and mammals) matters more to tourism within protected areas than country-wide. The spatially explicit model aimed at explaining tourism within protected areas, as measured by density of Flickr photographs-user-days (PUDs) per protected area shows a better fit (area under the curve [AUC] = 0.680; SI Appendix, Table S6 -Model 5) than the nation-wide model (AUC = 0.608; SI Appendix, Table S6 -Model 10). The protected area model shows that protected area size (50.8%), hotel density (27.1%), distance to roads (10.3%), and distance to water (11.1%) are the most important predictors in terms of permutation importance (SI Appendix, Table S6-Model 5). We find negative effects of protected area size, distance to roads, and distance to water, suggesting that larger parks and parks further from roads and water have lower densities of tourists compared to smaller parks with greater accessibility ( Fig. 1 and SI Appendix, Fig. S1 ).

The species richness of amphibians, reptiles, and mammals, is the least important predictor with the lowest values of permutation importance (0.6%) and shows a curved pattern, suggesting that parks with intermediate levels of vertebrate species richness have more tourists, though effects are muted (SI Appendix, Fig. S1 ). Importantly, the role of amphibians, reptiles and mammals is larger in explaining tourism inside protected areas compared to just bird species richness (0.4% permutation importance, SI Appendix, Table S6 -Model 3, AUC = 0.679), and much lower compared to the role of richness of threatened and endemic birds (10%, SI Appendix, Table S6 -Model 4, AUC = 0.690). Nonetheless, in the three models, the effect of biodiversity indicates that tourism is higher with intermediate levels of diversity (SI Appendix, Fig. S1 ). In all of the models, the main predictor is protected area size with a negative relationship (SI Appendix, Table S6 - .

Analyses at the country-wide scale reveal that accessibility is the main driver of tourism. Indeed, the main predictors explaining tourism patterns, as measured by Flickr PUDs, are distance to roads (∼50% of permutation importance in all three models, SI Appendix, Table S6 -Models 8-10), and hotel density (∼20% of permutation importance in all models, SI Appendix, Table S6 -Models 8-10). Distance to roads also has a negative effect, reaffirming that places closer to roads have more tourists. Hotel density has a quadratic effect, suggesting that people visit areas with intermediate levels of hotel density ( Fig. 2 and SI Appendix, Fig. S2 ). Distance to protected areas is the third most important predictor across models (∼20% of permutation importance in all three models) and it exhibits a negative linear effect, indicating that being close to a protected area yields more tourism.

In keeping with the protected area analysis, the biodiversity variable that better predicts overall tourism at the country-wide scale is the richness of threatened and endemic birds (5.3%; AUC = 0.612, SI Appendix, Table S6 -Model 9). The importance of amphibians, reptiles, and mammals (0.1%, SI Appendix, Table S6 -Model 10) and bird species richness (3.8%, SI Appendix, Table S6 -Model 8) are lower. Nonetheless, in all models, we find a positive effect of biodiversity in explaining tourism ( Fig. 2 and SI Appendix, Fig. S2 ). Last, distance to water has a permutation importance ranging from 1.9 to 3.1% and has a negative nonlinear effect, suggesting that people prefer being closer to water, within a few kilometers and generally not more.

Birdwatching tourism is patterned similarly to overall tourism but exhibits greater nonlinearities with regards to infrastructure and biodiversity variables. Models explaining birdwatching inside protected areas are the best performing models across our entire analysis, across the two different scales and two types of tourism (AUC = 0.776 for the model with all birds, and AUC = 0.780 for the model with threatened and endemic birds, SI Appendix, Table  S6 -Models 1, 2). The main predictors for birdwatching inside protected areas are distance to roads (61.1% and 54.2% for the model with bird species richness [Model 1] and for the model with richness of threatened and endemic birds [Model 2] respectively), and protected area size (13.1% and 14.8%), both showing negative nonlinear effects ( Fig. 1 and SI Appendix, Fig. S1 ). We find that the richness of threatened and endemic birds (7.1%) is a better predictor of birdwatching compared to all birds (2.1%). We find that birdwatching changes nonlinearly with respect to bird species richness, showing that birdwatching tourism peaks at 200 species and at~100 species of threatened and endemic birds (Fig. 1) .

Our model examining birdwatching tourism at the country level, as measured by eBird checklists with bird richness as the biodiversity predictor variable, has a similar fit compared to the model with threatened and endemic bird richness (AUC = 0.723, Model 6 vs. AUC = 0.721, Model 7). Distance to roads (∼58%) and hotel density (∼13%) are the main predictors of birdwatching tourism nationally (SI Appendix, Table S6 -Models 6 and 7). We find stronger negative relationships between predictor variables and birdwatching tourism than overall tourism (Fig. 2) . Specifically, we find negative quadratic effects for distance to roads (SI Appendix, Fig. S2 ) and higher negative slopes of distance to water (SI Appendix, Fig. S2 ) compared to birdwatching at the protected area level (SI Appendix, Fig. S1 ), suggesting that being close to a road or to water leads to more birdwatching up to certain distances (<50 m for roads and <50 m for water) before dropping significantly (Fig. 2) . Conversely, we find a positive saturating effect for hotel density, indicating that places with more hotels have more birdwatching but saturate at 0.75% of hotel density. Species richness of threatened and endemic birds is the third most important predictor in explaining birdwatching tourism at the countrywide scale (12.1%, SI Appendix, Table S6-Model 7) . Importantly, we find that the richness of threatened and endemic species is a better predictor of birdwatching patterns compared to all birds (10.3%, SI Appendix, Table S6-Model 6).

Contribution of Biodiversity to Spatial Patterns of Tourism. Applying our spatial regressions to national maps clearly demonstrates the importance of Costa Rica's major population centers and road networks to tourism patterns (areas in dark green in Fig. 3 A and D) . When the models are run with and without biodiversity, new patterns emerge. Patterns are stronger for Birdwatching and overall tourism decline when distance to roads increases (D), and with increased protected area size (E), as well as increased distance to water (F). Protected areas with intermediate levels of hotel density experience the highest tourism (G). These curves show the results of ten k-fold replicated MaxEnt runs using either Flickr PUDs (solid lines) or eBird CUDs (dotted lines) as the response variable and the corresponding environmental variable as the predictor. The line represents the mean of the ten replicates, while the shaded region is drawn between the minimum and maximum prediction of the 10 replicates.

birdwatching tourism (Fig. 3B ) than for overall tourism (Fig.  3E ), but in both cases, without considering biodiversity, the tourism value of the Cordillera de Talamanca (dark red areas in Fig. 3 B and E) is underestimated, while the tourism value of the Nicoya Peninsula is overestimated (blue areas in Fig. 3 B and E).

Models with and without infrastructure show that tourism in intact tracts of forest is overestimated, such as the Amistad International Park (dark blue area in the southeast), the Corcovado National Park in the Osa Peninsula (protected area in the southwest), and the Barra del Colorado Wildlife Refuge (northeast of the country), whereas areas along the Pacific coastline and near urban centers tourism are underestimated (Fig. 3 C  and F) .

Results show that tourism patterns are higher when both biodiversity and infrastructure are considered together, especially in areas near Monteverde, and the Volc an Po as National Park (Fig. 3 , center of the country). We find some spatial mismatches between tourism and birdwatching; in particular, birdwatching is more reliant on the road network than overall tourism, and certain protected areas are more important for tourism than for birdwatching, such as the Amistad International Park and the Barra del Colorado Wildlife Refuge (Fig. 3 A and D) .

Our results indicate that if Costa Rica invests in both infrastructure development and biodiversity conservation, they will continue deriving economic benefits from their natural assets (ceteris paribus). We show that while the role of infrastructure, such as hotels and roads, is larger in driving tourism compared to the role of biodiversity, tourism is highest in places where both biodiversity and infrastructure are present, particularly in mountainous areas. It is well-known that the tourism footprint is strongly associated with access and is particularly dependent on the presence of roads (31, 38, 39) . Our findings reaffirm that proximity to roads increases tourism inside protected areas and country-wide. Even for birdwatchers, being close to roads is a better predictor of where they go rather than the richness of threatened and endemic birds. Thus, our results echo prior studies, finding that accessibility is more important than the ecological value of a place in driving tourism, including naturebased tourism (40) . However, failing to adequately account for nature's contributions will overestimate the value of highly accessible but ecologically depauperate places, and will underestimate the importance of sustainable development to maintaining current levels of benefits.

Distance to water is less important in driving tourism patterns than we expected. Considering that ∼70% of all international visitors to Costa Rica list ecotourism, adventure, and beaches as their main motivation for visiting the country (34), we expected proximity to water to play a larger role in driving tourism country-wide than our analyses suggest. While beachgoers and surfers constitute a significant portion of tourism, many visitors are seeking other activities such as birdwatching that are not dependent on proximity to water. These results Fig. 2 . Across Costa Rica we find that tourism increases when species richness of birds (A), richness of threatened and endemic birds (B), and richness of amphibians, reptiles, and mammals increase (C). Tourism is predicted to decrease as distance to roads (D), distance to protected areas (E), and distance to water increases (F). The relationship between hotel density and tourism shows a nonlinear pattern, where intermediate levels of hotel density lead to highest levels of tourism (G). These curves show the results of ten k-fold replicated MaxEnt runs using either Flickr photographs (solid lines) or eBird CUDs (dotted lines) as the response variable and the corresponding environmental variable as the predictor. The line represents the mean of the ten replicates, while the shaded region is drawn between the minimum and maximum prediction of the ten replicates.

might indicate that tourist activities in Costa Rica are more diverse than previously reported (41, 42) .

Aligned with our predictions, we find that the richness of amphibians, reptiles, and mammals explains tourism inside protected areas better than country-wide. However, the contribution of this variable to tourism is smaller than we expected. The slightly bigger predictive power of these vertebrates, compared to birds, of tourism inside protected areas might be explained by visitors' preferences toward charismatic fauna, such as monkeys (e.g., White-faced capuchins, howlers) and spectacled caiman (Caiman crocodilus), which are some of Costa Rica's flagship species and main attractors of tourists (43, 44) .

Moreover, while the slope of the relationship is fairly shallow, we see that intermediate levels of diversity yield the highest tourism; this echoes results from Namibian conservancies, where more diverse ecological communities provide tourists with a greater range of viewing opportunities (45) . With respect to birdwatching, it is unsurprising that biodiversity is more important in explaining this form of nature-based tourism, compared to all tourism (22) . It is important to recognize that the nonlinear effect of species richness on birdwatching and overall tourism may indicate that the identity of the birds in communities is more important than total richness. These results are consistent with studies showing that species richness alone does not drive cultural value, as in the case of wildlife viewing and preferences for beautiful landscapes with flowering plants (11, 15) . Many birders seek out rare and endemic species, whose protection may depend on greater conservation effort that, in turn, may bring greater economic benefit.

Our spatial models identify some differences between where birdwatchers go and where general tourists go. While spatial congruence of species richness tends to be high between vertebrate taxa (e.g., places with more birds also have more mammals), it is not always congruent with the richness of rare and endemic vertebrate species (46) . Thus, in our study we show that the highlands of Costa Rica have more birdwatching tourism than any other place. Indeed, the highlands are where the endemic birds are more common (47) , and where the most preferred species by birdwatchers can be found (e.g., Resplendent Quetzal and Ornate Hawk-Eagle) (48) (Fig. 3) .

Moreover, we find that while birdwatching, like all tourism, increases with infrastructure, it is highest where both biodiversity and infrastructure co-occur. Thus, if biodiversity is not considered, places like Amistad International Park, and the Corcovado National Park in the Osa Peninsula would be predicted to have less tourism than currently observed (Fig. 3) . We conclude that a degradation of the country's ecosystems could negatively impact the country's economy. In typical years, Costa Rican tourism employs 160,000 people directly and 450,000 more people indirectly (33)-although we recognize that the COVID-19 pandemic has severely impacted the tourism sector and the local communities that rely on tourism in the past years (49) . Loss of biodiversity could mean further attrition in tourist's visits, with cascading impacts to revenues and jobs.

Funding from multilateral development banks, bilateral cooperation agencies, and the private sector could prioritize investing in infrastructure to support sustainable tourism in places where biodiversity is high (50) and to advance the international policy agendas of mainstreaming biodiversity in the tourism and infrastructure sectors (51, 52) . Nonetheless, it is critical to recognize that increased access to protected areas by expanding road networks in undisturbed areas (e.g., Amistad International Park), might cause biodiversity declines, as has been demonstrated in other tropical regions (53, 54) . Therefore, a crucial management challenge for Costa Rica, and other countries depending on tourism, is to plan infrastructure development in the least destructive ways. This means not only avoiding the most sensitive areas and maintaining the balance between an ecosystem's integrity and its accessibility to the public (55), but also policies, plans, investments, and cultural norms reinforcing the conservation of nature (56) .

The relationship between tourism, biodiversity, and infrastructure documented here may be sufficiently general to translate to other developing countries with similar species richness, but the relative importance of biodiversity and infrastructure for tourism might vary across spatial scales. For instance, for nature-based activities such as birdwatching and wildlife viewing, the presence of certain species and the diversity of species might determine where people go (15) . Birdwatchers, in particular, are known to be collectors who prefer endemic and rangerestricted species (57) , which might factor into the decision of which country to visit, whereas infrastructure such as hotels and roads might explain where people go birdwatching once they are in a country (5) .

We also note that Costa Rica is politically stable and has elected a series of governments that prioritize environmental management and protection (58) . Governance indeed is a key factor determining tourism success, even inside protected areas (59) . Thus, tourism in other places that have high biodiversity and high road density might not yield the same economic revenues due to the lack of political stability and security. Moreover, leakage can occur when some areas are developed for tourism while others with similar characteristics get degraded due to the lack of oversight and regulation (60) . Future studies could compare the role that biodiversity and infrastructure play in driving tourism relative to the role of governance in locations where tourism is a promising solution for sustainable development, but political conditions are less stable (e.g., Colombia and the Democratic Republic of Congo) (61, 62) .

Our research builds on a growing number of studies showing that social media are an effective source of information on tourism (19, 63, 64) . While MaxEnt was originally created to model the distributions of species based on environmental predictors such as climate and land use (65) , it has recently been used for mapping social values for ecosystem services (16, 17) . To date, MaxEnt applications involve either the mapping of cultural ecosystem services by identifying people's preferences for landscape attributes (16) , or species distribution modeling for conservation planning (66) . They are rarely done together, except in very small-scale studies (67) . Therefore, the methodological advancements presented in this paper link Earth observations (via land-use maps derived from satellite imagery) with vertebrate species distributions, and with social media data at national and subnational scales.

While this method is powerful for informing national planning decisions, particularly by highlighting the places with high tourism demand and the places where biodiversity has a disproportionate effect on tourism, it also has some limitations. One is the potential taxonomic bias in the Global Biodiversity Information Facility (GBIF) dataset used to generate our species distribution models. Specifically, many rare species may not have sufficient observations in GBIF to meet our stated requirements for inclusion (68) . For example, we were able to model many more birds than amphibians, reptiles, and mammals, and yet tourist's surveys indicate that these latter three taxa are their preferred animals to find (43) . Therefore, our model is likely underestimating the role of rare species, such as jaguars, in driving tourism (69) . Future studies are needed to evaluate the role of rare species on driving tourism, perhaps with complementary research methods (e.g., choice experiments).

In conclusion, by using Costa Rica as an example, we provide a pathway for future evaluations of the role of biodiversity in driving tourism and nature-based tourism at national scales, which is the most relevant scale for development planning and biodiversity conservation strategies (70) . Our proposed methodology and approach could inform UN SEEA accounts and other approaches aimed at mainstreaming biodiversity in the tourism sector and scaling it up to other countries (71) .

Study Area. Costa Rica's land surface covers 51,100 km 2 and harbors ∼5% of the world's macroscopic species, making it the world's number one country in terms of species density per 1,000 km 2 (72) . This small country encompasses 95 distinct climatic zones (73) , ranging from the submoist dry and very warm climate of the Nicoya Peninsula lowlands in the northwest to the very wet and cold climates in the Cordillera de Talamanca, which runs from the northwest to the southeast of the country and rises to an elevation of 3,819 m (73) . Costa Rican climate is influenced by both the Atlantic and the Pacific oceans. Its mean annual temperature varies from 26.7°C on the north Pacific coast, to 26°C on the Caribbean coast, to 6°C on the highest peak of the mountain range. Total annual precipitation varies from 1,300 mm to 7,467 mm, and average monthly relative humidity varies between 65 to 90% (73). The climatic variation results in a tremendous diversity of ecosystems. The most recent classification suggests 14 major ecosystem types that range from the p aramos (above treeline), to seasonal tropical dry forests, evergreen tropical rainforests, mangroves, rivers and estuaries, bogs, marshes, and swamps (74) .

Data Collection. To test the relationship between predictor variables and tourism and birdwatching we conducted analyses at two spatial scales. First, we evaluated such relationships at the country level, and second, we did so at the protected area level (n = 108 protected areas, SI Appendix, Table S1 ). We collected data for eight predictor variables (1-Distance to roads, 2-Hotel density, 3-Distance to water, 4-Distance to protected area, 5-Protected area size, 6-Richness of amphibians, reptiles and mammals, 7-Bird species richness, 8-Richness of threatened and endemic birds) and two response variables (1-Flickr photographs, 2-eBird checklists) (SI Appendix, Table S2 ). What follows presents the methods for calculating each predictor variable. Distance to roads. We extracted the roads in Costa Rica from the Global Roads Open Access Database in 2018. This layer was then rasterized at a 915 m resolution (in Universal Transverse Mercator [UTM] zone 16N), with all pixels intersecting roads being assigned a value of 1. We then calculated the distance to roads using the distance command in the R package raster. We used rasters of 915 × 915 m of resolution throughout the study. Hotel density. We first downloaded a list of hotels registered in Costa Rica by 2018 through the GeoNames database (https://www.geonames.org/). This list contains the latitude and longitude of hotels, nature lodges, and hostels throughout the country (SI Appendix, Table S3 ). We converted this layer to UTM zone 16N coordinates. We used the "kde2d" function from the MASS package in R (49) to perform a kernel density calculation. Using a Gaussian decay function with a SD of 25 km we generated a hotel density map for the entirety of Costa Rica matching the resolution of our other rasters (915 × 915 m). Distance to water. We created a composite layer of water bodies in Costa Rica by combining shapefiles of inland rivers and lakes with a global ocean shapefile, both downloaded from the Natural Earth database (https://www. naturalearthdata.com/). We converted this shapefile to a raster matching the resolution of the previous layers. All pixels that touched polyline features (rivers), and pixels that were covered by at least 50% of a polygon (lakes and oceans) were converted to values of 1 (indicating presence of water). Distance to water was then calculated for each pixel as the Euclidean distance to the nearest body of water using the distance command in the R package raster (48) . We extracted the mean value of each of the above variables for each protected area using zonal statistics. Distance to protected areas and protected area size. We rasterized the protected areas of Costa Rica (https://www.protectedplanet.net/) with those pixels being composed of at least 50% protected area at a 915 m resolution being assigned a value of 1. We then calculated the distance to the nearest protected area using the "distance" command in the R package raster (48) . The total surface area for the protected areas, which we downloaded from the world's protected area network database, was also used as a predictor variable (SI Appendix, Table S1 ). Biodiversity variables. We downloaded occurrence points for terrestrial vertebrate species (mammals, reptiles, amphibians, and birds) in Costa Rica from the Global Biodiversity Information Facility (GBIF) in July 2021 (GBIF.org, https:// doi.org/10.15468/dl.4a3zpk). We removed records that had geospatial issues as noted by GBIF (outside of assumed datum or rounded coordinates). We eliminated all occurrence records not resolved to species, with spatial resolution less than 1 km. We also excluded any species for which we had fewer than 25 occurrence records. We only included observations from 2005 to 2015 to match the temporal extent of the Flickr and eBird data (see below).

With these data, we estimated spatial terrestrial vertebrate species richness in Costa Rica using MaxEnt spatial modeling software version 3.1 (50) , which generates probabilistic range maps using the climatic and land cover conditions present where species are found. Models were fit using four climatic variables from WorldClim (annual temperature bio1, temperature seasonality bio4, annual precipitation bio12, and precipitation seasonality bio15). Climate data were acquired at a 30 arc-second resolution across the entire country (46) . These files were reprojected to UTM zone 16N from their native projection system (Wideband Global SATCOM, 1984). All spatial products were also reprojected to the UTM zone 16N system using bilinear interpretation.

We gathered Landsat satellite images (United States Geological Survey, United States) to create a land-use/land-cover map. We used fully-constrained spectral mixture analysis to create subpixel fractional cover maps of soil (fco-ver1), photosynthetic vegetation (fcover2), and impervious surfaces (fcover3), also from Landsat (75) . Soil and impervious spectra were resampled from a global spectral library (76) , and vegetation spectra were simulated using PRO-SAIL (77) . We created percent tree cover maps using random forest regression, with high resolution tree cover maps from (78) as training data and Landsat/ fcover as covariates. Each layer was created at 30 m resolution and upscaled them to 915 m resolution (47) . Thus, models were also fit using the land cover maps. Open source code to create these maps is available in GitHub (https:// github.com/earth-chris/earthlib/).

To account for uneven spatial distribution of sampling effort, we created a bias raster using all GBIF occurrences within the extent of Costa Rica. After applying the same quality controls described above, we reprojected the coordinates to UTM zone 16N. We then used the "kde2d" function of the MASS package with a Gaussian decay function and a 5-km SD to create a raster of the relative density of GBIF points across the landscape. We used this map to randomly sample 10,000 pseudo absences from the background, with the likelihood of a pixel being chosen being weighted by the value of that pixel in the bias map.

We allowed the model to be fit with linear, quadratic, multiplicative, and hinge terms using standard penalization techniques for more complex terms. The models were fit with four bioclimatic variables (mean annual temperature, mean annual precipitation, temperature seasonality, and precipitation seasonality), percent soil cover, percent photosynthetic cover, and percent tree cover. The run was replicated five times for each species, withholding 20% of the data for testing the models fit using AUC. All other parameters were kept as the MaxEnt default options. We excluded any species where the AUC of the species distribution model when evaluating the test dataset was less than 0.75, indicating a poor model fit. A total of 765 vertebrate species were ultimately included in the analysis (birds = 699 spp, mammals= 23 spp, reptiles = 21 spp, amphibians= 22 spp, see SI Appendix, Table S4 ). All models were conducted using the dismo package in R (51) . To subset species that were classified as threatened and endemic birds, we took the list of birds (n = 699) and consulted the BirdLife International and International Union for Conservation of Nature (IUCN) red list status. We classified a species as threatened or endemic if the species or a subspecies was endemic to Costa Rica, endemic to Central America or listed as vulnerable and endangered by IUCN. In total, we found 279 species that belonged in this category (SI Appendix, Table S5 ). We used these three biodiversity variables (1-bird species richness, 2-richness of threatened and endemic birds, and 3-amphibian, mammals, reptiles species richness) in separate models to avoid collinearity given the high correlation between them (SI Appendix, Figs. S3 and S4) .

We extracted the mean prevalence for each species, by first making a binary presence/absence map of each species using a thresholding method common in the literature which sets the probability of a false-positive and a false-negative to be equal (79) . These were then summed to create species richness maps using species belonging to each biodiversity variable (either bird species richness, endemic and threatened birds, or richness of amphibians, reptiles, and mammals).

Tourism. We measured tourism via Flickr photographs. Flickr is a social media platform for sharing photographs, some of which are tagged with the date and location where the user captured the image. Visitation was estimated by summing unique Flickr users who shared at least one image from a location per day. This metric, known as photograph-user-days (PUD) (19) corrects for potential biases created by Flickr users who upload multiple photographs from one location per day (19) , and it is the most common measure of visitation for research on recreation ecosystem services (79) . We used the Integrated Valuation of Ecosystem Services and Tradeoffs (InVEST) visitation model to compute PUD across all the available years of data (2005-2017) (80) . We used two separate InVEST runs to generate PUDs for each protected area and for each 915 m × 915 m pixel country-wide. We also created a threshold of 4 km 2 as the minimum protected area size. We did so because we identified that protected areas smaller than 4 km 2 had zero or very few Flickr photograph uploads and were inducing skewness in the data via zeroinflated patterns.

To measure birdwatching tourism, we used all eBird observations made from 2005 to 2017 (81). eBird is an online platform where birdwatchers share geolocated checklists of birds observed at a specific time. As we did for Flickr PUDs, we only included unique eBird checklist-user-days (CUD), so that a given user contributed no more than a single count per day in an individual area. We summed CUD across all years in each protected area for the protected area analysis and per pixel (915 m × 915 m) for the country-wide analysis.

Data Analysis. We used MaxEnt software (82) to map how tourism can be explained by predictor variables. We used species distribution models and instead of species occurrence variables we used the two tourism responses: Flickr PUD and eBird CUD. We applied square root transformations to distance variables (distance to water and distance to protected areas), and applied a log transformation to hotel density. We did so to conform with model assumptions and to assess goodness of fit for linear regressions (83) (SI Appendix, Figs. S6 and S7). We ran ten models in total, out of which five were at the countrywide scale and five at the protected area scale (see SI Appendix, Table S6 ).

We fit models using distance to roads, hotel density, distance to water, protected area size (for models within protected areas) or distance to protected areas (for the country-wide analyses), and the three biodiversity variables independently (see SI Appendix, Table S6 for all ten models). For all models, the selected five layers used as predictor variables did not show collinearity (variance inflation factor [VIF] <2, SI Appendix, Table S7 ) as measured by the R package usdm. For the models that considered tourism within protected areas, we clipped all input layers using polygons provided by the Protected Areas of Costa Rica. Unlike our species distribution models for biodiversity, we constrained the model to only select linear and/or quadratic terms. Models were replicated ten times, withholding 10% of the data for validation in a k-fold pattern. Pseudo absences were chosen at random across the entirety of the extent (either the whole country for country-wide models or just within protected areas). Otherwise default options were used in the models. All models were conducted using the R dismo package (84) .

In the main text, we report results of all ten models. In the SI Appendix, we report summary statistics for each model (e.g., model AUC, percent contribution, permutation importance for each variable, see SI Appendix, Table S6 ) and show graphically how each variable affects tourism (Figs. 1 and 2). We show this by plotting on the x-axis each predictor variable and on the y-axis the expected tourism value for MaxEnt models constructed using only this variable.

We projected tourism as either Flickr PUD or eBird CUD across the entirety of the country using the predict function in the dismo package in R (Fig. 3 A  and D) . Lastly, we created maps of the difference between the models created with all variables and those which withheld either biodiversity (bird richness, threatened and endemic birds, or amphibians, reptiles, and mammals) (Fig. 3 B and E), or infrastructure variables (distance to roads and hotel density) (Fig.  3 C and F) to see how excluding these variables affected the predicted patterns of tourism.

Data Availability. All original data is publicly available and referenced in the text. Intermediate data layers and final data layers are available on Dryad at: https://doi.org/10.5061/dryad.1ns1rn8w7 (85) .

Scripts to run the analyses and to create the figures are available on Dryad and on GitHub at: https://github.com/jeffreysmith-jrs/natcapCR/blob/main/bin/ tourismModel.R.

Scripts for calculating the fractional cover maps are available at: https:// github.com/earth-chris/earthlib/.

We are thankful to colleagues at the Natural Capital Project and at the Center for Conservation Biology at Stanford University for invaluable comments. Particularly, we thank Drs. Morgan Kain, Lisa Mandle, Adrian Vogl, and Mary Ruckelshaus for their feedback. We also thank the Costa Rican Tourism Institute and the Central Bank of

Linking livelihoods and conservation: A conceptual framework and scale for assessing the integration of human needs and biodiversity

A global perspective on trends in nature-based tourism

Complementary benefits of tourism and hunting to communal conservancies in Namibia

The role of transport infrastructure in international tourism development: A gravity model approach

A spatial-temporal analysis of hotels in urban tourism destination

Estimating economic losses to tourism in Africa from the illegal killing of elephants

Biodiversity and nature-based tourism at forest reserves in Uganda

Ecotourism marketing alternative to charismatic megafauna can also support biodiversity conservation

Assessing global popularity and threats to Important Bird and Biodiversity Areas using social media data

Advancing sustainable development and protected area management with social media-based tourism data

Effects of bird community dynamics on the seasonal distribution of cultural ecosystem services

Ecotourism for conservation?

Mismatches between supply and demand in wildlife tourism: Insights for assessing cultural ecosystem services

Assessing, mapping, and quantifying cultural ecosystem services at community level

Species richness alone does not predict cultural ecosystem service value

Validating a method for transferring social values of ecosystem services between public lands in the Rocky Mountain region

Or€ uc€ u, MaxEnt modelling of the potential distribution areas of cultural ecosystem services using social media data and GIS

Mart ın-L opez, Spatial patterns of cultural ecosystem services provision in Southern Patagonia

Using social media to quantify nature-based tourism and recreation

Beltr an, Cultural ecosystem services provided by rivers across diverse social-ecological landscapes: A social media analysis

Differential reporting of biodiversity in two citizen science platforms during COVID-19 lockdown in Colombia

Economic and conservation potential of birdwatching tourism in postconflict Colombia

Large mammal diversity matters for wildlife tourism in Southern African Protected Areas: Insights for management

Priority list of biodiversity metrics to observe from space

Framing the concept of satellite remote sensing essential biodiversity variables: Challenges and future directions

Satellite remote sensing of ecosystem functions: Opportunities, challenges and way forward

The value of species rarity in biodiversity recreation: A birdwatching example

Managing ecosystem services: What do we need to know about their ecology?

Determining beach user knowledge of rip currents in Costa Rica

Surf localism in Costa Rica: Exploring territoriality among Costa Rican and foreign resident surfers

Geolocated social media as a rapid indicator of park visitation and equitable park access

The non-market value of birding sites and the marginal value of additional species: Biodiversity in a random utility model of site choice by eBird members

Plan nacional de desarrollo tur ıstico de Costa Rica

Encuesta a erea de no residentes, Aeropuerto Internacional Daniel Oduber Quir os

Can ecotourism deliver real economic, social, and environmental benefits? A study of the Osa Peninsula, Costa Rica

The effect of land use change and ecotourism on biodiversity: A case study of Manuel Antonio

Global relationships between biodiversity and naturebased tourism in protected areas

Quantifying tourism booms and the increasing footprint in the Arctic with social media data

Social media reveal ecoregional variation in how weather influences visitor behavior in U.S. National Park Service units

Quantifying wildlife watchers' preferences to investigate the overlap between recreational and conservation value of natural areas

Ecotourism at a crossroads: The case of Costa Rica

Leveraging support for conservation from ecotourists: Can relational values play a role?

Species-specific responses to tourist interactions by white-faced capuchins (Cebus imitator) and mantled howlers (Alouatta palliata) in a Costa Rican wildlife refuge

Measuring the economic value and social impact of crocodile Tourism in T arcoles (FIU

Effect of diversity of large wildlife species on financial benefits to local communities in northwest Namibia

Global distribution and conservation of rare and threatened vertebrates

Reproduction, habitat use, and natural history of the Black-and-yellow Silky-Flycatcher (Phainoptila Melanoxantha), an endemic bird of the western Panama-Costa Rican highlands

Percepci on social mediante redes sociales para definir las aves m as bellas de Costa Rica

A loss of tourism threatens Costa Rica's lush paradise

Ecotourism and economic incentives-An empirical approach

The Cancun declaration on mainstreaming the sustainable use and conservation of biodiversity for well-being

Decision adopted by the Conference of the Parties to the Convention on Biological Diversity. 14/3 Mainstreaming biodiversity in the energy and mining, infrastructure, manufacturing and processing sectors

The impacts of roads and other infrastructure on mammal and bird populations: A meta-analysis

Effects of roads on terrestrial vertebrate species in Latin America

A framework for managing the negative impacts of wildlife tourism on wildlife

Amplifying small solutions for systemwide change" in Green Growth That Works: Natural Capital Policy and Finance Mechanisms From Around the World

Can avian functional traits predict cultural ecosystem services? People and Nature

The role of ecotourism in conservation: Panacea or Pandora's box? Biodivers

Vulnerability to climate change of nature-based tourism in the Nepalese Himalayas

Why keep lions instead of livestock? Assessing wildlife tourismbased payment for ecosystem services involving herders in the Maasai Mara

Valuing beaches to develop payment for ecosystem services schemes in Colombia's Seaflower marine protected area

Can nature tourism help finance protected areas in the Congo Basin?

Assessing the usability of social media data for visitor monitoring in protected areas

The Geographic spread and preferences of tourists revealed by user-generated information on Jeju Island

A maximum entropy approach to species distribution modeling

MaxEnt modeling for predicting suitable habitat for threatened and endangered tree Canacomyrica monticola in New Caledonia

Spatially explicit modeling of conflict zones between wildlife and snow sports: Prioritizing areas for winter refuges

Geographical, temporal and taxonomic biases in insect GBIF data on biodiversity and extinction

When do ecosystem services depend on rare species?

As far as possible and as appropriate": Implementing the Aichi biodiversity targets

Mainstreaming biodiversity: A review of national strategies

Estado del conocimiento y gesti on

Climate of Costa Rica

Costa Rican ecosystems: Setting the stage

Fully constrained least-squares based linear unmixing [hyperspectral image classification

A global spectral library to characterize the world's soil

PROSPECT+SAIL models: A review of use for vegetation characterization

Improving tree cover estimates for fine-scale landscape ecology

Uses and limitations of social media to inform visitor use management in parks and protected areas: A systematic review

InVEST user's guide (The Natural Capital Project

EOD -eBird observation dataset. Occurrence dataset

A brief tutorial on MaxEnt. Biodiversity Informatics. American Museum of Natural History

The analysis of biological data (Roberts and Company Publishers

Package "dismo": Species distribution modeling

Data from: Biodiversity and infrastructure interact to drive tourism to and within Costa Rica

Costa Rica for providing reports and perspectives on the tourism sector in the country. In particular, we thank Johnny Villareal, Roxana Arguedas, and Rodolfo Lizano. We thank our reviewers Drs. Natalia Ocampo Peñuela and Monica Turner for their thorough reviews and for helping us improve this manuscript. Lastly, we thank Jim Zook for helping us identify the endemic and threatened birds of Costa Rica, and Dr. Daniel S. Karp for support throughout this project. We acknowledge funding from NASA (80NSSC18K0434) and the Gordon and Betty Moore Foundation.